C
C There are a total of 26 entries in the algebraic variable array.
C There are a total of 8 entries in each of the rate and state variable arrays.
C There are a total of 46 entries in the constant variable array.
C
C
C VOI is time in component environment (second).
C CONSTS(1) is volume_i0 in component volume_i (volume).
C CONSTS(44) is volume_i in component volume_i (volume).
C CONSTS(2) is gamma in component volume_er (dimensionless).
C CONSTS(45) is volume_er in component volume_er (volume).
C STATES(1) is volumeCa_i in component volumeCa_i (microMVolume).
C ALGBRC(14) is Jipt in component Jipt (microMvolumepersec).
C ALGBRC(16) is Jryr in component Jryr (microMvolumepersec).
C ALGBRC(18) is Jer in component Jer (microMvolumepersec).
C ALGBRC(19) is Jserca in component Jserca (microMvolumepersec).
C ALGBRC(25) is Jin in component Jin (microMvolumepersec).
C ALGBRC(26) is Jpm in component Jpm (microMvolumepersec).
C ALGBRC(1) is Ca_i in component Ca_i (microM).
C STATES(2) is volumeCa_er in component volumeCa_er (microMVolume).
C ALGBRC(4) is Ca_er in component Ca_er (microM).
C STATES(3) is volume_iIP3 in component volume_iIP3 (microMVolume).
C ALGBRC(13) is J_ip3P in component J_ip3P (microMpersec).
C ALGBRC(17) is J_ip3D in component J_ip3D (microMpersec).
C ALGBRC(2) is IP3 in component IP3 (microM).
C STATES(4) is PIP2 in component PIP2 (dimensionless).
C CONSTS(3) is PIP2_Total in component PIP2 (microM).
C CONSTS(46) is J_ip3R in component J_ip3R (per_sec).
C ALGBRC(5) is IPX in component IPX (dimensionless).
C CONSTS(4) is PIP2_Total in component IPX (microM).
C CONSTS(5) is J_ip3R0 in component J_ip3R (per_sec).
C STATES(5) is ATP_e in component receptor (microM).
C ALGBRC(10) is V_IP3 in component receptor (microMpersec).
C CONSTS(6) is MechanicalStimulation in component receptor (dimensionless).
C ALGBRC(7) is V_mech in component receptor (dimensionless).
C CONSTS(7) is V_mech0 in component receptor (dimensionless).
C CONSTS(8) is C_1 in component receptor (microMpersec).
C CONSTS(9) is C_2 in component receptor (microM).
C CONSTS(10) is C_3 in component receptor (microMpersec).
C CONSTS(11) is C_4 in component receptor (microM).
C CONSTS(12) is V_ATP in component receptor (microMpersec).
C CONSTS(13) is K_ATP in component receptor (microM).
C CONSTS(14) is Beta_1 in component J_ip3D (microMpersec).
C CONSTS(15) is Beta_2 in component J_ip3D (microMpersec).
C ALGBRC(15) is F_Ca in component J_ip3D (dimensionless).
C CONSTS(16) is K_rc in component J_ip3D (microM).
C CONSTS(17) is K_IP3 in component J_ip3D (microM).
C ALGBRC(11) is Pipt in component Jipt (dimensionless).
C STATES(6) is O in component O (dimensionless).
C CONSTS(18) is Jipt0 in component Jipt (per_sec).
C CONSTS(19) is Ka in component Jryr (microM).
C CONSTS(20) is Kb in component Jryr (microM).
C CONSTS(21) is Kc in component Jryr (dimensionless).
C CONSTS(22) is Jryr0 in component Jryr (per_sec).
C ALGBRC(6) is Pryr in component Jryr (dimensionless).
C ALGBRC(3) is W_inf in component Jryr (dimensionless).
C CONSTS(23) is Jer0 in component Jer (per_sec).
C CONSTS(24) is Vserca in component Jserca (microM2persec).
C CONSTS(25) is Kserca in component Jserca (microM).
C CONSTS(26) is Vm in component Jin (milliV).
C ALGBRC(20) is V_Ca in component Jin (milliV).
C CONSTS(27) is R in component Jin (millijoule_per_mole_kelvin).
C CONSTS(28) is T in component Jin (kelvin).
C CONSTS(29) is F in component Jin (coulomb_per_mole).
C CONSTS(30) is z_Ca in component Jin (dimensionless).
C CONSTS(31) is Ca_ext in component Jin (microM).
C CONSTS(32) is P_pm_1 in component Jin (conductancepervolume).
C CONSTS(33) is P_pm_2 in component Jin (conductancepervolume).
C ALGBRC(21) is Iin_1 in component Jin (current).
C ALGBRC(22) is Iin_2 in component Jin (current).
C ALGBRC(23) is Jin_1 in component Jin (microMvolumepersec).
C ALGBRC(24) is Jin_2 in component Jin (microMvolumepersec).
C CONSTS(34) is Vpm in component Jpm (microMpersec).
C CONSTS(35) is Kpm in component Jpm (microM).
C ALGBRC(8) is S in component O (dimensionless).
C STATES(7) is I1 in component I1 (dimensionless).
C STATES(8) is I2 in component I2 (dimensionless).
C ALGBRC(12) is k1 in component constants (per_secmicroM).
C CONSTS(36) is k1_ in component constants (first_order_rate_constant).
C CONSTS(37) is k2 in component constants (first_order_rate_constant).
C CONSTS(38) is k3 in component constants (first_order_rate_constant).
C ALGBRC(9) is k4 in component constants (first_order_rate_constant).
C CONSTS(39) is k5 in component constants (first_order_rate_constant).
C CONSTS(40) is alpha_1 in component constants (per_secmicroM).
C CONSTS(41) is beta_1 in component constants (microM).
C CONSTS(42) is alpha_4 in component constants (first_order_rate_constant).
C CONSTS(43) is beta_4 in component constants (microM).
C RATES(1) is d/dt volumeCa_i in component volumeCa_i (microMVolume).
C RATES(2) is d/dt volumeCa_er in component volumeCa_er (microMVolume).
C RATES(3) is d/dt volume_iIP3 in component volume_iIP3 (microMVolume).
C RATES(4) is d/dt PIP2 in component PIP2 (dimensionless).
C RATES(5) is d/dt ATP_e in component receptor (microM).
C RATES(6) is d/dt O in component O (dimensionless).
C RATES(7) is d/dt I1 in component I1 (dimensionless).
C RATES(8) is d/dt I2 in component I2 (dimensionless).
C
      SUBROUTINE initConsts(CONSTS, RATES, STATES)
      REAL CONSTS(*), RATES(*), STATES(*)
      CONSTS(1) = 1450
      CONSTS(2) = 4.9
      STATES(1) = 75.8335890358239
      STATES(2) = 1199.73747802144
      STATES(3) = -3.28287434380887e-93
      STATES(4) = 1
      CONSTS(3) = 5
      CONSTS(4) = 25
      CONSTS(5) = 10
      STATES(5) = 0
      CONSTS(6) = 0
      CONSTS(7) = 1
      CONSTS(8) = 2.6
      CONSTS(9) = 0.02
      CONSTS(10) = 2.3
      CONSTS(11) = 10
      CONSTS(12) = 0.001
      CONSTS(13) = 0.0001
      CONSTS(14) = 0
      CONSTS(15) = 36
      CONSTS(16) = 0.3
      CONSTS(17) = 0.01
      STATES(6) = -2.42013500345263e-89
      CONSTS(18) = 28
      CONSTS(19) = 0.3722419
      CONSTS(20) = 0.636005
      CONSTS(21) = 0.0571428
      CONSTS(22) = 0
      CONSTS(23) = 0.0035
      CONSTS(24) = 0.09
      CONSTS(25) = 0.04
      CONSTS(26) = -39
      CONSTS(27) = 8314
      CONSTS(28) = 310
      CONSTS(29) = 96485
      CONSTS(30) = 2
      CONSTS(31) = 1300
      CONSTS(32) = 0.6
      CONSTS(33) = 28
      CONSTS(34) = 0.072
      CONSTS(35) = 0.6
      STATES(7) = 2.24343309680416e-81
      STATES(8) = 1.01276768045278e-40
      CONSTS(36) = 0.88
      CONSTS(37) = 0.5
      CONSTS(38) = 0.5
      CONSTS(39) = 0.02
      CONSTS(40) = 40
      CONSTS(41) = 0.8
      CONSTS(42) = 0.06
      CONSTS(43) = 0.01
      CONSTS(44) = CONSTS(1)
      CONSTS(45) = CONSTS(44)/CONSTS(2)
      CONSTS(46) = CONSTS(5)
      RETURN
      END
      SUBROUTINE computeRates(VOI, CONSTS,  RATES, STATES, ALGBRC)
      REAL VOI, CONSTS(*), RATES(*), STATES(*), ALGBRC(*)
      RATES(5) = ( - CONSTS(12)*STATES(5))/(CONSTS(13)+STATES(5))
      ALGBRC(2) = STATES(3)/CONSTS(44)
      ALGBRC(9) = ( CONSTS(42)*ALGBRC(2))/(CONSTS(43)+ALGBRC(2))
      RATES(7) =  CONSTS(37)*STATES(6) -  (CONSTS(38)+ALGBRC(9))*STATES(7)
      RATES(8) =  ALGBRC(9)*STATES(7) -  CONSTS(39)*STATES(8)
      ALGBRC(8) = 1.00000 - (STATES(6)+STATES(7)+STATES(8))
      ALGBRC(1) = STATES(1)/CONSTS(44)
      ALGBRC(12) = ( CONSTS(40)*ALGBRC(1) ** 3.00000)/(CONSTS(41) ** 3.00000+ALGBRC(1) ** 3.00000)
      RATES(6) =  ALGBRC(12)*ALGBRC(2)*ALGBRC(8) - ( CONSTS(36)*STATES(6)+ CONSTS(37)*STATES(6))
      ALGBRC(10) = ( CONSTS(8)*STATES(5))/(CONSTS(9)+STATES(5))+( CONSTS(10)*STATES(5) ** 2.00000)/(CONSTS(11) ** 2.00000+STATES(5) ** 2.00000)
      ALGBRC(13) =  ALGBRC(10)*STATES(4)
      ALGBRC(5) = (1.00000 - ALGBRC(2)/CONSTS(4)) - STATES(4)
      RATES(4) =  CONSTS(46)*ALGBRC(5) - ALGBRC(13)/CONSTS(3)
      ALGBRC(15) = ALGBRC(1)/(CONSTS(16)+ALGBRC(1))
      ALGBRC(17) = ( (CONSTS(14)+ CONSTS(15)*ALGBRC(15) ** 2.00000)*ALGBRC(2))/(CONSTS(17)+ALGBRC(2))
      RATES(3) =  (ALGBRC(13) - ALGBRC(17))*CONSTS(44)
      ALGBRC(4) = STATES(2)/CONSTS(45)
      ALGBRC(11) = STATES(6) ** 4.00000
      ALGBRC(14) =  CONSTS(18)*CONSTS(1)*ALGBRC(11)*(ALGBRC(4) - ALGBRC(1))
      ALGBRC(16) =  CONSTS(22)*CONSTS(1)*(ALGBRC(4) - ALGBRC(1))
      ALGBRC(18) =  CONSTS(23)*CONSTS(1)*(ALGBRC(4) - ALGBRC(1))
      ALGBRC(19) = ( (( CONSTS(1)*CONSTS(24)*ALGBRC(1) ** 2.00000)/(CONSTS(25) ** 2.00000+ALGBRC(1) ** 2.00000))*1.00000)/ALGBRC(4)
      RATES(2) = - ((ALGBRC(14)+ALGBRC(16)+ALGBRC(18)) - ALGBRC(19))
      ALGBRC(20) =  (( CONSTS(27)*CONSTS(28))/( CONSTS(30)*CONSTS(29)))*log(CONSTS(31)/ALGBRC(1))
      ALGBRC(21) =  CONSTS(1)*CONSTS(32)*(ALGBRC(20) - CONSTS(26))
      ALGBRC(23) = ALGBRC(21)/( CONSTS(30)*CONSTS(29))
      ALGBRC(7) = TERNRY(VOI.GT.10.0000.AND.VOI.LT.25.0000.AND.CONSTS(6).GT.0.00000, CONSTS(7), 0.00000)
      ALGBRC(22) =  CONSTS(1)*CONSTS(33)*ALGBRC(7)*(ALGBRC(20) - CONSTS(26))
      ALGBRC(24) = ALGBRC(22)/( CONSTS(30)*CONSTS(29))
      ALGBRC(25) = ALGBRC(23)+ALGBRC(24)
      ALGBRC(26) = ( CONSTS(34)*CONSTS(1)*ALGBRC(1) ** 2.00000)/(CONSTS(35) ** 2.00000+ALGBRC(1) ** 2.00000)
      RATES(1) = ((ALGBRC(14)+ALGBRC(16)+ALGBRC(18)) - ALGBRC(19))+(ALGBRC(25) - ALGBRC(26))
      RETURN
      END
      SUBROUTINE computeVariables(VOI, CONSTS, RATES, STATES, ALGBRC)
      REAL VOI, CONSTS(*), RATES(*), STATES(*), ALGBRC(*)
      ALGBRC(2) = STATES(3)/CONSTS(44)
      ALGBRC(9) = ( CONSTS(42)*ALGBRC(2))/(CONSTS(43)+ALGBRC(2))
      ALGBRC(8) = 1.00000 - (STATES(6)+STATES(7)+STATES(8))
      ALGBRC(1) = STATES(1)/CONSTS(44)
      ALGBRC(12) = ( CONSTS(40)*ALGBRC(1) ** 3.00000)/(CONSTS(41) ** 3.00000+ALGBRC(1) ** 3.00000)
      ALGBRC(10) = ( CONSTS(8)*STATES(5))/(CONSTS(9)+STATES(5))+( CONSTS(10)*STATES(5) ** 2.00000)/(CONSTS(11) ** 2.00000+STATES(5) ** 2.00000)
      ALGBRC(13) =  ALGBRC(10)*STATES(4)
      ALGBRC(5) = (1.00000 - ALGBRC(2)/CONSTS(4)) - STATES(4)
      ALGBRC(15) = ALGBRC(1)/(CONSTS(16)+ALGBRC(1))
      ALGBRC(17) = ( (CONSTS(14)+ CONSTS(15)*ALGBRC(15) ** 2.00000)*ALGBRC(2))/(CONSTS(17)+ALGBRC(2))
      ALGBRC(4) = STATES(2)/CONSTS(45)
      ALGBRC(11) = STATES(6) ** 4.00000
      ALGBRC(14) =  CONSTS(18)*CONSTS(1)*ALGBRC(11)*(ALGBRC(4) - ALGBRC(1))
      ALGBRC(16) =  CONSTS(22)*CONSTS(1)*(ALGBRC(4) - ALGBRC(1))
      ALGBRC(18) =  CONSTS(23)*CONSTS(1)*(ALGBRC(4) - ALGBRC(1))
      ALGBRC(19) = ( (( CONSTS(1)*CONSTS(24)*ALGBRC(1) ** 2.00000)/(CONSTS(25) ** 2.00000+ALGBRC(1) ** 2.00000))*1.00000)/ALGBRC(4)
      ALGBRC(20) =  (( CONSTS(27)*CONSTS(28))/( CONSTS(30)*CONSTS(29)))*log(CONSTS(31)/ALGBRC(1))
      ALGBRC(21) =  CONSTS(1)*CONSTS(32)*(ALGBRC(20) - CONSTS(26))
      ALGBRC(23) = ALGBRC(21)/( CONSTS(30)*CONSTS(29))
      ALGBRC(7) = TERNRY(VOI.GT.10.0000.AND.VOI.LT.25.0000.AND.CONSTS(6).GT.0.00000, CONSTS(7), 0.00000)
      ALGBRC(22) =  CONSTS(1)*CONSTS(33)*ALGBRC(7)*(ALGBRC(20) - CONSTS(26))
      ALGBRC(24) = ALGBRC(22)/( CONSTS(30)*CONSTS(29))
      ALGBRC(25) = ALGBRC(23)+ALGBRC(24)
      ALGBRC(26) = ( CONSTS(34)*CONSTS(1)*ALGBRC(1) ** 2.00000)/(CONSTS(35) ** 2.00000+ALGBRC(1) ** 2.00000)
      ALGBRC(3) = (1.00000+CONSTS(19)/ALGBRC(1) ** 4.00000+ALGBRC(1)/CONSTS(20) ** 3.00000)/(1.00000+1.00000/CONSTS(21)+ALGBRC(1)/CONSTS(20) ** 3.00000+CONSTS(19)/ALGBRC(1) ** 4.00000)
      ALGBRC(6) = ( ALGBRC(3)*(1.00000+CONSTS(19)/ALGBRC(1) ** 4.00000+ALGBRC(1)/CONSTS(20) ** 3.00000))/(1.00000+1.00000/CONSTS(21)+ALGBRC(1)/CONSTS(20) ** 3.00000+CONSTS(19)/ALGBRC(1) ** 4.00000)
      RETURN
      END
      REAL FUNCTION TERNRY(TEST, VALA, VALB)
      LOGICAL TEST
      REAL VALA, VALB
      IF (TEST) THEN
        TERNRY = VALA
      ELSE
        TERNRY = VALB
      ENDIF
      RETURN
      END