Generated Code

The following is python code generated by the CellML API from this CellML file. (Back to language selection)

The raw code is available.

# Size of variable arrays:
sizeAlgebraic = 90
sizeStates = 31
sizeConstants = 133
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] = "isotonic in component parameters (dimensionless)"
    legend_constants[1] = "alpha_1 in component parameters (per_micrometre)"
    legend_constants[2] = "beta_1 in component parameters (millinewton)"
    legend_constants[3] = "alpha_2 in component parameters (per_micrometre)"
    legend_constants[4] = "beta_2 in component parameters (millinewton)"
    legend_constants[5] = "alpha_3 in component parameters (per_micrometre)"
    legend_constants[6] = "beta_3 in component parameters (millinewton)"
    legend_constants[7] = "lambda in component parameters (millinewton)"
    legend_constants[8] = "mu in component parameters (dimensionless)"
    legend_constants[9] = "k_mu in component parameters (dimensionless)"
    legend_constants[10] = "kappa in component parameters (dimensionless)"
    legend_constants[11] = "kappa_0 in component parameters (dimensionless)"
    legend_constants[12] = "m_0 in component parameters (dimensionless)"
    legend_constants[13] = "v_max in component parameters (micrometre_per_second)"
    legend_constants[14] = "a in component parameters (dimensionless)"
    legend_constants[15] = "d_h in component parameters (dimensionless)"
    legend_constants[16] = "alpha_P in component parameters (dimensionless)"
    legend_algebraic[37] = "l in component length (micrometre)"
    legend_algebraic[26] = "F_muscle in component force (millinewton)"
    legend_constants[17] = "F_afterload in component isotonic (millinewton)"
    legend_algebraic[40] = "isotonic_mode in component isotonic (dimensionless)"
    legend_constants[18] = "l_0 in component length (micrometre)"
    legend_constants[19] = "stim_period in component membrane (second)"
    legend_constants[20] = "S_0 in component parameters_izakov_et_al_1991 (micrometre)"
    legend_algebraic[2] = "q_v in component parameters_izakov_et_al_1991 (per_second)"
    legend_constants[21] = "q_1 in component parameters_izakov_et_al_1991 (per_second)"
    legend_constants[22] = "q_2 in component parameters_izakov_et_al_1991 (per_second)"
    legend_constants[23] = "q_3 in component parameters_izakov_et_al_1991 (per_second)"
    legend_constants[24] = "q_4 in component parameters_izakov_et_al_1991 (per_second)"
    legend_constants[25] = "x_st in component parameters_izakov_et_al_1991 (dimensionless)"
    legend_constants[123] = "v_st in component parameters_izakov_et_al_1991 (micrometre_per_second)"
    legend_constants[124] = "v_1 in component parameters_izakov_et_al_1991 (micrometre_per_second)"
    legend_constants[26] = "alpha_G in component parameters_izakov_et_al_1991 (dimensionless)"
    legend_constants[27] = "k_A in component parameters_izakov_et_al_1991 (per_millimolar)"
    legend_states[0] = "v in component CE_velocity (micrometre_per_second)"
    legend_constants[28] = "alpha_Q in component parameters_izakov_et_al_1991 (dimensionless)"
    legend_constants[29] = "beta_Q in component parameters_izakov_et_al_1991 (dimensionless)"
    legend_algebraic[71] = "F_CE in component force (millinewton)"
    legend_algebraic[15] = "F_XSE in component force (millinewton)"
    legend_algebraic[0] = "F_SE in component force (millinewton)"
    legend_algebraic[1] = "F_PE in component force (millinewton)"
    legend_states[1] = "N in component crossbridge_kinetics (dimensionless)"
    legend_algebraic[44] = "k_P_vis in component CE_velocity (millinewton_second_per_micrometre)"
    legend_algebraic[51] = "k_S_vis in component PE_velocity (millinewton_second_per_micrometre)"
    legend_states[2] = "w in component PE_velocity (micrometre_per_second)"
    legend_states[3] = "l_1 in component length (micrometre)"
    legend_states[4] = "l_2 in component length (micrometre)"
    legend_states[5] = "l_3 in component length (micrometre)"
    legend_algebraic[69] = "p_v in component average_crossbridge_force (dimensionless)"
    legend_algebraic[63] = "K_kappa in component crossbridge_kinetics (per_second)"
    legend_algebraic[31] = "M_A in component crossbridge_kinetics (dimensionless)"
    legend_algebraic[33] = "n_1 in component crossbridge_kinetics (dimensionless)"
    legend_algebraic[35] = "L_oz in component crossbridge_kinetics (dimensionless)"
    legend_algebraic[59] = "k_p_v in component crossbridge_kinetics (per_second)"
    legend_algebraic[61] = "k_m_v in component crossbridge_kinetics (per_second)"
    legend_states[6] = "A in component intracellular_calcium_concentration (millimolar)"
    legend_constants[30] = "A_tot in component intracellular_calcium_concentration (millimolar)"
    legend_algebraic[57] = "G_star in component average_crossbridge_force (dimensionless)"
    legend_algebraic[53] = "P_star in component average_crossbridge_force (dimensionless)"
    legend_constants[31] = "g_1 in component crossbridge_kinetics (per_micrometre)"
    legend_constants[32] = "g_2 in component crossbridge_kinetics (dimensionless)"
    legend_algebraic[3] = "dl_1_dt in component length (micrometre_per_second)"
    legend_algebraic[54] = "dl_2_dt in component length (micrometre_per_second)"
    legend_algebraic[55] = "dl_3_dt in component length (micrometre_per_second)"
    legend_algebraic[47] = "phi_chi_2 in component CE_velocity (micrometre_per_second)"
    legend_constants[33] = "stim_start in component membrane (second)"
    legend_constants[34] = "stim_end in component membrane (second)"
    legend_constants[35] = "stim_duration in component membrane (second)"
    legend_constants[36] = "stim_amplitude in component membrane (nanoA)"
    legend_algebraic[74] = "phi_chi in component CE_velocity (micrometre_per_second2)"
    legend_algebraic[72] = "p_prime_v in component average_crossbridge_force (second_per_micrometre)"
    legend_constants[37] = "alpha_P_lengthening in component CE_velocity (per_micrometre)"
    legend_constants[38] = "beta_P_lengthening in component CE_velocity (millinewton_second_per_micrometre)"
    legend_constants[39] = "alpha_P_shortening in component CE_velocity (per_micrometre)"
    legend_constants[40] = "beta_P_shortening in component CE_velocity (millinewton_second_per_micrometre)"
    legend_algebraic[42] = "alp_p in component CE_velocity (per_micrometre)"
    legend_constants[41] = "alpha_S_lengthening in component PE_velocity (per_micrometre)"
    legend_constants[42] = "beta_S_lengthening in component PE_velocity (millinewton_second_per_micrometre)"
    legend_constants[43] = "alpha_S_shortening in component PE_velocity (per_micrometre)"
    legend_constants[44] = "beta_S_shortening in component PE_velocity (millinewton_second_per_micrometre)"
    legend_algebraic[49] = "alp_s in component PE_velocity (per_micrometre)"
    legend_constants[130] = "gamma in component average_crossbridge_force (dimensionless)"
    legend_constants[125] = "case_1 in component average_crossbridge_force (second_per_micrometre)"
    legend_algebraic[65] = "case_2 in component average_crossbridge_force (second_per_micrometre)"
    legend_constants[126] = "case_3 in component average_crossbridge_force (second_per_micrometre)"
    legend_algebraic[67] = "case_4 in component average_crossbridge_force (second_per_micrometre)"
    legend_states[7] = "V in component membrane (millivolt)"
    legend_constants[45] = "R in component membrane (joule_per_kilomole_kelvin)"
    legend_constants[46] = "T in component membrane (kelvin)"
    legend_constants[47] = "F in component membrane (coulomb_per_mole)"
    legend_constants[48] = "Cm in component membrane (microF)"
    legend_algebraic[52] = "i_K1 in component time_independent_potassium_current (nanoA)"
    legend_algebraic[79] = "i_to in component transient_outward_current (nanoA)"
    legend_algebraic[56] = "i_Kr in component rapid_delayed_rectifier_potassium_current (nanoA)"
    legend_algebraic[58] = "i_Ks in component slow_delayed_rectifier_potassium_current (nanoA)"
    legend_algebraic[68] = "i_Ca_L_K_cyt in component L_type_Ca_channel (nanoA)"
    legend_algebraic[75] = "i_Ca_L_K_ds in component L_type_Ca_channel (nanoA)"
    legend_algebraic[80] = "i_NaK in component sodium_potassium_pump (nanoA)"
    legend_algebraic[60] = "i_Na in component fast_sodium_current (nanoA)"
    legend_algebraic[64] = "i_b_Na in component sodium_background_current (nanoA)"
    legend_algebraic[62] = "i_p_Na in component persistent_sodium_current (nanoA)"
    legend_algebraic[70] = "i_Ca_L_Na_cyt in component L_type_Ca_channel (nanoA)"
    legend_algebraic[76] = "i_Ca_L_Na_ds in component L_type_Ca_channel (nanoA)"
    legend_algebraic[81] = "i_NaCa_cyt in component sodium_calcium_exchanger (nanoA)"
    legend_algebraic[82] = "i_NaCa_ds in component sodium_calcium_exchanger (nanoA)"
    legend_algebraic[66] = "i_Ca_L_Ca_cyt in component L_type_Ca_channel (nanoA)"
    legend_algebraic[73] = "i_Ca_L_Ca_ds in component L_type_Ca_channel (nanoA)"
    legend_algebraic[78] = "i_b_Ca in component calcium_background_current (nanoA)"
    legend_algebraic[39] = "i_Stim in component membrane (nanoA)"
    legend_algebraic[41] = "E_Na in component reversal_potentials (millivolt)"
    legend_algebraic[43] = "E_K in component reversal_potentials (millivolt)"
    legend_algebraic[46] = "E_Ks in component reversal_potentials (millivolt)"
    legend_algebraic[48] = "E_Ca in component reversal_potentials (millivolt)"
    legend_algebraic[50] = "E_mh in component reversal_potentials (millivolt)"
    legend_constants[49] = "P_kna in component reversal_potentials (dimensionless)"
    legend_states[8] = "K_o in component extracellular_potassium_concentration (millimolar)"
    legend_constants[50] = "Na_o in component extracellular_sodium_concentration (millimolar)"
    legend_states[9] = "K_i in component intracellular_potassium_concentration (millimolar)"
    legend_states[10] = "Na_i in component intracellular_sodium_concentration (millimolar)"
    legend_constants[51] = "Ca_o in component extracellular_calcium_concentration (millimolar)"
    legend_states[11] = "Ca_i in component intracellular_calcium_concentration (millimolar)"
    legend_constants[52] = "K_mk1 in component time_independent_potassium_current (millimolar)"
    legend_constants[53] = "g_K1 in component time_independent_potassium_current (microS)"
    legend_constants[54] = "g_Kr1 in component rapid_delayed_rectifier_potassium_current (microS)"
    legend_constants[55] = "g_Kr2 in component rapid_delayed_rectifier_potassium_current (microS)"
    legend_states[12] = "xr1 in component rapid_delayed_rectifier_potassium_current_xr1_gate (dimensionless)"
    legend_states[13] = "xr2 in component rapid_delayed_rectifier_potassium_current_xr2_gate (dimensionless)"
    legend_algebraic[4] = "alpha_xr1 in component rapid_delayed_rectifier_potassium_current_xr1_gate (per_second)"
    legend_algebraic[16] = "beta_xr1 in component rapid_delayed_rectifier_potassium_current_xr1_gate (per_second)"
    legend_algebraic[5] = "alpha_xr2 in component rapid_delayed_rectifier_potassium_current_xr2_gate (per_second)"
    legend_algebraic[17] = "beta_xr2 in component rapid_delayed_rectifier_potassium_current_xr2_gate (per_second)"
    legend_constants[56] = "g_Ks in component slow_delayed_rectifier_potassium_current (microS)"
    legend_states[14] = "xs in component slow_delayed_rectifier_potassium_current_xs_gate (dimensionless)"
    legend_algebraic[6] = "alpha_xs in component slow_delayed_rectifier_potassium_current_xs_gate (per_second)"
    legend_algebraic[18] = "beta_xs in component slow_delayed_rectifier_potassium_current_xs_gate (per_second)"
    legend_algebraic[45] = "i_KNa in component sodium_activated_potassium_current (nanoA)"
    legend_constants[57] = "g_K_Na in component sodium_activated_potassium_current (microS)"
    legend_constants[58] = "K_kna in component sodium_activated_potassium_current (millimolar)"
    legend_constants[59] = "g_Na in component fast_sodium_current (microS)"
    legend_states[15] = "m in component fast_sodium_current_m_gate (dimensionless)"
    legend_states[16] = "h in component fast_sodium_current_h_gate (dimensionless)"
    legend_algebraic[19] = "alpha_m in component fast_sodium_current_m_gate (per_second)"
    legend_algebraic[27] = "beta_m in component fast_sodium_current_m_gate (per_second)"
    legend_constants[60] = "delta_m in component fast_sodium_current_m_gate (millivolt)"
    legend_algebraic[7] = "E0_m in component fast_sodium_current_m_gate (millivolt)"
    legend_algebraic[8] = "alpha_h in component fast_sodium_current_h_gate (per_second)"
    legend_algebraic[20] = "beta_h in component fast_sodium_current_h_gate (per_second)"
    legend_constants[61] = "g_pna in component persistent_sodium_current (microS)"
    legend_constants[62] = "g_bna in component sodium_background_current (microS)"
    legend_algebraic[77] = "i_Ca_L in component L_type_Ca_channel (nanoA)"
    legend_constants[63] = "P_Ca_L in component L_type_Ca_channel (nanoA_per_millimolar)"
    legend_constants[64] = "P_CaK in component L_type_Ca_channel (dimensionless)"
    legend_constants[65] = "P_CaNa in component L_type_Ca_channel (dimensionless)"
    legend_states[17] = "Ca_ds in component intracellular_calcium_concentration (millimolar)"
    legend_states[18] = "d in component L_type_Ca_channel_d_gate (dimensionless)"
    legend_states[19] = "f in component L_type_Ca_channel_f_gate (dimensionless)"
    legend_states[20] = "f2 in component L_type_Ca_channel_f2_gate (dimensionless)"
    legend_states[21] = "f2ds in component L_type_Ca_channel_f2ds_gate (dimensionless)"
    legend_constants[66] = "Km_f2 in component L_type_Ca_channel (millimolar)"
    legend_constants[67] = "Km_f2ds in component L_type_Ca_channel (millimolar)"
    legend_constants[68] = "R_decay in component L_type_Ca_channel (per_second)"
    legend_constants[69] = "FrICa in component L_type_Ca_channel (dimensionless)"
    legend_algebraic[21] = "alpha_d in component L_type_Ca_channel_d_gate (per_second)"
    legend_algebraic[28] = "beta_d in component L_type_Ca_channel_d_gate (per_second)"
    legend_algebraic[9] = "E0_d in component L_type_Ca_channel_d_gate (millivolt)"
    legend_constants[70] = "speed_d in component L_type_Ca_channel_d_gate (dimensionless)"
    legend_algebraic[22] = "alpha_f in component L_type_Ca_channel_f_gate (per_second)"
    legend_algebraic[29] = "beta_f in component L_type_Ca_channel_f_gate (per_second)"
    legend_constants[71] = "speed_f in component L_type_Ca_channel_f_gate (dimensionless)"
    legend_constants[72] = "delta_f in component L_type_Ca_channel_f_gate (millivolt)"
    legend_algebraic[10] = "E0_f in component L_type_Ca_channel_f_gate (millivolt)"
    legend_constants[73] = "g_bca in component calcium_background_current (microS)"
    legend_constants[74] = "g_to in component transient_outward_current (microS)"
    legend_constants[75] = "g_tos in component transient_outward_current (dimensionless)"
    legend_states[22] = "s in component transient_outward_current_s_gate (dimensionless)"
    legend_states[23] = "r in component transient_outward_current_r_gate (dimensionless)"
    legend_algebraic[11] = "alpha_s in component transient_outward_current_s_gate (per_second)"
    legend_algebraic[23] = "beta_s in component transient_outward_current_s_gate (per_second)"
    legend_constants[76] = "i_NaK_max in component sodium_potassium_pump (nanoA)"
    legend_constants[77] = "K_mK in component sodium_potassium_pump (millimolar)"
    legend_constants[78] = "K_mNa in component sodium_potassium_pump (millimolar)"
    legend_algebraic[84] = "i_NaCa in component sodium_calcium_exchanger (nanoA)"
    legend_constants[79] = "k_NaCa in component sodium_calcium_exchanger (nanoA)"
    legend_constants[80] = "n_NaCa in component sodium_calcium_exchanger (dimensionless)"
    legend_constants[81] = "d_NaCa in component sodium_calcium_exchanger (dimensionless)"
    legend_constants[82] = "gamma in component sodium_calcium_exchanger (dimensionless)"
    legend_constants[83] = "FRiNaCa in component sodium_calcium_exchanger (dimensionless)"
    legend_algebraic[85] = "i_up in component sarcoplasmic_reticulum_calcium_pump (millimolar_per_second)"
    legend_constants[128] = "K_1 in component sarcoplasmic_reticulum_calcium_pump (dimensionless)"
    legend_algebraic[83] = "K_2 in component sarcoplasmic_reticulum_calcium_pump (millimolar)"
    legend_constants[84] = "K_cyca in component sarcoplasmic_reticulum_calcium_pump (millimolar)"
    legend_constants[85] = "K_xcs in component sarcoplasmic_reticulum_calcium_pump (dimensionless)"
    legend_constants[86] = "K_srca in component sarcoplasmic_reticulum_calcium_pump (millimolar)"
    legend_constants[87] = "alpha_up in component sarcoplasmic_reticulum_calcium_pump (millimolar_per_second)"
    legend_constants[88] = "beta_up in component sarcoplasmic_reticulum_calcium_pump (millimolar_per_second)"
    legend_states[24] = "Ca_up in component intracellular_calcium_concentration (millimolar)"
    legend_constants[89] = "flag_ingib in component sarcoplasmic_reticulum_calcium_pump (dimensionless)"
    legend_constants[90] = "K_inh in component sarcoplasmic_reticulum_calcium_pump (millimolar)"
    legend_algebraic[86] = "i_trans in component calcium_translocation (millimolar_per_second)"
    legend_states[25] = "Ca_rel in component intracellular_calcium_concentration (millimolar)"
    legend_constants[91] = "a_tr in component calcium_translocation (per_second)"
    legend_constants[92] = "alpha_CaS in component calcium_translocation (per_millimolar_second)"
    legend_constants[93] = "beta_CaS in component calcium_translocation (per_second)"
    legend_constants[94] = "CaS_tot in component calcium_translocation (millimolar)"
    legend_constants[129] = "beta in component calcium_translocation (millimolar)"
    legend_algebraic[87] = "i_rel in component calcium_release (millimolar_per_second)"
    legend_algebraic[12] = "VoltDep in component calcium_release (dimensionless)"
    legend_algebraic[30] = "RegBindSite in component calcium_release (dimensionless)"
    legend_algebraic[13] = "CaiReg in component calcium_release (dimensionless)"
    legend_algebraic[24] = "CadsReg in component calcium_release (dimensionless)"
    legend_algebraic[32] = "ActRate in component calcium_release (per_second)"
    legend_algebraic[34] = "InactRate in component calcium_release (per_second)"
    legend_constants[95] = "SRLeak in component calcium_release (per_second)"
    legend_constants[96] = "K_m_rel in component calcium_release (per_second)"
    legend_constants[97] = "K_m_Ca_cyt in component calcium_release (millimolar)"
    legend_constants[98] = "K_m_Ca_ds in component calcium_release (millimolar)"
    legend_algebraic[38] = "PrecFrac in component calcium_release (dimensionless)"
    legend_states[26] = "ActFrac in component calcium_release (dimensionless)"
    legend_states[27] = "ProdFrac in component calcium_release (dimensionless)"
    legend_algebraic[36] = "SpeedRel in component calcium_release (dimensionless)"
    legend_constants[131] = "V_i in component intracellular_calcium_concentration (micrometre3)"
    legend_constants[99] = "n_NaK in component intracellular_sodium_concentration (dimensionless)"
    legend_constants[100] = "K_b in component extracellular_potassium_concentration (millimolar)"
    legend_constants[132] = "V_e in component intracellular_calcium_concentration (micrometre3)"
    legend_constants[101] = "radius in component intracellular_calcium_concentration (micrometre)"
    legend_constants[102] = "length in component intracellular_calcium_concentration (micrometre)"
    legend_constants[127] = "V_Cell in component intracellular_calcium_concentration (micrometre3)"
    legend_constants[103] = "V_i_ratio in component intracellular_calcium_concentration (dimensionless)"
    legend_constants[104] = "V_ds_ratio in component intracellular_calcium_concentration (dimensionless)"
    legend_constants[105] = "V_rel_ratio in component intracellular_calcium_concentration (dimensionless)"
    legend_constants[106] = "V_e_ratio in component intracellular_calcium_concentration (dimensionless)"
    legend_constants[107] = "V_up_ratio in component intracellular_calcium_concentration (dimensionless)"
    legend_constants[108] = "Kdecay in component intracellular_calcium_concentration (per_second)"
    legend_algebraic[88] = "N_A in component intracellular_calcium_concentration (dimensionless)"
    legend_algebraic[89] = "pi_N_A in component intracellular_calcium_concentration (dimensionless)"
    legend_constants[109] = "pi_min in component intracellular_calcium_concentration (dimensionless)"
    legend_constants[110] = "a_on in component intracellular_calcium_concentration (per_millimolar_second)"
    legend_constants[111] = "a_off in component intracellular_calcium_concentration (per_second)"
    legend_states[28] = "B_1 in component intracellular_calcium_concentration (millimolar)"
    legend_constants[112] = "B_1_tot in component intracellular_calcium_concentration (millimolar)"
    legend_constants[113] = "b_1_on in component intracellular_calcium_concentration (per_millimolar_second)"
    legend_constants[114] = "b_1_off in component intracellular_calcium_concentration (per_second)"
    legend_states[29] = "B_2 in component intracellular_calcium_concentration (millimolar)"
    legend_constants[115] = "B_2_tot in component intracellular_calcium_concentration (millimolar)"
    legend_constants[116] = "b_2_on in component intracellular_calcium_concentration (per_millimolar_second)"
    legend_constants[117] = "b_2_off in component intracellular_calcium_concentration (per_second)"
    legend_constants[118] = "g_fibro_junct in component fibroblast (microS)"
    legend_constants[119] = "g_fibro in component fibroblast (microS)"
    legend_constants[120] = "c_fibro in component fibroblast (microF)"
    legend_constants[121] = "g_fibro_stretch in component fibroblast (microS)"
    legend_constants[122] = "E_fibro_stretch in component fibroblast (millivolt)"
    legend_states[30] = "V_fibro in component fibroblast (millivolt)"
    legend_algebraic[14] = "i_fibro in component fibroblast (nanoA)"
    legend_algebraic[25] = "i_fibro_junct in component fibroblast (nanoA)"
    legend_rates[1] = "d/dt N in component crossbridge_kinetics (dimensionless)"
    legend_rates[3] = "d/dt l_1 in component length (micrometre)"
    legend_rates[4] = "d/dt l_2 in component length (micrometre)"
    legend_rates[5] = "d/dt l_3 in component length (micrometre)"
    legend_rates[0] = "d/dt v in component CE_velocity (micrometre_per_second)"
    legend_rates[2] = "d/dt w in component PE_velocity (micrometre_per_second)"
    legend_rates[7] = "d/dt V in component membrane (millivolt)"
    legend_rates[12] = "d/dt xr1 in component rapid_delayed_rectifier_potassium_current_xr1_gate (dimensionless)"
    legend_rates[13] = "d/dt xr2 in component rapid_delayed_rectifier_potassium_current_xr2_gate (dimensionless)"
    legend_rates[14] = "d/dt xs in component slow_delayed_rectifier_potassium_current_xs_gate (dimensionless)"
    legend_rates[15] = "d/dt m in component fast_sodium_current_m_gate (dimensionless)"
    legend_rates[16] = "d/dt h in component fast_sodium_current_h_gate (dimensionless)"
    legend_rates[18] = "d/dt d in component L_type_Ca_channel_d_gate (dimensionless)"
    legend_rates[19] = "d/dt f in component L_type_Ca_channel_f_gate (dimensionless)"
    legend_rates[20] = "d/dt f2 in component L_type_Ca_channel_f2_gate (dimensionless)"
    legend_rates[21] = "d/dt f2ds in component L_type_Ca_channel_f2ds_gate (dimensionless)"
    legend_rates[22] = "d/dt s in component transient_outward_current_s_gate (dimensionless)"
    legend_rates[23] = "d/dt r in component transient_outward_current_r_gate (dimensionless)"
    legend_rates[26] = "d/dt ActFrac in component calcium_release (dimensionless)"
    legend_rates[27] = "d/dt ProdFrac in component calcium_release (dimensionless)"
    legend_rates[10] = "d/dt Na_i in component intracellular_sodium_concentration (millimolar)"
    legend_rates[8] = "d/dt K_o in component extracellular_potassium_concentration (millimolar)"
    legend_rates[9] = "d/dt K_i in component intracellular_potassium_concentration (millimolar)"
    legend_rates[6] = "d/dt A in component intracellular_calcium_concentration (millimolar)"
    legend_rates[28] = "d/dt B_1 in component intracellular_calcium_concentration (millimolar)"
    legend_rates[29] = "d/dt B_2 in component intracellular_calcium_concentration (millimolar)"
    legend_rates[11] = "d/dt Ca_i in component intracellular_calcium_concentration (millimolar)"
    legend_rates[17] = "d/dt Ca_ds in component intracellular_calcium_concentration (millimolar)"
    legend_rates[24] = "d/dt Ca_up in component intracellular_calcium_concentration (millimolar)"
    legend_rates[25] = "d/dt Ca_rel in component intracellular_calcium_concentration (millimolar)"
    legend_rates[30] = "d/dt V_fibro in component fibroblast (millivolt)"
    return (legend_states, legend_algebraic, legend_voi, legend_constants)

def initConsts():
    constants = [0.0] * sizeConstants; states = [0.0] * sizeStates;
    constants[0] = 0
    constants[1] = 14.6
    constants[2] = 0.84
    constants[3] = 14.6
    constants[4] = 0.0018
    constants[5] = 48
    constants[6] = 0.015
    constants[7] = 30
    constants[8] = 3
    constants[9] = 0.6
    constants[10] = 0.705
    constants[11] = 3
    constants[12] = 0.9
    constants[13] = 5.5
    constants[14] = 0.25
    constants[15] = 0.5
    constants[16] = 4
    constants[17] = 2
    constants[18] = 0.525139356105856
    constants[19] = 1
    constants[20] = 1.14
    constants[21] = 17.3
    constants[22] = 259
    constants[23] = 17.3
    constants[24] = 15
    constants[25] = 0.964285
    constants[26] = 1
    constants[27] = 40
    states[0] = 0
    constants[28] = 10
    constants[29] = 5
    states[1] = 7.455e-8
    states[2] = 0
    states[3] = 0.436333342969918
    states[4] = 0.436333525334166
    states[5] = 0.088805830771694
    states[6] = 0.00015
    constants[30] = 0.07
    constants[31] = 0.6
    constants[32] = 0.52
    constants[33] = 0.06
    constants[34] = 10000
    constants[35] = 0.0025
    constants[36] = -3
    constants[37] = 16
    constants[38] = 0.0015
    constants[39] = 16
    constants[40] = 0.0015
    constants[41] = 46
    constants[42] = 0
    constants[43] = 39
    constants[44] = 0
    states[7] = -93.658148
    constants[45] = 8314.472
    constants[46] = 310
    constants[47] = 96485.3415
    constants[48] = 9.5e-5
    constants[49] = 0.03
    states[8] = 3.988
    constants[50] = 140
    states[9] = 139.054
    states[10] = 5.18787513289509
    constants[51] = 2
    states[11] = 6.15e-6
    constants[52] = 10
    constants[53] = 0.5
    constants[54] = 0.0021
    constants[55] = 0.0013
    states[12] = 8.88859784542779e-6
    states[13] = 1.53745791069154e-7
    constants[56] = 0.0026
    states[14] = 0.001
    constants[57] = 0
    constants[58] = 20
    constants[59] = 2.5
    states[15] = 0.0015
    states[16] = 0.995
    constants[60] = 1e-5
    constants[61] = 0.004
    constants[62] = 0.0006
    constants[63] = 0.1
    constants[64] = 0.002
    constants[65] = 0.01
    states[17] = 2.55e-6
    states[18] = 0
    states[19] = 1
    states[20] = 1
    states[21] = 0.997
    constants[66] = 100000
    constants[67] = 0.001
    constants[68] = 20
    constants[69] = 1
    constants[70] = 3
    constants[71] = 0.3
    constants[72] = 0.0001
    constants[73] = 0.00025
    constants[74] = 0.006
    constants[75] = 0
    states[22] = 0.997044616031121
    states[23] = 1.63117173173398e-8
    constants[76] = 0.7
    constants[77] = 1
    constants[78] = 24.2
    constants[79] = 0.0005
    constants[80] = 3
    constants[81] = 0
    constants[82] = 0.5
    constants[83] = 0.001
    constants[84] = 0.00015
    constants[85] = 0.4
    constants[86] = 0.5
    constants[87] = 1
    constants[88] = 0.03
    states[24] = 0.994579
    constants[89] = 0
    constants[90] = 4
    states[25] = 0.989665
    constants[91] = 15
    constants[92] = 50000
    constants[93] = 32500
    constants[94] = 40
    constants[95] = 0.05
    constants[96] = 10000
    constants[97] = 0.0005
    constants[98] = 0.01
    states[26] = 0.001914
    states[27] = 0.2854569
    constants[99] = 1.5
    constants[100] = 4
    constants[101] = 12
    constants[102] = 74
    constants[103] = 0.49
    constants[104] = 0.1
    constants[105] = 0.003
    constants[106] = 0.4
    constants[107] = 0.03
    constants[108] = 10
    constants[109] = 0.03
    constants[110] = 70000
    constants[111] = 200
    states[28] = 0
    constants[112] = 0.08
    constants[113] = 100000
    constants[114] = 182
    states[29] = 0
    constants[115] = 0.1
    constants[116] = 1000
    constants[117] = 3
    constants[118] = 2.9e-4
    constants[119] = 2e-4
    constants[120] = 1e-5
    constants[121] = 0
    constants[122] = 0
    states[30] = -20
    constants[123] = constants[25]*constants[13]
    constants[124] = constants[13]/10.0000
    constants[125] = (constants[14]*(0.400000+0.400000*constants[14]))/(constants[13]*(power((constants[14]+1.00000)*0.400000, 2.00000)))
    constants[126] = (0.400000*constants[14]+1.00000)/(constants[14]*constants[13])
    constants[127] = (3.14159*(power(constants[101]/1000.00, 2.00000))*constants[102])/1000.00
    constants[128] = (constants[84]*constants[85])/constants[86]
    constants[129] = constants[93]/constants[92]
    constants[130] = (constants[14]*constants[15]*(power(constants[124]/constants[13], 2.00000)))/(3.00000*constants[14]*constants[15]-((constants[14]+1.00000)*constants[124])/constants[13])
    constants[131] = constants[127]*constants[103]
    constants[132] = constants[127]*constants[106]
    return (states, constants)

def computeRates(voi, states, constants):
    rates = [0.0] * sizeStates; algebraic = [0.0] * sizeAlgebraic
    rates[20] = 1.00000-1.00000*(states[11]/(constants[66]+states[11])+states[20])
    rates[21] = constants[68]*(1.00000-(states[17]/(constants[67]+states[17])+states[21]))
    rates[23] = 333.000*(1.00000/(1.00000+exp(-(states[7]+4.00000)/5.00000))-states[23])
    algebraic[3] = states[0]
    rates[3] = algebraic[3]
    algebraic[4] = 50.0000/(1.00000+exp(-(states[7]-5.00000)/9.00000))
    algebraic[16] = 0.0500000*exp(-(states[7]-20.0000)/15.0000)
    rates[12] = algebraic[4]*(1.00000-states[12])-algebraic[16]*states[12]
    algebraic[5] = 50.0000/(1.00000+exp(-(states[7]-5.00000)/9.00000))
    algebraic[17] = 0.400000*exp(-(power((states[7]+30.0000)/30.0000, 3.00000)))
    rates[13] = algebraic[5]*(1.00000-states[13])-algebraic[17]*states[13]
    algebraic[6] = 14.0000/(1.00000+exp(-(states[7]-40.0000)/9.00000))
    algebraic[18] = 1.00000*exp(-states[7]/45.0000)
    rates[14] = algebraic[6]*(1.00000-states[14])-algebraic[18]*states[14]
    algebraic[8] = 20.0000*exp(-0.125000*(states[7]+75.0000))
    algebraic[20] = 2000.00/(1.00000+320.000*exp(-0.100000*(states[7]+75.0000)))
    rates[16] = algebraic[8]*(1.00000-states[16])-algebraic[20]*states[16]
    algebraic[11] = 0.0330000*exp(-states[7]/17.0000)
    algebraic[23] = 33.0000/(1.00000+exp(-0.125000*(states[7]+10.0000)))
    rates[22] = algebraic[11]*(1.00000-states[22])-algebraic[23]*states[22]
    algebraic[14] = constants[119]*(states[30]+20.0000)+constants[121]*(states[30]-constants[122])
    algebraic[25] = -constants[118]*(states[7]-states[30])
    rates[30] = -(algebraic[14]+algebraic[25])/constants[120]
    algebraic[7] = states[7]+41.0000
    algebraic[19] = custom_piecewise([less(fabs(algebraic[7]) , constants[60]), 2000.00 , True, (200.000*algebraic[7])/(1.00000-exp(-0.100000*algebraic[7]))])
    algebraic[27] = 8000.00*exp(-0.0560000*(states[7]+66.0000))
    rates[15] = algebraic[19]*(1.00000-states[15])-algebraic[27]*states[15]
    algebraic[9] = (states[7]+24.0000)-5.00000
    algebraic[21] = custom_piecewise([less(fabs(algebraic[9]) , 0.000100000), 120.000 , True, (30.0000*algebraic[9])/(1.00000-exp(-algebraic[9]/4.00000))])
    algebraic[28] = custom_piecewise([less(fabs(algebraic[9]) , 0.000100000), 120.000 , True, (12.0000*algebraic[9])/(exp(algebraic[9]/10.0000)-1.00000)])
    rates[18] = constants[70]*(algebraic[21]*(1.00000-states[18])-algebraic[28]*states[18])
    algebraic[10] = states[7]+34.0000
    algebraic[22] = custom_piecewise([less(fabs(algebraic[10]) , constants[72]), 25.0000 , True, (6.25000*algebraic[10])/(exp(algebraic[10]/4.00000)-1.00000)])
    algebraic[29] = 12.0000/(1.00000+exp((-1.00000*(states[7]+34.0000))/4.00000))
    rates[19] = constants[71]*(algebraic[22]*(1.00000-states[19])-algebraic[29]*states[19])
    algebraic[13] = states[11]/(states[11]+constants[97])
    algebraic[24] = states[17]/(states[17]+constants[98])
    algebraic[30] = algebraic[13]+(1.00000-algebraic[13])*algebraic[24]
    algebraic[34] = 60.0000+500.000*(power(algebraic[30], 2.00000))
    algebraic[36] = custom_piecewise([less(states[7] , -50.0000), 5.00000 , True, 1.00000])
    rates[27] = states[26]*algebraic[36]*algebraic[34]-algebraic[36]*0.600000*states[27]
    algebraic[32] = 500.000*(power(algebraic[30], 2.00000))
    algebraic[38] = (1.00000-states[26])-states[27]
    rates[26] = algebraic[38]*algebraic[36]*algebraic[32]-states[26]*algebraic[36]*algebraic[34]
    algebraic[51] = custom_piecewise([less_equal(states[2] , states[0]), constants[42]*exp(constants[41]*(states[4]-states[3])) , True, constants[44]*exp(constants[43]*(states[4]-states[3]))])
    algebraic[37] = states[4]+states[5]
    algebraic[15] = constants[6]*(exp(constants[5]*states[5])-1.00000)
    algebraic[26] = algebraic[15]
    algebraic[40] = custom_piecewise([equal(constants[0] , 1.00000) & greater(algebraic[26] , constants[17]) & less_equal(algebraic[37] , constants[18]*(1.00000+0.000100000)), 1.00000 , True, 0.00000])
    algebraic[47] = custom_piecewise([equal(algebraic[40] , 1.00000), (constants[1]*constants[2]*exp(constants[1]*(states[4]-states[3]))*states[0])/(constants[1]*constants[2]*exp(constants[1]*(states[4]-states[3]))+constants[3]*constants[4]*exp(constants[3]*states[4])) , True, (constants[1]*constants[2]*exp(constants[1]*(states[4]-states[3]))*states[0])/(constants[1]*constants[2]*exp(constants[1]*(states[4]-states[3]))+constants[3]*constants[4]*exp(constants[3]*states[4])+constants[5]*constants[6]*exp(constants[5]*states[5]))])
    algebraic[54] = custom_piecewise([equal(algebraic[51] , 0.00000), algebraic[47] , True, states[2]])
    rates[4] = algebraic[54]
    algebraic[55] = custom_piecewise([equal(algebraic[40] , 1.00000), 0.00000 , equal(algebraic[40] , 0.00000) & equal(algebraic[51] , 0.00000), -algebraic[47] , equal(algebraic[40] , 0.00000) & algebraic[51] != 0.00000, -states[2] , True, float('nan')])
    rates[5] = custom_piecewise([greater_equal(voi , constants[33]) & less_equal(voi , constants[34]) & less_equal((voi-constants[33])-floor((voi-constants[33])/constants[19])*constants[19] , constants[35]), -((states[4]+states[5])-constants[18])/constants[35] , True, algebraic[55]])
    algebraic[31] = ((power(states[6]/constants[30], constants[8]))*(1.00000+power(constants[9], constants[8])))/(power(states[6]/constants[30], constants[8])+power(constants[9], constants[8]))
    algebraic[33] = custom_piecewise([less(constants[31]*states[3]+constants[32] , 0.00000), 0.00000 , less(constants[31]*states[3]+constants[32] , 1.00000), constants[31]*states[3]+constants[32] , True, 1.00000])
    algebraic[35] = custom_piecewise([greater(states[3] , 0.550000), (states[3]+constants[20])/(0.460000+constants[20]) , True, (constants[20]+0.550000)*1.00000])
    algebraic[2] = custom_piecewise([less_equal(states[0] , 0.00000), constants[21]-(constants[22]*states[0])/constants[13] , less_equal(states[0] , constants[123]) & less(0.00000 , states[0]), ((constants[24]-constants[23])*states[0])/constants[123]+constants[23] , True, constants[24]/(power(1.00000+(constants[29]*(states[0]-constants[123]))/constants[13], constants[28]))])
    algebraic[53] = custom_piecewise([less_equal(states[0] , 0.00000), (constants[14]*(1.00000+states[0]/constants[13]))/(constants[14]-states[0]/constants[13]) , True, (1.00000+constants[15])-((power(constants[15], 2.00000))*constants[14])/(((constants[14]*constants[15])/constants[130])*(power(states[0]/constants[13], 2.00000))+((constants[14]+1.00000)*states[0])/constants[13]+constants[14]*constants[15])])
    algebraic[57] = custom_piecewise([less_equal(-constants[13] , states[0]) & less_equal(states[0] , 0.00000), 1.00000+(0.600000*states[0])/constants[13] , less(0.00000 , states[0]) & less_equal(states[0] , constants[124]), algebraic[53]/(((0.400000*constants[14]+1.00000)*states[0])/(constants[14]*constants[13])+1.00000) , True, (algebraic[53]*exp(-constants[26]*(power((states[0]-constants[124])/constants[13], constants[16]))))/(((0.400000*constants[14]+1.00000)*states[0])/(constants[14]*constants[13])+1.00000)])
    algebraic[59] = constants[10]*constants[11]*algebraic[2]*constants[12]*algebraic[57]
    algebraic[61] = constants[11]*algebraic[2]*(1.00000-constants[10]*constants[12]*algebraic[57])
    algebraic[63] = algebraic[59]*algebraic[31]*algebraic[33]*algebraic[35]*(1.00000-states[1])-algebraic[61]*states[1]
    rates[1] = algebraic[63]
    algebraic[44] = custom_piecewise([less_equal(states[0] , 0.00000), constants[38]*exp(constants[37]*states[3]) , True, constants[40]*exp(constants[39]*states[3])])
    algebraic[69] = algebraic[53]/algebraic[57]
    algebraic[65] = (constants[14]*1.00000*(1.00000+0.400000*constants[14]+(1.20000*states[0])/constants[13]+0.600000*(power(states[0]/constants[13], 2.00000))))/(constants[13]*(power((constants[14]-states[0]/constants[13])*(1.00000+(0.600000*states[0])/constants[13]), 2.00000)))
    algebraic[67] = (1.00000/constants[13])*exp(-constants[26]*(power(states[0]/constants[13]-constants[124]/constants[13], constants[16])))*((0.400000*constants[14]+1.00000)/constants[14]+constants[26]*constants[16]*(1.00000+((0.400000*constants[14]+1.00000)*states[0])/(constants[14]*constants[13]))*(power(states[0]/constants[13]-constants[124]/constants[13], constants[16]-1.00000)))
    algebraic[72] = custom_piecewise([less_equal(states[0] , -constants[13]), constants[125] , less(-constants[13] , states[0]) & less_equal(states[0] , 0.00000), algebraic[65] , less(0.00000 , states[0]) & less_equal(states[0] , constants[124]), constants[126] , True, algebraic[67]])
    algebraic[42] = custom_piecewise([less_equal(states[0] , 0.00000), constants[37] , True, constants[39]])
    algebraic[74] = custom_piecewise([equal(algebraic[40] , 1.00000), -(constants[7]*algebraic[63]*algebraic[69]+algebraic[42]*algebraic[44]*(power(states[0], 2.00000))+constants[3]*constants[4]*exp(constants[3]*states[4])*states[2])/(constants[7]*states[1]*algebraic[72]+algebraic[44]) , True, -(constants[7]*algebraic[63]*algebraic[69]+algebraic[42]*algebraic[44]*(power(states[0], 2.00000))+(constants[3]*constants[4]*exp(constants[3]*states[4])+constants[5]*constants[6]*exp(constants[5]*states[5]))*states[2])/(constants[7]*states[1]*algebraic[72]+algebraic[44])])
    rates[0] = custom_piecewise([equal(algebraic[51] , 0.00000), (constants[1]*constants[2]*exp(constants[1]*(states[4]-states[3]))*(algebraic[47]-states[0])-(constants[7]*algebraic[63]*algebraic[69]+algebraic[42]*algebraic[44]*(power(states[0], 2.00000))))/(constants[7]*states[1]*algebraic[72]+algebraic[44]) , True, algebraic[74]])
    algebraic[49] = custom_piecewise([less_equal(states[2] , states[0]), constants[41] , True, constants[43]])
    rates[2] = custom_piecewise([equal(algebraic[40] , 1.00000) & algebraic[51] != 0.00000, ((algebraic[51]*(algebraic[74]-algebraic[49]*(power(states[2]-states[0], 2.00000)))-constants[1]*constants[2]*exp(constants[1]*(states[4]-states[3]))*(states[2]-states[0]))-constants[3]*constants[4]*exp(constants[3]*states[4])*states[2])/algebraic[51] , equal(algebraic[40] , 0.00000) & algebraic[51] != 0.00000, (algebraic[74]-algebraic[49]*(power(states[2]-states[0], 2.00000)))-(constants[1]*constants[2]*exp(constants[1]*(states[4]-states[3]))*(states[2]-states[0])+(constants[3]*constants[4]*exp(constants[3]*states[4])+constants[5]*constants[6]*exp(constants[5]*states[5]))*states[2])/algebraic[51] , equal(algebraic[51] , 0.00000), 0.00000 , True, float('nan')])
    algebraic[43] = ((constants[45]*constants[46])/constants[47])*log(states[8]/states[9])
    algebraic[52] = (((constants[53]*states[8])/(states[8]+constants[52]))*(states[7]-algebraic[43]))/(1.00000+exp((((states[7]-algebraic[43])-10.0000)*constants[47]*1.25000)/(constants[45]*constants[46])))
    algebraic[79] = constants[74]*(constants[75]+states[22]*(1.00000-constants[75]))*states[23]*(states[7]-algebraic[43])
    algebraic[56] = (((constants[54]*states[12]+constants[55]*states[13])*1.00000)/(1.00000+exp((states[7]+9.00000)/22.4000)))*(states[7]-algebraic[43])
    algebraic[46] = ((constants[45]*constants[46])/constants[47])*log((states[8]+constants[49]*constants[50])/(states[9]+constants[49]*states[10]))
    algebraic[58] = constants[56]*(power(states[14], 2.00000))*(states[7]-algebraic[46])
    algebraic[68] = ((((1.00000-constants[69])*constants[64]*constants[63]*states[18]*states[19]*states[20]*(states[7]-50.0000)*constants[47])/(constants[45]*constants[46]))/(1.00000-exp((-(states[7]-50.0000)*constants[47])/(constants[45]*constants[46]))))*(states[9]*exp((50.0000*constants[47])/(constants[45]*constants[46]))-states[8]*exp((-(states[7]-50.0000)*constants[47])/(constants[45]*constants[46])))
    algebraic[75] = (((constants[69]*constants[64]*constants[63]*states[18]*states[19]*states[21]*(states[7]-50.0000)*constants[47])/(constants[45]*constants[46]))/(1.00000-exp((-(states[7]-50.0000)*constants[47])/(constants[45]*constants[46]))))*(states[9]*exp((50.0000*constants[47])/(constants[45]*constants[46]))-states[8]*exp((-(states[7]-50.0000)*constants[47])/(constants[45]*constants[46])))
    algebraic[80] = (((constants[76]*states[8])/(constants[77]+states[8]))*states[10])/(constants[78]+states[10])
    rates[8] = (1.00000*(((algebraic[56]+algebraic[58]+algebraic[52]+algebraic[79])-(1.00000/(constants[99]-1.00000))*algebraic[80])+algebraic[68]+algebraic[75]))/(1.00000*constants[132]*constants[47])-0.700000*(states[8]-constants[100])
    rates[9] = (-1.00000/(1.00000*constants[131]*constants[47]))*((algebraic[52]+algebraic[56]+algebraic[58]+algebraic[68]+algebraic[75]+algebraic[79])-(1.00000/(constants[99]-1.00000))*algebraic[80])
    algebraic[50] = ((constants[45]*constants[46])/constants[47])*log((constants[50]+0.120000*states[8])/(states[10]+0.120000*states[9]))
    algebraic[60] = constants[59]*(power(states[15], 3.00000))*states[16]*(states[7]-algebraic[50])
    algebraic[41] = ((constants[45]*constants[46])/constants[47])*log(constants[50]/states[10])
    algebraic[64] = constants[62]*(states[7]-algebraic[41])
    algebraic[62] = ((constants[61]*1.00000)/(1.00000+exp(-(states[7]+52.0000)/8.00000)))*(states[7]-algebraic[41])
    algebraic[70] = ((((1.00000-constants[69])*constants[65]*constants[63]*states[18]*states[19]*states[20]*(states[7]-50.0000)*constants[47])/(constants[45]*constants[46]))/(1.00000-exp((-(states[7]-50.0000)*constants[47])/(constants[45]*constants[46]))))*(states[10]*exp((50.0000*constants[47])/(constants[45]*constants[46]))-constants[50]*exp((-(states[7]-50.0000)*constants[47])/(constants[45]*constants[46])))
    algebraic[76] = (((constants[69]*constants[65]*constants[63]*states[18]*states[19]*states[21]*(states[7]-50.0000)*constants[47])/(constants[45]*constants[46]))/(1.00000-exp((-(states[7]-50.0000)*constants[47])/(constants[45]*constants[46]))))*(states[10]*exp((50.0000*constants[47])/(constants[45]*constants[46]))-constants[50]*exp((-(states[7]-50.0000)*constants[47])/(constants[45]*constants[46])))
    algebraic[81] = ((1.00000-constants[83])*constants[79]*(exp((constants[82]*(constants[80]-2.00000)*states[7]*constants[47])/(constants[45]*constants[46]))*(power(states[10], constants[80]))*constants[51]-exp(((constants[82]-1.00000)*(constants[80]-2.00000)*states[7]*constants[47])/(constants[45]*constants[46]))*(power(constants[50], constants[80]))*states[11]))/((1.00000+constants[81]*(states[11]*(power(constants[50], constants[80]))+constants[51]*(power(states[10], constants[80]))))*(1.00000+states[11]/0.00690000))
    algebraic[82] = (constants[83]*constants[79]*(exp((constants[82]*(constants[80]-2.00000)*states[7]*constants[47])/(constants[45]*constants[46]))*(power(states[10], constants[80]))*constants[51]-exp(((constants[82]-1.00000)*(constants[80]-2.00000)*states[7]*constants[47])/(constants[45]*constants[46]))*(power(constants[50], constants[80]))*states[17]))/((1.00000+constants[81]*(states[17]*(power(constants[50], constants[80]))+constants[51]*(power(states[10], constants[80]))))*(1.00000+states[17]/0.00690000))
    algebraic[66] = ((((1.00000-constants[69])*4.00000*constants[63]*states[18]*states[19]*states[20]*(states[7]-50.0000)*constants[47])/(constants[45]*constants[46]))/(1.00000-exp((-(states[7]-50.0000)*constants[47]*2.00000)/(constants[45]*constants[46]))))*(states[11]*exp((100.000*constants[47])/(constants[45]*constants[46]))-constants[51]*exp((-(states[7]-50.0000)*constants[47]*2.00000)/(constants[45]*constants[46])))
    algebraic[73] = (((constants[69]*4.00000*constants[63]*states[18]*states[19]*states[21]*(states[7]-50.0000)*constants[47])/(constants[45]*constants[46]))/(1.00000-exp((-(states[7]-50.0000)*constants[47]*2.00000)/(constants[45]*constants[46]))))*(states[11]*exp((100.000*constants[47])/(constants[45]*constants[46]))-constants[51]*exp((-(states[7]-50.0000)*constants[47]*2.00000)/(constants[45]*constants[46])))
    algebraic[48] = ((0.500000*constants[45]*constants[46])/constants[47])*log(constants[51]/states[11])
    algebraic[78] = constants[73]*(states[7]-algebraic[48])
    algebraic[39] = custom_piecewise([greater_equal(voi , constants[33]) & less_equal(voi , constants[34]) & less_equal((voi-constants[33])-floor((voi-constants[33])/constants[19])*constants[19] , constants[35]), constants[36] , True, 0.00000])
    rates[7] = (-1.00000/constants[48])*(algebraic[39]+algebraic[52]+algebraic[79]+algebraic[56]+algebraic[58]+algebraic[80]+algebraic[60]+algebraic[64]+algebraic[62]+algebraic[70]+algebraic[76]+algebraic[81]+algebraic[82]+algebraic[66]+algebraic[73]+algebraic[68]+algebraic[75]+algebraic[78])
    rates[17] = (-1.00000*algebraic[73]+(2.00000*algebraic[82])/(constants[80]-2.00000))/(2.00000*1.00000*constants[104]*constants[131]*constants[47])-states[17]*constants[108]
    algebraic[84] = algebraic[81]+algebraic[82]
    rates[10] = (-1.00000/(1.00000*constants[131]*constants[47]))*(algebraic[60]+algebraic[62]+algebraic[64]+algebraic[70]+algebraic[76]+(constants[99]/(constants[99]-1.00000))*algebraic[80]+(constants[80]/(constants[80]-2.00000))*algebraic[84])
    algebraic[83] = states[11]+states[24]*constants[128]+constants[84]*constants[85]+constants[84]
    algebraic[85] = custom_piecewise([equal(constants[89] , 0.00000), (states[11]/algebraic[83])*constants[87]-((states[24]*constants[128])/algebraic[83])*constants[88] , True, ((states[11]/algebraic[83])*constants[87])/(1.00000+states[24]/constants[90])-((states[24]*constants[128])/algebraic[83])*constants[88]])
    algebraic[86] = constants[91]*(states[24]-states[25])
    rates[24] = (constants[103]/constants[107])*algebraic[85]-algebraic[86]
    algebraic[87] = ((power(states[26]/(states[26]+0.250000), 2.00000))*constants[96]+constants[95])*states[25]
    rates[25] = ((constants[107]/constants[105])*algebraic[86]-algebraic[87])/(1.00000+(constants[129]*constants[94])/(power(states[25]+constants[129], 2.00000)))
    algebraic[88] = (constants[30]*states[1])/(algebraic[35]*states[6])
    algebraic[89] = custom_piecewise([less_equal(algebraic[88] , 0.00000), constants[109] , less_equal(algebraic[88] , 1.00000), power(constants[109], algebraic[88]) , True, 1.00000])
    rates[6] = constants[110]*(constants[30]-states[6])*states[11]-constants[111]*exp(-constants[27]*states[6])*algebraic[89]*states[6]
    rates[28] = constants[113]*(constants[112]-states[28])*states[11]-constants[114]*states[28]
    rates[29] = constants[116]*(constants[115]-states[29])*states[11]-constants[117]*states[29]
    rates[11] = (((((-1.00000/(2.00000*1.00000*constants[131]*constants[47]))*((algebraic[66]+algebraic[78])-(2.00000/(constants[80]-2.00000))*algebraic[81])+states[17]*constants[104]*constants[108]+(algebraic[87]*constants[105])/constants[103])-rates[6])-rates[28])-rates[29])-algebraic[85]
    return(rates)

def computeAlgebraic(constants, states, voi):
    algebraic = array([[0.0] * len(voi)] * sizeAlgebraic)
    states = array(states)
    voi = array(voi)
    algebraic[3] = states[0]
    algebraic[4] = 50.0000/(1.00000+exp(-(states[7]-5.00000)/9.00000))
    algebraic[16] = 0.0500000*exp(-(states[7]-20.0000)/15.0000)
    algebraic[5] = 50.0000/(1.00000+exp(-(states[7]-5.00000)/9.00000))
    algebraic[17] = 0.400000*exp(-(power((states[7]+30.0000)/30.0000, 3.00000)))
    algebraic[6] = 14.0000/(1.00000+exp(-(states[7]-40.0000)/9.00000))
    algebraic[18] = 1.00000*exp(-states[7]/45.0000)
    algebraic[8] = 20.0000*exp(-0.125000*(states[7]+75.0000))
    algebraic[20] = 2000.00/(1.00000+320.000*exp(-0.100000*(states[7]+75.0000)))
    algebraic[11] = 0.0330000*exp(-states[7]/17.0000)
    algebraic[23] = 33.0000/(1.00000+exp(-0.125000*(states[7]+10.0000)))
    algebraic[14] = constants[119]*(states[30]+20.0000)+constants[121]*(states[30]-constants[122])
    algebraic[25] = -constants[118]*(states[7]-states[30])
    algebraic[7] = states[7]+41.0000
    algebraic[19] = custom_piecewise([less(fabs(algebraic[7]) , constants[60]), 2000.00 , True, (200.000*algebraic[7])/(1.00000-exp(-0.100000*algebraic[7]))])
    algebraic[27] = 8000.00*exp(-0.0560000*(states[7]+66.0000))
    algebraic[9] = (states[7]+24.0000)-5.00000
    algebraic[21] = custom_piecewise([less(fabs(algebraic[9]) , 0.000100000), 120.000 , True, (30.0000*algebraic[9])/(1.00000-exp(-algebraic[9]/4.00000))])
    algebraic[28] = custom_piecewise([less(fabs(algebraic[9]) , 0.000100000), 120.000 , True, (12.0000*algebraic[9])/(exp(algebraic[9]/10.0000)-1.00000)])
    algebraic[10] = states[7]+34.0000
    algebraic[22] = custom_piecewise([less(fabs(algebraic[10]) , constants[72]), 25.0000 , True, (6.25000*algebraic[10])/(exp(algebraic[10]/4.00000)-1.00000)])
    algebraic[29] = 12.0000/(1.00000+exp((-1.00000*(states[7]+34.0000))/4.00000))
    algebraic[13] = states[11]/(states[11]+constants[97])
    algebraic[24] = states[17]/(states[17]+constants[98])
    algebraic[30] = algebraic[13]+(1.00000-algebraic[13])*algebraic[24]
    algebraic[34] = 60.0000+500.000*(power(algebraic[30], 2.00000))
    algebraic[36] = custom_piecewise([less(states[7] , -50.0000), 5.00000 , True, 1.00000])
    algebraic[32] = 500.000*(power(algebraic[30], 2.00000))
    algebraic[38] = (1.00000-states[26])-states[27]
    algebraic[51] = custom_piecewise([less_equal(states[2] , states[0]), constants[42]*exp(constants[41]*(states[4]-states[3])) , True, constants[44]*exp(constants[43]*(states[4]-states[3]))])
    algebraic[37] = states[4]+states[5]
    algebraic[15] = constants[6]*(exp(constants[5]*states[5])-1.00000)
    algebraic[26] = algebraic[15]
    algebraic[40] = custom_piecewise([equal(constants[0] , 1.00000) & greater(algebraic[26] , constants[17]) & less_equal(algebraic[37] , constants[18]*(1.00000+0.000100000)), 1.00000 , True, 0.00000])
    algebraic[47] = custom_piecewise([equal(algebraic[40] , 1.00000), (constants[1]*constants[2]*exp(constants[1]*(states[4]-states[3]))*states[0])/(constants[1]*constants[2]*exp(constants[1]*(states[4]-states[3]))+constants[3]*constants[4]*exp(constants[3]*states[4])) , True, (constants[1]*constants[2]*exp(constants[1]*(states[4]-states[3]))*states[0])/(constants[1]*constants[2]*exp(constants[1]*(states[4]-states[3]))+constants[3]*constants[4]*exp(constants[3]*states[4])+constants[5]*constants[6]*exp(constants[5]*states[5]))])
    algebraic[54] = custom_piecewise([equal(algebraic[51] , 0.00000), algebraic[47] , True, states[2]])
    algebraic[55] = custom_piecewise([equal(algebraic[40] , 1.00000), 0.00000 , equal(algebraic[40] , 0.00000) & equal(algebraic[51] , 0.00000), -algebraic[47] , equal(algebraic[40] , 0.00000) & algebraic[51] != 0.00000, -states[2] , True, float('nan')])
    algebraic[31] = ((power(states[6]/constants[30], constants[8]))*(1.00000+power(constants[9], constants[8])))/(power(states[6]/constants[30], constants[8])+power(constants[9], constants[8]))
    algebraic[33] = custom_piecewise([less(constants[31]*states[3]+constants[32] , 0.00000), 0.00000 , less(constants[31]*states[3]+constants[32] , 1.00000), constants[31]*states[3]+constants[32] , True, 1.00000])
    algebraic[35] = custom_piecewise([greater(states[3] , 0.550000), (states[3]+constants[20])/(0.460000+constants[20]) , True, (constants[20]+0.550000)*1.00000])
    algebraic[2] = custom_piecewise([less_equal(states[0] , 0.00000), constants[21]-(constants[22]*states[0])/constants[13] , less_equal(states[0] , constants[123]) & less(0.00000 , states[0]), ((constants[24]-constants[23])*states[0])/constants[123]+constants[23] , True, constants[24]/(power(1.00000+(constants[29]*(states[0]-constants[123]))/constants[13], constants[28]))])
    algebraic[53] = custom_piecewise([less_equal(states[0] , 0.00000), (constants[14]*(1.00000+states[0]/constants[13]))/(constants[14]-states[0]/constants[13]) , True, (1.00000+constants[15])-((power(constants[15], 2.00000))*constants[14])/(((constants[14]*constants[15])/constants[130])*(power(states[0]/constants[13], 2.00000))+((constants[14]+1.00000)*states[0])/constants[13]+constants[14]*constants[15])])
    algebraic[57] = custom_piecewise([less_equal(-constants[13] , states[0]) & less_equal(states[0] , 0.00000), 1.00000+(0.600000*states[0])/constants[13] , less(0.00000 , states[0]) & less_equal(states[0] , constants[124]), algebraic[53]/(((0.400000*constants[14]+1.00000)*states[0])/(constants[14]*constants[13])+1.00000) , True, (algebraic[53]*exp(-constants[26]*(power((states[0]-constants[124])/constants[13], constants[16]))))/(((0.400000*constants[14]+1.00000)*states[0])/(constants[14]*constants[13])+1.00000)])
    algebraic[59] = constants[10]*constants[11]*algebraic[2]*constants[12]*algebraic[57]
    algebraic[61] = constants[11]*algebraic[2]*(1.00000-constants[10]*constants[12]*algebraic[57])
    algebraic[63] = algebraic[59]*algebraic[31]*algebraic[33]*algebraic[35]*(1.00000-states[1])-algebraic[61]*states[1]
    algebraic[44] = custom_piecewise([less_equal(states[0] , 0.00000), constants[38]*exp(constants[37]*states[3]) , True, constants[40]*exp(constants[39]*states[3])])
    algebraic[69] = algebraic[53]/algebraic[57]
    algebraic[65] = (constants[14]*1.00000*(1.00000+0.400000*constants[14]+(1.20000*states[0])/constants[13]+0.600000*(power(states[0]/constants[13], 2.00000))))/(constants[13]*(power((constants[14]-states[0]/constants[13])*(1.00000+(0.600000*states[0])/constants[13]), 2.00000)))
    algebraic[67] = (1.00000/constants[13])*exp(-constants[26]*(power(states[0]/constants[13]-constants[124]/constants[13], constants[16])))*((0.400000*constants[14]+1.00000)/constants[14]+constants[26]*constants[16]*(1.00000+((0.400000*constants[14]+1.00000)*states[0])/(constants[14]*constants[13]))*(power(states[0]/constants[13]-constants[124]/constants[13], constants[16]-1.00000)))
    algebraic[72] = custom_piecewise([less_equal(states[0] , -constants[13]), constants[125] , less(-constants[13] , states[0]) & less_equal(states[0] , 0.00000), algebraic[65] , less(0.00000 , states[0]) & less_equal(states[0] , constants[124]), constants[126] , True, algebraic[67]])
    algebraic[42] = custom_piecewise([less_equal(states[0] , 0.00000), constants[37] , True, constants[39]])
    algebraic[74] = custom_piecewise([equal(algebraic[40] , 1.00000), -(constants[7]*algebraic[63]*algebraic[69]+algebraic[42]*algebraic[44]*(power(states[0], 2.00000))+constants[3]*constants[4]*exp(constants[3]*states[4])*states[2])/(constants[7]*states[1]*algebraic[72]+algebraic[44]) , True, -(constants[7]*algebraic[63]*algebraic[69]+algebraic[42]*algebraic[44]*(power(states[0], 2.00000))+(constants[3]*constants[4]*exp(constants[3]*states[4])+constants[5]*constants[6]*exp(constants[5]*states[5]))*states[2])/(constants[7]*states[1]*algebraic[72]+algebraic[44])])
    algebraic[49] = custom_piecewise([less_equal(states[2] , states[0]), constants[41] , True, constants[43]])
    algebraic[43] = ((constants[45]*constants[46])/constants[47])*log(states[8]/states[9])
    algebraic[52] = (((constants[53]*states[8])/(states[8]+constants[52]))*(states[7]-algebraic[43]))/(1.00000+exp((((states[7]-algebraic[43])-10.0000)*constants[47]*1.25000)/(constants[45]*constants[46])))
    algebraic[79] = constants[74]*(constants[75]+states[22]*(1.00000-constants[75]))*states[23]*(states[7]-algebraic[43])
    algebraic[56] = (((constants[54]*states[12]+constants[55]*states[13])*1.00000)/(1.00000+exp((states[7]+9.00000)/22.4000)))*(states[7]-algebraic[43])
    algebraic[46] = ((constants[45]*constants[46])/constants[47])*log((states[8]+constants[49]*constants[50])/(states[9]+constants[49]*states[10]))
    algebraic[58] = constants[56]*(power(states[14], 2.00000))*(states[7]-algebraic[46])
    algebraic[68] = ((((1.00000-constants[69])*constants[64]*constants[63]*states[18]*states[19]*states[20]*(states[7]-50.0000)*constants[47])/(constants[45]*constants[46]))/(1.00000-exp((-(states[7]-50.0000)*constants[47])/(constants[45]*constants[46]))))*(states[9]*exp((50.0000*constants[47])/(constants[45]*constants[46]))-states[8]*exp((-(states[7]-50.0000)*constants[47])/(constants[45]*constants[46])))
    algebraic[75] = (((constants[69]*constants[64]*constants[63]*states[18]*states[19]*states[21]*(states[7]-50.0000)*constants[47])/(constants[45]*constants[46]))/(1.00000-exp((-(states[7]-50.0000)*constants[47])/(constants[45]*constants[46]))))*(states[9]*exp((50.0000*constants[47])/(constants[45]*constants[46]))-states[8]*exp((-(states[7]-50.0000)*constants[47])/(constants[45]*constants[46])))
    algebraic[80] = (((constants[76]*states[8])/(constants[77]+states[8]))*states[10])/(constants[78]+states[10])
    algebraic[50] = ((constants[45]*constants[46])/constants[47])*log((constants[50]+0.120000*states[8])/(states[10]+0.120000*states[9]))
    algebraic[60] = constants[59]*(power(states[15], 3.00000))*states[16]*(states[7]-algebraic[50])
    algebraic[41] = ((constants[45]*constants[46])/constants[47])*log(constants[50]/states[10])
    algebraic[64] = constants[62]*(states[7]-algebraic[41])
    algebraic[62] = ((constants[61]*1.00000)/(1.00000+exp(-(states[7]+52.0000)/8.00000)))*(states[7]-algebraic[41])
    algebraic[70] = ((((1.00000-constants[69])*constants[65]*constants[63]*states[18]*states[19]*states[20]*(states[7]-50.0000)*constants[47])/(constants[45]*constants[46]))/(1.00000-exp((-(states[7]-50.0000)*constants[47])/(constants[45]*constants[46]))))*(states[10]*exp((50.0000*constants[47])/(constants[45]*constants[46]))-constants[50]*exp((-(states[7]-50.0000)*constants[47])/(constants[45]*constants[46])))
    algebraic[76] = (((constants[69]*constants[65]*constants[63]*states[18]*states[19]*states[21]*(states[7]-50.0000)*constants[47])/(constants[45]*constants[46]))/(1.00000-exp((-(states[7]-50.0000)*constants[47])/(constants[45]*constants[46]))))*(states[10]*exp((50.0000*constants[47])/(constants[45]*constants[46]))-constants[50]*exp((-(states[7]-50.0000)*constants[47])/(constants[45]*constants[46])))
    algebraic[81] = ((1.00000-constants[83])*constants[79]*(exp((constants[82]*(constants[80]-2.00000)*states[7]*constants[47])/(constants[45]*constants[46]))*(power(states[10], constants[80]))*constants[51]-exp(((constants[82]-1.00000)*(constants[80]-2.00000)*states[7]*constants[47])/(constants[45]*constants[46]))*(power(constants[50], constants[80]))*states[11]))/((1.00000+constants[81]*(states[11]*(power(constants[50], constants[80]))+constants[51]*(power(states[10], constants[80]))))*(1.00000+states[11]/0.00690000))
    algebraic[82] = (constants[83]*constants[79]*(exp((constants[82]*(constants[80]-2.00000)*states[7]*constants[47])/(constants[45]*constants[46]))*(power(states[10], constants[80]))*constants[51]-exp(((constants[82]-1.00000)*(constants[80]-2.00000)*states[7]*constants[47])/(constants[45]*constants[46]))*(power(constants[50], constants[80]))*states[17]))/((1.00000+constants[81]*(states[17]*(power(constants[50], constants[80]))+constants[51]*(power(states[10], constants[80]))))*(1.00000+states[17]/0.00690000))
    algebraic[66] = ((((1.00000-constants[69])*4.00000*constants[63]*states[18]*states[19]*states[20]*(states[7]-50.0000)*constants[47])/(constants[45]*constants[46]))/(1.00000-exp((-(states[7]-50.0000)*constants[47]*2.00000)/(constants[45]*constants[46]))))*(states[11]*exp((100.000*constants[47])/(constants[45]*constants[46]))-constants[51]*exp((-(states[7]-50.0000)*constants[47]*2.00000)/(constants[45]*constants[46])))
    algebraic[73] = (((constants[69]*4.00000*constants[63]*states[18]*states[19]*states[21]*(states[7]-50.0000)*constants[47])/(constants[45]*constants[46]))/(1.00000-exp((-(states[7]-50.0000)*constants[47]*2.00000)/(constants[45]*constants[46]))))*(states[11]*exp((100.000*constants[47])/(constants[45]*constants[46]))-constants[51]*exp((-(states[7]-50.0000)*constants[47]*2.00000)/(constants[45]*constants[46])))
    algebraic[48] = ((0.500000*constants[45]*constants[46])/constants[47])*log(constants[51]/states[11])
    algebraic[78] = constants[73]*(states[7]-algebraic[48])
    algebraic[39] = custom_piecewise([greater_equal(voi , constants[33]) & less_equal(voi , constants[34]) & less_equal((voi-constants[33])-floor((voi-constants[33])/constants[19])*constants[19] , constants[35]), constants[36] , True, 0.00000])
    algebraic[84] = algebraic[81]+algebraic[82]
    algebraic[83] = states[11]+states[24]*constants[128]+constants[84]*constants[85]+constants[84]
    algebraic[85] = custom_piecewise([equal(constants[89] , 0.00000), (states[11]/algebraic[83])*constants[87]-((states[24]*constants[128])/algebraic[83])*constants[88] , True, ((states[11]/algebraic[83])*constants[87])/(1.00000+states[24]/constants[90])-((states[24]*constants[128])/algebraic[83])*constants[88]])
    algebraic[86] = constants[91]*(states[24]-states[25])
    algebraic[87] = ((power(states[26]/(states[26]+0.250000), 2.00000))*constants[96]+constants[95])*states[25]
    algebraic[88] = (constants[30]*states[1])/(algebraic[35]*states[6])
    algebraic[89] = custom_piecewise([less_equal(algebraic[88] , 0.00000), constants[109] , less_equal(algebraic[88] , 1.00000), power(constants[109], algebraic[88]) , True, 1.00000])
    algebraic[0] = constants[2]*(exp(constants[1]*(states[4]-states[3]))-1.00000)
    algebraic[1] = constants[4]*(exp(constants[3]*states[4])-1.00000)
    algebraic[12] = exp(0.0800000*(states[7]-40.0000))
    algebraic[45] = ((constants[57]*states[10])/(states[10]+constants[58]))*(states[7]-algebraic[43])
    algebraic[71] = constants[7]*algebraic[69]*states[1]
    algebraic[77] = algebraic[66]+algebraic[68]+algebraic[70]+algebraic[73]+algebraic[75]+algebraic[76]
    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)