# Size of variable arrays:
sizeAlgebraic = 107
sizeStates = 43
sizeConstants = 85
from math import *
from numpy import *

def createLegends():
    legend_states = [""] * sizeStates
    legend_rates = [""] * sizeStates
    legend_algebraic = [""] * sizeAlgebraic
    legend_voi = ""
    legend_constants = [""] * sizeConstants
    legend_constants[58] = "Aj_nj in component cell (um2)"
    legend_constants[0] = "BCa in component calcium (mM)"
    legend_constants[1] = "CSQN in component calcium (mM)"
    legend_states[0] = "CaSR1 in component calcium (mM)"
    legend_states[1] = "CaSR2 in component calcium (mM)"
    legend_states[2] = "CaSR3 in component calcium (mM)"
    legend_states[3] = "CaSR4 in component calcium (mM)"
    legend_states[4] = "Cai1 in component calcium (mM)"
    legend_states[5] = "Cai2 in component calcium (mM)"
    legend_states[6] = "Cai3 in component calcium (mM)"
    legend_states[7] = "Cai4 in component calcium (mM)"
    legend_states[8] = "Cass in component calcium (mM)"
    legend_constants[2] = "DCa in component calcium (m2_per_s_times_1e_minus_12)"
    legend_constants[3] = "DCaBm in component calcium (m2_per_s_times_1e_minus_12)"
    legend_constants[4] = "DCaSR in component calcium (m2_per_s_times_1e_minus_12)"
    legend_constants[5] = "F in component phys (C_per_mol)"
    legend_algebraic[58] = "ICaL in component ical (pA)"
    legend_algebraic[63] = "ICaP in component icap (pA)"
    legend_algebraic[81] = "ICab in component icab (pA)"
    legend_algebraic[71] = "INaCa in component inaca (pA)"
    legend_algebraic[75] = "JCa1 in component calcium (kat_times_1e_minus_12)"
    legend_algebraic[76] = "JCa2 in component calcium (kat_times_1e_minus_12)"
    legend_algebraic[77] = "JCa3 in component calcium (kat_times_1e_minus_12)"
    legend_algebraic[48] = "JCa4 in component calcium (kat_times_1e_minus_12)"
    legend_algebraic[93] = "JCass in component calcium (kat_times_1e_minus_12)"
    legend_algebraic[72] = "JSRCa1 in component calcium (kat_times_1e_minus_12)"
    legend_algebraic[73] = "JSRCa2 in component calcium (kat_times_1e_minus_12)"
    legend_algebraic[74] = "JSRCa3 in component calcium (kat_times_1e_minus_12)"
    legend_algebraic[91] = "JSRCa4 in component calcium (kat_times_1e_minus_12)"
    legend_algebraic[0] = "JSRCaleak1 in component calcium (kat_times_1e_minus_12)"
    legend_algebraic[1] = "JSRCaleak2 in component calcium (kat_times_1e_minus_12)"
    legend_algebraic[2] = "JSRCaleak3 in component calcium (kat_times_1e_minus_12)"
    legend_algebraic[3] = "JSRCaleakss in component calcium (kat_times_1e_minus_12)"
    legend_algebraic[64] = "J_SERCASR1 in component serca (kat_times_1e_minus_12)"
    legend_algebraic[65] = "J_SERCASR2 in component serca (kat_times_1e_minus_12)"
    legend_algebraic[66] = "J_SERCASR3 in component serca (kat_times_1e_minus_12)"
    legend_algebraic[88] = "J_SERCASRss in component serca (kat_times_1e_minus_12)"
    legend_algebraic[68] = "J_bulkSERCA1 in component serca (kat_times_1e_minus_12)"
    legend_algebraic[69] = "J_bulkSERCA2 in component serca (kat_times_1e_minus_12)"
    legend_algebraic[70] = "J_bulkSERCA3 in component serca (kat_times_1e_minus_12)"
    legend_algebraic[90] = "J_bulkSERCAss in component serca (kat_times_1e_minus_12)"
    legend_algebraic[31] = "Jj_nj in component calcium (kat_times_1e_minus_12)"
    legend_algebraic[60] = "Jrel1 in component ryr (kat_times_1e_minus_12)"
    legend_algebraic[61] = "Jrel2 in component ryr (kat_times_1e_minus_12)"
    legend_algebraic[62] = "Jrel3 in component ryr (kat_times_1e_minus_12)"
    legend_algebraic[86] = "Jrelss in component ryr (kat_times_1e_minus_12)"
    legend_constants[6] = "KdBCa in component calcium (mM)"
    legend_constants[7] = "KdCSQN in component calcium (mM)"
    legend_constants[8] = "KdSLhigh in component calcium (mM)"
    legend_constants[9] = "KdSLlow in component calcium (mM)"
    legend_constants[10] = "SLhigh in component calcium (mM)"
    legend_constants[11] = "SLlow in component calcium (mM)"
    legend_constants[72] = "VSR1 in component cell (nL)"
    legend_constants[77] = "VSR2 in component cell (nL)"
    legend_constants[80] = "VSR3 in component cell (nL)"
    legend_constants[83] = "VSR4 in component cell (nL)"
    legend_constants[59] = "Vnonjunct1 in component cell (nL)"
    legend_constants[71] = "Vnonjunct2 in component cell (nL)"
    legend_constants[76] = "Vnonjunct3 in component cell (nL)"
    legend_constants[79] = "Vnonjunct4 in component cell (nL)"
    legend_constants[12] = "Vss in component cell (nL)"
    legend_algebraic[28] = "calcium_CaSR1_beta in component calcium (dimensionless)"
    legend_algebraic[29] = "calcium_CaSR2_beta in component calcium (dimensionless)"
    legend_algebraic[30] = "calcium_CaSR3_beta in component calcium (dimensionless)"
    legend_algebraic[49] = "calcium_CaSR4_beta in component calcium (dimensionless)"
    legend_algebraic[50] = "calcium_Cai1_beta in component calcium (dimensionless)"
    legend_algebraic[45] = "calcium_Cai1_gamma in component calcium (dimensionless)"
    legend_algebraic[51] = "calcium_Cai2_beta in component calcium (dimensionless)"
    legend_algebraic[46] = "calcium_Cai2_gamma in component calcium (dimensionless)"
    legend_algebraic[52] = "calcium_Cai3_beta in component calcium (dimensionless)"
    legend_algebraic[47] = "calcium_Cai3_gamma in component calcium (dimensionless)"
    legend_algebraic[59] = "calcium_Cai4_beta in component calcium (dimensionless)"
    legend_algebraic[53] = "calcium_Cai4_gamma in component calcium (dimensionless)"
    legend_algebraic[54] = "calcium_Cass_beta in component calcium (dimensionless)"
    legend_algebraic[82] = "calcium_Cass_i_tot in component calcium (pA)"
    legend_constants[13] = "dx in component cell (um)"
    legend_constants[14] = "kSRleak in component calcium (hertz)"
    legend_voi = "time in component engine (second)"
    legend_constants[60] = "xj_nj in component cell (um)"
    legend_constants[15] = "Cm in component cell (nF)"
    legend_constants[84] = "Vcytosol in component cell (nL)"
    legend_constants[82] = "Vnonjunct_Nai in component cell (nL)"
    legend_constants[16] = "lcell in component cell (um)"
    legend_constants[17] = "pi in component cell (dimensionless)"
    legend_constants[18] = "rjunct in component cell (um)"
    legend_constants[61] = "xj_nj_Nai in component cell (um)"
    legend_constants[19] = "Cao in component extra (mM)"
    legend_constants[20] = "Ko in component extra (mM)"
    legend_constants[21] = "Nao in component extra (mM)"
    legend_algebraic[80] = "ECa in component nernst (mV)"
    legend_states[9] = "V in component membrane (mV)"
    legend_constants[22] = "gCab in component icab (nS)"
    legend_constants[23] = "ECa_app in component ical (mV)"
    legend_states[10] = "d in component ical (hertz)"
    legend_states[11] = "f1 in component ical (hertz)"
    legend_states[12] = "f2 in component ical (hertz)"
    legend_algebraic[4] = "f_inf in component ical (dimensionless)"
    legend_states[13] = "fca in component ical (hertz)"
    legend_constants[24] = "gCaL in component ical (nS)"
    legend_algebraic[5] = "ical_d_inf in component ical (dimensionless)"
    legend_algebraic[32] = "ical_d_tau in component ical (second)"
    legend_algebraic[33] = "ical_f1_tau in component ical (second)"
    legend_algebraic[34] = "ical_f2_tau in component ical (second)"
    legend_algebraic[6] = "ical_fca_inf in component ical (dimensionless)"
    legend_constants[25] = "ical_fca_tau in component ical (second)"
    legend_constants[26] = "kCa in component ical (mM)"
    legend_constants[27] = "kCan in component ical (dimensionless)"
    legend_constants[28] = "ICaPmax in component icap (pA)"
    legend_constants[29] = "kCaP in component icap (mM)"
    legend_algebraic[83] = "EK in component nernst (mV)"
    legend_algebraic[96] = "ENa in component nernst (mV)"
    legend_algebraic[98] = "If in component if (pA)"
    legend_algebraic[85] = "IfK in component if (pA)"
    legend_algebraic[97] = "IfNa in component if (pA)"
    legend_constants[30] = "gIf in component if (nS)"
    legend_algebraic[7] = "if_y_inf in component if (dimensionless)"
    legend_algebraic[35] = "if_y_tau in component if (second)"
    legend_states[14] = "y in component if (hertz)"
    legend_constants[73] = "FRT in component phys (per_mV)"
    legend_algebraic[87] = "IK1 in component ik1 (pA)"
    legend_constants[62] = "gK1 in component ik1 (nS)"
    legend_algebraic[89] = "IKr in component ikr (pA)"
    legend_constants[31] = "gKr in component ikr (nS)"
    legend_algebraic[8] = "ikr_pa_inf in component ikr (dimensionless)"
    legend_algebraic[36] = "ikr_pa_tau in component ikr (second)"
    legend_states[15] = "pa in component ikr (hertz)"
    legend_algebraic[67] = "pi in component ikr (dimensionless)"
    legend_algebraic[92] = "IKs in component iks (pA)"
    legend_constants[32] = "gKs in component iks (nS)"
    legend_algebraic[9] = "iks_n_inf in component iks (dimensionless)"
    legend_algebraic[37] = "iks_n_tau in component iks (second)"
    legend_states[16] = "n in component iks (hertz)"
    legend_algebraic[94] = "IKur in component ikur (pA)"
    legend_constants[63] = "gKur in component ikur (nS)"
    legend_algebraic[10] = "ikur_r_inf in component ikur (dimensionless)"
    legend_algebraic[38] = "ikur_r_tau in component ikur (second)"
    legend_algebraic[11] = "ikur_s_inf in component ikur (dimensionless)"
    legend_algebraic[39] = "ikur_s_tau in component ikur (second)"
    legend_states[17] = "r in component ikur (hertz)"
    legend_states[18] = "s in component ikur (hertz)"
    legend_algebraic[99] = "INa in component ina (pA)"
    legend_constants[33] = "PNa in component ina (m3_per_s_times_1e_minus_12)"
    legend_states[19] = "h1 in component ina (hertz)"
    legend_states[20] = "h2 in component ina (hertz)"
    legend_algebraic[12] = "h_inf in component ina (dimensionless)"
    legend_algebraic[40] = "ina_h1_tau in component ina (second)"
    legend_algebraic[41] = "ina_h2_tau in component ina (second)"
    legend_algebraic[13] = "ina_m_inf in component ina (dimensionless)"
    legend_algebraic[42] = "ina_m_tau in component ina (second)"
    legend_states[21] = "m in component ina (hertz)"
    legend_algebraic[100] = "INab in component inab (pA)"
    legend_constants[34] = "gNab in component inab (nS)"
    legend_states[22] = "Nass in component sodium (mM)"
    legend_constants[35] = "dNaCa in component inaca (m12_per_mol4)"
    legend_constants[36] = "fCaNCX in component inaca (dimensionless)"
    legend_constants[37] = "gam in component inaca (dimensionless)"
    legend_constants[38] = "kNaCa in component inaca (m12_A_per_mol4_times_1e_minus_12)"
    legend_algebraic[79] = "INaK in component inak (pA)"
    legend_constants[39] = "INaKmax in component inak (pA)"
    legend_algebraic[78] = "Nass15 in component inak (dimensionless)"
    legend_constants[40] = "kNaKK in component inak (mM)"
    legend_constants[41] = "kNaKNa in component inak (dimensionless)"
    legend_algebraic[95] = "It in component it (pA)"
    legend_constants[64] = "gt in component it (nS)"
    legend_algebraic[14] = "it_r_inf in component it (dimensionless)"
    legend_algebraic[43] = "it_r_tau in component it (second)"
    legend_algebraic[15] = "it_s_inf in component it (dimensionless)"
    legend_algebraic[44] = "it_s_tau in component it (second)"
    legend_states[23] = "r in component it (hertz)"
    legend_states[24] = "s in component it (hertz)"
    legend_algebraic[101] = "i_ion in component membrane (pA)"
    legend_algebraic[103] = "i_stim in component stimulus (pA)"
    legend_states[25] = "Ki in component potassium (mM)"
    legend_constants[65] = "RTF in component phys (mV)"
    legend_constants[42] = "R in component phys (mJ_per_mol_per_K)"
    legend_constants[43] = "T in component phys (kelvin)"
    legend_algebraic[105] = "i_tot in component potassium (pA)"
    legend_algebraic[55] = "SRCa1 in component ryr (dimensionless)"
    legend_algebraic[56] = "SRCa2 in component ryr (dimensionless)"
    legend_algebraic[57] = "SRCa3 in component ryr (dimensionless)"
    legend_algebraic[84] = "SRCass in component ryr (dimensionless)"
    legend_states[26] = "a1 in component ryr (hertz)"
    legend_states[27] = "a2 in component ryr (hertz)"
    legend_states[28] = "a3 in component ryr (hertz)"
    legend_algebraic[16] = "ainf1 in component ryr (dimensionless)"
    legend_algebraic[17] = "ainf2 in component ryr (dimensionless)"
    legend_algebraic[18] = "ainf3 in component ryr (dimensionless)"
    legend_algebraic[19] = "ainfss in component ryr (dimensionless)"
    legend_states[29] = "ass in component ryr (hertz)"
    legend_states[30] = "c1 in component ryr (hertz)"
    legend_states[31] = "c2 in component ryr (hertz)"
    legend_states[32] = "c3 in component ryr (hertz)"
    legend_algebraic[20] = "cinf1 in component ryr (dimensionless)"
    legend_algebraic[21] = "cinf2 in component ryr (dimensionless)"
    legend_algebraic[22] = "cinf3 in component ryr (dimensionless)"
    legend_algebraic[23] = "cinfss in component ryr (dimensionless)"
    legend_states[33] = "css in component ryr (hertz)"
    legend_constants[74] = "nu1 in component ryr (m3_per_s_times_1e_minus_12)"
    legend_constants[78] = "nu2 in component ryr (m3_per_s_times_1e_minus_12)"
    legend_constants[81] = "nu3 in component ryr (m3_per_s_times_1e_minus_12)"
    legend_constants[66] = "nuss in component ryr (m3_per_s_times_1e_minus_12)"
    legend_states[34] = "o1 in component ryr (hertz)"
    legend_states[35] = "o2 in component ryr (hertz)"
    legend_states[36] = "o3 in component ryr (hertz)"
    legend_algebraic[24] = "oinf1 in component ryr (dimensionless)"
    legend_algebraic[25] = "oinf2 in component ryr (dimensionless)"
    legend_algebraic[26] = "oinf3 in component ryr (dimensionless)"
    legend_algebraic[27] = "oinfss in component ryr (dimensionless)"
    legend_states[37] = "oss in component ryr (hertz)"
    legend_constants[44] = "tau_act in component ryr (second)"
    legend_constants[45] = "tau_actss in component ryr (second)"
    legend_constants[46] = "tau_adapt in component ryr (second)"
    legend_constants[47] = "tau_inact in component ryr (second)"
    legend_constants[48] = "tau_inactss in component ryr (second)"
    legend_constants[49] = "SERCAKmf in component serca (mM)"
    legend_constants[50] = "SERCAKmr in component serca (mM)"
    legend_states[38] = "a1 in component serca (mM)"
    legend_states[39] = "a2 in component serca (mM)"
    legend_states[40] = "a3 in component serca (mM)"
    legend_states[41] = "ass in component serca (mM)"
    legend_constants[51] = "cpumps in component serca (mM)"
    legend_constants[67] = "k1 in component serca (m6_per_s_per_mol2)"
    legend_constants[75] = "k2 in component serca (hertz)"
    legend_constants[68] = "k3 in component serca (m6_per_s_per_mol2)"
    legend_constants[52] = "k4 in component serca (hertz)"
    legend_constants[69] = "BNa in component sodium (mM)"
    legend_constants[53] = "DNa in component sodium (m2_per_s_times_1e_minus_12)"
    legend_algebraic[102] = "JNa in component sodium (kat_times_1e_minus_12)"
    legend_constants[54] = "KdBNa in component sodium (mM)"
    legend_states[42] = "Nai in component sodium (mM)"
    legend_algebraic[104] = "betaNass in component sodium (dimensionless)"
    legend_algebraic[106] = "i_ss in component sodium (pA)"
    legend_constants[70] = "amplitude in component stimulus (pA)"
    legend_constants[55] = "duration in component stimulus (second)"
    legend_constants[56] = "offset in component stimulus (second)"
    legend_constants[57] = "period in component stimulus (second)"
    legend_rates[0] = "d/dt CaSR1 in component calcium (mM)"
    legend_rates[1] = "d/dt CaSR2 in component calcium (mM)"
    legend_rates[2] = "d/dt CaSR3 in component calcium (mM)"
    legend_rates[3] = "d/dt CaSR4 in component calcium (mM)"
    legend_rates[4] = "d/dt Cai1 in component calcium (mM)"
    legend_rates[5] = "d/dt Cai2 in component calcium (mM)"
    legend_rates[6] = "d/dt Cai3 in component calcium (mM)"
    legend_rates[7] = "d/dt Cai4 in component calcium (mM)"
    legend_rates[8] = "d/dt Cass in component calcium (mM)"
    legend_rates[10] = "d/dt d in component ical (hertz)"
    legend_rates[11] = "d/dt f1 in component ical (hertz)"
    legend_rates[12] = "d/dt f2 in component ical (hertz)"
    legend_rates[13] = "d/dt fca in component ical (hertz)"
    legend_rates[14] = "d/dt y in component if (hertz)"
    legend_rates[15] = "d/dt pa in component ikr (hertz)"
    legend_rates[16] = "d/dt n in component iks (hertz)"
    legend_rates[17] = "d/dt r in component ikur (hertz)"
    legend_rates[18] = "d/dt s in component ikur (hertz)"
    legend_rates[19] = "d/dt h1 in component ina (hertz)"
    legend_rates[20] = "d/dt h2 in component ina (hertz)"
    legend_rates[21] = "d/dt m in component ina (hertz)"
    legend_rates[23] = "d/dt r in component it (hertz)"
    legend_rates[24] = "d/dt s in component it (hertz)"
    legend_rates[9] = "d/dt V in component membrane (mV)"
    legend_rates[25] = "d/dt Ki in component potassium (mM)"
    legend_rates[26] = "d/dt a1 in component ryr (hertz)"
    legend_rates[27] = "d/dt a2 in component ryr (hertz)"
    legend_rates[28] = "d/dt a3 in component ryr (hertz)"
    legend_rates[29] = "d/dt ass in component ryr (hertz)"
    legend_rates[30] = "d/dt c1 in component ryr (hertz)"
    legend_rates[31] = "d/dt c2 in component ryr (hertz)"
    legend_rates[32] = "d/dt c3 in component ryr (hertz)"
    legend_rates[33] = "d/dt css in component ryr (hertz)"
    legend_rates[34] = "d/dt o1 in component ryr (hertz)"
    legend_rates[35] = "d/dt o2 in component ryr (hertz)"
    legend_rates[36] = "d/dt o3 in component ryr (hertz)"
    legend_rates[37] = "d/dt oss in component ryr (hertz)"
    legend_rates[38] = "d/dt a1 in component serca (mM)"
    legend_rates[39] = "d/dt a2 in component serca (mM)"
    legend_rates[40] = "d/dt a3 in component serca (mM)"
    legend_rates[41] = "d/dt ass in component serca (mM)"
    legend_rates[42] = "d/dt Nai in component sodium (mM)"
    legend_rates[22] = "d/dt Nass in component sodium (mM)"
    return (legend_states, legend_algebraic, legend_voi, legend_constants)

def initConsts():
    constants = [0.0] * sizeConstants; states = [0.0] * sizeStates;
    constants[0] = 0.024
    constants[1] = 6.7
    states[0] = 0.6189225
    states[1] = 0.6076289
    states[2] = 0.5905266
    states[3] = 0.5738108
    states[4] = 1.35496500000000013e-04
    states[5] = 1.38142100000000014e-04
    states[6] = 1.44208699999999994e-04
    states[7] = 1.56184399999999995e-04
    states[8] = 1.61937700000000013e-04
    constants[2] = 780.0
    constants[3] = 25.0
    constants[4] = 44.0
    constants[5] = 96487.0
    constants[6] = 0.00238
    constants[7] = 0.8
    constants[8] = 0.013
    constants[9] = 1.1
    constants[10] = 13.0
    constants[11] = 165.0
    constants[12] = 4.99231999999999966e-05
    constants[13] = 1.625
    constants[14] = 0.006
    constants[15] = 0.05
    constants[16] = 122.051
    constants[17] = 3.14159265358979312e+00
    constants[18] = 6.5
    constants[19] = 1.8
    constants[20] = 5.4
    constants[21] = 130.0
    states[9] = -75.42786
    constants[22] = 0.0952
    constants[23] = 60.0
    states[10] = 1.06091699999999996e-05
    states[11] = 0.9988566
    states[12] = 0.9988624
    states[13] = 0.9744374
    constants[24] = 25.3125
    constants[25] = 0.002
    constants[26] = 0.001
    constants[27] = 2.0
    constants[28] = 2.0
    constants[29] = 0.0005
    constants[30] = 1.0
    states[14] = 5.62066499999999969e-02
    constants[31] = 0.5
    states[15] = 4.18941700000000008e-05
    constants[32] = 1.0
    states[16] = 4.10975100000000003e-03
    states[17] = 3.11170299999999984e-04
    states[18] = 0.9751094
    constants[33] = 0.0018
    states[19] = 0.90391
    states[20] = 0.9039673
    states[21] = 2.77581199999999990e-03
    constants[34] = 0.060599
    states[22] = 8.691504
    constants[35] = 0.0003
    constants[36] = 1.0
    constants[37] = 0.45
    constants[38] = 0.0084
    constants[39] = 70.8253
    constants[40] = 1.0
    constants[41] = 11.0
    states[23] = 9.59425800000000026e-04
    states[24] = 0.954338
    states[25] = 134.6313
    constants[42] = 8314.0
    constants[43] = 306.15
    states[26] = 0.1925362
    states[27] = 0.2010345
    states[28] = 0.2163122
    states[29] = 0.2455297
    states[30] = 0.9993722
    states[31] = 0.9995086
    states[32] = 0.9995604
    states[33] = 0.9999717
    states[34] = 9.47851400000000044e-05
    states[35] = 7.76550300000000031e-05
    states[36] = 5.67494700000000006e-05
    states[37] = 3.97509699999999973e-05
    constants[44] = 0.01875
    constants[45] = 0.005
    constants[46] = 1.0
    constants[47] = 0.0875
    constants[48] = 0.015
    constants[49] = 0.00025
    constants[50] = 1.8
    states[38] = 4.63856499999999988e-03
    states[39] = 4.51207800000000010e-03
    states[40] = 4.32640899999999981e-03
    states[41] = 4.25044500000000026e-03
    constants[51] = 0.04
    constants[52] = 7.5
    constants[53] = 0.12
    constants[54] = 10.0
    states[42] = 9.28686
    constants[55] = 0.001
    constants[56] = 0.01
    constants[57] = 1.0
    constants[58] = (((constants[17]*constants[18])*2.00000)*constants[16])*0.500000
    constants[59] = ((((power(1.00000*constants[13], 2.00000)-power(0.00000*constants[13], 2.00000))*constants[17])*constants[16])*0.500000)*1.00000e-06
    constants[60] = 0.0200000/2.00000+constants[13]/2.00000
    constants[61] = 0.0200000/2.00000+2.00000*constants[13]
    constants[62] = 3.82500*0.900000
    constants[63] = 0.890000*2.75000
    constants[64] = 1.09000*7.50000
    constants[65] = (constants[42]*constants[43])/constants[5]
    constants[66] = 625.000*constants[12]
    constants[67] = (power(1000.00, 2.00000))*constants[52]
    constants[68] = constants[52]/(power(constants[50], 2.00000))
    constants[69] = 0.490000*2.31000
    constants[70] = -2500.00
    constants[71] = ((((power(2.00000*constants[13], 2.00000)-power(1.00000*constants[13], 2.00000))*constants[17])*constants[16])*0.500000)*1.00000e-06
    constants[72] = ((0.0500000*constants[59])/2.00000)*0.900000
    constants[73] = 1.00000/constants[65]
    constants[74] = 1.00000*constants[59]
    constants[75] = constants[67]*(power(constants[49], 2.00000))
    constants[76] = ((((power(3.00000*constants[13], 2.00000)-power(2.00000*constants[13], 2.00000))*constants[17])*constants[16])*0.500000)*1.00000e-06
    constants[77] = ((0.0500000*constants[71])/2.00000)*0.900000
    constants[78] = 1.00000*constants[71]
    constants[79] = ((((power(4.00000*constants[13], 2.00000)-power(3.00000*constants[13], 2.00000))*constants[17])*constants[16])*0.500000)*1.00000e-06
    constants[80] = ((0.0500000*constants[76])/2.00000)*0.900000
    constants[81] = 1.00000*constants[76]
    constants[82] = ((constants[59]+constants[71])+constants[76])+constants[79]
    constants[83] = ((0.0500000*constants[79])/2.00000)*0.900000
    constants[84] = constants[82]+constants[12]
    return (states, constants)

def computeRates(voi, states, constants):
    rates = [0.0] * sizeStates; algebraic = [0.0] * sizeAlgebraic
    algebraic[6] = 1.00000-1.00000/(1.00000+power(constants[26]/states[8], constants[27]))
    rates[13] = (algebraic[6]-states[13])/constants[25]
    algebraic[16] = 0.505000-0.427000/(1.00000+exp((states[4]*1000.00-0.290000)/0.0820000))
    rates[26] = (algebraic[16]-states[26])/constants[46]
    algebraic[17] = 0.505000-0.427000/(1.00000+exp((states[5]*1000.00-0.290000)/0.0820000))
    rates[27] = (algebraic[17]-states[27])/constants[46]
    algebraic[18] = 0.505000-0.427000/(1.00000+exp((states[6]*1000.00-0.290000)/0.0820000))
    rates[28] = (algebraic[18]-states[28])/constants[46]
    algebraic[19] = 0.505000-0.427000/(1.00000+exp((states[8]*1000.00-0.290000)/0.0820000))
    rates[29] = (algebraic[19]-states[29])/constants[46]
    algebraic[20] = 1.00000/(1.00000+exp((states[4]*1000.00-(states[26]+0.0200000))/0.0100000))
    rates[30] = (algebraic[20]-states[30])/constants[47]
    algebraic[21] = 1.00000/(1.00000+exp((states[5]*1000.00-(states[27]+0.0200000))/0.0100000))
    rates[31] = (algebraic[21]-states[31])/constants[47]
    algebraic[22] = 1.00000/(1.00000+exp((states[6]*1000.00-(states[28]+0.0200000))/0.0100000))
    rates[32] = (algebraic[22]-states[32])/constants[47]
    algebraic[23] = 1.00000/(1.00000+exp((states[8]*1000.00-(states[29]+0.0200000))/0.0100000))
    rates[33] = (algebraic[23]-states[33])/constants[48]
    algebraic[24] = 1.00000-1.00000/(1.00000+exp((states[4]*1000.00-(states[26]+0.220000))/0.0300000))
    rates[34] = (algebraic[24]-states[34])/constants[44]
    algebraic[25] = 1.00000-1.00000/(1.00000+exp((states[5]*1000.00-(states[27]+0.220000))/0.0300000))
    rates[35] = (algebraic[25]-states[35])/constants[44]
    algebraic[26] = 1.00000-1.00000/(1.00000+exp((states[6]*1000.00-(states[28]+0.220000))/0.0300000))
    rates[36] = (algebraic[26]-states[36])/constants[44]
    algebraic[27] = 1.00000-1.00000/(1.00000+exp((states[8]*1000.00-(states[29]+0.220000))/0.0300000))
    rates[37] = (algebraic[27]-states[37])/constants[45]
    algebraic[5] = 1.00000/(1.00000+exp((states[9]+9.00000)/-5.80000))
    algebraic[32] = 0.00270000*exp(-(power((states[9]+35.0000)/30.0000, 2.00000)))+0.00200000
    rates[10] = (algebraic[5]-states[10])/algebraic[32]
    algebraic[4] = 1.00000/(1.00000+exp((states[9]+27.4000)/7.10000))
    algebraic[33] = ((0.986980*exp(-(power((states[9]+30.1605)/7.09396, 2.00000)))+0.0427500/(1.00000+exp((states[9]-51.6155)/-80.6133)))+0.0357600/(1.00000+exp((states[9]+29.5727)/13.2176)))-0.00821000
    rates[11] = (algebraic[4]-states[11])/algebraic[33]
    algebraic[34] = 1.33230*exp(-(power((states[9]+40.0000)/14.2000, 2.00000)))+0.0626000
    rates[12] = (algebraic[4]-states[12])/algebraic[34]
    algebraic[7] = 1.00000/(1.00000+exp((states[9]+97.8287)/12.4802))
    algebraic[35] = 1.00000/(0.00332000*exp(-states[9]/16.5410)+23.7184*exp(states[9]/16.5410))
    rates[14] = (algebraic[7]-states[14])/algebraic[35]
    algebraic[8] = 1.00000/(1.00000+exp((states[9]+15.0000)/-6.00000))
    algebraic[36] = 0.217180*exp(-(power((states[9]+20.1376)/22.1996, 2.00000)))+0.0311800
    rates[15] = (algebraic[8]-states[15])/algebraic[36]
    algebraic[9] = 1.00000/(1.00000+exp((states[9]-19.9000)/-12.7000))
    algebraic[37] = 0.400000*exp(-(power((states[9]-20.0000)/20.0000, 2.00000)))+0.700000
    rates[16] = (algebraic[9]-states[16])/algebraic[37]
    algebraic[10] = 1.00000/(1.00000+exp((states[9]+6.00000)/-8.60000))
    algebraic[38] = 0.00900000/(1.00000+exp((states[9]+5.00000)/12.0000))+0.000500000
    rates[17] = (algebraic[10]-states[17])/algebraic[38]
    algebraic[11] = 1.00000/(1.00000+exp((states[9]+7.50000)/10.0000))
    algebraic[39] = 0.590000/(1.00000+exp((states[9]+60.0000)/10.0000))+3.05000
    rates[18] = (algebraic[11]-states[18])/algebraic[39]
    algebraic[12] = 1.00000/(1.00000+exp((states[9]+63.6000)/5.30000))
    algebraic[40] = 0.0300000/(1.00000+exp((states[9]+35.1000)/3.20000))+0.000300000
    rates[19] = (algebraic[12]-states[19])/algebraic[40]
    algebraic[41] = 0.120000/(1.00000+exp((states[9]+35.1000)/3.20000))+0.00300000
    rates[20] = (algebraic[12]-states[20])/algebraic[41]
    algebraic[13] = 1.00000/(1.00000+exp((states[9]+27.1200)/-8.21000))
    algebraic[42] = 4.20000e-05*exp(-(power((states[9]+25.5700)/28.8000, 2.00000)))+2.40000e-05
    rates[21] = (algebraic[13]-states[21])/algebraic[42]
    algebraic[14] = 1.00000/(1.00000+exp((states[9]-1.00000)/-11.0000))
    algebraic[43] = 0.00350000*exp(-(power((states[9]+0.00000)/30.0000, 2.00000)))+0.00150000
    rates[23] = (algebraic[14]-states[23])/algebraic[43]
    algebraic[15] = 1.00000/(1.00000+exp((states[9]+40.5000)/11.5000))
    algebraic[44] = 0.0256350*exp(-(power((states[9]+52.4500)/15.8827, 2.00000)))+0.0141400
    rates[24] = (algebraic[15]-states[24])/algebraic[44]
    algebraic[31] = (((constants[2]*constants[58])/constants[60])*(states[8]-states[7]))*1.00000e-06
    algebraic[48] = algebraic[31]
    algebraic[53] = (constants[0]*constants[6])/(power(states[7]+constants[6], 2.00000))
    algebraic[59] = 1.00000/(1.00000+algebraic[53])
    rates[7] = ((algebraic[59]*(constants[2]+algebraic[53]*constants[3]))*(((states[7]-2.00000*states[7])+states[6])/(power(constants[13], 2.00000))+(states[7]-states[6])/((2.00000*4.00000)*(power(constants[13], 2.00000))))-((((2.00000*algebraic[59])*algebraic[53])*constants[3])/(constants[6]+states[7]))*(power((states[7]-states[6])/(2.00000*constants[13]), 2.00000)))+(algebraic[48]/constants[79])*algebraic[59]
    algebraic[0] = (constants[14]*(states[0]-states[4]))*constants[59]
    algebraic[64] = (((-constants[68]*(power(states[0], 2.00000)))*(constants[51]-states[38])+constants[52]*states[38])*constants[59])*2.00000
    algebraic[55] = 1.00000-1.00000/(1.00000+exp((states[0]-0.300000)/0.100000))
    algebraic[60] = (((constants[74]*states[34])*states[30])*algebraic[55])*(states[0]-states[4])
    algebraic[72] = (algebraic[64]-algebraic[0])-algebraic[60]
    algebraic[28] = 1.00000/(1.00000+(constants[1]*constants[7])/(power(states[0]+constants[7], 2.00000)))
    rates[0] = (algebraic[28]*constants[4])*(((states[1]-2.00000*states[0])+states[0])/(power(constants[13], 2.00000))+(states[1]-states[0])/((2.00000*1.00000)*(power(constants[13], 2.00000))))+(algebraic[72]/constants[72])*algebraic[28]
    algebraic[1] = (constants[14]*(states[1]-states[5]))*constants[71]
    algebraic[65] = (((-constants[68]*(power(states[1], 2.00000)))*(constants[51]-states[39])+constants[52]*states[39])*constants[71])*2.00000
    algebraic[56] = 1.00000-1.00000/(1.00000+exp((states[1]-0.300000)/0.100000))
    algebraic[61] = (((constants[78]*states[35])*states[31])*algebraic[56])*(states[1]-states[5])
    algebraic[73] = (algebraic[65]-algebraic[1])-algebraic[61]
    algebraic[29] = 1.00000/(1.00000+(constants[1]*constants[7])/(power(states[1]+constants[7], 2.00000)))
    rates[1] = (algebraic[29]*constants[4])*(((states[2]-2.00000*states[1])+states[0])/(power(constants[13], 2.00000))+(states[2]-states[0])/((2.00000*2.00000)*(power(constants[13], 2.00000))))+(algebraic[73]/constants[77])*algebraic[29]
    algebraic[2] = (constants[14]*(states[2]-states[6]))*constants[76]
    algebraic[66] = (((-constants[68]*(power(states[2], 2.00000)))*(constants[51]-states[40])+constants[52]*states[40])*constants[76])*2.00000
    algebraic[57] = 1.00000-1.00000/(1.00000+exp((states[2]-0.300000)/0.100000))
    algebraic[62] = (((constants[81]*states[36])*states[32])*algebraic[57])*(states[2]-states[6])
    algebraic[74] = (algebraic[66]-algebraic[2])-algebraic[62]
    algebraic[30] = 1.00000/(1.00000+(constants[1]*constants[7])/(power(states[2]+constants[7], 2.00000)))
    rates[2] = (algebraic[30]*constants[4])*(((states[3]-2.00000*states[2])+states[1])/(power(constants[13], 2.00000))+(states[3]-states[1])/((2.00000*3.00000)*(power(constants[13], 2.00000))))+(algebraic[74]/constants[80])*algebraic[30]
    algebraic[68] = (((constants[67]*(power(states[4], 2.00000)))*(constants[51]-states[38])-constants[75]*states[38])*constants[59])*2.00000
    rates[38] = (0.500000*(-algebraic[64]+algebraic[68]))/constants[59]
    algebraic[69] = (((constants[67]*(power(states[5], 2.00000)))*(constants[51]-states[39])-constants[75]*states[39])*constants[71])*2.00000
    rates[39] = (0.500000*(-algebraic[65]+algebraic[69]))/constants[71]
    algebraic[70] = (((constants[67]*(power(states[6], 2.00000)))*(constants[51]-states[40])-constants[75]*states[40])*constants[76])*2.00000
    rates[40] = (0.500000*(-algebraic[66]+algebraic[70]))/constants[76]
    algebraic[75] = (-algebraic[68]+algebraic[0])+algebraic[60]
    algebraic[45] = (constants[0]*constants[6])/(power(states[4]+constants[6], 2.00000))
    algebraic[50] = 1.00000/(1.00000+algebraic[45])
    rates[4] = ((algebraic[50]*(constants[2]+algebraic[45]*constants[3]))*(((states[5]-2.00000*states[4])+states[4])/(power(constants[13], 2.00000))+(states[5]-states[4])/((2.00000*1.00000)*(power(constants[13], 2.00000))))-((((2.00000*algebraic[50])*algebraic[45])*constants[3])/(constants[6]+states[4]))*(power((states[5]-states[4])/(2.00000*constants[13]), 2.00000)))+(algebraic[75]/constants[59])*algebraic[50]
    algebraic[76] = (-algebraic[69]+algebraic[1])+algebraic[61]
    algebraic[46] = (constants[0]*constants[6])/(power(states[5]+constants[6], 2.00000))
    algebraic[51] = 1.00000/(1.00000+algebraic[46])
    rates[5] = ((algebraic[51]*(constants[2]+algebraic[46]*constants[3]))*(((states[6]-2.00000*states[5])+states[4])/(power(constants[13], 2.00000))+(states[6]-states[4])/((2.00000*2.00000)*(power(constants[13], 2.00000))))-((((2.00000*algebraic[51])*algebraic[46])*constants[3])/(constants[6]+states[5]))*(power((states[6]-states[4])/(2.00000*constants[13]), 2.00000)))+(algebraic[76]/constants[71])*algebraic[51]
    algebraic[77] = (-algebraic[70]+algebraic[2])+algebraic[62]
    algebraic[47] = (constants[0]*constants[6])/(power(states[6]+constants[6], 2.00000))
    algebraic[52] = 1.00000/(1.00000+algebraic[47])
    rates[6] = ((algebraic[52]*(constants[2]+algebraic[47]*constants[3]))*(((states[7]-2.00000*states[6])+states[5])/(power(constants[13], 2.00000))+(states[7]-states[5])/((2.00000*3.00000)*(power(constants[13], 2.00000))))-((((2.00000*algebraic[52])*algebraic[47])*constants[3])/(constants[6]+states[6]))*(power((states[7]-states[5])/(2.00000*constants[13]), 2.00000)))+(algebraic[77]/constants[76])*algebraic[52]
    algebraic[3] = (constants[14]*(states[3]-states[8]))*constants[12]
    algebraic[88] = (((-constants[68]*(power(states[3], 2.00000)))*(constants[51]-states[41])+constants[52]*states[41])*constants[12])*2.00000
    algebraic[84] = 1.00000-1.00000/(1.00000+exp((states[3]-0.300000)/0.100000))
    algebraic[86] = (((constants[66]*states[37])*states[33])*algebraic[84])*(states[3]-states[8])
    algebraic[91] = (algebraic[88]-algebraic[3])-algebraic[86]
    algebraic[49] = 1.00000/(1.00000+(constants[1]*constants[7])/(power(states[3]+constants[7], 2.00000)))
    rates[3] = (algebraic[49]*constants[4])*(((states[3]-2.00000*states[3])+states[2])/(power(constants[13], 2.00000))+(states[3]-states[2])/((2.00000*4.00000)*(power(constants[13], 2.00000))))+(algebraic[91]/constants[83])*algebraic[49]
    algebraic[90] = (((constants[67]*(power(states[8], 2.00000)))*(constants[51]-states[41])-constants[75]*states[41])*constants[12])*2.00000
    rates[41] = (0.500000*(-algebraic[88]+algebraic[90]))/constants[12]
    algebraic[93] = ((-algebraic[31]+algebraic[3])-algebraic[90])+algebraic[86]
    algebraic[54] = 1.00000/(((1.00000+(constants[11]*constants[9])/(power(states[8]+constants[9], 2.00000)))+(constants[10]*constants[8])/(power(states[8]+constants[8], 2.00000)))+(constants[0]*constants[6])/(power(states[8]+constants[6], 2.00000)))
    algebraic[58] = ((((constants[24]*states[10])*states[13])*states[11])*states[12])*(states[9]-constants[23])
    algebraic[63] = (constants[28]*states[8])/(constants[29]+states[8])
    algebraic[80] = (constants[65]*log(constants[19]/states[8]))/2.00000
    algebraic[81] = constants[22]*(states[9]-algebraic[80])
    algebraic[71] = (constants[38]*((exp((constants[37]*states[9])*constants[73])*(power(states[22], 3.00000)))*constants[19]-((exp(((constants[37]-1.00000)*states[9])*constants[73])*(power(constants[21], 3.00000)))*states[8])*constants[36]))/(1.00000+constants[35]*(((power(constants[21], 3.00000))*states[8])*constants[36]+(power(states[22], 3.00000))*constants[19]))
    algebraic[82] = ((-algebraic[58]-algebraic[81])-algebraic[63])+2.00000*algebraic[71]
    rates[8] = algebraic[54]*(algebraic[93]/constants[12]+algebraic[82]/((2.00000*constants[12])*constants[5]))
    algebraic[102] = (((constants[53]*constants[58])/constants[61])*(states[22]-states[42]))*1.00000e-06
    rates[42] = algebraic[102]/constants[82]
    algebraic[83] = constants[65]*log(constants[20]/states[25])
    algebraic[85] = (constants[30]*states[14])*((1.00000-0.267700)*(states[9]-algebraic[83]))
    algebraic[96] = constants[65]*log(constants[21]/states[22])
    algebraic[97] = (constants[30]*states[14])*(0.267700*(states[9]-algebraic[96]))
    algebraic[98] = algebraic[85]+algebraic[97]
    algebraic[87] = ((constants[62]*(power(constants[20]*1.00000, 0.445700)))*(states[9]-algebraic[83]))/(1.00000+exp((1.50000*((states[9]-algebraic[83])+3.60000))*constants[73]))
    algebraic[67] = 1.00000/(1.00000+exp((states[9]+55.0000)/24.0000))
    algebraic[89] = ((constants[31]*states[15])*algebraic[67])*(states[9]-algebraic[83])
    algebraic[92] = (constants[32]*states[16])*(states[9]-algebraic[83])
    algebraic[94] = ((constants[63]*states[17])*states[18])*(states[9]-algebraic[83])
    algebraic[99] = (((((((constants[33]*(power(states[21], 3.00000)))*(0.900000*states[19]+0.100000*states[20]))*constants[21])*states[9])*constants[5])*constants[73])*(exp((states[9]-algebraic[96])*constants[73])-1.00000))/(exp(states[9]*constants[73])-1.00000)
    algebraic[100] = constants[34]*(states[9]-algebraic[96])
    algebraic[78] = power(states[22]*1.00000, 1.50000)
    algebraic[79] = (((((constants[39]*constants[20])/(constants[20]+constants[40]))*algebraic[78])/(algebraic[78]+power(constants[41], 1.50000)))*(states[9]+150.000))/(states[9]+200.000)
    algebraic[95] = ((constants[64]*states[23])*states[24])*(states[9]-algebraic[83])
    algebraic[101] = ((((((((((((algebraic[99])+algebraic[58])+algebraic[95])+algebraic[94])+algebraic[87])+algebraic[89])+algebraic[92])+algebraic[100])+algebraic[81])+algebraic[79])+algebraic[63])+algebraic[71])+algebraic[98]
    algebraic[103] = (custom_piecewise([less((voi-constants[56])-constants[57]*floor((voi-constants[56])/constants[57]) , constants[55]), 1.00000 , True, 0.00000]))*constants[70]
    rates[9] = -(algebraic[101]+algebraic[103])/constants[15]
    algebraic[105] = (((((((algebraic[95])+algebraic[94])+algebraic[87])+algebraic[89])+algebraic[92])-2.00000*algebraic[79])+algebraic[85])+algebraic[103]
    rates[25] = -algebraic[105]/(constants[84]*constants[5])
    algebraic[104] = 1.00000/(1.00000+(constants[69]*constants[54])/(power(states[22]+constants[54], 2.00000)))
    algebraic[106] = ((((algebraic[99])+algebraic[100])+3.00000*algebraic[79])+3.00000*algebraic[71])+algebraic[97]
    rates[22] = algebraic[104]*(-algebraic[102]/constants[12]-algebraic[106]/(constants[12]*constants[5]))
    return(rates)

def computeAlgebraic(constants, states, voi):
    algebraic = array([[0.0] * len(voi)] * sizeAlgebraic)
    states = array(states)
    voi = array(voi)
    algebraic[6] = 1.00000-1.00000/(1.00000+power(constants[26]/states[8], constants[27]))
    algebraic[16] = 0.505000-0.427000/(1.00000+exp((states[4]*1000.00-0.290000)/0.0820000))
    algebraic[17] = 0.505000-0.427000/(1.00000+exp((states[5]*1000.00-0.290000)/0.0820000))
    algebraic[18] = 0.505000-0.427000/(1.00000+exp((states[6]*1000.00-0.290000)/0.0820000))
    algebraic[19] = 0.505000-0.427000/(1.00000+exp((states[8]*1000.00-0.290000)/0.0820000))
    algebraic[20] = 1.00000/(1.00000+exp((states[4]*1000.00-(states[26]+0.0200000))/0.0100000))
    algebraic[21] = 1.00000/(1.00000+exp((states[5]*1000.00-(states[27]+0.0200000))/0.0100000))
    algebraic[22] = 1.00000/(1.00000+exp((states[6]*1000.00-(states[28]+0.0200000))/0.0100000))
    algebraic[23] = 1.00000/(1.00000+exp((states[8]*1000.00-(states[29]+0.0200000))/0.0100000))
    algebraic[24] = 1.00000-1.00000/(1.00000+exp((states[4]*1000.00-(states[26]+0.220000))/0.0300000))
    algebraic[25] = 1.00000-1.00000/(1.00000+exp((states[5]*1000.00-(states[27]+0.220000))/0.0300000))
    algebraic[26] = 1.00000-1.00000/(1.00000+exp((states[6]*1000.00-(states[28]+0.220000))/0.0300000))
    algebraic[27] = 1.00000-1.00000/(1.00000+exp((states[8]*1000.00-(states[29]+0.220000))/0.0300000))
    algebraic[5] = 1.00000/(1.00000+exp((states[9]+9.00000)/-5.80000))
    algebraic[32] = 0.00270000*exp(-(power((states[9]+35.0000)/30.0000, 2.00000)))+0.00200000
    algebraic[4] = 1.00000/(1.00000+exp((states[9]+27.4000)/7.10000))
    algebraic[33] = ((0.986980*exp(-(power((states[9]+30.1605)/7.09396, 2.00000)))+0.0427500/(1.00000+exp((states[9]-51.6155)/-80.6133)))+0.0357600/(1.00000+exp((states[9]+29.5727)/13.2176)))-0.00821000
    algebraic[34] = 1.33230*exp(-(power((states[9]+40.0000)/14.2000, 2.00000)))+0.0626000
    algebraic[7] = 1.00000/(1.00000+exp((states[9]+97.8287)/12.4802))
    algebraic[35] = 1.00000/(0.00332000*exp(-states[9]/16.5410)+23.7184*exp(states[9]/16.5410))
    algebraic[8] = 1.00000/(1.00000+exp((states[9]+15.0000)/-6.00000))
    algebraic[36] = 0.217180*exp(-(power((states[9]+20.1376)/22.1996, 2.00000)))+0.0311800
    algebraic[9] = 1.00000/(1.00000+exp((states[9]-19.9000)/-12.7000))
    algebraic[37] = 0.400000*exp(-(power((states[9]-20.0000)/20.0000, 2.00000)))+0.700000
    algebraic[10] = 1.00000/(1.00000+exp((states[9]+6.00000)/-8.60000))
    algebraic[38] = 0.00900000/(1.00000+exp((states[9]+5.00000)/12.0000))+0.000500000
    algebraic[11] = 1.00000/(1.00000+exp((states[9]+7.50000)/10.0000))
    algebraic[39] = 0.590000/(1.00000+exp((states[9]+60.0000)/10.0000))+3.05000
    algebraic[12] = 1.00000/(1.00000+exp((states[9]+63.6000)/5.30000))
    algebraic[40] = 0.0300000/(1.00000+exp((states[9]+35.1000)/3.20000))+0.000300000
    algebraic[41] = 0.120000/(1.00000+exp((states[9]+35.1000)/3.20000))+0.00300000
    algebraic[13] = 1.00000/(1.00000+exp((states[9]+27.1200)/-8.21000))
    algebraic[42] = 4.20000e-05*exp(-(power((states[9]+25.5700)/28.8000, 2.00000)))+2.40000e-05
    algebraic[14] = 1.00000/(1.00000+exp((states[9]-1.00000)/-11.0000))
    algebraic[43] = 0.00350000*exp(-(power((states[9]+0.00000)/30.0000, 2.00000)))+0.00150000
    algebraic[15] = 1.00000/(1.00000+exp((states[9]+40.5000)/11.5000))
    algebraic[44] = 0.0256350*exp(-(power((states[9]+52.4500)/15.8827, 2.00000)))+0.0141400
    algebraic[31] = (((constants[2]*constants[58])/constants[60])*(states[8]-states[7]))*1.00000e-06
    algebraic[48] = algebraic[31]
    algebraic[53] = (constants[0]*constants[6])/(power(states[7]+constants[6], 2.00000))
    algebraic[59] = 1.00000/(1.00000+algebraic[53])
    algebraic[0] = (constants[14]*(states[0]-states[4]))*constants[59]
    algebraic[64] = (((-constants[68]*(power(states[0], 2.00000)))*(constants[51]-states[38])+constants[52]*states[38])*constants[59])*2.00000
    algebraic[55] = 1.00000-1.00000/(1.00000+exp((states[0]-0.300000)/0.100000))
    algebraic[60] = (((constants[74]*states[34])*states[30])*algebraic[55])*(states[0]-states[4])
    algebraic[72] = (algebraic[64]-algebraic[0])-algebraic[60]
    algebraic[28] = 1.00000/(1.00000+(constants[1]*constants[7])/(power(states[0]+constants[7], 2.00000)))
    algebraic[1] = (constants[14]*(states[1]-states[5]))*constants[71]
    algebraic[65] = (((-constants[68]*(power(states[1], 2.00000)))*(constants[51]-states[39])+constants[52]*states[39])*constants[71])*2.00000
    algebraic[56] = 1.00000-1.00000/(1.00000+exp((states[1]-0.300000)/0.100000))
    algebraic[61] = (((constants[78]*states[35])*states[31])*algebraic[56])*(states[1]-states[5])
    algebraic[73] = (algebraic[65]-algebraic[1])-algebraic[61]
    algebraic[29] = 1.00000/(1.00000+(constants[1]*constants[7])/(power(states[1]+constants[7], 2.00000)))
    algebraic[2] = (constants[14]*(states[2]-states[6]))*constants[76]
    algebraic[66] = (((-constants[68]*(power(states[2], 2.00000)))*(constants[51]-states[40])+constants[52]*states[40])*constants[76])*2.00000
    algebraic[57] = 1.00000-1.00000/(1.00000+exp((states[2]-0.300000)/0.100000))
    algebraic[62] = (((constants[81]*states[36])*states[32])*algebraic[57])*(states[2]-states[6])
    algebraic[74] = (algebraic[66]-algebraic[2])-algebraic[62]
    algebraic[30] = 1.00000/(1.00000+(constants[1]*constants[7])/(power(states[2]+constants[7], 2.00000)))
    algebraic[68] = (((constants[67]*(power(states[4], 2.00000)))*(constants[51]-states[38])-constants[75]*states[38])*constants[59])*2.00000
    algebraic[69] = (((constants[67]*(power(states[5], 2.00000)))*(constants[51]-states[39])-constants[75]*states[39])*constants[71])*2.00000
    algebraic[70] = (((constants[67]*(power(states[6], 2.00000)))*(constants[51]-states[40])-constants[75]*states[40])*constants[76])*2.00000
    algebraic[75] = (-algebraic[68]+algebraic[0])+algebraic[60]
    algebraic[45] = (constants[0]*constants[6])/(power(states[4]+constants[6], 2.00000))
    algebraic[50] = 1.00000/(1.00000+algebraic[45])
    algebraic[76] = (-algebraic[69]+algebraic[1])+algebraic[61]
    algebraic[46] = (constants[0]*constants[6])/(power(states[5]+constants[6], 2.00000))
    algebraic[51] = 1.00000/(1.00000+algebraic[46])
    algebraic[77] = (-algebraic[70]+algebraic[2])+algebraic[62]
    algebraic[47] = (constants[0]*constants[6])/(power(states[6]+constants[6], 2.00000))
    algebraic[52] = 1.00000/(1.00000+algebraic[47])
    algebraic[3] = (constants[14]*(states[3]-states[8]))*constants[12]
    algebraic[88] = (((-constants[68]*(power(states[3], 2.00000)))*(constants[51]-states[41])+constants[52]*states[41])*constants[12])*2.00000
    algebraic[84] = 1.00000-1.00000/(1.00000+exp((states[3]-0.300000)/0.100000))
    algebraic[86] = (((constants[66]*states[37])*states[33])*algebraic[84])*(states[3]-states[8])
    algebraic[91] = (algebraic[88]-algebraic[3])-algebraic[86]
    algebraic[49] = 1.00000/(1.00000+(constants[1]*constants[7])/(power(states[3]+constants[7], 2.00000)))
    algebraic[90] = (((constants[67]*(power(states[8], 2.00000)))*(constants[51]-states[41])-constants[75]*states[41])*constants[12])*2.00000
    algebraic[93] = ((-algebraic[31]+algebraic[3])-algebraic[90])+algebraic[86]
    algebraic[54] = 1.00000/(((1.00000+(constants[11]*constants[9])/(power(states[8]+constants[9], 2.00000)))+(constants[10]*constants[8])/(power(states[8]+constants[8], 2.00000)))+(constants[0]*constants[6])/(power(states[8]+constants[6], 2.00000)))
    algebraic[58] = ((((constants[24]*states[10])*states[13])*states[11])*states[12])*(states[9]-constants[23])
    algebraic[63] = (constants[28]*states[8])/(constants[29]+states[8])
    algebraic[80] = (constants[65]*log(constants[19]/states[8]))/2.00000
    algebraic[81] = constants[22]*(states[9]-algebraic[80])
    algebraic[71] = (constants[38]*((exp((constants[37]*states[9])*constants[73])*(power(states[22], 3.00000)))*constants[19]-((exp(((constants[37]-1.00000)*states[9])*constants[73])*(power(constants[21], 3.00000)))*states[8])*constants[36]))/(1.00000+constants[35]*(((power(constants[21], 3.00000))*states[8])*constants[36]+(power(states[22], 3.00000))*constants[19]))
    algebraic[82] = ((-algebraic[58]-algebraic[81])-algebraic[63])+2.00000*algebraic[71]
    algebraic[102] = (((constants[53]*constants[58])/constants[61])*(states[22]-states[42]))*1.00000e-06
    algebraic[83] = constants[65]*log(constants[20]/states[25])
    algebraic[85] = (constants[30]*states[14])*((1.00000-0.267700)*(states[9]-algebraic[83]))
    algebraic[96] = constants[65]*log(constants[21]/states[22])
    algebraic[97] = (constants[30]*states[14])*(0.267700*(states[9]-algebraic[96]))
    algebraic[98] = algebraic[85]+algebraic[97]
    algebraic[87] = ((constants[62]*(power(constants[20]*1.00000, 0.445700)))*(states[9]-algebraic[83]))/(1.00000+exp((1.50000*((states[9]-algebraic[83])+3.60000))*constants[73]))
    algebraic[67] = 1.00000/(1.00000+exp((states[9]+55.0000)/24.0000))
    algebraic[89] = ((constants[31]*states[15])*algebraic[67])*(states[9]-algebraic[83])
    algebraic[92] = (constants[32]*states[16])*(states[9]-algebraic[83])
    algebraic[94] = ((constants[63]*states[17])*states[18])*(states[9]-algebraic[83])
    algebraic[99] = (((((((constants[33]*(power(states[21], 3.00000)))*(0.900000*states[19]+0.100000*states[20]))*constants[21])*states[9])*constants[5])*constants[73])*(exp((states[9]-algebraic[96])*constants[73])-1.00000))/(exp(states[9]*constants[73])-1.00000)
    algebraic[100] = constants[34]*(states[9]-algebraic[96])
    algebraic[78] = power(states[22]*1.00000, 1.50000)
    algebraic[79] = (((((constants[39]*constants[20])/(constants[20]+constants[40]))*algebraic[78])/(algebraic[78]+power(constants[41], 1.50000)))*(states[9]+150.000))/(states[9]+200.000)
    algebraic[95] = ((constants[64]*states[23])*states[24])*(states[9]-algebraic[83])
    algebraic[101] = ((((((((((((algebraic[99])+algebraic[58])+algebraic[95])+algebraic[94])+algebraic[87])+algebraic[89])+algebraic[92])+algebraic[100])+algebraic[81])+algebraic[79])+algebraic[63])+algebraic[71])+algebraic[98]
    algebraic[103] = (custom_piecewise([less((voi-constants[56])-constants[57]*floor((voi-constants[56])/constants[57]) , constants[55]), 1.00000 , True, 0.00000]))*constants[70]
    algebraic[105] = (((((((algebraic[95])+algebraic[94])+algebraic[87])+algebraic[89])+algebraic[92])-2.00000*algebraic[79])+algebraic[85])+algebraic[103]
    algebraic[104] = 1.00000/(1.00000+(constants[69]*constants[54])/(power(states[22]+constants[54], 2.00000)))
    algebraic[106] = ((((algebraic[99])+algebraic[100])+3.00000*algebraic[79])+3.00000*algebraic[71])+algebraic[97]
    return algebraic

def custom_piecewise(cases):
    """Compute result of a piecewise function"""
    return select(cases[0::2],cases[1::2])

def solve_model():
    """Solve model with ODE solver"""
    from scipy.integrate import ode
    # Initialise constants and state variables
    (init_states, constants) = initConsts()

    # Set timespan to solve over
    voi = linspace(0, 10, 500)

    # Construct ODE object to solve
    r = ode(computeRates)
    r.set_integrator('vode', method='bdf', atol=1e-06, rtol=1e-06, max_step=1)
    r.set_initial_value(init_states, voi[0])
    r.set_f_params(constants)

    # Solve model
    states = array([[0.0] * len(voi)] * sizeStates)
    states[:,0] = init_states
    for (i,t) in enumerate(voi[1:]):
        if r.successful():
            r.integrate(t)
            states[:,i+1] = r.y
        else:
            break

    # Compute algebraic variables
    algebraic = computeAlgebraic(constants, states, voi)
    return (voi, states, algebraic)

def plot_model(voi, states, algebraic):
    """Plot variables against variable of integration"""
    import pylab
    (legend_states, legend_algebraic, legend_voi, legend_constants) = createLegends()
    pylab.figure(1)
    pylab.plot(voi,vstack((states,algebraic)).T)
    pylab.xlabel(legend_voi)
    pylab.legend(legend_states + legend_algebraic, loc='best')
    pylab.show()

if __name__ == "__main__":
    (voi, states, algebraic) = solve_model()
    plot_model(voi, states, algebraic)