C
C There are a total of 11 entries in the algebraic variable array.
C There are a total of 4 entries in each of the rate and state variable arrays.
C There are a total of 21 entries in the constant variable array.
C
C
C VOI is time in component environment (hour).
C ALGBRC(1) is Metin in component environment (flux).
C STATES(1) is Met in component Met (micromolar).
C ALGBRC(8) is V_MS in component V_MS (flux).
C ALGBRC(11) is V_BHMT in component V_BHMT (flux).
C ALGBRC(2) is V_MATI in component V_MATI (flux).
C ALGBRC(4) is V_MATIII in component V_MATIII (flux).
C STATES(2) is AdoMet in component AdoMet (micromolar).
C ALGBRC(7) is V_METH in component V_METH (flux).
C ALGBRC(5) is V_GNMT in component V_GNMT (flux).
C STATES(3) is AdoHcy in component AdoHcy (micromolar).
C ALGBRC(9) is V_AH in component V_AH (flux).
C STATES(4) is Hcy in component Hcy (micromolar).
C ALGBRC(10) is V_CBS in component V_CBS (flux).
C CONSTS(1) is V_MATImax in component V_MATI (flux).
C CONSTS(2) is Km_MATI in component V_MATI (micromolar).
C CONSTS(3) is Ki_MATI in component V_MATI (micromolar).
C CONSTS(4) is V_MATIIImax in component V_MATIII (flux).
C ALGBRC(3) is Km1_MATIII in component V_MATIII (micromolar).
C CONSTS(5) is Km2_MATIII in component V_MATIII (micromolar).
C CONSTS(6) is V_GNMTmax in component V_GNMT (flux).
C CONSTS(7) is Km_GNMT in component V_GNMT (micromolar).
C CONSTS(8) is Ki_GNMT in component V_GNMT (micromolar).
C CONSTS(9) is V_METHmax in component V_METH (flux).
C ALGBRC(6) is Km1_METH in component V_METH (micromolar).
C CONSTS(10) is Km2_METH_A in component V_METH (dimensionless).
C CONSTS(11) is five_mTHF in component V_MS (micromolar).
C CONSTS(12) is V_MSmax in component V_MS (flux).
C CONSTS(13) is Kd_MS in component V_MS (micromolar).
C CONSTS(14) is Km_Hcy_MS in component V_MS (micromolar).
C CONSTS(15) is Km_five_mTHF_MS in component V_MS (micromolar).
C CONSTS(16) is alpha1 in component V_AH (first_order_rate_constant).
C CONSTS(17) is alpha2 in component V_AH (dimensionless).
C CONSTS(18) is beta1 in component V_CBS (second_order_rate_constant).
C CONSTS(19) is beta2 in component V_CBS (first_order_rate_constant).
C CONSTS(20) is V_BHMTmax in component V_BHMT (flux).
C CONSTS(21) is Km_BHMT in component V_BHMT (micromolar).
C RATES(1) is d/dt Met in component Met (micromolar).
C RATES(2) is d/dt AdoMet in component AdoMet (micromolar).
C RATES(3) is d/dt AdoHcy in component AdoHcy (micromolar).
C RATES(4) is d/dt Hcy in component Hcy (micromolar).
C
      SUBROUTINE initConsts(CONSTS, RATES, STATES)
      REAL CONSTS(*), RATES(*), STATES(*)
      STATES(1) = 53.5
      STATES(2) = 137.6
      STATES(3) = 13.2
      STATES(4) = 0.88
      CONSTS(1) = 561
      CONSTS(2) = 41
      CONSTS(3) = 50
      CONSTS(4) = 22870
      CONSTS(5) = 21.1
      CONSTS(6) = 10600
      CONSTS(7) = 4500
      CONSTS(8) = 20
      CONSTS(9) = 4521
      CONSTS(10) = 10
      CONSTS(11) = 5.2
      CONSTS(12) = 500
      CONSTS(13) = 1
      CONSTS(14) = 0.1
      CONSTS(15) = 25
      CONSTS(16) = 100
      CONSTS(17) = 10
      CONSTS(18) = 1.7
      CONSTS(19) = 30
      CONSTS(20) = 2500
      CONSTS(21) = 12
      RETURN
      END
      SUBROUTINE computeRates(VOI, CONSTS,  RATES, STATES, ALGBRC)
      REAL VOI, CONSTS(*), RATES(*), STATES(*), ALGBRC(*)
      ALGBRC(2) = CONSTS(1)/(1.00000+ (CONSTS(2)/STATES(1))*(1.00000+STATES(2)/CONSTS(3)))
      ALGBRC(3) = 20000.0/(1.00000+ 5.70000*STATES(2)/(STATES(2)+600.000) ** 2.00000)
      ALGBRC(4) = CONSTS(4)/(1.00000+( ALGBRC(3)*CONSTS(5))/(STATES(1) ** 2.00000+ STATES(1)*CONSTS(5)))
      ALGBRC(6) =  1.00000*(1.00000+STATES(3)/4.00000)
      ALGBRC(7) = CONSTS(9)/(1.00000+ALGBRC(6)/STATES(2)+CONSTS(10)+( CONSTS(10)*ALGBRC(6))/STATES(2))
      ALGBRC(5) = ( (CONSTS(6)/(1.00000+CONSTS(7)/STATES(2) ** 2.30000))*1.00000)/(1.00000+STATES(3)/CONSTS(8))
      RATES(2) = (ALGBRC(2)+ALGBRC(4)) - (ALGBRC(7)+ALGBRC(5))
      ALGBRC(9) =  CONSTS(16)*(STATES(3) -  CONSTS(17)*STATES(4))
      RATES(3) = (ALGBRC(7)+ALGBRC(5)) - ALGBRC(9)
      ALGBRC(1) = TERNRY(VOI.LT.2.00000.OR.VOI.GE.8.00000, 200.000, TERNRY(VOI.GE.2.00000.AND.VOI.LT.5.00000, 300.000, TERNRY(VOI.GE.5.00000.AND.VOI.LT.8.00000, 100.000, 200.000)
      ALGBRC(8) = ( CONSTS(12)*CONSTS(11)*STATES(4))/( CONSTS(13)*CONSTS(14)+ CONSTS(14)*CONSTS(11)+ CONSTS(15)*STATES(4)+ CONSTS(11)*STATES(4))
      ALGBRC(11) = ( (0.700000 -  0.0250000*((STATES(2)+STATES(3)) - 150.000))*CONSTS(20)*STATES(4))/(CONSTS(21)+STATES(4))
      RATES(1) = (ALGBRC(8)+ALGBRC(11)+ALGBRC(1)) - (ALGBRC(2)+ALGBRC(4))
      ALGBRC(10) =  ( CONSTS(18)*(STATES(2)+STATES(3)) - CONSTS(19))*STATES(4)
      RATES(4) = ALGBRC(9) - (ALGBRC(10)+ALGBRC(8)+ALGBRC(11))
      RETURN
      END
      SUBROUTINE computeVariables(VOI, CONSTS, RATES, STATES, ALGBRC)
      REAL VOI, CONSTS(*), RATES(*), STATES(*), ALGBRC(*)
      ALGBRC(2) = CONSTS(1)/(1.00000+ (CONSTS(2)/STATES(1))*(1.00000+STATES(2)/CONSTS(3)))
      ALGBRC(3) = 20000.0/(1.00000+ 5.70000*STATES(2)/(STATES(2)+600.000) ** 2.00000)
      ALGBRC(4) = CONSTS(4)/(1.00000+( ALGBRC(3)*CONSTS(5))/(STATES(1) ** 2.00000+ STATES(1)*CONSTS(5)))
      ALGBRC(6) =  1.00000*(1.00000+STATES(3)/4.00000)
      ALGBRC(7) = CONSTS(9)/(1.00000+ALGBRC(6)/STATES(2)+CONSTS(10)+( CONSTS(10)*ALGBRC(6))/STATES(2))
      ALGBRC(5) = ( (CONSTS(6)/(1.00000+CONSTS(7)/STATES(2) ** 2.30000))*1.00000)/(1.00000+STATES(3)/CONSTS(8))
      ALGBRC(9) =  CONSTS(16)*(STATES(3) -  CONSTS(17)*STATES(4))
      ALGBRC(1) = TERNRY(VOI.LT.2.00000.OR.VOI.GE.8.00000, 200.000, TERNRY(VOI.GE.2.00000.AND.VOI.LT.5.00000, 300.000, TERNRY(VOI.GE.5.00000.AND.VOI.LT.8.00000, 100.000, 200.000)
      ALGBRC(8) = ( CONSTS(12)*CONSTS(11)*STATES(4))/( CONSTS(13)*CONSTS(14)+ CONSTS(14)*CONSTS(11)+ CONSTS(15)*STATES(4)+ CONSTS(11)*STATES(4))
      ALGBRC(11) = ( (0.700000 -  0.0250000*((STATES(2)+STATES(3)) - 150.000))*CONSTS(20)*STATES(4))/(CONSTS(21)+STATES(4))
      ALGBRC(10) =  ( CONSTS(18)*(STATES(2)+STATES(3)) - CONSTS(19))*STATES(4)
      RETURN
      END
      REAL FUNCTION TERNRY(TEST, VALA, VALB)
      LOGICAL TEST
      REAL VALA, VALB
      IF (TEST) THEN
        TERNRY = VALA
      ELSE
        TERNRY = VALB
      ENDIF
      RETURN
      END