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 = 86
sizeStates = 44
sizeConstants = 96
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 Circ_Environment (ms)"
    legend_constants[0] = "HR in component Circ_Environment (ratepm)"
    legend_constants[91] = "hrf in component Circ_Environment (Hz)"
    legend_algebraic[0] = "beattime in component LATiming (ms)"
    legend_constants[1] = "Esys in component LAElastanceFunction (elastance)"
    legend_constants[2] = "Edia in component LAElastanceFunction (elastance)"
    legend_constants[3] = "TsK in component LAElastanceFunction (ms)"
    legend_constants[92] = "Ts in component LAElastanceFunction (ms)"
    legend_algebraic[1] = "E_LA in component LAElastanceFunction (elastance)"
    legend_algebraic[82] = "P_PV in component PulmonaryVein (mmHg)"
    legend_algebraic[9] = "F_LV in component LeftVentricle (flow)"
    legend_constants[4] = "P_LA_ext in component LeftAtrium (mmHg)"
    legend_states[0] = "V_LA in component LeftAtrium (mL)"
    legend_constants[5] = "V_LA_0 in component LeftAtrium (mL)"
    legend_algebraic[3] = "P_LA in component LeftAtrium (mmHg)"
    legend_algebraic[83] = "F_LA in component LeftAtrium (flow)"
    legend_constants[6] = "R_LA in component LeftAtrium (resistance)"
    legend_algebraic[2] = "conc_X_LA in component LeftAtrium (concentration)"
    legend_algebraic[81] = "conc_X_PV in component PulmonaryVein (concentration)"
    legend_states[1] = "sub_X_LA in component LeftAtrium (mL)"
    legend_constants[7] = "PRint in component LVTiming (ms)"
    legend_algebraic[4] = "beattime in component LVTiming (ms)"
    legend_constants[8] = "Esys in component LVElastanceFunction (elastance)"
    legend_constants[9] = "Edia in component LVElastanceFunction (elastance)"
    legend_constants[10] = "TsK in component LVElastanceFunction (second)"
    legend_constants[93] = "Ts in component LVElastanceFunction (second)"
    legend_algebraic[5] = "E_LV in component LVElastanceFunction (elastance)"
    legend_algebraic[12] = "F_AO in component Aorta (flow)"
    legend_algebraic[15] = "F_CO in component CoronaryCirc (flow)"
    legend_constants[11] = "pLVPext in component LeftVentricle (mmHg)"
    legend_states[2] = "V_LV in component LeftVentricle (mL)"
    legend_constants[12] = "V_LV_0 in component LeftVentricle (mL)"
    legend_algebraic[7] = "P_LV in component LeftVentricle (mmHg)"
    legend_constants[13] = "R_LV in component LeftVentricle (resistance)"
    legend_algebraic[6] = "conc_X_LV in component LeftVentricle (concentration)"
    legend_states[3] = "sub_X_LV in component LeftVentricle (mL)"
    legend_algebraic[8] = "fiber_stress in component LeftVentricle (mmHg)"
    legend_algebraic[19] = "F_CR in component CarotidCirc (flow)"
    legend_algebraic[22] = "F_AR in component Arteries (flow)"
    legend_constants[14] = "C_AO in component Aorta (capacitance)"
    legend_constants[15] = "L_AO in component Aorta (inductance)"
    legend_constants[16] = "R_AO in component Aorta (resistance)"
    legend_constants[17] = "V_AO_0 in component Aorta (mL)"
    legend_states[4] = "V_AO in component Aorta (mL)"
    legend_states[5] = "F_AOalways in component Aorta (flow)"
    legend_algebraic[11] = "P_AO in component Aorta (mmHg)"
    legend_algebraic[10] = "conc_X_AO in component Aorta (concentration)"
    legend_states[6] = "sub_X_AO in component Aorta (mL)"
    legend_algebraic[62] = "P_RA in component RightAtrium (mmHg)"
    legend_constants[18] = "R_CO in component CoronaryCirc (resistance)"
    legend_constants[19] = "R_COV in component CoronaryCirc (resistance)"
    legend_constants[20] = "V_CO_0 in component CoronaryCirc (mL)"
    legend_constants[21] = "C_CO in component CoronaryCirc (capacitance)"
    legend_algebraic[63] = "F_COV in component CoronaryCirc (flow)"
    legend_states[7] = "V_CO in component CoronaryCirc (mL)"
    legend_algebraic[14] = "P_CO in component CoronaryCirc (mmHg)"
    legend_algebraic[13] = "conc_X_CO in component CoronaryCirc (concentration)"
    legend_states[8] = "sub_X_CO in component CoronaryCirc (mL)"
    legend_algebraic[55] = "P_VC in component VenaCava (mmHg)"
    legend_constants[22] = "R_CR in component CarotidCirc (resistance)"
    legend_constants[23] = "R_CRV in component CarotidCirc (resistance)"
    legend_constants[24] = "V_CR_0 in component CarotidCirc (mL)"
    legend_constants[25] = "C_CR in component CarotidCirc (capacitance)"
    legend_algebraic[56] = "F_CRV in component CarotidCirc (flow)"
    legend_states[9] = "V_CR in component CarotidCirc (mL)"
    legend_algebraic[18] = "P_CR in component CarotidCirc (mmHg)"
    legend_algebraic[17] = "conc_X_CR in component CarotidCirc (concentration)"
    legend_states[10] = "sub_X_CR in component CarotidCirc (mL)"
    legend_algebraic[16] = "conc_X_CRtissue in component CarotidCirc (concentration)"
    legend_constants[26] = "D_CRtissue in component CarotidCirc (flow)"
    legend_constants[27] = "V_CRtissue in component CarotidCirc (mL)"
    legend_states[11] = "sub_X_CRtissue in component CarotidCirc (mL)"
    legend_algebraic[25] = "F_AD in component Adipose (flow)"
    legend_algebraic[28] = "F_MU in component Muscle (flow)"
    legend_algebraic[31] = "F_GI in component GutIntestine (flow)"
    legend_algebraic[35] = "F_LI in component Liver (flow)"
    legend_algebraic[38] = "F_KI in component Kidney (flow)"
    legend_algebraic[44] = "F_OT in component OtherTissue (flow)"
    legend_algebraic[41] = "F_SK in component Skin (flow)"
    legend_constants[28] = "R_AR in component Arteries (resistance)"
    legend_constants[29] = "C_AR in component Arteries (capacitance)"
    legend_constants[30] = "V_AR_0 in component Arteries (mL)"
    legend_states[12] = "V_AR in component Arteries (mL)"
    legend_algebraic[21] = "P_AR in component Arteries (mmHg)"
    legend_algebraic[20] = "conc_X_AR in component Arteries (concentration)"
    legend_states[13] = "sub_X_AR in component Arteries (mL)"
    legend_algebraic[46] = "P_VE in component Veins (mmHg)"
    legend_constants[31] = "C_AD in component Adipose (capacitance)"
    legend_constants[32] = "R_AD in component Adipose (resistance)"
    legend_constants[33] = "R_ADV in component Adipose (resistance)"
    legend_constants[34] = "V_AD_0 in component Adipose (mL)"
    legend_states[14] = "V_AD in component Adipose (mL)"
    legend_algebraic[47] = "F_ADV in component Adipose (flow)"
    legend_algebraic[24] = "P_AD in component Adipose (mmHg)"
    legend_algebraic[23] = "conc_X_AD in component Adipose (concentration)"
    legend_states[15] = "sub_X_AD in component Adipose (mL)"
    legend_constants[35] = "C_MU in component Muscle (capacitance)"
    legend_constants[36] = "R_MU in component Muscle (resistance)"
    legend_constants[37] = "R_MUV in component Muscle (resistance)"
    legend_constants[38] = "V_MU_0 in component Muscle (mL)"
    legend_states[16] = "V_MU in component Muscle (mL)"
    legend_algebraic[48] = "F_MUV in component Muscle (flow)"
    legend_algebraic[27] = "P_MU in component Muscle (mmHg)"
    legend_algebraic[26] = "conc_X_MU in component Muscle (concentration)"
    legend_states[17] = "sub_X_MU in component Muscle (mL)"
    legend_algebraic[33] = "P_LI in component Liver (mmHg)"
    legend_constants[39] = "C_GI in component GutIntestine (capacitance)"
    legend_constants[40] = "R_GI in component GutIntestine (resistance)"
    legend_constants[41] = "R_GIV in component GutIntestine (resistance)"
    legend_constants[42] = "V_GI_0 in component GutIntestine (mL)"
    legend_states[18] = "V_GI in component GutIntestine (mL)"
    legend_algebraic[34] = "F_GIV in component GutIntestine (flow)"
    legend_algebraic[30] = "P_GI in component GutIntestine (mmHg)"
    legend_algebraic[29] = "conc_X_GI in component GutIntestine (concentration)"
    legend_states[19] = "sub_X_GI in component GutIntestine (mL)"
    legend_constants[43] = "C_LI in component Liver (capacitance)"
    legend_constants[44] = "R_LI in component Liver (resistance)"
    legend_constants[45] = "R_LIV in component Liver (resistance)"
    legend_constants[46] = "V_LI_0 in component Liver (mL)"
    legend_states[20] = "V_LI in component Liver (mL)"
    legend_algebraic[49] = "F_LIV in component Liver (flow)"
    legend_algebraic[32] = "conc_X_LI in component Liver (concentration)"
    legend_states[21] = "sub_X_LI in component Liver (mL)"
    legend_constants[47] = "metabolism in component Liver (flow)"
    legend_constants[48] = "C_KI in component Kidney (capacitance)"
    legend_constants[49] = "R_KI in component Kidney (resistance)"
    legend_constants[50] = "V_KI_0 in component Kidney (mL)"
    legend_constants[51] = "R_KIV in component Kidney (resistance)"
    legend_states[22] = "V_KI in component Kidney (mL)"
    legend_algebraic[50] = "F_KIV in component Kidney (flow)"
    legend_algebraic[37] = "P_KI in component Kidney (mmHg)"
    legend_algebraic[36] = "conc_X_KI in component Kidney (concentration)"
    legend_states[23] = "sub_X_KI in component Kidney (mL)"
    legend_constants[52] = "C_SK in component Skin (capacitance)"
    legend_constants[53] = "R_SK in component Skin (resistance)"
    legend_constants[54] = "R_SKV in component Skin (resistance)"
    legend_constants[55] = "V_SK_0 in component Skin (mL)"
    legend_states[24] = "V_SK in component Skin (mL)"
    legend_algebraic[51] = "F_SKV in component Skin (flow)"
    legend_algebraic[40] = "P_SK in component Skin (mmHg)"
    legend_algebraic[39] = "conc_X_SK in component Skin (concentration)"
    legend_states[25] = "sub_X_SK in component Skin (mL)"
    legend_constants[56] = "C_OT in component OtherTissue (capacitance)"
    legend_constants[57] = "R_OT in component OtherTissue (resistance)"
    legend_constants[58] = "R_OTV in component OtherTissue (resistance)"
    legend_constants[59] = "V_OT_0 in component OtherTissue (mL)"
    legend_states[26] = "V_OT in component OtherTissue (mL)"
    legend_algebraic[52] = "F_OTV in component OtherTissue (flow)"
    legend_algebraic[43] = "P_OT in component OtherTissue (mmHg)"
    legend_algebraic[42] = "conc_X_OT in component OtherTissue (concentration)"
    legend_states[27] = "sub_X_OT in component OtherTissue (mL)"
    legend_constants[60] = "R_VEV in component Veins (resistance)"
    legend_states[28] = "V_VE in component Veins (mL)"
    legend_constants[61] = "C_VE in component Veins (capacitance)"
    legend_constants[62] = "V_VE_0 in component Veins (mL)"
    legend_algebraic[53] = "F_VE in component Veins (flow)"
    legend_algebraic[57] = "F_VEV in component Veins (flow)"
    legend_algebraic[45] = "conc_X_VE in component Veins (concentration)"
    legend_states[29] = "sub_X_VE in component Veins (mL)"
    legend_constants[63] = "R_VCV in component VenaCava (resistance)"
    legend_states[30] = "V_VC in component VenaCava (mL)"
    legend_constants[64] = "V_VC_0 in component VenaCava (mL)"
    legend_constants[65] = "C_VC in component VenaCava (capacitance)"
    legend_algebraic[58] = "F_VC in component VenaCava (flow)"
    legend_algebraic[64] = "F_VCV in component VenaCava (flow)"
    legend_algebraic[54] = "conc_X_VC in component VenaCava (concentration)"
    legend_states[31] = "sub_X_VC in component VenaCava (mL)"
    legend_algebraic[59] = "beattime in component RATiming (ms)"
    legend_constants[66] = "Esys in component RAElastanceFunction (elastance)"
    legend_constants[67] = "Edia in component RAElastanceFunction (elastance)"
    legend_constants[68] = "TsK in component RAElastanceFunction (second)"
    legend_constants[94] = "Ts in component RAElastanceFunction (second)"
    legend_algebraic[60] = "E_RA in component RAElastanceFunction (elastance)"
    legend_algebraic[70] = "F_RV in component RightVentricle (flow)"
    legend_constants[69] = "P_RA_ext in component RightAtrium (mmHg)"
    legend_states[32] = "V_RA in component RightAtrium (mL)"
    legend_constants[70] = "V_RA_0 in component RightAtrium (mL)"
    legend_algebraic[65] = "F_RA in component RightAtrium (flow)"
    legend_algebraic[61] = "conc_X_RA in component RightAtrium (concentration)"
    legend_states[33] = "sub_X_RA in component RightAtrium (mL)"
    legend_constants[71] = "PRinterval in component RVTiming (second)"
    legend_algebraic[66] = "beattime in component RVTiming (ms)"
    legend_constants[72] = "Esys in component RVElastanceFunction (elastance)"
    legend_constants[73] = "Edia in component RVElastanceFunction (elastance)"
    legend_constants[74] = "TsK in component RVElastanceFunction (second)"
    legend_constants[95] = "Ts in component RVElastanceFunction (second)"
    legend_algebraic[67] = "E_RV in component RVElastanceFunction (elastance)"
    legend_algebraic[73] = "F_PA in component PulmonaryArtery (flow)"
    legend_constants[75] = "P_RV_ext in component RightVentricle (mmHg)"
    legend_states[34] = "V_RV in component RightVentricle (mL)"
    legend_constants[76] = "V_RV_0 in component RightVentricle (mL)"
    legend_algebraic[69] = "P_RV in component RightVentricle (mmHg)"
    legend_constants[77] = "R_RV in component RightVentricle (resistance)"
    legend_algebraic[68] = "conc_X_RV in component RightVentricle (concentration)"
    legend_states[35] = "sub_X_RV in component RightVentricle (mL)"
    legend_algebraic[76] = "F_Pa in component PulmonaryArteries (flow)"
    legend_constants[78] = "R_PA in component PulmonaryArtery (resistance)"
    legend_constants[79] = "C_PA in component PulmonaryArtery (capacitance)"
    legend_constants[80] = "V_PA_0 in component PulmonaryArtery (mL)"
    legend_states[36] = "V_PA in component PulmonaryArtery (mL)"
    legend_algebraic[72] = "P_PA in component PulmonaryArtery (mmHg)"
    legend_algebraic[71] = "conc_X_PA in component PulmonaryArtery (concentration)"
    legend_states[37] = "sub_X_PA in component PulmonaryArtery (mL)"
    legend_algebraic[80] = "F_PC in component PulmonaryCapillaries (flow)"
    legend_algebraic[79] = "F_SH in component Shunt (flow)"
    legend_states[38] = "V_Pa in component PulmonaryArteries (mL)"
    legend_constants[81] = "C_Pa in component PulmonaryArteries (capacitance)"
    legend_constants[82] = "R_Pa in component PulmonaryArteries (resistance)"
    legend_constants[83] = "V_Pa_0 in component PulmonaryArteries (mL)"
    legend_algebraic[75] = "P_Pa in component PulmonaryArteries (mmHg)"
    legend_algebraic[74] = "conc_X_Pa in component PulmonaryArteries (concentration)"
    legend_states[39] = "sub_X_Pa in component PulmonaryArteries (mL)"
    legend_constants[84] = "R_SH in component Shunt (resistance)"
    legend_algebraic[78] = "P_PC in component PulmonaryCapillaries (mmHg)"
    legend_constants[85] = "C_PC in component PulmonaryCapillaries (capacitance)"
    legend_constants[86] = "R_PC in component PulmonaryCapillaries (resistance)"
    legend_constants[87] = "R_PCV in component PulmonaryCapillaries (resistance)"
    legend_constants[88] = "V_PC_0 in component PulmonaryCapillaries (mL)"
    legend_states[40] = "V_PC in component PulmonaryCapillaries (mL)"
    legend_algebraic[84] = "F_PCV in component PulmonaryCapillaries (flow)"
    legend_algebraic[77] = "conc_X_PC in component PulmonaryCapillaries (concentration)"
    legend_states[41] = "sub_X_PC in component PulmonaryCapillaries (mL)"
    legend_constants[89] = "C_PV in component PulmonaryVein (capacitance)"
    legend_constants[90] = "V_PV_0 in component PulmonaryVein (mL)"
    legend_states[42] = "V_PV in component PulmonaryVein (mL)"
    legend_algebraic[85] = "F_PV in component PulmonaryVein (flow)"
    legend_states[43] = "sub_X_PV in component PulmonaryVein (mL)"
    legend_rates[1] = "d/dt sub_X_LA in component LeftAtrium (mL)"
    legend_rates[0] = "d/dt V_LA in component LeftAtrium (mL)"
    legend_rates[3] = "d/dt sub_X_LV in component LeftVentricle (mL)"
    legend_rates[2] = "d/dt V_LV in component LeftVentricle (mL)"
    legend_rates[6] = "d/dt sub_X_AO in component Aorta (mL)"
    legend_rates[5] = "d/dt F_AOalways in component Aorta (flow)"
    legend_rates[4] = "d/dt V_AO in component Aorta (mL)"
    legend_rates[8] = "d/dt sub_X_CO in component CoronaryCirc (mL)"
    legend_rates[7] = "d/dt V_CO in component CoronaryCirc (mL)"
    legend_rates[11] = "d/dt sub_X_CRtissue in component CarotidCirc (mL)"
    legend_rates[10] = "d/dt sub_X_CR in component CarotidCirc (mL)"
    legend_rates[9] = "d/dt V_CR in component CarotidCirc (mL)"
    legend_rates[13] = "d/dt sub_X_AR in component Arteries (mL)"
    legend_rates[12] = "d/dt V_AR in component Arteries (mL)"
    legend_rates[15] = "d/dt sub_X_AD in component Adipose (mL)"
    legend_rates[14] = "d/dt V_AD in component Adipose (mL)"
    legend_rates[17] = "d/dt sub_X_MU in component Muscle (mL)"
    legend_rates[16] = "d/dt V_MU in component Muscle (mL)"
    legend_rates[19] = "d/dt sub_X_GI in component GutIntestine (mL)"
    legend_rates[18] = "d/dt V_GI in component GutIntestine (mL)"
    legend_rates[21] = "d/dt sub_X_LI in component Liver (mL)"
    legend_rates[20] = "d/dt V_LI in component Liver (mL)"
    legend_rates[23] = "d/dt sub_X_KI in component Kidney (mL)"
    legend_rates[22] = "d/dt V_KI in component Kidney (mL)"
    legend_rates[25] = "d/dt sub_X_SK in component Skin (mL)"
    legend_rates[24] = "d/dt V_SK in component Skin (mL)"
    legend_rates[27] = "d/dt sub_X_OT in component OtherTissue (mL)"
    legend_rates[26] = "d/dt V_OT in component OtherTissue (mL)"
    legend_rates[29] = "d/dt sub_X_VE in component Veins (mL)"
    legend_rates[28] = "d/dt V_VE in component Veins (mL)"
    legend_rates[31] = "d/dt sub_X_VC in component VenaCava (mL)"
    legend_rates[30] = "d/dt V_VC in component VenaCava (mL)"
    legend_rates[33] = "d/dt sub_X_RA in component RightAtrium (mL)"
    legend_rates[32] = "d/dt V_RA in component RightAtrium (mL)"
    legend_rates[35] = "d/dt sub_X_RV in component RightVentricle (mL)"
    legend_rates[34] = "d/dt V_RV in component RightVentricle (mL)"
    legend_rates[37] = "d/dt sub_X_PA in component PulmonaryArtery (mL)"
    legend_rates[36] = "d/dt V_PA in component PulmonaryArtery (mL)"
    legend_rates[39] = "d/dt sub_X_Pa in component PulmonaryArteries (mL)"
    legend_rates[38] = "d/dt V_Pa in component PulmonaryArteries (mL)"
    legend_rates[41] = "d/dt sub_X_PC in component PulmonaryCapillaries (mL)"
    legend_rates[40] = "d/dt V_PC in component PulmonaryCapillaries (mL)"
    legend_rates[43] = "d/dt sub_X_PV in component PulmonaryVein (mL)"
    legend_rates[42] = "d/dt V_PV in component PulmonaryVein (mL)"
    return (legend_states, legend_algebraic, legend_voi, legend_constants)

def initConsts():
    constants = [0.0] * sizeConstants; states = [0.0] * sizeStates;
    constants[0] = 77
    constants[1] = 0.145
    constants[2] = 0.131
    constants[3] = 0.2
    constants[4] = 0
    states[0] = 144.038898121062
    constants[5] = 15
    constants[6] = 0.01
    states[1] = 0
    constants[7] = 0.12
    constants[8] = 5.6
    constants[9] = 0.19
    constants[10] = 0.35
    constants[11] = 0
    states[2] = 118.75315693072
    constants[12] = 30
    constants[13] = 0.001
    states[3] = 0
    constants[14] = 1.2
    constants[15] = 5e-5
    constants[16] = 0.01
    constants[17] = 200
    states[4] = 302.714548324312
    states[5] = -6898.97886842467
    states[6] = 0
    constants[18] = 1
    constants[19] = 0.2
    constants[20] = 10
    constants[21] = 0.1
    states[7] = 11.4012164759732
    states[8] = 0
    constants[22] = 5
    constants[23] = 0.02
    constants[24] = 10
    constants[25] = 1
    states[9] = 27.2147903053333
    states[10] = 0
    constants[26] = 1
    constants[27] = 1
    states[11] = 0
    constants[28] = 1.5
    constants[29] = 0.5
    constants[30] = 370
    states[12] = 383.557125288364
    states[13] = 0
    constants[31] = 1.81
    constants[32] = 1
    constants[33] = 0.33
    constants[34] = 75
    states[14] = 111.639088723426
    states[15] = 0
    constants[35] = 1.81
    constants[36] = 1
    constants[37] = 0.33
    constants[38] = 75
    states[16] = 111.639088723426
    states[17] = 0
    constants[39] = 1.81
    constants[40] = 1
    constants[41] = 0.33
    constants[42] = 75
    states[18] = 111.639088723426
    states[19] = 0
    constants[43] = 1.81
    constants[44] = 1
    constants[45] = 0.33
    constants[46] = 75
    states[20] = 111.639088723426
    states[21] = 0
    constants[47] = 0
    constants[48] = 1.81
    constants[49] = 1
    constants[50] = 75
    constants[51] = 0.33
    states[22] = 111.639088723426
    states[23] = 0
    constants[52] = 1.81
    constants[53] = 1
    constants[54] = 0.33
    constants[55] = 75
    states[24] = 111.639088723426
    states[25] = 0
    constants[56] = 1.81
    constants[57] = 1
    constants[58] = 0.33
    constants[59] = 75
    states[26] = 111.639088723426
    states[27] = 0
    constants[60] = 0.0223
    states[28] = 833.291684587657
    constants[61] = 13.24
    constants[62] = 596
    states[29] = 5
    constants[63] = 0.0267
    states[30] = 1849.93788284599
    constants[64] = 600
    constants[65] = 73.88
    states[31] = 0
    constants[66] = 0.11
    constants[67] = 0.099
    constants[68] = 0.2
    constants[69] = 0
    states[32] = 152.150765432928
    constants[70] = 15
    states[33] = 0
    constants[71] = 0.12
    constants[72] = 0.67
    constants[73] = 0.104
    constants[74] = 0.35
    constants[75] = 0
    states[34] = 154.337916417774
    constants[76] = 25
    constants[77] = 0.002
    states[35] = 0
    constants[78] = 0.0227
    constants[79] = 2.222
    constants[80] = 50
    states[36] = 112.059036699839
    states[37] = 0
    states[38] = 66.9910537203609
    constants[81] = 1.481
    constants[82] = 0.053
    constants[83] = 30
    states[39] = 0
    constants[84] = 2
    constants[85] = 1.778
    constants[86] = 0.0379
    constants[87] = 0.0379
    constants[88] = 53
    states[40] = 91.7784213531534
    states[41] = 0
    constants[89] = 5
    constants[90] = 150
    states[42] = 239.25897115602
    states[43] = 0
    constants[91] = constants[0]/0.0600000
    constants[92] = constants[3]*(power(constants[91]/1000.00, 1.0/2))
    constants[93] = constants[10]*(power(constants[91]/1000.00, 1.0/2))
    constants[94] = constants[68]*(power(constants[91]/1000.00, 1.0/2))
    constants[95] = constants[74]*(power(constants[91]/1000.00, 1.0/2))
    return (states, constants)

def computeRates(voi, states, constants):
    rates = [0.0] * sizeStates; algebraic = [0.0] * sizeAlgebraic
    algebraic[4] = (voi-floor(voi/constants[91])*constants[91])/1000.00-constants[7]
    algebraic[5] = custom_piecewise([greater_equal(algebraic[4]*0.00100000 , 0.00000) & less_equal(algebraic[4]*0.00100000 , constants[93]), constants[9]+((constants[8]-constants[9])*(1.00000-cos(( pi*(algebraic[4]*0.00100000))/constants[93])))/2.00000 , less(algebraic[4]*0.00100000 , 1.50000*constants[93]) & greater_equal(algebraic[4]*0.00100000 , constants[93]), constants[9]+((constants[8]-constants[9])*(1.00000+cos((2.00000* pi*(algebraic[4]*0.00100000-constants[93]))/constants[93])))/2.00000 , True, constants[9]])
    algebraic[7] = algebraic[5]*(states[2]-constants[12])+constants[11]
    algebraic[11] = (1.00000/constants[14])*(states[4]-constants[17])
    rates[5] = ((1.00000/constants[15])*(algebraic[7]-algebraic[11])-(constants[16]/constants[15])*states[5])/1000.00
    algebraic[0] = (voi-floor(voi/constants[91])*constants[91])/1000.00
    algebraic[1] = custom_piecewise([greater_equal(algebraic[0] , 0.00000) & less_equal(algebraic[0] , constants[92]), constants[2]+((constants[1]-constants[2])*(1.00000-cos(( pi*algebraic[0])/constants[92])))/2.00000 , less(algebraic[0] , 1.50000*constants[92]) & greater_equal(algebraic[0] , constants[92]), constants[2]+((constants[1]-constants[2])*(1.00000+cos((2.00000* pi*(algebraic[0]-constants[92]))/constants[92])))/2.00000 , True, constants[2]])
    algebraic[3] = algebraic[1]*(states[0]-constants[5])+constants[4]
    algebraic[9] = custom_piecewise([greater(algebraic[3] , algebraic[7]), (algebraic[3]-algebraic[7])/constants[13] , True, 0.00000])
    algebraic[2] = states[1]/states[0]
    algebraic[12] = custom_piecewise([greater(algebraic[7] , algebraic[11]), states[5] , True, 0.00000])
    algebraic[14] = (1.00000/constants[21])*(states[7]-constants[20])
    algebraic[15] = custom_piecewise([greater(algebraic[7] , algebraic[14]), (algebraic[7]-algebraic[14])/constants[18] , True, 0.00000])
    algebraic[6] = states[3]/states[2]
    rates[3] = (algebraic[2]*algebraic[9]-algebraic[6]*(algebraic[12]+algebraic[15]))/1000.00
    rates[2] = ((algebraic[9]-algebraic[12])-algebraic[15])/1000.00
    algebraic[17] = states[10]/states[9]
    algebraic[16] = states[11]/constants[27]
    rates[11] = ((algebraic[17]-algebraic[16])*constants[26])/1000.00
    algebraic[18] = (1.00000/constants[25])*(states[9]-constants[24])
    algebraic[19] = (algebraic[11]-algebraic[18])/constants[22]
    algebraic[21] = (1.00000/constants[29])*(states[12]-constants[30])
    algebraic[22] = (algebraic[11]-algebraic[21])/constants[28]
    algebraic[10] = states[6]/states[4]
    rates[6] = (algebraic[6]*algebraic[12]-algebraic[10]*(algebraic[19]+algebraic[22]))/1000.00
    rates[4] = ((algebraic[12]-algebraic[22])-algebraic[19])/1000.00
    algebraic[30] = (1.00000/constants[39])*(states[18]-constants[42])
    algebraic[31] = (algebraic[21]-algebraic[30])/constants[40]
    algebraic[20] = states[13]/states[12]
    algebraic[33] = (1.00000/constants[43])*(states[20]-constants[46])
    algebraic[34] = (algebraic[30]-algebraic[33])/constants[41]
    algebraic[29] = states[19]/states[18]
    rates[19] = ((algebraic[20]*algebraic[31]-algebraic[29]*algebraic[34])/1000.00)*0.00100000
    rates[18] = ((algebraic[31]-algebraic[34])/1000.00)*0.00100000
    algebraic[24] = (1.00000/constants[31])*(states[14]-constants[34])
    algebraic[25] = (algebraic[21]-algebraic[24])/constants[32]
    algebraic[27] = (1.00000/constants[35])*(states[16]-constants[38])
    algebraic[28] = (algebraic[21]-algebraic[27])/constants[36]
    algebraic[35] = (algebraic[21]-algebraic[33])/constants[44]
    algebraic[37] = (1.00000/constants[48])*(states[22]-constants[50])
    algebraic[38] = (algebraic[21]-algebraic[37])/constants[49]
    algebraic[43] = (1.00000/constants[56])*(states[26]-constants[59])
    algebraic[44] = (algebraic[21]-algebraic[43])/constants[57]
    algebraic[40] = (1.00000/constants[52])*(states[24]-constants[55])
    algebraic[41] = (algebraic[21]-algebraic[40])/constants[53]
    rates[13] = (algebraic[10]*algebraic[22]-algebraic[20]*(algebraic[41]+algebraic[25]+algebraic[28]+algebraic[31]+algebraic[35]+algebraic[38]+algebraic[44]))/1000.00
    rates[12] = (((((((algebraic[22]-algebraic[25])-algebraic[28])-algebraic[31])-algebraic[35])-algebraic[38])-algebraic[44])-algebraic[41])/1000.00
    algebraic[46] = (1.00000/constants[61])*(states[28]-constants[62])
    algebraic[47] = (algebraic[24]-algebraic[46])/constants[33]
    algebraic[23] = states[15]/states[14]
    rates[15] = ((algebraic[20]*algebraic[25]-algebraic[23]*algebraic[47])/1000.00)*0.00100000
    rates[14] = ((algebraic[25]-algebraic[47])/1000.00)*0.00100000
    algebraic[48] = (algebraic[27]-algebraic[46])/constants[37]
    algebraic[26] = states[17]/states[16]
    rates[17] = ((algebraic[20]*algebraic[28]-algebraic[26]*algebraic[48])/1000.00)*0.00100000
    rates[16] = ((algebraic[28]-algebraic[48])/1000.00)*0.00100000
    algebraic[49] = (algebraic[33]-algebraic[46])/constants[45]
    algebraic[32] = states[21]/states[20]
    rates[21] = ((((algebraic[29]*algebraic[34]+algebraic[20]*algebraic[35])-algebraic[32]*algebraic[49])+algebraic[32]*constants[47])/1000.00)*0.00100000
    rates[20] = (((algebraic[34]+algebraic[35])-algebraic[49])/1000.00)*0.00100000
    algebraic[50] = (algebraic[37]-algebraic[46])/constants[51]
    algebraic[36] = states[23]/states[22]
    rates[23] = ((algebraic[20]*algebraic[38]-algebraic[36]*algebraic[50])/1000.00)*0.00100000
    rates[22] = ((algebraic[38]-algebraic[50])/1000.00)*0.00100000
    algebraic[51] = (algebraic[40]-algebraic[46])/constants[54]
    algebraic[39] = states[25]/states[24]
    rates[25] = ((algebraic[20]*algebraic[41]-algebraic[39]*algebraic[51])/1000.00)*0.00100000
    rates[24] = ((algebraic[41]-algebraic[51])/1000.00)*0.00100000
    algebraic[52] = (algebraic[43]-algebraic[46])/constants[58]
    algebraic[42] = states[27]/states[26]
    rates[27] = ((algebraic[20]*algebraic[44]-algebraic[42]*algebraic[52])/1000.00)*0.00100000
    rates[26] = ((algebraic[44]-algebraic[52])/1000.00)*0.00100000
    algebraic[55] = (1.00000/constants[65])*(states[30]-constants[64])
    algebraic[56] = (algebraic[18]-algebraic[55])/constants[23]
    rates[10] = ((algebraic[10]*algebraic[19]-algebraic[17]*algebraic[56])-(algebraic[17]-algebraic[16])*constants[26])/1000.00
    rates[9] = (algebraic[19]-algebraic[56])/1000.00
    algebraic[57] = (algebraic[46]-algebraic[55])/constants[60]
    algebraic[45] = states[29]/states[28]
    rates[29] = ((algebraic[39]*algebraic[51]+algebraic[23]*algebraic[47]+algebraic[26]*algebraic[48]+algebraic[32]*algebraic[49]+algebraic[36]*algebraic[50]+algebraic[42]*algebraic[52])-algebraic[45]*algebraic[57])/1000.00
    algebraic[53] = algebraic[51]+algebraic[47]+algebraic[48]+algebraic[49]+algebraic[50]+algebraic[52]
    rates[28] = (algebraic[53]-algebraic[57])/1000.00
    algebraic[59] = (voi-floor(voi/constants[91])*constants[91])/1000.00
    algebraic[60] = custom_piecewise([greater_equal(algebraic[59]*0.00100000 , 0.00000) & less_equal(algebraic[59]*0.00100000 , constants[94]), constants[67]+((constants[66]-constants[67])*(1.00000-cos(( pi*(algebraic[59]*0.00100000))/constants[94])))/2.00000 , less(algebraic[59]*0.00100000 , 1.50000*constants[94]) & greater_equal(algebraic[59]*0.00100000 , constants[94]), constants[67]+((constants[66]-constants[67])*(1.00000+cos((2.00000* pi*(algebraic[59]*0.00100000-constants[94]))/constants[94])))/2.00000 , True, constants[67]])
    algebraic[62] = algebraic[60]*(states[32]-constants[70])+constants[69]
    algebraic[63] = (algebraic[14]-algebraic[62])/constants[19]
    algebraic[13] = states[8]/states[7]
    rates[8] = (algebraic[6]*algebraic[15]-algebraic[13]*algebraic[63])/1000.00
    rates[7] = (algebraic[15]-algebraic[63])/1000.00
    algebraic[64] = (algebraic[55]-algebraic[62])/constants[63]
    algebraic[54] = states[31]/states[30]
    rates[31] = ((algebraic[45]*algebraic[57]+algebraic[17]*algebraic[56])-algebraic[54]*algebraic[64])/1000.00
    algebraic[58] = algebraic[56]+algebraic[57]
    rates[30] = (algebraic[58]-algebraic[64])/1000.00
    algebraic[66] = (voi-floor(voi/constants[91])*constants[91])/1000.00-constants[71]
    algebraic[67] = custom_piecewise([greater_equal(algebraic[66]*0.00100000 , 0.00000) & less_equal(algebraic[66]*0.00100000 , constants[95]), constants[73]+((constants[72]-constants[73])*(1.00000-cos(( pi*(algebraic[66]*0.00100000))/constants[95])))/2.00000 , less(algebraic[66]*0.00100000 , 1.50000*constants[95]) & greater_equal(algebraic[66]*0.00100000 , constants[95]), constants[73]+((constants[72]-constants[73])*(1.00000+cos((2.00000* pi*(algebraic[66]*0.00100000-constants[95]))/constants[95])))/2.00000 , True, constants[73]])
    algebraic[69] = algebraic[67]*(states[34]-constants[76])+constants[75]
    algebraic[70] = custom_piecewise([greater(algebraic[62] , algebraic[69]), (algebraic[62]-algebraic[69])/constants[77] , True, 0.00000])
    algebraic[61] = states[33]/states[32]
    rates[33] = ((algebraic[54]*algebraic[64]+algebraic[13]*algebraic[63])-algebraic[61]*algebraic[70])/1000.00
    algebraic[65] = algebraic[64]+algebraic[63]
    rates[32] = (algebraic[65]-algebraic[70])/1000.00
    algebraic[72] = (1.00000/constants[79])*(states[36]-constants[80])
    algebraic[73] = custom_piecewise([greater(algebraic[69] , algebraic[72]), (algebraic[69]-algebraic[72])/constants[78] , True, 0.00000])
    algebraic[68] = states[35]/states[34]
    rates[35] = (algebraic[61]*algebraic[70]-algebraic[68]*algebraic[73])/1000.00
    rates[34] = (algebraic[70]-algebraic[73])/1000.00
    algebraic[75] = (1.00000/constants[81])*(states[38]-constants[83])
    algebraic[76] = (algebraic[72]-algebraic[75])/constants[82]
    algebraic[71] = states[37]/states[36]
    rates[37] = (algebraic[68]*algebraic[73]-algebraic[71]*algebraic[76])/1000.00
    rates[36] = (algebraic[73]-algebraic[76])/1000.00
    algebraic[78] = (1.00000/constants[85])*(states[40]-constants[88])
    algebraic[80] = (algebraic[75]-algebraic[78])/constants[86]
    algebraic[79] = (algebraic[75]-algebraic[78])/constants[84]
    algebraic[74] = states[39]/states[38]
    rates[39] = (algebraic[71]*algebraic[76]-algebraic[74]*(algebraic[80]+algebraic[79]))/1000.00
    rates[38] = ((algebraic[76]-algebraic[80])-algebraic[79])/1000.00
    algebraic[82] = (1.00000/constants[89])*(states[42]-constants[90])
    algebraic[83] = (algebraic[82]-algebraic[3])/constants[6]
    algebraic[81] = states[43]/states[42]
    rates[1] = (algebraic[81]*algebraic[83]-algebraic[2]*algebraic[9])/1000.00
    rates[0] = (algebraic[83]-algebraic[9])/1000.00
    algebraic[84] = (algebraic[78]-algebraic[82])/constants[87]
    algebraic[77] = states[41]/states[40]
    rates[41] = (algebraic[74]*algebraic[80]-algebraic[77]*algebraic[84])/1000.00
    rates[40] = (algebraic[80]-algebraic[84])/1000.00
    rates[43] = ((algebraic[77]*algebraic[84]+algebraic[74]*algebraic[79])-algebraic[81]*algebraic[83])/1000.00
    algebraic[85] = algebraic[84]+algebraic[79]
    rates[42] = (algebraic[85]-algebraic[83])/1000.00
    return(rates)

def computeAlgebraic(constants, states, voi):
    algebraic = array([[0.0] * len(voi)] * sizeAlgebraic)
    states = array(states)
    voi = array(voi)
    algebraic[4] = (voi-floor(voi/constants[91])*constants[91])/1000.00-constants[7]
    algebraic[5] = custom_piecewise([greater_equal(algebraic[4]*0.00100000 , 0.00000) & less_equal(algebraic[4]*0.00100000 , constants[93]), constants[9]+((constants[8]-constants[9])*(1.00000-cos(( pi*(algebraic[4]*0.00100000))/constants[93])))/2.00000 , less(algebraic[4]*0.00100000 , 1.50000*constants[93]) & greater_equal(algebraic[4]*0.00100000 , constants[93]), constants[9]+((constants[8]-constants[9])*(1.00000+cos((2.00000* pi*(algebraic[4]*0.00100000-constants[93]))/constants[93])))/2.00000 , True, constants[9]])
    algebraic[7] = algebraic[5]*(states[2]-constants[12])+constants[11]
    algebraic[11] = (1.00000/constants[14])*(states[4]-constants[17])
    algebraic[0] = (voi-floor(voi/constants[91])*constants[91])/1000.00
    algebraic[1] = custom_piecewise([greater_equal(algebraic[0] , 0.00000) & less_equal(algebraic[0] , constants[92]), constants[2]+((constants[1]-constants[2])*(1.00000-cos(( pi*algebraic[0])/constants[92])))/2.00000 , less(algebraic[0] , 1.50000*constants[92]) & greater_equal(algebraic[0] , constants[92]), constants[2]+((constants[1]-constants[2])*(1.00000+cos((2.00000* pi*(algebraic[0]-constants[92]))/constants[92])))/2.00000 , True, constants[2]])
    algebraic[3] = algebraic[1]*(states[0]-constants[5])+constants[4]
    algebraic[9] = custom_piecewise([greater(algebraic[3] , algebraic[7]), (algebraic[3]-algebraic[7])/constants[13] , True, 0.00000])
    algebraic[2] = states[1]/states[0]
    algebraic[12] = custom_piecewise([greater(algebraic[7] , algebraic[11]), states[5] , True, 0.00000])
    algebraic[14] = (1.00000/constants[21])*(states[7]-constants[20])
    algebraic[15] = custom_piecewise([greater(algebraic[7] , algebraic[14]), (algebraic[7]-algebraic[14])/constants[18] , True, 0.00000])
    algebraic[6] = states[3]/states[2]
    algebraic[17] = states[10]/states[9]
    algebraic[16] = states[11]/constants[27]
    algebraic[18] = (1.00000/constants[25])*(states[9]-constants[24])
    algebraic[19] = (algebraic[11]-algebraic[18])/constants[22]
    algebraic[21] = (1.00000/constants[29])*(states[12]-constants[30])
    algebraic[22] = (algebraic[11]-algebraic[21])/constants[28]
    algebraic[10] = states[6]/states[4]
    algebraic[30] = (1.00000/constants[39])*(states[18]-constants[42])
    algebraic[31] = (algebraic[21]-algebraic[30])/constants[40]
    algebraic[20] = states[13]/states[12]
    algebraic[33] = (1.00000/constants[43])*(states[20]-constants[46])
    algebraic[34] = (algebraic[30]-algebraic[33])/constants[41]
    algebraic[29] = states[19]/states[18]
    algebraic[24] = (1.00000/constants[31])*(states[14]-constants[34])
    algebraic[25] = (algebraic[21]-algebraic[24])/constants[32]
    algebraic[27] = (1.00000/constants[35])*(states[16]-constants[38])
    algebraic[28] = (algebraic[21]-algebraic[27])/constants[36]
    algebraic[35] = (algebraic[21]-algebraic[33])/constants[44]
    algebraic[37] = (1.00000/constants[48])*(states[22]-constants[50])
    algebraic[38] = (algebraic[21]-algebraic[37])/constants[49]
    algebraic[43] = (1.00000/constants[56])*(states[26]-constants[59])
    algebraic[44] = (algebraic[21]-algebraic[43])/constants[57]
    algebraic[40] = (1.00000/constants[52])*(states[24]-constants[55])
    algebraic[41] = (algebraic[21]-algebraic[40])/constants[53]
    algebraic[46] = (1.00000/constants[61])*(states[28]-constants[62])
    algebraic[47] = (algebraic[24]-algebraic[46])/constants[33]
    algebraic[23] = states[15]/states[14]
    algebraic[48] = (algebraic[27]-algebraic[46])/constants[37]
    algebraic[26] = states[17]/states[16]
    algebraic[49] = (algebraic[33]-algebraic[46])/constants[45]
    algebraic[32] = states[21]/states[20]
    algebraic[50] = (algebraic[37]-algebraic[46])/constants[51]
    algebraic[36] = states[23]/states[22]
    algebraic[51] = (algebraic[40]-algebraic[46])/constants[54]
    algebraic[39] = states[25]/states[24]
    algebraic[52] = (algebraic[43]-algebraic[46])/constants[58]
    algebraic[42] = states[27]/states[26]
    algebraic[55] = (1.00000/constants[65])*(states[30]-constants[64])
    algebraic[56] = (algebraic[18]-algebraic[55])/constants[23]
    algebraic[57] = (algebraic[46]-algebraic[55])/constants[60]
    algebraic[45] = states[29]/states[28]
    algebraic[53] = algebraic[51]+algebraic[47]+algebraic[48]+algebraic[49]+algebraic[50]+algebraic[52]
    algebraic[59] = (voi-floor(voi/constants[91])*constants[91])/1000.00
    algebraic[60] = custom_piecewise([greater_equal(algebraic[59]*0.00100000 , 0.00000) & less_equal(algebraic[59]*0.00100000 , constants[94]), constants[67]+((constants[66]-constants[67])*(1.00000-cos(( pi*(algebraic[59]*0.00100000))/constants[94])))/2.00000 , less(algebraic[59]*0.00100000 , 1.50000*constants[94]) & greater_equal(algebraic[59]*0.00100000 , constants[94]), constants[67]+((constants[66]-constants[67])*(1.00000+cos((2.00000* pi*(algebraic[59]*0.00100000-constants[94]))/constants[94])))/2.00000 , True, constants[67]])
    algebraic[62] = algebraic[60]*(states[32]-constants[70])+constants[69]
    algebraic[63] = (algebraic[14]-algebraic[62])/constants[19]
    algebraic[13] = states[8]/states[7]
    algebraic[64] = (algebraic[55]-algebraic[62])/constants[63]
    algebraic[54] = states[31]/states[30]
    algebraic[58] = algebraic[56]+algebraic[57]
    algebraic[66] = (voi-floor(voi/constants[91])*constants[91])/1000.00-constants[71]
    algebraic[67] = custom_piecewise([greater_equal(algebraic[66]*0.00100000 , 0.00000) & less_equal(algebraic[66]*0.00100000 , constants[95]), constants[73]+((constants[72]-constants[73])*(1.00000-cos(( pi*(algebraic[66]*0.00100000))/constants[95])))/2.00000 , less(algebraic[66]*0.00100000 , 1.50000*constants[95]) & greater_equal(algebraic[66]*0.00100000 , constants[95]), constants[73]+((constants[72]-constants[73])*(1.00000+cos((2.00000* pi*(algebraic[66]*0.00100000-constants[95]))/constants[95])))/2.00000 , True, constants[73]])
    algebraic[69] = algebraic[67]*(states[34]-constants[76])+constants[75]
    algebraic[70] = custom_piecewise([greater(algebraic[62] , algebraic[69]), (algebraic[62]-algebraic[69])/constants[77] , True, 0.00000])
    algebraic[61] = states[33]/states[32]
    algebraic[65] = algebraic[64]+algebraic[63]
    algebraic[72] = (1.00000/constants[79])*(states[36]-constants[80])
    algebraic[73] = custom_piecewise([greater(algebraic[69] , algebraic[72]), (algebraic[69]-algebraic[72])/constants[78] , True, 0.00000])
    algebraic[68] = states[35]/states[34]
    algebraic[75] = (1.00000/constants[81])*(states[38]-constants[83])
    algebraic[76] = (algebraic[72]-algebraic[75])/constants[82]
    algebraic[71] = states[37]/states[36]
    algebraic[78] = (1.00000/constants[85])*(states[40]-constants[88])
    algebraic[80] = (algebraic[75]-algebraic[78])/constants[86]
    algebraic[79] = (algebraic[75]-algebraic[78])/constants[84]
    algebraic[74] = states[39]/states[38]
    algebraic[82] = (1.00000/constants[89])*(states[42]-constants[90])
    algebraic[83] = (algebraic[82]-algebraic[3])/constants[6]
    algebraic[81] = states[43]/states[42]
    algebraic[84] = (algebraic[78]-algebraic[82])/constants[87]
    algebraic[77] = states[41]/states[40]
    algebraic[85] = algebraic[84]+algebraic[79]
    algebraic[8] = algebraic[7]*(1.00000+(3.00000*states[2])/200.000)
    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)