C
C There are a total of 12 entries in the algebraic variable array.
C There are a total of 5 entries in each of the rate and state variable arrays.
C There are a total of 32 entries in the constant variable array.
C
C
C VOI is time in component environment (millisecond).
C STATES(1) is V in component membrane (millivolt).
C CONSTS(1) is Cm in component membrane (femtofarad).
C ALGBRC(6) is Ica in component Ica (femtoampere).
C ALGBRC(8) is Is1 in component Is1 (femtoampere).
C ALGBRC(12) is Is2 in component Is2 (femtoampere).
C ALGBRC(11) is Il in component Il (femtoampere).
C ALGBRC(9) is Ik in component Ik (femtoampere).
C CONSTS(2) is gCa in component Ica (picosiemens).
C CONSTS(3) is VCa in component Ica (millivolt).
C ALGBRC(1) is m_infinity in component m (dimensionless).
C CONSTS(4) is vm in component m (millivolt).
C CONSTS(5) is sm in component m (millivolt).
C CONSTS(6) is gs1 in component Is1 (picosiemens).
C CONSTS(7) is VK in component Ik (millivolt).
C STATES(2) is s1 in component s1 (dimensionless).
C ALGBRC(2) is s1_infinity in component s1 (dimensionless).
C CONSTS(8) is autos1 in component s1 (dimensionless).
C CONSTS(9) is s1knot in component s1 (dimensionless).
C CONSTS(10) is tau_s1 in component s1 (millisecond).
C CONSTS(11) is vs1 in component s1 (millivolt).
C CONSTS(12) is ss1 in component s1 (millivolt).
C CONSTS(13) is gK in component Ik (picosiemens).
C STATES(3) is n in component n (dimensionless).
C ALGBRC(3) is n_infinity in component n (dimensionless).
C CONSTS(14) is tau_n_bar in component n (millisecond).
C ALGBRC(7) is tau_n in component n (millisecond).
C CONSTS(15) is vn in component n (millivolt).
C CONSTS(16) is sn in component n (millivolt).
C CONSTS(17) is gl in component Il (picosiemens).
C CONSTS(18) is Vl in component Il (millivolt).
C ALGBRC(10) is q in component Il (dimensionless).
C STATES(4) is p in component Il (dimensionless).
C CONSTS(19) is alpha_p in component Il (dimensionless).
C CONSTS(20) is tau_p in component Il (millisecond).
C CONSTS(31) is beta_p in component Il (dimensionless).
C CONSTS(21) is p0 in component Il (dimensionless).
C CONSTS(22) is noise in component Il (dimensionless).
C ALGBRC(4) is sigma in component Il (dimensionless).
C CONSTS(23) is w in component Il (dimensionless).
C CONSTS(32) is nstoc in component Il (dimensionless).
C CONSTS(24) is delNoise in component Il (dimensionless).
C CONSTS(25) is gs2 in component Is2 (picosiemens).
C STATES(5) is s2 in component s2 (dimensionless).
C ALGBRC(5) is s2_infinity in component s2 (dimensionless).
C CONSTS(26) is autos2 in component s2 (dimensionless).
C CONSTS(27) is s2knot in component s2 (dimensionless).
C CONSTS(28) is tau_s2 in component s2 (millisecond).
C CONSTS(29) is vs2 in component s2 (millivolt).
C CONSTS(30) is ss2 in component s2 (millivolt).
C RATES(1) is d/dt V in component membrane (millivolt).
C RATES(2) is d/dt s1 in component s1 (dimensionless).
C RATES(3) is d/dt n in component n (dimensionless).
C RATES(4) is d/dt p in component Il (dimensionless).
C RATES(5) is d/dt s2 in component s2 (dimensionless).
C
      SUBROUTINE initConsts(CONSTS, RATES, STATES)
      REAL CONSTS(*), RATES(*), STATES(*)
      STATES(1) = -40.0
      CONSTS(1) = 4525.0
      CONSTS(2) = 280.0
      CONSTS(3) = 100.0
      CONSTS(4) = -22.0
      CONSTS(5) = 7.5
      CONSTS(6) = 22.0
      CONSTS(7) = -80.0
      STATES(2) = 0.9
      CONSTS(8) = 1
      CONSTS(9) = 1
      CONSTS(10) = 1000.0
      CONSTS(11) = -50.0
      CONSTS(12) = 5
      CONSTS(13) = 1300.0
      STATES(3) = 0.0
      CONSTS(14) = 8.25
      CONSTS(15) = -9.0
      CONSTS(16) = 10.0
      CONSTS(17) = 41.0
      CONSTS(18) = -40.0
      STATES(4) = 0.14
      CONSTS(19) = 1.0
      CONSTS(20) = 100.0
      CONSTS(21) = 0.2
      CONSTS(22) = 1
      CONSTS(23) = 1
      CONSTS(24) = 3
      CONSTS(25) = 16
      STATES(5) = 0.5
      CONSTS(26) = 1
      CONSTS(27) = 0.47
      CONSTS(28) = 30000.0
      CONSTS(29) = -40.0
      CONSTS(30) = 15
      CONSTS(31) =  CONSTS(19)*(1.00000/CONSTS(21) - 1.00000)
      CONSTS(32) = 1000.00/CONSTS(24) ** 2.00000
      RETURN
      END
      SUBROUTINE computeRates(VOI, CONSTS,  RATES, STATES, ALGBRC)
      REAL VOI, CONSTS(*), RATES(*), STATES(*), ALGBRC(*)
      ALGBRC(2) = 1.00000/(1.00000+EXP((CONSTS(11) - STATES(1))/CONSTS(12)))
      RATES(2) =  CONSTS(8)*((ALGBRC(2) - STATES(2))/CONSTS(10))+ (1.00000 - CONSTS(8))*(CONSTS(9) - STATES(2))
      ALGBRC(4) =  (( CONSTS(19)*(1.00000 - STATES(4))+ CONSTS(31)*STATES(4))/( CONSTS(20)*CONSTS(32))) ** (1.0 / 2)
      RATES(4) = ( CONSTS(19)*(1.00000 - STATES(4)) -  CONSTS(31)*STATES(4))/CONSTS(20)+ CONSTS(22)*CONSTS(23)*ALGBRC(4)
      ALGBRC(5) = 1.00000/(1.00000+EXP((CONSTS(29) - STATES(1))/CONSTS(30)))
      RATES(5) =  CONSTS(26)*((ALGBRC(5) - STATES(5))/CONSTS(28))+ (1.00000 - CONSTS(26))*(CONSTS(27) - STATES(5))
      ALGBRC(3) = 1.00000/(1.00000+EXP((CONSTS(15) - STATES(1))/CONSTS(16)))
      ALGBRC(7) = CONSTS(14)/(1.00000+EXP((STATES(1) - CONSTS(15))/CONSTS(16)))
      RATES(3) = (ALGBRC(3) - STATES(3))/ALGBRC(7)
      ALGBRC(1) = 1.00000/(1.00000+EXP((CONSTS(4) - STATES(1))/CONSTS(5)))
      ALGBRC(6) =  CONSTS(2)*ALGBRC(1)*(STATES(1) - CONSTS(3))
      ALGBRC(8) =  CONSTS(6)*STATES(2)*(STATES(1) - CONSTS(7))
      ALGBRC(12) =  CONSTS(25)*STATES(5)*(STATES(1) - CONSTS(7))
      ALGBRC(10) = (1.00000+STATES(4))/2.00000
      ALGBRC(11) =  CONSTS(17)*ALGBRC(10)*(STATES(1) - CONSTS(18))
      ALGBRC(9) =  CONSTS(13)*STATES(3)*(STATES(1) - CONSTS(7))
      RATES(1) = - (ALGBRC(6)+ALGBRC(8)+ALGBRC(12)+ALGBRC(11)+ALGBRC(9))/CONSTS(1)
      RETURN
      END
      SUBROUTINE computeVariables(VOI, CONSTS, RATES, STATES, ALGBRC)
      REAL VOI, CONSTS(*), RATES(*), STATES(*), ALGBRC(*)
      ALGBRC(2) = 1.00000/(1.00000+EXP((CONSTS(11) - STATES(1))/CONSTS(12)))
      ALGBRC(4) =  (( CONSTS(19)*(1.00000 - STATES(4))+ CONSTS(31)*STATES(4))/( CONSTS(20)*CONSTS(32))) ** (1.0 / 2)
      ALGBRC(5) = 1.00000/(1.00000+EXP((CONSTS(29) - STATES(1))/CONSTS(30)))
      ALGBRC(3) = 1.00000/(1.00000+EXP((CONSTS(15) - STATES(1))/CONSTS(16)))
      ALGBRC(7) = CONSTS(14)/(1.00000+EXP((STATES(1) - CONSTS(15))/CONSTS(16)))
      ALGBRC(1) = 1.00000/(1.00000+EXP((CONSTS(4) - STATES(1))/CONSTS(5)))
      ALGBRC(6) =  CONSTS(2)*ALGBRC(1)*(STATES(1) - CONSTS(3))
      ALGBRC(8) =  CONSTS(6)*STATES(2)*(STATES(1) - CONSTS(7))
      ALGBRC(12) =  CONSTS(25)*STATES(5)*(STATES(1) - CONSTS(7))
      ALGBRC(10) = (1.00000+STATES(4))/2.00000
      ALGBRC(11) =  CONSTS(17)*ALGBRC(10)*(STATES(1) - CONSTS(18))
      ALGBRC(9) =  CONSTS(13)*STATES(3)*(STATES(1) - CONSTS(7))
      RETURN
      END