C
C There are a total of 9 entries in the algebraic variable array.
C There are a total of 1 entries in each of the rate and state variable arrays.
C There are a total of 29 entries in the constant variable array.
C
C
C VOI is time in component environment (second).
C ALGBRC(9) is F_CE in component F_CE (newton).
C ALGBRC(1) is f_L_CE in component f_L_CE (newton).
C ALGBRC(6) is g_V_CE in component g_V_CE (dimensionless).
C CONSTS(1) is a in component user_defined_constants (dimensionless).
C CONSTS(2) is F_min in component f_L_CE (newton).
C CONSTS(3) is F_max in component user_defined_constants (newton).
C STATES(1) is L_CE in component L_CE (metre).
C CONSTS(4) is L_CE_opt in component user_defined_constants (metre).
C CONSTS(5) is W in component f_L_CE (dimensionless).
C CONSTS(29) is lambda_a in component lambda_a (second_per_metre).
C CONSTS(6) is V_max in component g_V_CE (metre_per_second).
C ALGBRC(7) is V_CE in component V_CE (dimensionless).
C CONSTS(7) is A in component g_V_CE (dimensionless).
C CONSTS(8) is g_max in component g_V_CE (dimensionless).
C CONSTS(25) is d1 in component d1 (dimensionless).
C CONSTS(27) is d2 in component d2 (dimensionless).
C CONSTS(28) is d3 in component d3 (dimensionless).
C CONSTS(9) is gamma in component g_V_CE (dimensionless).
C CONSTS(10) is V_max in component d1 (metre_per_second).
C CONSTS(11) is A in component d1 (dimensionless).
C CONSTS(12) is g_max in component d1 (dimensionless).
C CONSTS(13) is S in component d1 (metre_per_second).
C CONSTS(14) is S in component d2 (metre_per_second).
C CONSTS(15) is A in component d2 (dimensionless).
C CONSTS(16) is V_max in component d2 (metre_per_second).
C CONSTS(17) is gamma in component d2 (dimensionless).
C CONSTS(18) is g_max in component d3 (dimensionless).
C CONSTS(19) is gamma in component d3 (dimensionless).
C ALGBRC(5) is F_SEE in component F_SEE (newton).
C CONSTS(20) is k_SEE in component F_SEE (newton_per_metre2).
C ALGBRC(4) is L_SEE in component L_SEE (metre).
C CONSTS(21) is L_slack in component F_SEE (metre).
C ALGBRC(2) is F_PEE in component F_PEE (newton).
C CONSTS(26) is k_PEE in component k_PEE (newton_per_metre2).
C CONSTS(22) is L_slack in component F_PEE (metre).
C CONSTS(23) is W in component k_PEE (dimensionless).
C CONSTS(24) is L_CE_opt in component k_PEE (metre).
C ALGBRC(3) is L_m in component L_m (metre).
C ALGBRC(8) is F_m in component F_m (newton).
C RATES(1) is d/dt L_CE in component L_CE (metre).
C
      SUBROUTINE initConsts(CONSTS, RATES, STATES)
      REAL CONSTS(*), RATES(*), STATES(*)
      CONSTS(1) = 0.8
      CONSTS(2) = 10
      CONSTS(3) = 7000
      STATES(1) = 0.038
      CONSTS(4) = 0.093
      CONSTS(5) = 0.63
      CONSTS(6) = 0.93
      CONSTS(7) = 0.25
      CONSTS(8) = 1.5
      CONSTS(9) = 5.67
      CONSTS(10) = 0.93
      CONSTS(11) = 0.25
      CONSTS(12) = 1.5
      CONSTS(13) = 2
      CONSTS(14) = 2
      CONSTS(15) = 0.25
      CONSTS(16) = 0.93
      CONSTS(17) = 5.67
      CONSTS(18) = 1.5
      CONSTS(19) = 5.67
      CONSTS(20) = 1000000
      CONSTS(21) = 0.0025
      CONSTS(22) = 0.0025
      CONSTS(23) = 0.63
      CONSTS(24) = 0.01
      CONSTS(25) = ( CONSTS(10)*CONSTS(11)*(CONSTS(12) - 1.00000))/( CONSTS(13)*(CONSTS(11)+1.00000))
      CONSTS(26) = CONSTS(3)/ CONSTS(23)*CONSTS(24) ** 2.00000
      CONSTS(27) = ( CONSTS(14)*(CONSTS(15)+1.00000))/( CONSTS(16)*CONSTS(15)*CONSTS(17)+1.00000 ** 2.00000)
      CONSTS(28) = ( (CONSTS(18) - 1.00000)*CONSTS(19) ** 2.00000)/CONSTS(19)+1.00000 ** 2.00000+1.00000
      CONSTS(29) =  1.00000*((1.00000 - EXP( - 3.82000*CONSTS(1)))+ CONSTS(1)*EXP(- 3.82000))
      RETURN
      END
      SUBROUTINE computeRates(VOI, CONSTS,  RATES, STATES, ALGBRC)
      REAL VOI, CONSTS(*), RATES(*), STATES(*), ALGBRC(*)
      ALGBRC(1) = ( CONSTS(3)*( 1.00000*(1.00000 - STATES(1)) - CONSTS(4) ** 2.00000))/( CONSTS(5) ** 2.00000*CONSTS(4) ** 2.00000)
      ALGBRC(3) = TERNRY(VOI.LE.1.00000, 0.0380000, TERNRY(VOI.GT.1.00000.AND.VOI.LT.2.00000, 0.0380000+ 0.00200000*(VOI - 1.00000), 0.0400000)
      ALGBRC(4) = ALGBRC(3) - STATES(1)
      ALGBRC(5) = TERNRY(ALGBRC(4).LE.CONSTS(21), 0.00000,  CONSTS(20)*ALGBRC(4) - CONSTS(21) ** 2.00000)
      ALGBRC(2) = TERNRY(STATES(1).LE.CONSTS(22), 0.00000,  CONSTS(26)*STATES(1) - CONSTS(22) ** 2.00000)
CALL minimize(minfunc_0, CONSTS, VARIABLES, <VAR>)
      RATES(1) =  1.00000*ALGBRC(7)
      RETURN
      END
      SUBROUTINE computeVariables(VOI, CONSTS, RATES, STATES, ALGBRC)
      REAL VOI, CONSTS(*), RATES(*), STATES(*), ALGBRC(*)
      ALGBRC(1) = ( CONSTS(3)*( 1.00000*(1.00000 - STATES(1)) - CONSTS(4) ** 2.00000))/( CONSTS(5) ** 2.00000*CONSTS(4) ** 2.00000)
      ALGBRC(3) = TERNRY(VOI.LE.1.00000, 0.0380000, TERNRY(VOI.GT.1.00000.AND.VOI.LT.2.00000, 0.0380000+ 0.00200000*(VOI - 1.00000), 0.0400000)
      ALGBRC(4) = ALGBRC(3) - STATES(1)
      ALGBRC(5) = TERNRY(ALGBRC(4).LE.CONSTS(21), 0.00000,  CONSTS(20)*ALGBRC(4) - CONSTS(21) ** 2.00000)
      ALGBRC(2) = TERNRY(STATES(1).LE.CONSTS(22), 0.00000,  CONSTS(26)*STATES(1) - CONSTS(22) ** 2.00000)
      ALGBRC(8) = ALGBRC(5)
      ALGBRC(9) =  ALGBRC(1)*ALGBRC(6)*CONSTS(1)
      RETURN
      END
REAL FUNCTION minfunc_0(CONSTS, VARIABLES)
            REAL CONSTS(*), VARIABLES(*) 
      minfunc_0[1] = abs(ALGBRC(6) - (TERNRY(ALGBRC(7).LE.0.00000, ( CONSTS(29)*CONSTS(6)+ALGBRC(7))/( CONSTS(29)*CONSTS(6) - ALGBRC(7)/CONSTS(7)), TERNRY(0.00000.LT.ALGBRC(7).AND.ALGBRC(7).LE. CONSTS(9)*CONSTS(25), ( CONSTS(8)*ALGBRC(7)+CONSTS(25))/(ALGBRC(7)+CONSTS(25)), TERNRY(ALGBRC(7).GT. CONSTS(9)*CONSTS(25), CONSTS(28)+ CONSTS(27)*ALGBRC(7), 0.0/0.0)))
      minfunc_0[2] = abs(ALGBRC(7) -  1.00000*(( (1.00000/ALGBRC(6))*( ALGBRC(5)*(ALGBRC(3) - STATES(1)) -  ALGBRC(2)*STATES(1)))/( CONSTS(1)*ALGBRC(1))))
      RETURN
      END
      REAL FUNCTION TERNRY(TEST, VALA, VALB)
      LOGICAL TEST
      REAL VALA, VALB
      IF (TEST) THEN
        TERNRY = VALA
      ELSE
        TERNRY = VALB
      ENDIF
      RETURN
      END