C
C There are a total of 50 entries in the algebraic variable array.
C There are a total of 11 entries in each of the rate and state variable arrays.
C There are a total of 71 entries in the constant variable array.
C
C
C VOI is time in component environment (millisecond).
C ALGBRC(4) is SOVFThick in component sarcomere_geometry (dimensionless).
C ALGBRC(5) is SOVFThin in component sarcomere_geometry (dimensionless).
C ALGBRC(1) is sovr_ze in component sarcomere_geometry (micrometre).
C ALGBRC(2) is sovr_cle in component sarcomere_geometry (micrometre).
C ALGBRC(3) is len_sovr in component sarcomere_geometry (micrometre).
C CONSTS(1) is len_thin in component model_parameters (micrometre).
C CONSTS(2) is len_thick in component model_parameters (micrometre).
C CONSTS(3) is len_hbare in component model_parameters (micrometre).
C STATES(1) is SL in component normalised_active_and_passive_force (micrometre).
C STATES(2) is TRPNCaL in component Ca_binding_to_troponin_to_thin_filament_regulation (dimensionless).
C STATES(3) is TRPNCaH in component Ca_binding_to_troponin_to_thin_filament_regulation (dimensionless).
C ALGBRC(8) is dTRPNCaL in component Ca_binding_to_troponin_to_thin_filament_regulation (first_order_rate_constant).
C ALGBRC(11) is dTRPNCaH in component Ca_binding_to_troponin_to_thin_filament_regulation (first_order_rate_constant).
C ALGBRC(13) is kn_pT in component Ca_binding_to_troponin_to_thin_filament_regulation (first_order_rate_constant).
C ALGBRC(19) is kp_nT in component Ca_binding_to_troponin_to_thin_filament_regulation (first_order_rate_constant).
C CONSTS(61) is H in component Ca_binding_to_troponin_to_thin_filament_regulation (micromolar).
C CONSTS(63) is H_cons in component Ca_binding_to_troponin_to_thin_filament_regulation (micromolar).
C CONSTS(64) is konT in component Ca_binding_to_troponin_to_thin_filament_regulation (second_order_rate_constant).
C CONSTS(59) is koffLT in component Ca_binding_to_troponin_to_thin_filament_regulation (first_order_rate_constant).
C CONSTS(60) is koffHT in component Ca_binding_to_troponin_to_thin_filament_regulation (first_order_rate_constant).
C CONSTS(4) is Qkon in component Ca_binding_to_troponin_to_thin_filament_regulation (dimensionless).
C CONSTS(5) is Qkoff in component Ca_binding_to_troponin_to_thin_filament_regulation (dimensionless).
C CONSTS(6) is Qkn_p in component Ca_binding_to_troponin_to_thin_filament_regulation (dimensionless).
C CONSTS(7) is Qkp_n in component Ca_binding_to_troponin_to_thin_filament_regulation (dimensionless).
C CONSTS(8) is kon in component Ca_binding_to_troponin_to_thin_filament_regulation (second_order_rate_constant).
C CONSTS(9) is koffL in component Ca_binding_to_troponin_to_thin_filament_regulation (first_order_rate_constant).
C CONSTS(10) is koffH in component Ca_binding_to_troponin_to_thin_filament_regulation (first_order_rate_constant).
C CONSTS(11) is perm50 in component Ca_binding_to_troponin_to_thin_filament_regulation (dimensionless).
C CONSTS(12) is nperm in component Ca_binding_to_troponin_to_thin_filament_regulation (dimensionless).
C CONSTS(13) is kn_p in component Ca_binding_to_troponin_to_thin_filament_regulation (first_order_rate_constant).
C CONSTS(14) is kp_n in component Ca_binding_to_troponin_to_thin_filament_regulation (first_order_rate_constant).
C CONSTS(15) is koffmod in component Ca_binding_to_troponin_to_thin_filament_regulation (dimensionless).
C ALGBRC(7) is Tropreg in component Ca_binding_to_troponin_to_thin_filament_regulation (dimensionless).
C ALGBRC(10) is permtot in component Ca_binding_to_troponin_to_thin_filament_regulation (dimensionless).
C ALGBRC(16) is inprmt in component Ca_binding_to_troponin_to_thin_filament_regulation (dimensionless).
C CONSTS(16) is pH in component Ca_binding_to_troponin_to_thin_filament_regulation (dimensionless).
C CONSTS(17) is m in component Ca_binding_to_troponin_to_thin_filament_regulation (dimensionless).
C CONSTS(18) is kdHCa in component Ca_binding_to_troponin_to_thin_filament_regulation (micromolar).
C CONSTS(19) is TmpC in component model_parameters (celsius).
C CONSTS(20) is Cai in component model_parameters (micromolar).
C CONSTS(65) is fappT in component thin_filament_regulation_and_crossbridge_cycling_rates (first_order_rate_constant).
C ALGBRC(17) is gappT in component thin_filament_regulation_and_crossbridge_cycling_rates (first_order_rate_constant).
C ALGBRC(24) is hfT in component thin_filament_regulation_and_crossbridge_cycling_rates (first_order_rate_constant).
C ALGBRC(25) is hbT in component thin_filament_regulation_and_crossbridge_cycling_rates (first_order_rate_constant).
C ALGBRC(27) is gxbT in component thin_filament_regulation_and_crossbridge_cycling_rates (first_order_rate_constant).
C CONSTS(21) is fapp in component thin_filament_regulation_and_crossbridge_cycling_rates (first_order_rate_constant).
C CONSTS(22) is gapp in component thin_filament_regulation_and_crossbridge_cycling_rates (first_order_rate_constant).
C CONSTS(23) is hf in component thin_filament_regulation_and_crossbridge_cycling_rates (first_order_rate_constant).
C CONSTS(24) is hb in component thin_filament_regulation_and_crossbridge_cycling_rates (first_order_rate_constant).
C CONSTS(25) is gxb in component thin_filament_regulation_and_crossbridge_cycling_rates (first_order_rate_constant).
C CONSTS(26) is gslmod in component thin_filament_regulation_and_crossbridge_cycling_rates (dimensionless).
C ALGBRC(20) is hfmd in component thin_filament_regulation_and_crossbridge_cycling_rates (dimensionless).
C ALGBRC(22) is hbmd in component thin_filament_regulation_and_crossbridge_cycling_rates (dimensionless).
C CONSTS(27) is hfmdc in component thin_filament_regulation_and_crossbridge_cycling_rates (dimensionless).
C CONSTS(28) is hbmdc in component thin_filament_regulation_and_crossbridge_cycling_rates (dimensionless).
C CONSTS(29) is sigmap in component thin_filament_regulation_and_crossbridge_cycling_rates (dimensionless).
C CONSTS(30) is sigman in component thin_filament_regulation_and_crossbridge_cycling_rates (dimensionless).
C CONSTS(31) is xbmodsp in component thin_filament_regulation_and_crossbridge_cycling_rates (dimensionless).
C CONSTS(32) is Qfapp in component thin_filament_regulation_and_crossbridge_cycling_rates (dimensionless).
C CONSTS(33) is Qgapp in component thin_filament_regulation_and_crossbridge_cycling_rates (dimensionless).
C CONSTS(34) is Qhf in component thin_filament_regulation_and_crossbridge_cycling_rates (dimensionless).
C CONSTS(35) is Qhb in component thin_filament_regulation_and_crossbridge_cycling_rates (dimensionless).
C CONSTS(36) is Qgxb in component thin_filament_regulation_and_crossbridge_cycling_rates (dimensionless).
C ALGBRC(26) is gxbmd in component thin_filament_regulation_and_crossbridge_cycling_rates (dimensionless).
C ALGBRC(14) is gapslmd in component thin_filament_regulation_and_crossbridge_cycling_rates (dimensionless).
C CONSTS(37) is x_0 in component model_parameters (micrometre).
C STATES(4) is xXBpostr in component mean_strain_of_strongly_bound_states (micrometre).
C STATES(5) is xXBprer in component mean_strain_of_strongly_bound_states (micrometre).
C STATES(6) is XBpostr in component regulation_and_crossbridge_cycling_state_equations (dimensionless).
C STATES(7) is XBprer in component regulation_and_crossbridge_cycling_state_equations (dimensionless).
C ALGBRC(35) is dXBpostr in component regulation_and_crossbridge_cycling_state_equations (first_order_rate_constant).
C ALGBRC(32) is dXBprer in component regulation_and_crossbridge_cycling_state_equations (first_order_rate_constant).
C CONSTS(66) is alpha1_plus in component regulation_and_crossbridge_cycling_state_equations (first_order_rate_constant).
C ALGBRC(28) is alpha2_plus in component regulation_and_crossbridge_cycling_state_equations (first_order_rate_constant).
C ALGBRC(29) is alpha3_plus in component regulation_and_crossbridge_cycling_state_equations (first_order_rate_constant).
C ALGBRC(30) is alpha1_minus in component regulation_and_crossbridge_cycling_state_equations (first_order_rate_constant).
C ALGBRC(31) is alpha2_minus in component regulation_and_crossbridge_cycling_state_equations (first_order_rate_constant).
C ALGBRC(34) is alpha3_minus in component regulation_and_crossbridge_cycling_state_equations (first_order_rate_constant).
C CONSTS(38) is kMgATP in component regulation_and_crossbridge_cycling_state_equations (micromolar2).
C CONSTS(39) is kdADP in component regulation_and_crossbridge_cycling_state_equations (micromolar).
C CONSTS(40) is xPi_cons in component regulation_and_crossbridge_cycling_state_equations (micromolar).
C CONSTS(41) is MgATP_cons in component regulation_and_crossbridge_cycling_state_equations (micromolar).
C ALGBRC(33) is fxbT in component regulation_and_crossbridge_cycling_state_equations (first_order_rate_constant).
C STATES(8) is N_NoXB in component regulation_and_crossbridge_cycling_state_equations (dimensionless).
C STATES(9) is P_NoXB in component regulation_and_crossbridge_cycling_state_equations (dimensionless).
C STATES(10) is P in component regulation_and_crossbridge_cycling_state_equations (dimensionless).
C ALGBRC(23) is N in component regulation_and_crossbridge_cycling_state_equations (dimensionless).
C CONSTS(42) is MgADP_cons in component regulation_and_crossbridge_cycling_state_equations (micromolar).
C CONSTS(43) is xPi in component regulation_and_crossbridge_cycling_state_equations (micromolar).
C CONSTS(44) is MgATP in component regulation_and_crossbridge_cycling_state_equations (micromolar).
C CONSTS(45) is MgADP in component regulation_and_crossbridge_cycling_state_equations (micromolar).
C ALGBRC(39) is dxXBpostr in component mean_strain_of_strongly_bound_states (micrometre_per_millisecond).
C ALGBRC(37) is dxXBprer in component mean_strain_of_strongly_bound_states (micrometre_per_millisecond).
C CONSTS(46) is xPsi in component mean_strain_of_strongly_bound_states (dimensionless).
C ALGBRC(36) is dutyprer in component mean_strain_of_strongly_bound_states (dimensionless).
C ALGBRC(38) is dutypostr in component mean_strain_of_strongly_bound_states (dimensionless).
C CONSTS(62) is dSL in component normalised_active_and_passive_force (micrometre_per_millisecond).
C CONSTS(68) is SSXBpostr in component normalised_active_and_passive_force (dimensionless).
C ALGBRC(40) is SSXBprer in component normalised_active_and_passive_force (dimensionless).
C CONSTS(47) is kxb in component normalised_active_and_passive_force (millinewton_per_millimetre2).
C CONSTS(69) is Fnordv in component normalised_active_and_passive_force (millinewton_micrometre_per_millimetre2).
C ALGBRC(6) is force in component normalised_active_and_passive_force (millinewton_micrometre_per_millimetre2).
C ALGBRC(9) is active in component normalised_active_and_passive_force (unit_normalised_force).
C ALGBRC(18) is ppforce in component normalised_active_and_passive_force (unit_normalised_force).
C ALGBRC(12) is ppforce_t in component normalised_active_and_passive_force (unit_normalised_force).
C ALGBRC(15) is ppforce_c in component normalised_active_and_passive_force (unit_normalised_force).
C CONSTS(70) is preload in component normalised_active_and_passive_force (unit_normalised_force).
C ALGBRC(21) is afterload in component normalised_active_and_passive_force (unit_normalised_force).
C STATES(11) is intf in component normalised_active_and_passive_force (unit_normalised_force_millisecond).
C CONSTS(48) is SL_c in component normalised_active_and_passive_force (micrometre).
C CONSTS(49) is SLrest in component normalised_active_and_passive_force (micrometre).
C CONSTS(50) is SLset in component normalised_active_and_passive_force (micrometre).
C CONSTS(51) is PCon_t in component normalised_active_and_passive_force (unit_normalised_force).
C CONSTS(52) is PExp_t in component normalised_active_and_passive_force (per_micrometre).
C CONSTS(53) is PCon_c in component normalised_active_and_passive_force (unit_normalised_force).
C CONSTS(54) is PExp_c in component normalised_active_and_passive_force (per_micrometre).
C CONSTS(55) is KSE in component normalised_active_and_passive_force (unit_normalised_force_per_micrometre).
C CONSTS(67) is fxb in component normalised_active_and_passive_force (first_order_rate_constant).
C CONSTS(56) is SEon in component normalised_active_and_passive_force (dimensionless).
C ALGBRC(41) is FrSBXB in component calculation_of_micromolar_per_millisecondes_of_Ca_for_apparent_Ca_binding (dimensionless).
C ALGBRC(42) is dFrSBXB in component calculation_of_micromolar_per_millisecondes_of_Ca_for_apparent_Ca_binding (first_order_rate_constant).
C ALGBRC(44) is dsovr_ze in component calculation_of_micromolar_per_millisecondes_of_Ca_for_apparent_Ca_binding (micrometre_per_millisecond).
C ALGBRC(45) is dsovr_cle in component calculation_of_micromolar_per_millisecondes_of_Ca_for_apparent_Ca_binding (micrometre_per_millisecond).
C ALGBRC(46) is dlen_sovr in component calculation_of_micromolar_per_millisecondes_of_Ca_for_apparent_Ca_binding (micrometre_per_millisecond).
C ALGBRC(47) is dSOVFThick in component calculation_of_micromolar_per_millisecondes_of_Ca_for_apparent_Ca_binding (first_order_rate_constant).
C ALGBRC(48) is dSOVFThin in component calculation_of_micromolar_per_millisecondes_of_Ca_for_apparent_Ca_binding (first_order_rate_constant).
C CONSTS(57) is kxb in component calculation_of_micromolar_per_millisecondes_of_Ca_for_apparent_Ca_binding (millinewton_per_millimetre2).
C ALGBRC(49) is dforce in component calculation_of_micromolar_per_millisecondes_of_Ca_for_apparent_Ca_binding (millinewton_micrometre_per_millimetre2_per_millisecond).
C CONSTS(58) is Trop_conc in component calculation_of_micromolar_per_millisecondes_of_Ca_for_apparent_Ca_binding (micromolar).
C ALGBRC(43) is TropTot in component calculation_of_micromolar_per_millisecondes_of_Ca_for_apparent_Ca_binding (micromolar).
C ALGBRC(50) is dTropTot in component calculation_of_micromolar_per_millisecondes_of_Ca_for_apparent_Ca_binding (micromolar_per_millisecond).
C RATES(2) is d/dt TRPNCaL in component Ca_binding_to_troponin_to_thin_filament_regulation (dimensionless).
C RATES(3) is d/dt TRPNCaH in component Ca_binding_to_troponin_to_thin_filament_regulation (dimensionless).
C RATES(8) is d/dt N_NoXB in component regulation_and_crossbridge_cycling_state_equations (dimensionless).
C RATES(9) is d/dt P_NoXB in component regulation_and_crossbridge_cycling_state_equations (dimensionless).
C RATES(10) is d/dt P in component regulation_and_crossbridge_cycling_state_equations (dimensionless).
C RATES(7) is d/dt XBprer in component regulation_and_crossbridge_cycling_state_equations (dimensionless).
C RATES(6) is d/dt XBpostr in component regulation_and_crossbridge_cycling_state_equations (dimensionless).
C RATES(5) is d/dt xXBprer in component mean_strain_of_strongly_bound_states (micrometre).
C RATES(4) is d/dt xXBpostr in component mean_strain_of_strongly_bound_states (micrometre).
C RATES(1) is d/dt SL in component normalised_active_and_passive_force (micrometre).
C RATES(11) is d/dt intf in component normalised_active_and_passive_force (unit_normalised_force_millisecond).
C
      SUBROUTINE initConsts(CONSTS, RATES, STATES)
      REAL CONSTS(*), RATES(*), STATES(*)
      CONSTS(1) = 1.2
      CONSTS(2) = 1.65
      CONSTS(3) = 0.1
      STATES(1) = 2.2
      STATES(2) = 0.0147730085063734
      STATES(3) = 0.13066096561522
      CONSTS(4) = 1.5
      CONSTS(5) = 1.3
      CONSTS(6) = 1.6
      CONSTS(7) = 1.6
      CONSTS(8) = 0.05
      CONSTS(9) = 0.25
      CONSTS(10) = 0.025
      CONSTS(11) = 0.5
      CONSTS(12) = 15
      CONSTS(13) = 0.5
      CONSTS(14) = 0.05
      CONSTS(15) = 1
      CONSTS(16) = 7.15
      CONSTS(17) = 1
      CONSTS(18) = 2e-2
      CONSTS(19) = 22
      CONSTS(20) = 200.0
      CONSTS(21) = 0.5
      CONSTS(22) = 0.07
      CONSTS(23) = 2
      CONSTS(24) = 0.4
      CONSTS(25) = 0.07
      CONSTS(26) = 6
      CONSTS(27) = 5
      CONSTS(28) = 0
      CONSTS(29) = 8
      CONSTS(30) = 1
      CONSTS(31) = 0.2
      CONSTS(32) = 6.25
      CONSTS(33) = 2.5
      CONSTS(34) = 6.25
      CONSTS(35) = 6.25
      CONSTS(36) = 6.25
      CONSTS(37) = 0.007
      STATES(4) = 0.00700005394873882
      STATES(5) = 3.41212828972468e-8
      STATES(6) = 1.81017564383744e-6
      STATES(7) = 3.0494964880038e-7
      CONSTS(38) = 15400e6
      CONSTS(39) = 4
      CONSTS(40) = 2e3
      CONSTS(41) = 5e3
      STATES(8) = 0.999999959256274
      STATES(9) = 4.07437173988636e-8
      STATES(10) = 0.999997834540066
      CONSTS(42) = 36
      CONSTS(43) = 2e3
      CONSTS(44) = 5e3
      CONSTS(45) = 36.3
      CONSTS(46) = 2
      CONSTS(47) = 120
      STATES(11) = -4.5113452510363e-6
      CONSTS(48) = 2.25
      CONSTS(49) = 1.85
      CONSTS(50) = 1.9
      CONSTS(51) = 0.002
      CONSTS(52) = 10
      CONSTS(53) = 0.02
      CONSTS(54) = 70
      CONSTS(55) = 1
      CONSTS(56) = 1
      CONSTS(57) = 120
      CONSTS(58) = 70
      CONSTS(59) =  CONSTS(9)*CONSTS(15)*CONSTS(5) ** (CONSTS(19) - 37.0000)/10.0000
      CONSTS(60) =  CONSTS(10)*CONSTS(15)*CONSTS(5) ** (CONSTS(19) - 37.0000)/10.0000
      CONSTS(61) =  1.00000e+07*10.0000 ** - CONSTS(16)
      CONSTS(62) = 0.00000
      CONSTS(63) =  1.00000e+07*10.0000 ** -7.15000
      CONSTS(70) = CONSTS(62)
      CONSTS(64) =  ((CONSTS(18) ** CONSTS(17)+CONSTS(63) ** CONSTS(17))/(CONSTS(18) ** CONSTS(17)+CONSTS(61) ** CONSTS(17)))*( CONSTS(8)*CONSTS(4) ** (CONSTS(19) - 37.0000)/10.0000)
      CONSTS(65) =  CONSTS(21)*CONSTS(31)*CONSTS(32) ** (CONSTS(19) - 37.0000)/10.0000
      CONSTS(66) = CONSTS(65)
      CONSTS(67) = ( CONSTS(39)*CONSTS(21)*CONSTS(23)*(CONSTS(25)/CONSTS(41)))/( (CONSTS(22)/CONSTS(40))*(CONSTS(24)/CONSTS(63))*CONSTS(38))
      CONSTS(68) = ( CONSTS(21)*CONSTS(23)+ CONSTS(67)*CONSTS(22)+ CONSTS(67)*CONSTS(24))/( CONSTS(23)*CONSTS(25)+ CONSTS(24)*CONSTS(22)+ CONSTS(25)*CONSTS(22)+ CONSTS(67)*CONSTS(24)+ CONSTS(25)*CONSTS(21)+ CONSTS(24)*CONSTS(21)+ CONSTS(67)*CONSTS(22)+ CONSTS(21)*CONSTS(23)+ CONSTS(67)*CONSTS(24))
      CONSTS(69) =  CONSTS(47)*CONSTS(37)*CONSTS(68)
      CONSTS(70) =  (ABS(CONSTS(50) - CONSTS(49))/(CONSTS(50) - CONSTS(49)))*CONSTS(51)*(EXP( CONSTS(52)*ABS(CONSTS(50) - CONSTS(49))) - 1.00000)
      RETURN
      END
      SUBROUTINE computeRates(VOI, CONSTS,  RATES, STATES, ALGBRC)
      REAL VOI, CONSTS(*), RATES(*), STATES(*), ALGBRC(*)
      RATES(1) = CONSTS(70)
      ALGBRC(8) =  CONSTS(64)*CONSTS(20)*(1.00000 - STATES(2)) -  CONSTS(59)*STATES(2)
      RATES(2) = ALGBRC(8)
      ALGBRC(11) =  CONSTS(64)*CONSTS(20)*(1.00000 - STATES(3)) -  CONSTS(60)*STATES(3)
      RATES(3) = ALGBRC(11)
      ALGBRC(1) = TERNRY(CONSTS(2)/2.00000.LT.STATES(1)/2.00000, CONSTS(2)/2.00000, STATES(1)/2.00000)
      ALGBRC(2) = TERNRY(STATES(1)/2.00000 - (STATES(1) - CONSTS(1)).GT.CONSTS(3)/2.00000, STATES(1)/2.00000 - (STATES(1) - CONSTS(1)), CONSTS(3)/2.00000)
      ALGBRC(3) = ALGBRC(1) - ALGBRC(2)
      ALGBRC(5) = ALGBRC(3)/CONSTS(1)
      ALGBRC(7) =  (1.00000 - ALGBRC(5))*STATES(2)+ ALGBRC(5)*STATES(3)
      ALGBRC(10) =  ABS(1.00000/(1.00000+CONSTS(11)/ALGBRC(7) ** CONSTS(12))) ** (1.0 / 2)
      ALGBRC(13) =  CONSTS(13)*ALGBRC(10)*CONSTS(6) ** (CONSTS(19) - 37.0000)/10.0000
      ALGBRC(16) = TERNRY(1.00000/ALGBRC(10).LT.100.000, 1.00000/ALGBRC(10), 100.000)
      ALGBRC(19) =  CONSTS(14)*ALGBRC(16)*CONSTS(7) ** (CONSTS(19) - 37.0000)/10.0000
      RATES(8) =  ALGBRC(19)*STATES(9) -  ALGBRC(13)*STATES(8)
      RATES(9) =  ALGBRC(13)*STATES(8) -  ALGBRC(19)*STATES(9)
      ALGBRC(4) = ( ALGBRC(3)*2.00000)/(CONSTS(2) - CONSTS(3))
      ALGBRC(6) =  CONSTS(47)*ALGBRC(4)*( STATES(4)*STATES(6)+ STATES(5)*STATES(7))
      ALGBRC(9) = ( 1.00000*ALGBRC(6))/CONSTS(69)
      ALGBRC(12) =  ((STATES(1) - CONSTS(49))/ABS(STATES(1) - CONSTS(49)))*CONSTS(51)*(EXP( CONSTS(52)*ABS(STATES(1) - CONSTS(49))) - 1.00000)
      ALGBRC(15) = TERNRY(STATES(1).GT.CONSTS(48),  CONSTS(53)*(EXP( CONSTS(54)*ABS(STATES(1) - CONSTS(48))) - 1.00000), 0.00000)
      ALGBRC(18) = ALGBRC(12)+ALGBRC(15)
      ALGBRC(21) = TERNRY(CONSTS(56).EQ.1.00000,  CONSTS(55)*(CONSTS(50) - STATES(1)), 0.00000)
      RATES(11) = (CONSTS(70)+ALGBRC(21)) - (ALGBRC(18)+ALGBRC(9))
      ALGBRC(23) = 1.00000 - (STATES(10)+STATES(7)+STATES(6))
      RATES(10) =  ALGBRC(13)*ALGBRC(23) -  (ALGBRC(19)+CONSTS(66))*STATES(10)
      ALGBRC(20) = EXP( (- STATES(5)/ABS(STATES(5)))*CONSTS(27)*STATES(5)/CONSTS(37) ** 2.00000)
      ALGBRC(24) =  CONSTS(23)*ALGBRC(20)*CONSTS(31)*CONSTS(34) ** (CONSTS(19) - 37.0000)/10.0000
      ALGBRC(28) = ALGBRC(24)
      ALGBRC(14) = 1.00000+ (1.00000 - ALGBRC(4))*CONSTS(26)
      ALGBRC(17) =  CONSTS(22)*ALGBRC(14)*CONSTS(31)*CONSTS(33) ** (CONSTS(19) - 37.0000)/10.0000
      ALGBRC(30) =  CONSTS(43)*(ALGBRC(17)/CONSTS(40))
      ALGBRC(22) = EXP( ((STATES(4) - CONSTS(37))/ABS(STATES(4) - CONSTS(37)))*CONSTS(28)*(STATES(4) - CONSTS(37))/CONSTS(37) ** 2.00000)
      ALGBRC(25) =  CONSTS(24)*ALGBRC(22)*CONSTS(31)*CONSTS(35) ** (CONSTS(19) - 37.0000)/10.0000
      ALGBRC(31) =  CONSTS(61)*(ALGBRC(25)/CONSTS(63))*((CONSTS(39)+CONSTS(42))/CONSTS(42))*(CONSTS(45)/(CONSTS(39)+CONSTS(45)))
      ALGBRC(32) = ( CONSTS(66)*STATES(10)+ ALGBRC(31)*STATES(6)) -  (ALGBRC(30)+ALGBRC(28))*STATES(7)
      RATES(7) = ALGBRC(32)
      ALGBRC(26) = TERNRY(STATES(4).LT.CONSTS(37), EXP( CONSTS(29)*(CONSTS(37) - STATES(4))/CONSTS(37) ** 2.00000), EXP( CONSTS(30)*(STATES(4) - CONSTS(37))/CONSTS(37) ** 2.00000))
      ALGBRC(27) =  CONSTS(25)*ALGBRC(26)*CONSTS(31)*CONSTS(36) ** (CONSTS(19) - 37.0000)/10.0000
      ALGBRC(29) =  CONSTS(44)*(ALGBRC(27)/CONSTS(41))*((CONSTS(39)+CONSTS(42))/(CONSTS(39)+CONSTS(45)))
      ALGBRC(33) = ( CONSTS(39)*CONSTS(65)*ALGBRC(24)*(ALGBRC(27)/CONSTS(41)))/( (ALGBRC(17)/CONSTS(40))*(ALGBRC(25)/CONSTS(63))*CONSTS(38))
      ALGBRC(34) = ALGBRC(33)
      ALGBRC(35) = ( ALGBRC(34)*STATES(10)+ ALGBRC(28)*STATES(7)) -  (ALGBRC(31)+ALGBRC(29))*STATES(6)
      RATES(6) = ALGBRC(35)
      ALGBRC(36) = ( ALGBRC(34)*ALGBRC(31)+ ALGBRC(29)*CONSTS(66)+ ALGBRC(31)*CONSTS(66))/( CONSTS(66)*ALGBRC(28)+ ALGBRC(34)*ALGBRC(30)+ ALGBRC(34)*ALGBRC(28)+ ALGBRC(34)*ALGBRC(31)+ ALGBRC(29)*CONSTS(66)+ ALGBRC(31)*CONSTS(66)+ ALGBRC(28)*ALGBRC(29)+ ALGBRC(34)*ALGBRC(30)+ ALGBRC(29)*ALGBRC(30))
      ALGBRC(37) = CONSTS(62)/2.00000+ (CONSTS(46)/ALGBRC(36))*(- ( CONSTS(66)*STATES(5))+ ALGBRC(31)*(STATES(4) - (CONSTS(37)+STATES(5))))
      RATES(5) = ALGBRC(37)
      ALGBRC(38) = ( CONSTS(66)*ALGBRC(28)+ ALGBRC(34)*ALGBRC(30)+ ALGBRC(34)*ALGBRC(28))/( CONSTS(66)*ALGBRC(28)+ ALGBRC(34)*ALGBRC(30)+ ALGBRC(34)*ALGBRC(28)+ ALGBRC(34)*ALGBRC(31)+ ALGBRC(29)*CONSTS(66)+ ALGBRC(31)*CONSTS(66)+ ALGBRC(28)*ALGBRC(29)+ ALGBRC(34)*ALGBRC(30)+ ALGBRC(29)*ALGBRC(30))
      ALGBRC(39) = CONSTS(62)/2.00000+ (CONSTS(46)/ALGBRC(38))*( ALGBRC(28)*((STATES(5)+CONSTS(37)) - STATES(4)))
      RATES(4) = ALGBRC(39)
      RETURN
      END
      SUBROUTINE computeVariables(VOI, CONSTS, RATES, STATES, ALGBRC)
      REAL VOI, CONSTS(*), RATES(*), STATES(*), ALGBRC(*)
      ALGBRC(8) =  CONSTS(64)*CONSTS(20)*(1.00000 - STATES(2)) -  CONSTS(59)*STATES(2)
      ALGBRC(11) =  CONSTS(64)*CONSTS(20)*(1.00000 - STATES(3)) -  CONSTS(60)*STATES(3)
      ALGBRC(1) = TERNRY(CONSTS(2)/2.00000.LT.STATES(1)/2.00000, CONSTS(2)/2.00000, STATES(1)/2.00000)
      ALGBRC(2) = TERNRY(STATES(1)/2.00000 - (STATES(1) - CONSTS(1)).GT.CONSTS(3)/2.00000, STATES(1)/2.00000 - (STATES(1) - CONSTS(1)), CONSTS(3)/2.00000)
      ALGBRC(3) = ALGBRC(1) - ALGBRC(2)
      ALGBRC(5) = ALGBRC(3)/CONSTS(1)
      ALGBRC(7) =  (1.00000 - ALGBRC(5))*STATES(2)+ ALGBRC(5)*STATES(3)
      ALGBRC(10) =  ABS(1.00000/(1.00000+CONSTS(11)/ALGBRC(7) ** CONSTS(12))) ** (1.0 / 2)
      ALGBRC(13) =  CONSTS(13)*ALGBRC(10)*CONSTS(6) ** (CONSTS(19) - 37.0000)/10.0000
      ALGBRC(16) = TERNRY(1.00000/ALGBRC(10).LT.100.000, 1.00000/ALGBRC(10), 100.000)
      ALGBRC(19) =  CONSTS(14)*ALGBRC(16)*CONSTS(7) ** (CONSTS(19) - 37.0000)/10.0000
      ALGBRC(4) = ( ALGBRC(3)*2.00000)/(CONSTS(2) - CONSTS(3))
      ALGBRC(6) =  CONSTS(47)*ALGBRC(4)*( STATES(4)*STATES(6)+ STATES(5)*STATES(7))
      ALGBRC(9) = ( 1.00000*ALGBRC(6))/CONSTS(69)
      ALGBRC(12) =  ((STATES(1) - CONSTS(49))/ABS(STATES(1) - CONSTS(49)))*CONSTS(51)*(EXP( CONSTS(52)*ABS(STATES(1) - CONSTS(49))) - 1.00000)
      ALGBRC(15) = TERNRY(STATES(1).GT.CONSTS(48),  CONSTS(53)*(EXP( CONSTS(54)*ABS(STATES(1) - CONSTS(48))) - 1.00000), 0.00000)
      ALGBRC(18) = ALGBRC(12)+ALGBRC(15)
      ALGBRC(21) = TERNRY(CONSTS(56).EQ.1.00000,  CONSTS(55)*(CONSTS(50) - STATES(1)), 0.00000)
      ALGBRC(23) = 1.00000 - (STATES(10)+STATES(7)+STATES(6))
      ALGBRC(20) = EXP( (- STATES(5)/ABS(STATES(5)))*CONSTS(27)*STATES(5)/CONSTS(37) ** 2.00000)
      ALGBRC(24) =  CONSTS(23)*ALGBRC(20)*CONSTS(31)*CONSTS(34) ** (CONSTS(19) - 37.0000)/10.0000
      ALGBRC(28) = ALGBRC(24)
      ALGBRC(14) = 1.00000+ (1.00000 - ALGBRC(4))*CONSTS(26)
      ALGBRC(17) =  CONSTS(22)*ALGBRC(14)*CONSTS(31)*CONSTS(33) ** (CONSTS(19) - 37.0000)/10.0000
      ALGBRC(30) =  CONSTS(43)*(ALGBRC(17)/CONSTS(40))
      ALGBRC(22) = EXP( ((STATES(4) - CONSTS(37))/ABS(STATES(4) - CONSTS(37)))*CONSTS(28)*(STATES(4) - CONSTS(37))/CONSTS(37) ** 2.00000)
      ALGBRC(25) =  CONSTS(24)*ALGBRC(22)*CONSTS(31)*CONSTS(35) ** (CONSTS(19) - 37.0000)/10.0000
      ALGBRC(31) =  CONSTS(61)*(ALGBRC(25)/CONSTS(63))*((CONSTS(39)+CONSTS(42))/CONSTS(42))*(CONSTS(45)/(CONSTS(39)+CONSTS(45)))
      ALGBRC(32) = ( CONSTS(66)*STATES(10)+ ALGBRC(31)*STATES(6)) -  (ALGBRC(30)+ALGBRC(28))*STATES(7)
      ALGBRC(26) = TERNRY(STATES(4).LT.CONSTS(37), EXP( CONSTS(29)*(CONSTS(37) - STATES(4))/CONSTS(37) ** 2.00000), EXP( CONSTS(30)*(STATES(4) - CONSTS(37))/CONSTS(37) ** 2.00000))
      ALGBRC(27) =  CONSTS(25)*ALGBRC(26)*CONSTS(31)*CONSTS(36) ** (CONSTS(19) - 37.0000)/10.0000
      ALGBRC(29) =  CONSTS(44)*(ALGBRC(27)/CONSTS(41))*((CONSTS(39)+CONSTS(42))/(CONSTS(39)+CONSTS(45)))
      ALGBRC(33) = ( CONSTS(39)*CONSTS(65)*ALGBRC(24)*(ALGBRC(27)/CONSTS(41)))/( (ALGBRC(17)/CONSTS(40))*(ALGBRC(25)/CONSTS(63))*CONSTS(38))
      ALGBRC(34) = ALGBRC(33)
      ALGBRC(35) = ( ALGBRC(34)*STATES(10)+ ALGBRC(28)*STATES(7)) -  (ALGBRC(31)+ALGBRC(29))*STATES(6)
      ALGBRC(36) = ( ALGBRC(34)*ALGBRC(31)+ ALGBRC(29)*CONSTS(66)+ ALGBRC(31)*CONSTS(66))/( CONSTS(66)*ALGBRC(28)+ ALGBRC(34)*ALGBRC(30)+ ALGBRC(34)*ALGBRC(28)+ ALGBRC(34)*ALGBRC(31)+ ALGBRC(29)*CONSTS(66)+ ALGBRC(31)*CONSTS(66)+ ALGBRC(28)*ALGBRC(29)+ ALGBRC(34)*ALGBRC(30)+ ALGBRC(29)*ALGBRC(30))
      ALGBRC(37) = CONSTS(62)/2.00000+ (CONSTS(46)/ALGBRC(36))*(- ( CONSTS(66)*STATES(5))+ ALGBRC(31)*(STATES(4) - (CONSTS(37)+STATES(5))))
      ALGBRC(38) = ( CONSTS(66)*ALGBRC(28)+ ALGBRC(34)*ALGBRC(30)+ ALGBRC(34)*ALGBRC(28))/( CONSTS(66)*ALGBRC(28)+ ALGBRC(34)*ALGBRC(30)+ ALGBRC(34)*ALGBRC(28)+ ALGBRC(34)*ALGBRC(31)+ ALGBRC(29)*CONSTS(66)+ ALGBRC(31)*CONSTS(66)+ ALGBRC(28)*ALGBRC(29)+ ALGBRC(34)*ALGBRC(30)+ ALGBRC(29)*ALGBRC(30))
      ALGBRC(39) = CONSTS(62)/2.00000+ (CONSTS(46)/ALGBRC(38))*( ALGBRC(28)*((STATES(5)+CONSTS(37)) - STATES(4)))
      ALGBRC(40) = ( CONSTS(67)*CONSTS(24)+ CONSTS(25)*CONSTS(21)+ ALGBRC(25)*CONSTS(21))/( CONSTS(23)*CONSTS(25)+ CONSTS(24)*CONSTS(22)+ CONSTS(25)*CONSTS(22)+ CONSTS(67)*CONSTS(24)+ CONSTS(25)*CONSTS(21)+ CONSTS(24)*CONSTS(21)+ CONSTS(67)*CONSTS(22)+ CONSTS(21)*CONSTS(23)+ CONSTS(67)*CONSTS(24))
      ALGBRC(41) = (STATES(6)+STATES(7))/(CONSTS(68)+ALGBRC(40))
      ALGBRC(42) = (ALGBRC(35)+ALGBRC(32))/(CONSTS(68)+ALGBRC(40))
      ALGBRC(43) =  CONSTS(58)*( (1.00000 - ALGBRC(5))*STATES(2)+ ALGBRC(5)*( ALGBRC(41)*STATES(3)+ (1.00000 - ALGBRC(41))*STATES(2)))
      ALGBRC(44) = TERNRY(STATES(1).LT.CONSTS(2),  - 0.500000*CONSTS(62), 0.00000)
      ALGBRC(45) = TERNRY( 2.00000*CONSTS(1) - STATES(1).GT.CONSTS(3),  - 0.500000*CONSTS(62), 0.00000)
      ALGBRC(46) = ALGBRC(44) - ALGBRC(45)
      ALGBRC(47) = ( 2.00000*ALGBRC(46))/(CONSTS(2) - CONSTS(3))
      ALGBRC(48) = ALGBRC(46)/CONSTS(1)
      ALGBRC(49) =  CONSTS(57)*ALGBRC(47)*( STATES(4)*STATES(6)+ STATES(5)*STATES(7))+ CONSTS(57)*ALGBRC(4)*( ALGBRC(39)*STATES(6)+ STATES(4)*ALGBRC(35)+ ALGBRC(37)*STATES(7)+ STATES(5)*ALGBRC(32))
      ALGBRC(50) =  CONSTS(58)*( - ALGBRC(48)*STATES(2)+ (1.00000 - ALGBRC(5))*ALGBRC(8)+ ALGBRC(48)*( ALGBRC(41)*STATES(3)+ (1.00000 - ALGBRC(41))*STATES(2))+ ALGBRC(5)*(( ALGBRC(42)*STATES(3)+ ALGBRC(41)*ALGBRC(11)+ (1.00000 - ALGBRC(41))*ALGBRC(8)) -  ALGBRC(42)*STATES(2)))
      RETURN
      END
      REAL FUNCTION TERNRY(TEST, VALA, VALB)
      LOGICAL TEST
      REAL VALA, VALB
      IF (TEST) THEN
        TERNRY = VALA
      ELSE
        TERNRY = VALB
      ENDIF
      RETURN
      END