Functions and Subroutines in FORTRAN

Functions and Subroutines: It is difficult to write the lengthy program corresponding to a complex problem/project as one single unit. Employing functions and subroutines, we can divide the complex problem in several small units, and can then write small subprograms, such as (newly defined) functions and/or subroutines, corresponding individually to these units. The main program then is written to tackle the main problem, with the subsidiary or minor aspects being taken care of by the subprograms, as and when they are called by the main program. Such programming method is called modular programming. The subprograms become particularly useful if some tasks (such as determination of roots of a quadratic equation - as in the following example) are to be done repeatedly. A function calculates the result of mathematical operations on some quantities called its arguments, which must be defined when it is called. A subroutine acts as a part of a program, and generally also have some arguments (for a subroutine, the output variables, if any, must also be mentioned as its arguments).

Problem: For chemical-reaction equilibria of the form H2(g) + I2(g) = 2HI(g) with products not initially added, the concentration of the product x is related to the initial reactant concentrations a & b by the obvious relation 4x2/{(a–x)(b–x)} = KC, i.e., the quadratic equation  (KC –4)x2 – (a+b)KC.x + abKC = 0, with only a real positive root, also less than both a & b, only being meaningful. Now, given the value of KC (the equilibrium constant), let us calculate this x for values of a & b each ranging from 0.01 to 0.04 mol dm–3, in steps of 0.01 mol dm–3. To do this main problem using FORTRAN, let us take the help of a subroutine giving the coefficients P, Q, R of the quadratic equation Px2 + Qx + R = 0 from the current values of a, b & KC, two functions to get the two real roots (under condition that they exist), and another subroutine to judge which of the two roots to accept, as exemplified below. [Let the user run the program by entering the value 49.8, say, for the constant KC].

C     Calculation of amount x of product for reaction between H2 & I2
*     (for various concentrations of the two reactants)
      WRITE(*,*) 'Enter the value of the equilibrium constant Kc:'
      WRITE(*,*) ' '
      READ(*,*) EqCnKc
      WRITE(*,*) ' '
      WRITE(*,*) 'H2-Conc.   I2-Conc.   HI-Conc.'
      WRITE(*,*) ' '
      ConcH2 = 0.00
      DO 11 I = 1, 4, 1
        ConcH2 = ConcH2 + 0.01
        ConcI2 = 0.00
        DO 12 J = 1, 4, 1
          ConcI2 = ConcI2 + 0.01
          CALL COEFFS(ConcH2,ConcI2,EqCnKc,P,Q,R)
          QQM4PR = Q*Q - 4.0*P*R
          IF (QQM4PR .GE. 0.0) THEN
            X1G = ROOTGR(P,Q,R)
            X2L = ROOTLR(P,Q,R)
            CALL JudgeR(X1G,X2L,ConcH2,ConcI2,Xgood)
            WRITE(*,3) ConcH2,ConcI2,Xgood
  3         FORMAT(F9.6,2X,F9.6,2X,F9.6)
          ELSE
            WRITE(*,*) 'Both roots are found to be complex!'
          ENDIF
 12     CONTINUE
 11   CONTINUE
      STOP
      END
      SUBROUTINE COEFFS(A,B,EC,P2,Q1,R0)
      P2 = EC - 4.0
      Q1 = -(A+B)*EC
      R0 = A*B*EC
      RETURN
      END
      FUNCTION ROOTGR(P2,Q1,R0)
      ROOTGR = (-Q1 + SQRT(Q1**2 - 4.0*P2*R0))/(2.0*P2)
      RETURN
      END
      FUNCTION ROOTLR(P2,Q1,R0)
      ROOTLR = (-Q1 - SQRT(Q1**2 - 4.0*P2*R0))/(2.0*P2)
      RETURN
      END
      SUBROUTINE JudgeR(RT1G,RT2L,C1REAC,C2REAC,RT0)
      IF (RT2L .LE. 0.0) THEN
        RT0 = RT1G
      ELSE
        CLess = MIN(C1REAC,C2REAC)
C       Note: Under this conditional state, both roots are positive
        IF (RT1G .GE. CLess) THEN
          RT0 = RT2L
        ELSE
C         Note: Under this conditional state, even greater-root < CLess
          RT0 = RT1G
        ENDIF
      ENDIF
      RETURN
      END

Now, let us note the following rules about writing and using functions and subroutines in FORTRAN:
        (1) All the functions and subroutines must be written after the END statement of the main program. However, no blank line(s) should separate the main program and the subprograms, or the various subprograms themselves.
        (2) Each subroutine must be written starting with a statement line SUBROUTINE SRName(Arguments); each function must be written starting with a statement line FUNCTION FNName(Arguments). The rules for naming FORTRAN functions and subroutines are same as those for the variables. Each function or subroutine must be concluded with a RETURN statement followed by an END statement (now, you probably understand why there were two concluding statements STOP and END in the main programs!). The arguments for a subroutine need include the input variables as well as the output ones, while those for a function includes only the input ones.
        (3) Within the main program, a subroutine can be called (i.e., invited to perform its specialized job) by the CALL SRName(Arguments) statement, while a newly defined function can be used just as one uses a pre-defined function, e.g., the SQRT(Var1) function or the MIN(Var1,Var2) function (used above). Note that a function may have either one argument or more than one arguments, but the value of the function is only one (e.g., the value of MIN(6,77) is 6 only even though it has two arguments, the value of SQRT(4.0) is only 2.0 but not –2.0 simultaneously).
        (4) The oddest thing you may perceive about functions and subroutines is that their arguments in the written function/subroutine are not, in general, the same as their arguments when they are called within the main program! Thus the arguments of the ROOTGR function (that gives the greater real root of a quadratic equation) in the written function is P2, Q1 and R0 (check that above), while within the main program they are P, Q and R (check that). Similarly, the arguments of the COEFFS subroutine in the written subroutine are A, B, EC, P2, Q1 and R0 (check that) while when called by the main program they are ConcH2, ConcI2, EqCnKc, P, Q and R. What matters then, regarding the arguments of functions and subroutines? Yes, it is only their order that matters. Thus, regarding the COEFFS subroutine, it is its sixth argument variable that will be allotted a value equal to the product of its first three arguments, as per the statement R0 = A*B*EC within its definition segment (check that).