# Size of variable arrays:
sizeAlgebraic = 93
sizeStates = 41
sizeConstants = 148
from math import *
from numpy import *

def createLegends():
    legend_states = [""] * sizeStates
    legend_rates = [""] * sizeStates
    legend_algebraic = [""] * sizeAlgebraic
    legend_voi = ""
    legend_constants = [""] * sizeConstants
    legend_voi = "time in component environment (second)"
    legend_constants[0] = "R in component model_parameters (joule_per_kilomole_kelvin)"
    legend_constants[1] = "T in component model_parameters (kelvin)"
    legend_constants[2] = "F in component model_parameters (coulomb_per_mole)"
    legend_constants[3] = "Na_e in component model_parameters (millimolar)"
    legend_constants[4] = "Ca_e in component model_parameters (millimolar)"
    legend_constants[5] = "K_e in component model_parameters (millimolar)"
    legend_constants[145] = "Vt in component model_parameters (cm3)"
    legend_constants[112] = "Vd in component model_parameters (cm3)"
    legend_constants[113] = "Vmyo in component model_parameters (cm3)"
    legend_constants[116] = "Sms in component model_parameters (cm2)"
    legend_constants[117] = "Smt in component model_parameters (cm2)"
    legend_constants[147] = "Cms in component model_parameters (microF)"
    legend_constants[146] = "Cmt in component model_parameters (microF)"
    legend_constants[106] = "Vc in component model_parameters (cm3)"
    legend_constants[114] = "VSRup in component model_parameters (cm3)"
    legend_constants[115] = "VSRrel in component model_parameters (cm3)"
    legend_constants[118] = "pt in component model_parameters (dimensionless)"
    legend_constants[143] = "Lt in component model_parameters (cm)"
    legend_constants[144] = "Rst in component model_parameters (ohm)"
    legend_constants[6] = "Rot in component model_parameters (ohm_per_cm2)"
    legend_constants[7] = "radiust in component model_parameters (cm)"
    legend_constants[8] = "ptcm in component model_parameters (per_cm)"
    legend_constants[9] = "Smtot in component model_parameters (cm2)"
    legend_algebraic[1] = "i_circ in component common_current (microA)"
    legend_constants[10] = "R_st in component common_current (kilo_ohm)"
    legend_states[0] = "Vm_s in component membrane (millivolt)"
    legend_states[1] = "Vm_t in component membrane (millivolt)"
    legend_algebraic[26] = "J_CaSRrel in component SR_Ca_release_channel (millimolar_per_second)"
    legend_constants[11] = "kap in component SR_Ca_release_channel (per_millimolar4_per_second)"
    legend_constants[12] = "kam in component SR_Ca_release_channel (per_second)"
    legend_constants[13] = "kbp in component SR_Ca_release_channel (per_millimolar3_per_second)"
    legend_constants[14] = "kbm in component SR_Ca_release_channel (per_second)"
    legend_constants[15] = "kcp in component SR_Ca_release_channel (per_second)"
    legend_constants[16] = "kcm in component SR_Ca_release_channel (per_second)"
    legend_states[2] = "F2 in component SR_Ca_release_channel (dimensionless)"
    legend_states[3] = "F3 in component SR_Ca_release_channel (dimensionless)"
    legend_states[4] = "F1 in component SR_Ca_release_channel (dimensionless)"
    legend_states[5] = "F4 in component SR_Ca_release_channel (dimensionless)"
    legend_states[6] = "Ca_ss in component ion_concentrations (millimolar)"
    legend_states[7] = "CaSRrel in component ion_concentrations (millimolar)"
    legend_algebraic[49] = "J_CaSRup in component SERCA2a_pump (millimolar_per_second)"
    legend_constants[17] = "K_fb in component SERCA2a_pump (millimolar)"
    legend_constants[18] = "K_rb in component SERCA2a_pump (millimolar)"
    legend_algebraic[47] = "fb in component SERCA2a_pump (dimensionless)"
    legend_algebraic[48] = "rb in component SERCA2a_pump (dimensionless)"
    legend_constants[19] = "Vmaxf in component SERCA2a_pump (millimolar_per_second)"
    legend_constants[20] = "Vmaxr in component SERCA2a_pump (millimolar_per_second)"
    legend_constants[21] = "K_SR in component SERCA2a_pump (dimensionless)"
    legend_constants[22] = "N_fb in component SERCA2a_pump (dimensionless)"
    legend_constants[23] = "N_rb in component SERCA2a_pump (dimensionless)"
    legend_states[8] = "Ca_i in component ion_concentrations (millimolar)"
    legend_states[9] = "CaSRup in component ion_concentrations (millimolar)"
    legend_algebraic[50] = "J_tr in component intracellular_and_SR_Ca_fluxes (millimolar_per_second)"
    legend_algebraic[0] = "J_xfer in component intracellular_and_SR_Ca_fluxes (millimolar_per_second)"
    legend_algebraic[53] = "J_trpn in component intracellular_and_SR_Ca_fluxes (millimolar_per_second)"
    legend_algebraic[51] = "J_BTRH in component intracellular_and_SR_Ca_fluxes (millimolar_per_second)"
    legend_algebraic[52] = "J_BTRL in component intracellular_and_SR_Ca_fluxes (millimolar_per_second)"
    legend_constants[24] = "BTRHmax in component intracellular_and_SR_Ca_fluxes (millimolar)"
    legend_constants[25] = "BTRLmax in component intracellular_and_SR_Ca_fluxes (millimolar)"
    legend_constants[26] = "tau_tr in component intracellular_and_SR_Ca_fluxes (second)"
    legend_constants[27] = "tau_xfer in component intracellular_and_SR_Ca_fluxes (second)"
    legend_states[10] = "BTRH in component intracellular_and_SR_Ca_fluxes (millimolar)"
    legend_states[11] = "BTRL in component intracellular_and_SR_Ca_fluxes (millimolar)"
    legend_constants[28] = "k_htrpn_plus in component intracellular_and_SR_Ca_fluxes (per_millimolar_per_second)"
    legend_constants[29] = "k_htrpn_minus in component intracellular_and_SR_Ca_fluxes (per_second)"
    legend_constants[30] = "k_ltrpn_plus in component intracellular_and_SR_Ca_fluxes (per_millimolar_per_second)"
    legend_constants[31] = "k_ltrpn_minus in component intracellular_and_SR_Ca_fluxes (per_second)"
    legend_states[12] = "Na_i in component ion_concentrations (millimolar)"
    legend_states[13] = "K_i in component ion_concentrations (millimolar)"
    legend_states[14] = "Na_t in component ion_concentrations (millimolar)"
    legend_states[15] = "Ca_t in component ion_concentrations (millimolar)"
    legend_states[16] = "K_t in component ion_concentrations (millimolar)"
    legend_algebraic[54] = "i_Kext in component ion_concentrations (microA)"
    legend_constants[32] = "K_mCMDN in component ion_concentrations (millimolar)"
    legend_constants[33] = "K_mCSQN in component ion_concentrations (millimolar)"
    legend_constants[34] = "K_mEGTA in component ion_concentrations (millimolar)"
    legend_constants[35] = "CMDN_tot in component ion_concentrations (millimolar)"
    legend_constants[36] = "CSQN_tot in component ion_concentrations (millimolar)"
    legend_constants[37] = "EGTA_tot in component ion_concentrations (millimolar)"
    legend_constants[38] = "tau_d in component ion_concentrations (second)"
    legend_constants[39] = "STIMULATION in component ion_concentrations (dimensionless)"
    legend_constants[40] = "time0 in component ion_concentrations (second)"
    legend_algebraic[59] = "i_Na_s in component sodium_current (microA)"
    legend_algebraic[61] = "i_Ca_s in component calcium_current (microA)"
    legend_algebraic[68] = "i_BNa_s in component background_currents (microA)"
    legend_algebraic[74] = "i_NaCa_s in component Na_Ca_ion_exchanger_current (microA)"
    legend_algebraic[72] = "i_NaK_s in component sodium_potassium_pump (microA)"
    legend_algebraic[60] = "i_fNa_s in component hyperpolarisation_activated_current (microA)"
    legend_algebraic[66] = "i_fK_s in component hyperpolarisation_activated_current (microA)"
    legend_algebraic[70] = "i_BK_s in component background_currents (microA)"
    legend_algebraic[65] = "i_K1_s in component inward_rectifier (microA)"
    legend_algebraic[63] = "i_Kto_s in component Ca_independent_transient_outward_K_current (microA)"
    legend_algebraic[64] = "i_Kss_s in component steady_state_outward_K_current (microA)"
    legend_algebraic[73] = "i_CaP_s in component sarcolemmal_calcium_pump_current (microA)"
    legend_algebraic[69] = "i_BCa_s in component background_currents (microA)"
    legend_algebraic[76] = "i_Na_t in component sodium_current (microA)"
    legend_algebraic[78] = "i_Ca_t in component calcium_current (microA)"
    legend_algebraic[85] = "i_BNa_t in component background_currents (microA)"
    legend_algebraic[92] = "i_NaCa_t in component Na_Ca_ion_exchanger_current (microA)"
    legend_algebraic[90] = "i_NaK_t in component sodium_potassium_pump (microA)"
    legend_algebraic[77] = "i_fNa_t in component hyperpolarisation_activated_current (microA)"
    legend_algebraic[83] = "i_fK_t in component hyperpolarisation_activated_current (microA)"
    legend_algebraic[87] = "i_BK_t in component background_currents (microA)"
    legend_algebraic[82] = "i_K1_t in component inward_rectifier (microA)"
    legend_algebraic[80] = "i_Kto_t in component Ca_independent_transient_outward_K_current (microA)"
    legend_algebraic[81] = "i_Kss_t in component steady_state_outward_K_current (microA)"
    legend_algebraic[91] = "i_CaP_t in component sarcolemmal_calcium_pump_current (microA)"
    legend_algebraic[86] = "i_BCa_t in component background_currents (microA)"
    legend_algebraic[55] = "Jte_Na in component t_tubular_ion_fluxes (millimolar_per_second)"
    legend_algebraic[56] = "Jte_Ca in component t_tubular_ion_fluxes (millimolar_per_second)"
    legend_algebraic[57] = "Jte_K in component t_tubular_ion_fluxes (millimolar_per_second)"
    legend_constants[41] = "tau_Na in component t_tubular_ion_fluxes (second)"
    legend_constants[42] = "tau_Ca in component t_tubular_ion_fluxes (second)"
    legend_constants[43] = "tau_K in component t_tubular_ion_fluxes (second)"
    legend_algebraic[67] = "i_f_s in component hyperpolarisation_activated_current (microA)"
    legend_algebraic[71] = "i_B_s in component background_currents (microA)"
    legend_algebraic[58] = "E_Na in component sodium_current (millivolt)"
    legend_constants[44] = "g_Na in component sodium_current (milliS_per_cm2)"
    legend_constants[119] = "g_Nas in component sodium_current (milliS)"
    legend_constants[45] = "fNat in component sodium_current (dimensionless)"
    legend_states[17] = "m in component sodium_current_m_gate (dimensionless)"
    legend_states[18] = "h in component sodium_current_h_gate (dimensionless)"
    legend_algebraic[2] = "alpha_m in component sodium_current_m_gate (per_second)"
    legend_algebraic[27] = "beta_m in component sodium_current_m_gate (per_second)"
    legend_algebraic[3] = "alpha_h in component sodium_current_h_gate (per_second)"
    legend_algebraic[28] = "beta_h in component sodium_current_h_gate (per_second)"
    legend_states[19] = "Co in component calcium_current (dimensionless)"
    legend_states[20] = "C_cast in component calcium_current (dimensionless)"
    legend_states[21] = "C1 in component calcium_current (dimensionless)"
    legend_states[22] = "C2 in component calcium_current (dimensionless)"
    legend_algebraic[29] = "tauCo in component calcium_current (second)"
    legend_constants[46] = "tauC_cast in component calcium_current (second)"
    legend_algebraic[30] = "tauC1 in component calcium_current (second)"
    legend_algebraic[31] = "tauC2 in component calcium_current (second)"
    legend_algebraic[4] = "Co_inf in component calcium_current (dimensionless)"
    legend_algebraic[5] = "C_cast_inf in component calcium_current (dimensionless)"
    legend_algebraic[6] = "C1_inf in component calcium_current (dimensionless)"
    legend_algebraic[7] = "C2_inf in component calcium_current (dimensionless)"
    legend_constants[47] = "g_Ca in component calcium_current (cm_per_second)"
    legend_constants[120] = "g_Cas in component calcium_current (cm3_per_second)"
    legend_constants[48] = "fCat in component calcium_current (dimensionless)"
    legend_algebraic[62] = "E_K in component Ca_independent_transient_outward_K_current (millivolt)"
    legend_constants[49] = "a in component Ca_independent_transient_outward_K_current (dimensionless)"
    legend_constants[50] = "b in component Ca_independent_transient_outward_K_current (dimensionless)"
    legend_constants[51] = "g_Kto in component Ca_independent_transient_outward_K_current (milliS_per_cm2)"
    legend_constants[121] = "g_Ktos in component Ca_independent_transient_outward_K_current (milliS)"
    legend_constants[52] = "fKtot in component Ca_independent_transient_outward_K_current (dimensionless)"
    legend_states[23] = "r in component Ca_independent_transient_outward_K_current_r_gate (dimensionless)"
    legend_states[24] = "s in component Ca_independent_transient_outward_K_current_s_gate (dimensionless)"
    legend_states[25] = "s_slow in component Ca_independent_transient_outward_K_current_s_slow_gate (dimensionless)"
    legend_algebraic[32] = "tau_r in component Ca_independent_transient_outward_K_current_r_gate (second)"
    legend_algebraic[8] = "r_max in component Ca_independent_transient_outward_K_current_r_gate (dimensionless)"
    legend_algebraic[33] = "tau_s in component Ca_independent_transient_outward_K_current_s_gate (second)"
    legend_algebraic[9] = "s_max in component Ca_independent_transient_outward_K_current_s_gate (dimensionless)"
    legend_algebraic[34] = "tau_s_slow in component Ca_independent_transient_outward_K_current_s_slow_gate (second)"
    legend_algebraic[10] = "s_slow_max in component Ca_independent_transient_outward_K_current_s_slow_gate (dimensionless)"
    legend_constants[53] = "g_Kss in component steady_state_outward_K_current (milliS_per_cm2)"
    legend_constants[122] = "g_Ksss in component steady_state_outward_K_current (milliS)"
    legend_constants[54] = "fKsst in component steady_state_outward_K_current (dimensionless)"
    legend_states[26] = "r_ss in component steady_state_outward_K_current_r_ss_gate (dimensionless)"
    legend_states[27] = "s_ss in component steady_state_outward_K_current_s_ss_gate (dimensionless)"
    legend_algebraic[35] = "tau_r_ss in component steady_state_outward_K_current_r_ss_gate (second)"
    legend_algebraic[11] = "r_ss_max in component steady_state_outward_K_current_r_ss_gate (dimensionless)"
    legend_constants[107] = "tau_s_ss in component steady_state_outward_K_current_s_ss_gate (second)"
    legend_algebraic[12] = "s_ss_max in component steady_state_outward_K_current_s_ss_gate (dimensionless)"
    legend_constants[55] = "g_K1 in component inward_rectifier (milliS_per_cm2)"
    legend_constants[123] = "g_K1s in component inward_rectifier (milliS)"
    legend_constants[56] = "fK1t in component inward_rectifier (dimensionless)"
    legend_constants[57] = "g_f in component hyperpolarisation_activated_current (milliS_per_cm2)"
    legend_constants[124] = "g_fs in component hyperpolarisation_activated_current (milliS)"
    legend_constants[58] = "fft in component hyperpolarisation_activated_current (dimensionless)"
    legend_constants[59] = "f_Na in component hyperpolarisation_activated_current (dimensionless)"
    legend_constants[108] = "f_K in component hyperpolarisation_activated_current (dimensionless)"
    legend_states[28] = "y in component hyperpolarisation_activated_current_y_gate (dimensionless)"
    legend_algebraic[36] = "tau_y in component hyperpolarisation_activated_current_y_gate (second)"
    legend_algebraic[13] = "y_infinity in component hyperpolarisation_activated_current_y_gate (dimensionless)"
    legend_constants[60] = "g_B_Ca in component background_currents (milliS_per_cm2)"
    legend_constants[61] = "g_B_Na in component background_currents (milliS_per_cm2)"
    legend_constants[62] = "g_B_K in component background_currents (milliS_per_cm2)"
    legend_constants[125] = "g_B_Cas in component background_currents (milliS)"
    legend_constants[126] = "g_B_Nas in component background_currents (milliS)"
    legend_constants[127] = "g_B_Ks in component background_currents (milliS)"
    legend_constants[63] = "fCabt in component background_currents (dimensionless)"
    legend_constants[64] = "fNabt in component background_currents (dimensionless)"
    legend_constants[65] = "fKbt in component background_currents (dimensionless)"
    legend_constants[66] = "E_Ca in component background_currents (millivolt)"
    legend_constants[67] = "i_NaK_max in component sodium_potassium_pump (microA_per_cm2)"
    legend_constants[128] = "i_NaK_max_s in component sodium_potassium_pump (microA)"
    legend_constants[68] = "fNaKt in component sodium_potassium_pump (dimensionless)"
    legend_constants[69] = "K_m_K in component sodium_potassium_pump (millimolar)"
    legend_constants[70] = "K_m_Na in component sodium_potassium_pump (millimolar)"
    legend_constants[109] = "sigma in component sodium_potassium_pump (dimensionless)"
    legend_constants[71] = "i_CaP_max in component sarcolemmal_calcium_pump_current (microA_per_cm2)"
    legend_constants[129] = "i_CaP_max_s in component sarcolemmal_calcium_pump_current (microA)"
    legend_constants[72] = "fCaPt in component sarcolemmal_calcium_pump_current (dimensionless)"
    legend_constants[73] = "i_NaCa_max in component Na_Ca_ion_exchanger_current (microA_per_cm2)"
    legend_constants[130] = "i_NaCa_max_s in component Na_Ca_ion_exchanger_current (microA)"
    legend_constants[74] = "fNaCat in component Na_Ca_ion_exchanger_current (dimensionless)"
    legend_algebraic[84] = "i_f_t in component hyperpolarisation_activated_current (microA)"
    legend_algebraic[88] = "i_B_t in component background_currents (microA)"
    legend_algebraic[75] = "E_Na in component sodium_current (millivolt)"
    legend_constants[75] = "g_Na in component sodium_current (milliS_per_cm2)"
    legend_constants[131] = "g_Nat in component sodium_current (milliS)"
    legend_constants[76] = "fNat in component sodium_current (dimensionless)"
    legend_states[29] = "m in component sodium_current_m_gate (dimensionless)"
    legend_states[30] = "h in component sodium_current_h_gate (dimensionless)"
    legend_algebraic[14] = "alpha_m in component sodium_current_m_gate (per_second)"
    legend_algebraic[37] = "beta_m in component sodium_current_m_gate (per_second)"
    legend_algebraic[15] = "alpha_h in component sodium_current_h_gate (per_second)"
    legend_algebraic[38] = "beta_h in component sodium_current_h_gate (per_second)"
    legend_states[31] = "Co in component calcium_current (dimensionless)"
    legend_states[32] = "C_cast in component calcium_current (dimensionless)"
    legend_states[33] = "C1 in component calcium_current (dimensionless)"
    legend_states[34] = "C2 in component calcium_current (dimensionless)"
    legend_algebraic[39] = "tauCo in component calcium_current (second)"
    legend_constants[77] = "tauC_cast in component calcium_current (second)"
    legend_algebraic[40] = "tauC1 in component calcium_current (second)"
    legend_algebraic[41] = "tauC2 in component calcium_current (second)"
    legend_algebraic[16] = "Co_inf in component calcium_current (dimensionless)"
    legend_algebraic[17] = "C_cast_inf in component calcium_current (dimensionless)"
    legend_algebraic[18] = "C1_inf in component calcium_current (dimensionless)"
    legend_algebraic[19] = "C2_inf in component calcium_current (dimensionless)"
    legend_constants[78] = "g_Ca in component calcium_current (cm_per_second)"
    legend_constants[132] = "g_Cat in component calcium_current (cm3_per_second)"
    legend_constants[79] = "fCat in component calcium_current (dimensionless)"
    legend_algebraic[79] = "E_K in component Ca_independent_transient_outward_K_current (millivolt)"
    legend_constants[80] = "a in component Ca_independent_transient_outward_K_current (dimensionless)"
    legend_constants[81] = "b in component Ca_independent_transient_outward_K_current (dimensionless)"
    legend_constants[82] = "g_Kto in component Ca_independent_transient_outward_K_current (milliS_per_cm2)"
    legend_constants[133] = "g_Ktot in component Ca_independent_transient_outward_K_current (milliS)"
    legend_constants[83] = "fKtot in component Ca_independent_transient_outward_K_current (dimensionless)"
    legend_states[35] = "r in component Ca_independent_transient_outward_K_current_r_gate (dimensionless)"
    legend_states[36] = "s in component Ca_independent_transient_outward_K_current_s_gate (dimensionless)"
    legend_states[37] = "s_slow in component Ca_independent_transient_outward_K_current_s_slow_gate (dimensionless)"
    legend_algebraic[42] = "tau_r in component Ca_independent_transient_outward_K_current_r_gate (second)"
    legend_algebraic[20] = "r_max in component Ca_independent_transient_outward_K_current_r_gate (dimensionless)"
    legend_algebraic[43] = "tau_t in component Ca_independent_transient_outward_K_current_s_gate (second)"
    legend_algebraic[21] = "s_max in component Ca_independent_transient_outward_K_current_s_gate (dimensionless)"
    legend_algebraic[44] = "tau_s_slow in component Ca_independent_transient_outward_K_current_s_slow_gate (second)"
    legend_algebraic[22] = "s_slow_max in component Ca_independent_transient_outward_K_current_s_slow_gate (dimensionless)"
    legend_constants[84] = "g_Kss in component steady_state_outward_K_current (milliS_per_cm2)"
    legend_constants[134] = "g_Ksst in component steady_state_outward_K_current (milliS)"
    legend_constants[85] = "fKsst in component steady_state_outward_K_current (dimensionless)"
    legend_states[38] = "r_ss in component steady_state_outward_K_current_r_ss_gate (dimensionless)"
    legend_states[39] = "s_ss in component steady_state_outward_K_current_s_ss_gate (dimensionless)"
    legend_algebraic[45] = "tau_r_ss in component steady_state_outward_K_current_r_ss_gate (second)"
    legend_algebraic[23] = "r_ss_max in component steady_state_outward_K_current_r_ss_gate (dimensionless)"
    legend_constants[110] = "tau_s_ss in component steady_state_outward_K_current_s_ss_gate (second)"
    legend_algebraic[24] = "s_ss_max in component steady_state_outward_K_current_s_ss_gate (dimensionless)"
    legend_constants[86] = "g_K1 in component inward_rectifier (milliS_per_cm2)"
    legend_constants[135] = "g_K1t in component inward_rectifier (milliS)"
    legend_constants[87] = "fK1t in component inward_rectifier (dimensionless)"
    legend_constants[88] = "g_f in component hyperpolarisation_activated_current (milliS_per_cm2)"
    legend_constants[136] = "g_ft in component hyperpolarisation_activated_current (milliS)"
    legend_constants[89] = "fft in component hyperpolarisation_activated_current (dimensionless)"
    legend_constants[90] = "f_Na in component hyperpolarisation_activated_current (dimensionless)"
    legend_constants[111] = "f_K in component hyperpolarisation_activated_current (dimensionless)"
    legend_states[40] = "y in component hyperpolarisation_activated_current_y_gate (dimensionless)"
    legend_algebraic[46] = "tau_y in component hyperpolarisation_activated_current_y_gate (second)"
    legend_algebraic[25] = "y_infinity in component hyperpolarisation_activated_current_y_gate (dimensionless)"
    legend_constants[91] = "g_B_Ca in component background_currents (milliS_per_cm2)"
    legend_constants[92] = "g_B_Na in component background_currents (milliS_per_cm2)"
    legend_constants[93] = "g_B_K in component background_currents (milliS_per_cm2)"
    legend_constants[137] = "g_B_Cat in component background_currents (milliS)"
    legend_constants[138] = "g_B_Nat in component background_currents (milliS)"
    legend_constants[139] = "g_B_Kt in component background_currents (milliS)"
    legend_constants[94] = "fCabt in component background_currents (dimensionless)"
    legend_constants[95] = "fNabt in component background_currents (dimensionless)"
    legend_constants[96] = "fKbt in component background_currents (dimensionless)"
    legend_constants[97] = "E_Ca in component background_currents (millivolt)"
    legend_constants[98] = "i_NaK_max in component sodium_potassium_pump (microA_per_cm2)"
    legend_constants[140] = "i_NaK_max_t in component sodium_potassium_pump (microA)"
    legend_constants[99] = "fNaKt in component sodium_potassium_pump (dimensionless)"
    legend_constants[100] = "K_m_K in component sodium_potassium_pump (millimolar)"
    legend_constants[101] = "K_m_Na in component sodium_potassium_pump (millimolar)"
    legend_algebraic[89] = "sigma in component sodium_potassium_pump (dimensionless)"
    legend_constants[102] = "i_CaP_max in component sarcolemmal_calcium_pump_current (microA_per_cm2)"
    legend_constants[141] = "i_CaP_max_t in component sarcolemmal_calcium_pump_current (microA)"
    legend_constants[103] = "fCaPt in component sarcolemmal_calcium_pump_current (dimensionless)"
    legend_constants[104] = "i_NaCa_max in component Na_Ca_ion_exchanger_current (microA_per_cm2)"
    legend_constants[142] = "i_NaCa_max_t in component Na_Ca_ion_exchanger_current (microA)"
    legend_constants[105] = "fNaCat in component Na_Ca_ion_exchanger_current (dimensionless)"
    legend_rates[4] = "d/dt F1 in component SR_Ca_release_channel (dimensionless)"
    legend_rates[2] = "d/dt F2 in component SR_Ca_release_channel (dimensionless)"
    legend_rates[3] = "d/dt F3 in component SR_Ca_release_channel (dimensionless)"
    legend_rates[5] = "d/dt F4 in component SR_Ca_release_channel (dimensionless)"
    legend_rates[10] = "d/dt BTRH in component intracellular_and_SR_Ca_fluxes (millimolar)"
    legend_rates[11] = "d/dt BTRL in component intracellular_and_SR_Ca_fluxes (millimolar)"
    legend_rates[12] = "d/dt Na_i in component ion_concentrations (millimolar)"
    legend_rates[14] = "d/dt Na_t in component ion_concentrations (millimolar)"
    legend_rates[13] = "d/dt K_i in component ion_concentrations (millimolar)"
    legend_rates[16] = "d/dt K_t in component ion_concentrations (millimolar)"
    legend_rates[8] = "d/dt Ca_i in component ion_concentrations (millimolar)"
    legend_rates[15] = "d/dt Ca_t in component ion_concentrations (millimolar)"
    legend_rates[6] = "d/dt Ca_ss in component ion_concentrations (millimolar)"
    legend_rates[7] = "d/dt CaSRrel in component ion_concentrations (millimolar)"
    legend_rates[9] = "d/dt CaSRup in component ion_concentrations (millimolar)"
    legend_rates[0] = "d/dt Vm_s in component membrane (millivolt)"
    legend_rates[17] = "d/dt m in component sodium_current_m_gate (dimensionless)"
    legend_rates[18] = "d/dt h in component sodium_current_h_gate (dimensionless)"
    legend_rates[19] = "d/dt Co in component calcium_current (dimensionless)"
    legend_rates[21] = "d/dt C1 in component calcium_current (dimensionless)"
    legend_rates[22] = "d/dt C2 in component calcium_current (dimensionless)"
    legend_rates[20] = "d/dt C_cast in component calcium_current (dimensionless)"
    legend_rates[23] = "d/dt r in component Ca_independent_transient_outward_K_current_r_gate (dimensionless)"
    legend_rates[24] = "d/dt s in component Ca_independent_transient_outward_K_current_s_gate (dimensionless)"
    legend_rates[25] = "d/dt s_slow in component Ca_independent_transient_outward_K_current_s_slow_gate (dimensionless)"
    legend_rates[26] = "d/dt r_ss in component steady_state_outward_K_current_r_ss_gate (dimensionless)"
    legend_rates[27] = "d/dt s_ss in component steady_state_outward_K_current_s_ss_gate (dimensionless)"
    legend_rates[28] = "d/dt y in component hyperpolarisation_activated_current_y_gate (dimensionless)"
    legend_rates[1] = "d/dt Vm_t in component membrane (millivolt)"
    legend_rates[29] = "d/dt m in component sodium_current_m_gate (dimensionless)"
    legend_rates[30] = "d/dt h in component sodium_current_h_gate (dimensionless)"
    legend_rates[31] = "d/dt Co in component calcium_current (dimensionless)"
    legend_rates[33] = "d/dt C1 in component calcium_current (dimensionless)"
    legend_rates[34] = "d/dt C2 in component calcium_current (dimensionless)"
    legend_rates[32] = "d/dt C_cast in component calcium_current (dimensionless)"
    legend_rates[35] = "d/dt r in component Ca_independent_transient_outward_K_current_r_gate (dimensionless)"
    legend_rates[36] = "d/dt s in component Ca_independent_transient_outward_K_current_s_gate (dimensionless)"
    legend_rates[37] = "d/dt s_slow in component Ca_independent_transient_outward_K_current_s_slow_gate (dimensionless)"
    legend_rates[38] = "d/dt r_ss in component steady_state_outward_K_current_r_ss_gate (dimensionless)"
    legend_rates[39] = "d/dt s_ss in component steady_state_outward_K_current_s_ss_gate (dimensionless)"
    legend_rates[40] = "d/dt y in component hyperpolarisation_activated_current_y_gate (dimensionless)"
    return (legend_states, legend_algebraic, legend_voi, legend_constants)

def initConsts():
    constants = [0.0] * sizeConstants; states = [0.0] * sizeStates;
    constants[0] = 8310
    constants[1] = 295
    constants[2] = 96500
    constants[3] = 140
    constants[4] = 1.2
    constants[5] = 5.4
    constants[6] = 83.33
    constants[7] = 150e-7
    constants[8] = 30000000
    constants[9] = 0.0001
    constants[10] = 20.0
    states[0] = -7.9709547e1
    states[1] = -7.9709536e1
    constants[11] = 12.15e12
    constants[12] = 576
    constants[13] = 4.05e9
    constants[14] = 1930
    constants[15] = 100
    constants[16] = 0.8
    states[2] = 3.8333331e-4
    states[3] = 6.4516879e-11
    states[4] = 6.7606483e-1
    states[5] = 3.2355462e-1
    states[6] = 4.3112114e-5
    states[7] = 2.2638169e-1
    constants[17] = 0.000168
    constants[18] = 3.29
    constants[19] = 0.4
    constants[20] = 0.9
    constants[21] = 1
    constants[22] = 1.2
    constants[23] = 1
    states[8] = 3.5486498e-5
    states[9] = 2.2755058e-1
    constants[24] = 0.14
    constants[25] = 0.07
    constants[26] = 0.02
    constants[27] = 0.008
    states[10] = 9.2804169e-1
    states[11] = 3.4317829e-2
    constants[28] = 20000
    constants[29] = 0.07
    constants[30] = 40000
    constants[31] = 40
    states[12] = 6.9435814e0
    states[13] = 1.3753311e2
    states[14] = 1.4001386e2
    states[15] = 1.2271399e0
    states[16] = 5.3627676e0
    constants[32] = 0.00238
    constants[33] = 0.8
    constants[34] = 0.00015
    constants[35] = 0.05
    constants[36] = 15
    constants[37] = 10
    constants[38] = 0.008
    constants[39] = 1
    constants[40] = 0
    constants[41] = 0.15
    constants[42] = 0.5
    constants[43] = 0.15
    constants[44] = 10.0
    constants[45] = 0.56
    states[17] = 2.1832558e-6
    states[18] = 5.8232781e-1
    states[19] = 2.5436465e-6
    states[20] = 9.9570063e-1
    states[21] = 9.9994546e-1
    states[22] = 9.9994546e-1
    constants[46] = 0.009
    constants[47] = 0.00021
    constants[48] = 0.87
    constants[49] = 0.886
    constants[50] = 0.114
    constants[51] = 0.35
    constants[52] = 0.56
    states[23] = 2.3485072e-3
    states[24] = 9.8315843e-1
    states[25] = 6.7483370e-1
    constants[53] = 0.07
    constants[54] = 0.76
    states[26] = 3.1079213e-3
    states[27] = 3.0773393e-1
    constants[55] = 0.24
    constants[56] = 0.56
    constants[57] = 0.0145
    constants[58] = 0.56
    constants[59] = 0.2
    states[28] = 3.4357356e-3
    constants[60] = 0.000648
    constants[61] = 0.0008015
    constants[62] = 0.00138
    constants[63] = 0.56
    constants[64] = 0.56
    constants[65] = 0.56
    constants[66] = 65
    constants[67] = 1
    constants[68] = 0.59
    constants[69] = 1.5
    constants[70] = 10
    constants[71] = 0.85
    constants[72] = 0.56
    constants[73] = 1.8e-4
    constants[74] = 0.81
    constants[75] = 10.0
    constants[76] = 0.56
    states[29] = 2.1832694e-6
    states[30] = 5.8232717e-1
    states[31] = 2.5436520e-6
    states[32] = 9.9570063e-1
    states[33] = 9.9785393e-1
    states[34] = 9.9994546e-1
    constants[77] = 0.009
    constants[78] = 0.00021
    constants[79] = 0.87
    constants[80] = 0.886
    constants[81] = 0.114
    constants[82] = 0.35
    constants[83] = 0.56
    states[35] = 2.3485094e-3
    states[36] = 9.8315865e-1
    states[37] = 6.7478532e-1
    constants[84] = 0.07
    constants[85] = 0.76
    states[38] = 3.1079241e-3
    states[39] = 3.0773368e-1
    constants[86] = 0.24
    constants[87] = 0.56
    constants[88] = 0.0145
    constants[89] = 0.56
    constants[90] = 0.2
    states[40] = 3.4357175e-3
    constants[91] = 0.000648
    constants[92] = 0.0008015
    constants[93] = 0.00138
    constants[94] = 0.56
    constants[95] = 0.56
    constants[96] = 0.56
    constants[97] = 65
    constants[98] = 1
    constants[99] = 0.59
    constants[100] = 1.5
    constants[101] = 10
    constants[102] = 0.85
    constants[103] = 0.56
    constants[104] = 1.8e-4
    constants[105] = 0.81
    constants[106] = constants[9]/7850.00
    constants[107] = 2.10000
    constants[108] = 1.00000-constants[59]
    constants[109] = (exp(constants[3]/67.3000)-1.00000)/7.00000
    constants[110] = 2.10000
    constants[111] = 1.00000-constants[90]
    constants[112] = constants[106]*7.50000e-05
    constants[113] = constants[106]*0.585000
    constants[114] = constants[106]*0.0315000
    constants[115] = constants[106]*0.00350000
    constants[116] = constants[9]-4400.00*constants[106]
    constants[117] = constants[9]-constants[116]
    constants[118] = constants[8]*constants[116]
    constants[119] = constants[44]*(constants[116]+constants[117])*(1.00000-constants[45])
    constants[120] = constants[47]*(constants[116]+constants[117])*(1.00000-constants[48])
    constants[121] = constants[51]*(constants[116]+constants[117])*(1.00000-constants[52])
    constants[122] = constants[53]*(constants[116]+constants[117])*(1.00000-constants[54])
    constants[123] = constants[55]*(constants[116]+constants[117])*(1.00000-constants[56])
    constants[124] = constants[57]*(constants[116]+constants[117])*(1.00000-constants[58])
    constants[125] = constants[60]*(constants[116]+constants[117])*(1.00000-constants[63])
    constants[126] = constants[61]*(constants[116]+constants[117])*(1.00000-constants[64])
    constants[127] = constants[62]*(constants[116]+constants[117])*(1.00000-constants[65])
    constants[128] = constants[67]*(constants[116]+constants[117])*(1.00000-constants[68])
    constants[129] = constants[71]*(constants[116]+constants[117])*(1.00000-constants[72])
    constants[130] = constants[73]*(constants[116]+constants[117])*(1.00000-constants[74])
    constants[131] = constants[75]*(constants[116]+constants[117])*constants[76]
    constants[132] = constants[78]*(constants[116]+constants[117])*constants[79]
    constants[133] = constants[82]*(constants[116]+constants[117])*constants[83]
    constants[134] = constants[84]*(constants[116]+constants[117])*constants[85]
    constants[135] = constants[86]*(constants[116]+constants[117])*constants[87]
    constants[136] = constants[88]*(constants[116]+constants[117])*constants[89]
    constants[137] = constants[91]*(constants[116]+constants[117])*constants[94]
    constants[138] = constants[92]*(constants[116]+constants[117])*constants[95]
    constants[139] = constants[93]*(constants[116]+constants[117])*constants[96]
    constants[140] = constants[98]*(constants[116]+constants[117])*constants[99]
    constants[141] = constants[102]*(constants[116]+constants[117])*constants[103]
    constants[142] = constants[104]*(constants[116]+constants[117])*constants[105]
    constants[143] = constants[117]/(2.00000*3.14159*constants[7]*constants[118])
    constants[144] = (constants[6]*(constants[143]/2.00000))/(3.14159*(power(constants[7], 2.00000))*constants[118])
    constants[145] = 3.14159*(power(constants[7], 2.00000))*constants[143]*constants[118]
    constants[146] = constants[117]*1.00000
    constants[147] = constants[116]*1.00000
    return (states, constants)

def computeRates(voi, states, constants):
    rates = [0.0] * sizeStates; algebraic = [0.0] * sizeAlgebraic
    rates[4] = -constants[11]*(power(states[6], 4.00000))*states[4]+constants[12]*states[2]
    rates[2] = (constants[11]*(power(states[6], 4.00000))*states[4]-(constants[12]*states[2]+constants[13]*(power(states[6], 3.00000))*states[2]+constants[15]*states[2]))+constants[14]*states[3]+constants[16]*states[5]
    rates[3] = constants[13]*(power(states[6], 3.00000))*states[2]-constants[14]*states[3]
    rates[5] = constants[15]*states[2]-constants[16]*states[5]
    algebraic[5] = 1.00000/(1.00000+states[6]/0.0100000)
    rates[20] = (algebraic[5]-states[20])/constants[46]
    algebraic[12] = 1.00000/(1.00000+exp((states[0]+87.5000)/10.3000))
    rates[27] = (algebraic[12]-states[27])/constants[107]
    algebraic[17] = 1.00000/(1.00000+states[6]/0.0100000)
    rates[32] = (algebraic[17]-states[32])/constants[77]
    algebraic[24] = 1.00000/(1.00000+exp((states[1]+87.5000)/10.3000))
    rates[39] = (algebraic[24]-states[39])/constants[110]
    algebraic[2] = (117.260*(states[0]+59.3000))/(1.00000-exp(-0.550000*(states[0]+59.3000)))
    algebraic[27] = 3800.00*exp(-0.0720000*(states[0]+61.0000))
    rates[17] = algebraic[2]*(1.00000-states[17])-algebraic[27]*states[17]
    algebraic[3] = 284.400/(1.00000+exp(0.0812000*(states[0]+115.900)))
    algebraic[28] = (18.7700*(states[0]+64.4000))/(1.00000-exp(-0.220000*(states[0]+64.4000)))
    rates[18] = algebraic[3]*(1.00000-states[18])-algebraic[28]*states[18]
    algebraic[29] = 0.00305000*exp(-0.00450000*(power(states[0]+7.00000, 2.00000)))+0.00105000*exp(-0.00200000*(power(states[0]-18.0000, 2.00000)))+0.000250000
    algebraic[4] = 1.00000/(1.00000+exp((states[0]+15.3000)/-5.00000))
    rates[19] = (algebraic[4]-states[19])/algebraic[29]
    algebraic[30] = 0.105000*exp(-(power((states[0]+45.0000)/12.0000, 2.00000)))+0.0400000/(1.00000+exp((-states[0]+25.0000)/25.0000))+0.0150000/(1.00000+exp((states[0]+75.0000)/25.0000))+0.00170000
    algebraic[6] = 1.00000/(1.00000+exp((states[0]+26.7000)/5.40000))
    rates[21] = (algebraic[6]-states[21])/algebraic[30]
    algebraic[31] = 0.0410000*exp(-(power((states[0]+47.0000)/12.0000, 2.00000)))+0.0800000/(1.00000+exp((states[0]+55.0000)/-5.00000))+0.0150000/(1.00000+exp((states[0]+75.0000)/25.0000))+0.00170000
    algebraic[7] = 1.00000/(1.00000+exp((states[0]+26.7000)/5.40000))
    rates[22] = (algebraic[7]-states[22])/algebraic[31]
    algebraic[32] = 1.00000/(45.1600*exp(0.0357700*(states[0]+50.0000))+98.9000*exp(-0.100000*(states[0]+38.0000)))
    algebraic[8] = 1.00000/(1.00000+exp((states[0]+10.6000)/-11.4200))
    rates[23] = (algebraic[8]-states[23])/algebraic[32]
    algebraic[33] = 0.350000*exp(-(power(0.0666000*(states[0]+70.0000), 2.00000)))+0.0350000
    algebraic[9] = 1.00000/(1.00000+exp((states[0]+45.3000)/6.88410))
    rates[24] = (algebraic[9]-states[24])/algebraic[33]
    algebraic[34] = 3.70000*exp(-(power(0.0333000*(states[0]+70.0000), 2.00000)))+0.0350000
    algebraic[10] = 1.00000/(1.00000+exp((states[0]+45.3000)/6.88410))
    rates[25] = (algebraic[10]-states[25])/algebraic[34]
    algebraic[35] = 10.0000/(45.1600*exp(0.0357700*(states[0]+50.0000))+98.9000*exp(-0.100000*(states[0]+38.0000)))
    algebraic[11] = 1.00000/(1.00000+exp(-(states[0]+11.5000)/11.8200))
    rates[26] = (algebraic[11]-states[26])/algebraic[35]
    algebraic[36] = 1.00000/(0.118850*exp((states[0]+80.0000)/28.3700)+0.562300*exp((states[0]+80.0000)/-14.1900))
    algebraic[13] = 1.00000/(1.00000+exp((states[0]+138.600)/10.4800))
    rates[28] = (algebraic[13]-states[28])/algebraic[36]
    algebraic[14] = (117.260*(states[1]+59.3000))/(1.00000-exp(-0.550000*(states[1]+59.3000)))
    algebraic[37] = 3800.00*exp(-0.0720000*(states[1]+61.0000))
    rates[29] = algebraic[14]*(1.00000-states[29])-algebraic[37]*states[29]
    algebraic[15] = 284.400/(1.00000+exp(0.0812000*(states[1]+115.900)))
    algebraic[38] = (18.7700*(states[1]+64.4000))/(1.00000-exp(-0.220000*(states[1]+64.4000)))
    rates[30] = algebraic[15]*(1.00000-states[30])-algebraic[38]*states[30]
    algebraic[39] = 0.00305000*exp(-0.00450000*(power(states[1]+7.00000, 2.00000)))+0.00105000*exp(-0.00200000*(power(states[1]-18.0000, 2.00000)))+0.000250000
    algebraic[16] = 1.00000/(1.00000+exp((states[1]+15.3000)/-5.00000))
    rates[31] = (algebraic[16]-states[31])/algebraic[39]
    algebraic[40] = 0.105000*exp(-(power((states[1]+45.0000)/12.0000, 2.00000)))+0.0400000/(1.00000+exp((-states[1]+25.0000)/25.0000))+0.0150000/(1.00000+exp((states[1]+75.0000)/25.0000))+0.00170000
    algebraic[18] = 1.00000/(1.00000+exp((states[1]+26.7000)/5.40000))
    rates[33] = (algebraic[18]-states[33])/algebraic[40]
    algebraic[41] = 0.0410000*exp(-(power((states[1]+47.0000)/12.0000, 2.00000)))+0.0800000/(1.00000+exp((states[1]+55.0000)/-5.00000))+0.0150000/(1.00000+exp((states[1]+75.0000)/25.0000))+0.00170000
    algebraic[19] = 1.00000/(1.00000+exp((states[1]+26.7000)/5.40000))
    rates[34] = (algebraic[19]-states[34])/algebraic[41]
    algebraic[42] = 1.00000/(45.1600*exp(0.0357700*(states[1]+50.0000))+98.9000*exp(-0.100000*(states[1]+38.0000)))
    algebraic[20] = 1.00000/(1.00000+exp((states[1]+10.6000)/-11.4200))
    rates[35] = (algebraic[20]-states[35])/algebraic[42]
    algebraic[43] = 0.350000*exp(-(power(0.0666000*(states[1]+70.0000), 2.00000)))+0.0350000
    algebraic[21] = 1.00000/(1.00000+exp((states[1]+45.3000)/6.88410))
    rates[36] = (algebraic[21]-states[36])/algebraic[43]
    algebraic[44] = 3.70000*exp(-(power(0.0333000*(states[1]+70.0000), 2.00000)))+0.0350000
    algebraic[22] = 1.00000/(1.00000+exp((states[1]+45.3000)/6.88410))
    rates[37] = (algebraic[22]-states[37])/algebraic[44]
    algebraic[45] = 10.0000/(45.1600*exp(0.0357700*(states[1]+50.0000))+98.9000*exp(-0.100000*(states[1]+38.0000)))
    algebraic[23] = 1.00000/(1.00000+exp(-(states[1]+11.5000)/11.8200))
    rates[38] = (algebraic[23]-states[38])/algebraic[45]
    algebraic[46] = 1.00000/(0.118850*exp((states[1]+80.0000)/28.3700)+0.562300*exp((states[1]+80.0000)/-14.1900))
    algebraic[25] = 1.00000/(1.00000+exp((states[1]+138.600)/10.4800))
    rates[40] = (algebraic[25]-states[40])/algebraic[46]
    algebraic[51] = constants[28]*states[8]*(1.00000-states[10])-constants[29]*states[10]
    rates[10] = algebraic[51]
    algebraic[26] = 1800.00*(states[2]+states[3])*(states[7]-states[6])
    algebraic[50] = (states[9]-states[7])/constants[26]
    rates[7] = (1.00000/(1.00000+(constants[36]*constants[33])/(power(constants[33]+states[7], 2.00000))))*(algebraic[50]-algebraic[26])
    algebraic[47] = power(states[8]/constants[17], constants[22])
    algebraic[48] = power(states[9]/constants[18], constants[23])
    algebraic[49] = (constants[21]*(constants[19]*algebraic[47]-constants[20]*algebraic[48]))/(1.00000+algebraic[47]+algebraic[48])
    rates[9] = algebraic[49]-algebraic[50]
    algebraic[52] = constants[30]*states[8]*(1.00000-states[11])-constants[31]*states[11]
    rates[11] = algebraic[52]
    algebraic[1] = 1000.00*((states[0]-states[1])/constants[10])
    algebraic[54] = custom_piecewise([equal(constants[39] , 1.00000) & less(voi-constants[40] , 0.00100000), 52.0000*(constants[116]+constants[117]) , True, 0.00000])
    algebraic[59] = constants[119]*(power(states[17], 3.00000))*states[18]*(states[0]-log((constants[3]+0.120000*constants[5])/(states[12]+0.120000*states[13]), 10)*((constants[0]*constants[1])/constants[2]))
    algebraic[61] = constants[120]*4.00000*((states[0]*(power(constants[2], 2.00000)))/(constants[0]*constants[1]))*((states[6]*exp((2.00000*states[0]*(power(constants[2], 2.00000)))/(constants[0]*constants[1]))-0.341000*constants[4])/(exp((2.00000*states[0]*(power(constants[2], 2.00000)))/(constants[0]*constants[1]))-1.00000))*states[19]*((0.900000+states[20]/10.0000)*states[21]+(0.100000+states[20]/10.0000)*states[22])
    algebraic[74] = (constants[130]*(exp(0.0187000*states[0])*(power(states[12], 3.00000))*constants[4]-exp(-0.0187000*states[0])*(power(constants[3], 3.00000))*states[8]))/(1.00000+0.000100000*(states[8]*(power(constants[3], 3.00000))+constants[4]*(power(states[12], 3.00000))))
    algebraic[72] = (((constants[128]/(1.00000+0.124500*exp((-0.100000*states[0]*constants[2])/(constants[0]*constants[1]))+0.0365000*constants[109]*exp((-states[0]*constants[2])/(constants[0]*constants[1]))))*constants[5])/(constants[5]+constants[69]))/(1.00000+power(constants[70]/states[12], 1.50000))
    algebraic[62] = ((constants[0]*constants[1])/constants[2])*log(constants[5]/states[13])
    algebraic[65] = ((constants[123]*(states[0]-(algebraic[62]+1.73000)))/(1.00000+exp(((1.61300*constants[2])/(constants[0]*constants[1]))*(states[0]-(algebraic[62]+1.73000)))))*(1.00000+exp(-(constants[5]-0.998800)/0.124000))
    algebraic[63] = constants[121]*states[23]*(constants[49]*states[24]+constants[50]*states[25])*(states[0]-algebraic[62])
    algebraic[64] = constants[122]*states[26]*states[27]*(states[0]-algebraic[62])
    algebraic[73] = (constants[129]*states[8])/(states[8]+0.000400000)
    algebraic[58] = ((constants[0]*constants[1])/constants[2])*log(constants[3]/states[12])
    algebraic[67] = constants[124]*states[28]*(0.200000*(states[0]-algebraic[58])+0.800000*(states[0]-algebraic[62]))
    algebraic[68] = constants[126]*(states[0]-algebraic[58])
    algebraic[70] = constants[127]*(states[0]-algebraic[62])
    algebraic[69] = constants[125]*(states[0]-constants[66])
    algebraic[71] = algebraic[68]+algebraic[69]+algebraic[70]
    rates[0] = 1000.00*(((algebraic[54]+algebraic[1])-(algebraic[59]+algebraic[61]+algebraic[63]+algebraic[64]+algebraic[67]+algebraic[65]+algebraic[71]+algebraic[72]+algebraic[74]+algebraic[73]))/constants[147])
    algebraic[78] = constants[132]*4.00000*((states[1]*(power(constants[2], 2.00000)))/(constants[0]*constants[1]))*((states[6]*exp((2.00000*states[1]*(power(constants[2], 2.00000)))/(constants[0]*constants[1]))-0.341000*states[15])/(exp((2.00000*states[1]*(power(constants[2], 2.00000)))/(constants[0]*constants[1]))-1.00000))*states[31]*((0.900000+states[32]/10.0000)*states[33]+(0.100000+states[32]/10.0000)*states[34])
    rates[6] = (1.00000/(1.00000+constants[35]*(constants[32]/(power(constants[32]+states[6], 2.00000)))))*(-(algebraic[61]+algebraic[78])/((2.00000*constants[112]*constants[2]+algebraic[26]/constants[112])-(((states[6]-states[8])*constants[113])/constants[38])/constants[112]))
    algebraic[66] = constants[124]*states[28]*constants[108]*(states[0]-algebraic[62])
    algebraic[89] = (exp(states[14]/67.3000)-1.00000)/7.00000
    algebraic[90] = (((constants[140]/(1.00000+0.124500*exp((-0.100000*states[1]*constants[2])/(constants[0]*constants[1]))+0.0365000*algebraic[89]*exp((-states[1]*constants[2])/(constants[0]*constants[1]))))*states[16])/(states[16]+constants[100]))/(1.00000+power(constants[101]/states[12], 1.50000))
    algebraic[79] = ((constants[0]*constants[1])/constants[2])*log(states[16]/states[13])
    algebraic[83] = constants[136]*states[40]*constants[111]*(states[1]-algebraic[79])
    algebraic[87] = constants[139]*(states[1]-algebraic[79])
    algebraic[82] = ((constants[135]*(states[1]-(algebraic[79]+1.73000)))/(1.00000+exp(((1.61300*constants[2])/(constants[0]*constants[1]))*(states[1]-(algebraic[79]+1.73000)))))*(1.00000+exp(-(states[16]-0.998800)/0.124000))
    algebraic[80] = constants[133]*states[35]*(constants[80]*states[36]+constants[81]*states[37])*(states[1]-algebraic[79])
    algebraic[81] = constants[134]*states[38]*states[39]*(states[1]-algebraic[79])
    rates[13] = -(-algebraic[54]+algebraic[66]+algebraic[83]+algebraic[64]+algebraic[81]+algebraic[65]+algebraic[82]+algebraic[70]+algebraic[87]+algebraic[72]*-2.00000+algebraic[90]*-2.00000+algebraic[63]+algebraic[80])/(constants[113]*constants[2])
    algebraic[57] = (1.00000*constants[145])/(constants[43]*(states[16]-constants[5]))
    rates[16] = ((algebraic[83]+algebraic[81]+algebraic[82]+algebraic[87]+algebraic[90]*-2.00000+algebraic[80])/(constants[2]-algebraic[57]))/constants[145]
    algebraic[76] = constants[131]*(power(states[29], 3.00000))*states[30]*(states[1]-log((states[14]+0.120000*states[16])/(states[12]+0.120000*states[13]), 10)*((constants[0]*constants[1])/constants[2]))
    algebraic[75] = ((constants[0]*constants[1])/constants[2])*log(states[14]/states[12])
    algebraic[85] = constants[138]*(states[1]-algebraic[75])
    algebraic[92] = (constants[142]*(exp(0.0187000*states[1])*(power(states[12], 3.00000))*states[15]-exp(-0.0187000*states[1])*(power(states[14], 3.00000))*states[8]))/(1.00000+0.000100000*(states[8]*(power(states[14], 3.00000))+states[15]*(power(states[12], 3.00000))))
    rates[12] = -(algebraic[59]+algebraic[76]+algebraic[68]+algebraic[85]+algebraic[74]*3.00000+algebraic[92]*3.00000+algebraic[72]*3.00000+algebraic[90]*3.00000)/(constants[113]*constants[2])
    algebraic[55] = (1.00000*constants[145])/(constants[41]*(states[14]-constants[3]))
    rates[14] = ((algebraic[76]+algebraic[85]+algebraic[92]*3.00000+algebraic[90]*3.00000)/(constants[2]-algebraic[55]))/constants[145]
    algebraic[91] = (constants[141]*states[8])/(states[8]+0.000400000)
    algebraic[86] = constants[137]*(states[1]-constants[97])
    rates[8] = (1.00000/(1.00000+constants[35]*(constants[32]/(power(constants[32]+states[8], 2.00000)))))*(((-2.00000*(algebraic[74]+algebraic[92])+(algebraic[69]+algebraic[86]+algebraic[73]+algebraic[91]))/(-2.00000*constants[113]*constants[2])+((((states[6]-states[8])*constants[113])/constants[38]+algebraic[26])-algebraic[49])/constants[113])-(algebraic[51]*constants[24]+algebraic[52]*constants[25]))
    algebraic[56] = (1.00000*constants[145])/(constants[42]*(states[15]-constants[4]))
    rates[15] = ((algebraic[78]+algebraic[91]+algebraic[86]+algebraic[92]*-2.00000)/(2.00000*constants[2]-algebraic[56]))/constants[145]
    algebraic[84] = constants[136]*states[40]*(0.200000*(states[1]-algebraic[75])+0.800000*(states[1]-algebraic[79]))
    algebraic[88] = algebraic[85]+algebraic[86]+algebraic[87]
    rates[1] = 1000.00*(-(algebraic[1]+algebraic[76]+algebraic[78]+algebraic[80]+algebraic[81]+algebraic[84]+algebraic[82]+algebraic[88]+algebraic[90]+algebraic[92]+algebraic[91])/constants[146])
    return(rates)

def computeAlgebraic(constants, states, voi):
    algebraic = array([[0.0] * len(voi)] * sizeAlgebraic)
    states = array(states)
    voi = array(voi)
    algebraic[5] = 1.00000/(1.00000+states[6]/0.0100000)
    algebraic[12] = 1.00000/(1.00000+exp((states[0]+87.5000)/10.3000))
    algebraic[17] = 1.00000/(1.00000+states[6]/0.0100000)
    algebraic[24] = 1.00000/(1.00000+exp((states[1]+87.5000)/10.3000))
    algebraic[2] = (117.260*(states[0]+59.3000))/(1.00000-exp(-0.550000*(states[0]+59.3000)))
    algebraic[27] = 3800.00*exp(-0.0720000*(states[0]+61.0000))
    algebraic[3] = 284.400/(1.00000+exp(0.0812000*(states[0]+115.900)))
    algebraic[28] = (18.7700*(states[0]+64.4000))/(1.00000-exp(-0.220000*(states[0]+64.4000)))
    algebraic[29] = 0.00305000*exp(-0.00450000*(power(states[0]+7.00000, 2.00000)))+0.00105000*exp(-0.00200000*(power(states[0]-18.0000, 2.00000)))+0.000250000
    algebraic[4] = 1.00000/(1.00000+exp((states[0]+15.3000)/-5.00000))
    algebraic[30] = 0.105000*exp(-(power((states[0]+45.0000)/12.0000, 2.00000)))+0.0400000/(1.00000+exp((-states[0]+25.0000)/25.0000))+0.0150000/(1.00000+exp((states[0]+75.0000)/25.0000))+0.00170000
    algebraic[6] = 1.00000/(1.00000+exp((states[0]+26.7000)/5.40000))
    algebraic[31] = 0.0410000*exp(-(power((states[0]+47.0000)/12.0000, 2.00000)))+0.0800000/(1.00000+exp((states[0]+55.0000)/-5.00000))+0.0150000/(1.00000+exp((states[0]+75.0000)/25.0000))+0.00170000
    algebraic[7] = 1.00000/(1.00000+exp((states[0]+26.7000)/5.40000))
    algebraic[32] = 1.00000/(45.1600*exp(0.0357700*(states[0]+50.0000))+98.9000*exp(-0.100000*(states[0]+38.0000)))
    algebraic[8] = 1.00000/(1.00000+exp((states[0]+10.6000)/-11.4200))
    algebraic[33] = 0.350000*exp(-(power(0.0666000*(states[0]+70.0000), 2.00000)))+0.0350000
    algebraic[9] = 1.00000/(1.00000+exp((states[0]+45.3000)/6.88410))
    algebraic[34] = 3.70000*exp(-(power(0.0333000*(states[0]+70.0000), 2.00000)))+0.0350000
    algebraic[10] = 1.00000/(1.00000+exp((states[0]+45.3000)/6.88410))
    algebraic[35] = 10.0000/(45.1600*exp(0.0357700*(states[0]+50.0000))+98.9000*exp(-0.100000*(states[0]+38.0000)))
    algebraic[11] = 1.00000/(1.00000+exp(-(states[0]+11.5000)/11.8200))
    algebraic[36] = 1.00000/(0.118850*exp((states[0]+80.0000)/28.3700)+0.562300*exp((states[0]+80.0000)/-14.1900))
    algebraic[13] = 1.00000/(1.00000+exp((states[0]+138.600)/10.4800))
    algebraic[14] = (117.260*(states[1]+59.3000))/(1.00000-exp(-0.550000*(states[1]+59.3000)))
    algebraic[37] = 3800.00*exp(-0.0720000*(states[1]+61.0000))
    algebraic[15] = 284.400/(1.00000+exp(0.0812000*(states[1]+115.900)))
    algebraic[38] = (18.7700*(states[1]+64.4000))/(1.00000-exp(-0.220000*(states[1]+64.4000)))
    algebraic[39] = 0.00305000*exp(-0.00450000*(power(states[1]+7.00000, 2.00000)))+0.00105000*exp(-0.00200000*(power(states[1]-18.0000, 2.00000)))+0.000250000
    algebraic[16] = 1.00000/(1.00000+exp((states[1]+15.3000)/-5.00000))
    algebraic[40] = 0.105000*exp(-(power((states[1]+45.0000)/12.0000, 2.00000)))+0.0400000/(1.00000+exp((-states[1]+25.0000)/25.0000))+0.0150000/(1.00000+exp((states[1]+75.0000)/25.0000))+0.00170000
    algebraic[18] = 1.00000/(1.00000+exp((states[1]+26.7000)/5.40000))
    algebraic[41] = 0.0410000*exp(-(power((states[1]+47.0000)/12.0000, 2.00000)))+0.0800000/(1.00000+exp((states[1]+55.0000)/-5.00000))+0.0150000/(1.00000+exp((states[1]+75.0000)/25.0000))+0.00170000
    algebraic[19] = 1.00000/(1.00000+exp((states[1]+26.7000)/5.40000))
    algebraic[42] = 1.00000/(45.1600*exp(0.0357700*(states[1]+50.0000))+98.9000*exp(-0.100000*(states[1]+38.0000)))
    algebraic[20] = 1.00000/(1.00000+exp((states[1]+10.6000)/-11.4200))
    algebraic[43] = 0.350000*exp(-(power(0.0666000*(states[1]+70.0000), 2.00000)))+0.0350000
    algebraic[21] = 1.00000/(1.00000+exp((states[1]+45.3000)/6.88410))
    algebraic[44] = 3.70000*exp(-(power(0.0333000*(states[1]+70.0000), 2.00000)))+0.0350000
    algebraic[22] = 1.00000/(1.00000+exp((states[1]+45.3000)/6.88410))
    algebraic[45] = 10.0000/(45.1600*exp(0.0357700*(states[1]+50.0000))+98.9000*exp(-0.100000*(states[1]+38.0000)))
    algebraic[23] = 1.00000/(1.00000+exp(-(states[1]+11.5000)/11.8200))
    algebraic[46] = 1.00000/(0.118850*exp((states[1]+80.0000)/28.3700)+0.562300*exp((states[1]+80.0000)/-14.1900))
    algebraic[25] = 1.00000/(1.00000+exp((states[1]+138.600)/10.4800))
    algebraic[51] = constants[28]*states[8]*(1.00000-states[10])-constants[29]*states[10]
    algebraic[26] = 1800.00*(states[2]+states[3])*(states[7]-states[6])
    algebraic[50] = (states[9]-states[7])/constants[26]
    algebraic[47] = power(states[8]/constants[17], constants[22])
    algebraic[48] = power(states[9]/constants[18], constants[23])
    algebraic[49] = (constants[21]*(constants[19]*algebraic[47]-constants[20]*algebraic[48]))/(1.00000+algebraic[47]+algebraic[48])
    algebraic[52] = constants[30]*states[8]*(1.00000-states[11])-constants[31]*states[11]
    algebraic[1] = 1000.00*((states[0]-states[1])/constants[10])
    algebraic[54] = custom_piecewise([equal(constants[39] , 1.00000) & less(voi-constants[40] , 0.00100000), 52.0000*(constants[116]+constants[117]) , True, 0.00000])
    algebraic[59] = constants[119]*(power(states[17], 3.00000))*states[18]*(states[0]-log((constants[3]+0.120000*constants[5])/(states[12]+0.120000*states[13]), 10)*((constants[0]*constants[1])/constants[2]))
    algebraic[61] = constants[120]*4.00000*((states[0]*(power(constants[2], 2.00000)))/(constants[0]*constants[1]))*((states[6]*exp((2.00000*states[0]*(power(constants[2], 2.00000)))/(constants[0]*constants[1]))-0.341000*constants[4])/(exp((2.00000*states[0]*(power(constants[2], 2.00000)))/(constants[0]*constants[1]))-1.00000))*states[19]*((0.900000+states[20]/10.0000)*states[21]+(0.100000+states[20]/10.0000)*states[22])
    algebraic[74] = (constants[130]*(exp(0.0187000*states[0])*(power(states[12], 3.00000))*constants[4]-exp(-0.0187000*states[0])*(power(constants[3], 3.00000))*states[8]))/(1.00000+0.000100000*(states[8]*(power(constants[3], 3.00000))+constants[4]*(power(states[12], 3.00000))))
    algebraic[72] = (((constants[128]/(1.00000+0.124500*exp((-0.100000*states[0]*constants[2])/(constants[0]*constants[1]))+0.0365000*constants[109]*exp((-states[0]*constants[2])/(constants[0]*constants[1]))))*constants[5])/(constants[5]+constants[69]))/(1.00000+power(constants[70]/states[12], 1.50000))
    algebraic[62] = ((constants[0]*constants[1])/constants[2])*log(constants[5]/states[13])
    algebraic[65] = ((constants[123]*(states[0]-(algebraic[62]+1.73000)))/(1.00000+exp(((1.61300*constants[2])/(constants[0]*constants[1]))*(states[0]-(algebraic[62]+1.73000)))))*(1.00000+exp(-(constants[5]-0.998800)/0.124000))
    algebraic[63] = constants[121]*states[23]*(constants[49]*states[24]+constants[50]*states[25])*(states[0]-algebraic[62])
    algebraic[64] = constants[122]*states[26]*states[27]*(states[0]-algebraic[62])
    algebraic[73] = (constants[129]*states[8])/(states[8]+0.000400000)
    algebraic[58] = ((constants[0]*constants[1])/constants[2])*log(constants[3]/states[12])
    algebraic[67] = constants[124]*states[28]*(0.200000*(states[0]-algebraic[58])+0.800000*(states[0]-algebraic[62]))
    algebraic[68] = constants[126]*(states[0]-algebraic[58])
    algebraic[70] = constants[127]*(states[0]-algebraic[62])
    algebraic[69] = constants[125]*(states[0]-constants[66])
    algebraic[71] = algebraic[68]+algebraic[69]+algebraic[70]
    algebraic[78] = constants[132]*4.00000*((states[1]*(power(constants[2], 2.00000)))/(constants[0]*constants[1]))*((states[6]*exp((2.00000*states[1]*(power(constants[2], 2.00000)))/(constants[0]*constants[1]))-0.341000*states[15])/(exp((2.00000*states[1]*(power(constants[2], 2.00000)))/(constants[0]*constants[1]))-1.00000))*states[31]*((0.900000+states[32]/10.0000)*states[33]+(0.100000+states[32]/10.0000)*states[34])
    algebraic[66] = constants[124]*states[28]*constants[108]*(states[0]-algebraic[62])
    algebraic[89] = (exp(states[14]/67.3000)-1.00000)/7.00000
    algebraic[90] = (((constants[140]/(1.00000+0.124500*exp((-0.100000*states[1]*constants[2])/(constants[0]*constants[1]))+0.0365000*algebraic[89]*exp((-states[1]*constants[2])/(constants[0]*constants[1]))))*states[16])/(states[16]+constants[100]))/(1.00000+power(constants[101]/states[12], 1.50000))
    algebraic[79] = ((constants[0]*constants[1])/constants[2])*log(states[16]/states[13])
    algebraic[83] = constants[136]*states[40]*constants[111]*(states[1]-algebraic[79])
    algebraic[87] = constants[139]*(states[1]-algebraic[79])
    algebraic[82] = ((constants[135]*(states[1]-(algebraic[79]+1.73000)))/(1.00000+exp(((1.61300*constants[2])/(constants[0]*constants[1]))*(states[1]-(algebraic[79]+1.73000)))))*(1.00000+exp(-(states[16]-0.998800)/0.124000))
    algebraic[80] = constants[133]*states[35]*(constants[80]*states[36]+constants[81]*states[37])*(states[1]-algebraic[79])
    algebraic[81] = constants[134]*states[38]*states[39]*(states[1]-algebraic[79])
    algebraic[57] = (1.00000*constants[145])/(constants[43]*(states[16]-constants[5]))
    algebraic[76] = constants[131]*(power(states[29], 3.00000))*states[30]*(states[1]-log((states[14]+0.120000*states[16])/(states[12]+0.120000*states[13]), 10)*((constants[0]*constants[1])/constants[2]))
    algebraic[75] = ((constants[0]*constants[1])/constants[2])*log(states[14]/states[12])
    algebraic[85] = constants[138]*(states[1]-algebraic[75])
    algebraic[92] = (constants[142]*(exp(0.0187000*states[1])*(power(states[12], 3.00000))*states[15]-exp(-0.0187000*states[1])*(power(states[14], 3.00000))*states[8]))/(1.00000+0.000100000*(states[8]*(power(states[14], 3.00000))+states[15]*(power(states[12], 3.00000))))
    algebraic[55] = (1.00000*constants[145])/(constants[41]*(states[14]-constants[3]))
    algebraic[91] = (constants[141]*states[8])/(states[8]+0.000400000)
    algebraic[86] = constants[137]*(states[1]-constants[97])
    algebraic[56] = (1.00000*constants[145])/(constants[42]*(states[15]-constants[4]))
    algebraic[84] = constants[136]*states[40]*(0.200000*(states[1]-algebraic[75])+0.800000*(states[1]-algebraic[79]))
    algebraic[88] = algebraic[85]+algebraic[86]+algebraic[87]
    algebraic[0] = (states[6]-states[8])/constants[27]
    algebraic[53] = algebraic[51]+algebraic[52]
    algebraic[60] = constants[124]*states[28]*constants[59]*(states[0]-algebraic[58])
    algebraic[77] = constants[136]*states[40]*constants[90]*(states[1]-algebraic[75])
    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)