/* There are a total of 21 entries in the algebraic variable array. There are a total of 6 entries in each of the rate and state variable arrays. There are a total of 40 entries in the constant variable array. */ /* * VOI is time in component environment (second). * CONSTANTS[0] is R_mt in component heart_parameters (kPa_second_per_liter). * CONSTANTS[1] is R_av in component heart_parameters (kPa_second_per_liter). * CONSTANTS[2] is R_tc in component heart_parameters (kPa_second_per_liter). * CONSTANTS[3] is R_pv in component heart_parameters (kPa_second_per_liter). * CONSTANTS[4] is R_pul in component heart_parameters (kPa_second_per_liter). * CONSTANTS[5] is R_sys in component heart_parameters (kPa_second_per_liter). * CONSTANTS[6] is HR in component heart_parameters (dimensionless). * CONSTANTS[7] is V_tot in component heart_parameters (liter). * CONSTANTS[8] is P_pl in component heart_parameters (kPa). * ALGEBRAIC[1] is e_t in component driver_function (dimensionless). * CONSTANTS[9] is A in component driver_function (dimensionless). * CONSTANTS[10] is B in component driver_function (dimensionless). * CONSTANTS[11] is C in component driver_function (dimensionless). * ALGEBRAIC[0] is tau in component driver_function (second). * CONSTANTS[12] is period in component driver_function (dimensionless). * ALGEBRAIC[2] is V_pcd in component pericardium (liter). * ALGEBRAIC[3] is P_pcd in component pericardium (kPa). * ALGEBRAIC[4] is P_peri in component pericardium (kPa). * STATES[0] is V_lv in component left_ventricle (liter). * STATES[1] is V_rv in component right_ventricle (liter). * CONSTANTS[13] is P_0_pcd in component pericardium (kPa). * CONSTANTS[14] is V_0_pcd in component pericardium (liter). * CONSTANTS[15] is lambda_pcd in component pericardium (per_liter). * ALGEBRAIC[9] is V_lvf in component left_ventricle (liter). * ALGEBRAIC[10] is P_lvf in component left_ventricle (kPa). * ALGEBRAIC[18] is P_lv in component left_ventricle (kPa). * ALGEBRAIC[11] is V_spt in component septum (liter). * ALGEBRAIC[12] is P_es_lvf in component lvf_calculator (kPa). * ALGEBRAIC[13] is P_ed_lvf in component lvf_calculator (kPa). * ALGEBRAIC[6] is P_pu in component pulmonary_vein (kPa). * ALGEBRAIC[7] is P_ao in component aorta (kPa). * CONSTANTS[16] is E_es_lvf in component lvf_calculator (kPa_per_liter). * CONSTANTS[17] is V_d_lvf in component lvf_calculator (liter). * CONSTANTS[18] is P_0_lvf in component lvf_calculator (kPa). * CONSTANTS[19] is lambda_lvf in component lvf_calculator (per_liter). * CONSTANTS[20] is V_0_lvf in component lvf_calculator (liter). * ALGEBRAIC[14] is V_rvf in component right_ventricle (liter). * ALGEBRAIC[15] is P_rvf in component right_ventricle (kPa). * ALGEBRAIC[19] is P_rv in component right_ventricle (kPa). * ALGEBRAIC[16] is P_es_rvf in component rvf_calculator (kPa). * ALGEBRAIC[17] is P_ed_rvf in component rvf_calculator (kPa). * ALGEBRAIC[5] is P_pa in component pulmonary_artery (kPa). * ALGEBRAIC[8] is P_vc in component vena_cava (kPa). * CONSTANTS[21] is E_es_rvf in component rvf_calculator (kPa_per_liter). * CONSTANTS[22] is V_d_rvf in component rvf_calculator (liter). * CONSTANTS[23] is P_0_rvf in component rvf_calculator (kPa). * CONSTANTS[24] is lambda_rvf in component rvf_calculator (per_liter). * CONSTANTS[25] is V_0_rvf in component rvf_calculator (liter). * ALGEBRAIC[20] is P_sept in component septum (kPa). * CONSTANTS[26] is E_es_spt in component septum (kPa_per_liter). * CONSTANTS[27] is V_d_spt in component septum (liter). * CONSTANTS[28] is P_0_spt in component septum (kPa). * CONSTANTS[29] is lambda_spt in component septum (per_liter). * CONSTANTS[30] is V_0_spt in component septum (liter). * CONSTANTS[31] is one in component septum (dimensionless). * CONSTANTS[32] is E_es_pa in component pulmonary_artery (kPa_per_liter). * STATES[2] is V_pa in component pulmonary_artery (liter). * CONSTANTS[33] is V_d_pa in component pulmonary_artery (liter). * CONSTANTS[34] is E_es_pu in component pulmonary_vein (kPa_per_liter). * STATES[3] is V_pu in component pulmonary_vein (liter). * CONSTANTS[35] is V_d_pu in component pulmonary_vein (liter). * CONSTANTS[36] is E_es_ao in component aorta (kPa_per_liter). * STATES[4] is V_ao in component aorta (liter). * CONSTANTS[37] is V_d_ao in component aorta (liter). * CONSTANTS[38] is E_es_vc in component vena_cava (kPa_per_liter). * STATES[5] is V_vc in component vena_cava (liter). * CONSTANTS[39] is V_d_vc in component vena_cava (liter). * RATES[0] is d/dt V_lv in component left_ventricle (liter). * RATES[1] is d/dt V_rv in component right_ventricle (liter). * RATES[2] is d/dt V_pa in component pulmonary_artery (liter). * RATES[3] is d/dt V_pu in component pulmonary_vein (liter). * RATES[4] is d/dt V_ao in component aorta (liter). * RATES[5] is d/dt V_vc in component vena_cava (liter). */ void initConsts(double* CONSTANTS, double* RATES, double *STATES) { CONSTANTS[0] = 0.06; CONSTANTS[1] = 1.4; CONSTANTS[2] = 0.18; CONSTANTS[3] = 0.48; CONSTANTS[4] = 19; CONSTANTS[5] = 140; CONSTANTS[6] = 80; CONSTANTS[7] = 5.5; CONSTANTS[8] = -0.533289474; CONSTANTS[9] = 1; CONSTANTS[10] = 80; CONSTANTS[11] = 0.27; CONSTANTS[12] = 0.405; STATES[0] = 0.005; STATES[1] = 0.005; CONSTANTS[13] = 0.067; CONSTANTS[14] = 0.2; CONSTANTS[15] = 30; CONSTANTS[16] = 454; CONSTANTS[17] = 0.005; CONSTANTS[18] = 0.17; CONSTANTS[19] = 15; CONSTANTS[20] = 0.005; CONSTANTS[21] = 87; CONSTANTS[22] = 0.005; CONSTANTS[23] = 0.16; CONSTANTS[24] = 15; CONSTANTS[25] = 0.005; CONSTANTS[26] = 6500; CONSTANTS[27] = 0.002; CONSTANTS[28] = 0.148; CONSTANTS[29] = 435; CONSTANTS[30] = 0.002; CONSTANTS[31] = 1; CONSTANTS[32] = 45; STATES[2] = 0.16; CONSTANTS[33] = 0.16; CONSTANTS[34] = 0.8; STATES[3] = 0.2; CONSTANTS[35] = 0.2; CONSTANTS[36] = 94; STATES[4] = 0.8; CONSTANTS[37] = 0.8; CONSTANTS[38] = 1.5; STATES[5] = 2.83; CONSTANTS[39] = 2.83; } void computeRates(double VOI, double* CONSTANTS, double* RATES, double* STATES, double* ALGEBRAIC) { ALGEBRAIC[2] = STATES[0]+STATES[1]; ALGEBRAIC[3] = CONSTANTS[13]*(exp( CONSTANTS[15]*(ALGEBRAIC[2] - CONSTANTS[14])) - 1.00000); ALGEBRAIC[4] = ALGEBRAIC[3]+CONSTANTS[8]; ALGEBRAIC[0] = (int)(VOI) % (int)(CONSTANTS[12]); ALGEBRAIC[1] = CONSTANTS[9]*exp( - CONSTANTS[10]*pow(( ALGEBRAIC[0]*CONSTANTS[6])/60.0000 - CONSTANTS[11], 2.00000)); rootfind_0(VOI, CONSTANTS, RATES, STATES, ALGEBRAIC, pret); ALGEBRAIC[18] = ALGEBRAIC[10]+ALGEBRAIC[4]; ALGEBRAIC[6] = CONSTANTS[34]*(STATES[3] - CONSTANTS[35]); ALGEBRAIC[7] = CONSTANTS[36]*(STATES[4] - CONSTANTS[37]); RATES[0] = (ALGEBRAIC[6] - ALGEBRAIC[18]<0.00000&&ALGEBRAIC[18] - ALGEBRAIC[7]<0.00000 ? 0.00000 : ALGEBRAIC[6] - ALGEBRAIC[18]<0.00000 ? - (ALGEBRAIC[18] - ALGEBRAIC[7])/CONSTANTS[1] : ALGEBRAIC[18] - ALGEBRAIC[7]<0.00000 ? (ALGEBRAIC[6] - ALGEBRAIC[18])/CONSTANTS[0] : (ALGEBRAIC[6] - ALGEBRAIC[18])/CONSTANTS[0] - (ALGEBRAIC[18] - ALGEBRAIC[7])/CONSTANTS[1]); ALGEBRAIC[5] = CONSTANTS[32]*(STATES[2] - CONSTANTS[33]); RATES[3] = (ALGEBRAIC[5] - ALGEBRAIC[6]<0.00000&&ALGEBRAIC[6] - ALGEBRAIC[18]<0.00000 ? 0.00000 : ALGEBRAIC[5] - ALGEBRAIC[6]<0.00000 ? - (ALGEBRAIC[6] - ALGEBRAIC[18])/CONSTANTS[0] : ALGEBRAIC[6] - ALGEBRAIC[18]<0.00000 ? (ALGEBRAIC[5] - ALGEBRAIC[6])/CONSTANTS[4] : (ALGEBRAIC[5] - ALGEBRAIC[6])/CONSTANTS[4] - (ALGEBRAIC[6] - ALGEBRAIC[18])/CONSTANTS[0]); ALGEBRAIC[8] = CONSTANTS[38]*(STATES[5] - CONSTANTS[39]); RATES[4] = (ALGEBRAIC[18] - ALGEBRAIC[7]<0.00000&&ALGEBRAIC[7] - ALGEBRAIC[8]<0.00000 ? 0.00000 : ALGEBRAIC[18] - ALGEBRAIC[7]<0.00000 ? - (ALGEBRAIC[7] - ALGEBRAIC[8])/CONSTANTS[5] : ALGEBRAIC[7] - ALGEBRAIC[8]<0.00000 ? (ALGEBRAIC[18] - ALGEBRAIC[7])/CONSTANTS[1] : (ALGEBRAIC[18] - ALGEBRAIC[7])/CONSTANTS[1] - (ALGEBRAIC[7] - ALGEBRAIC[8])/CONSTANTS[5]); ALGEBRAIC[19] = ALGEBRAIC[15]+ALGEBRAIC[4]; RATES[1] = (ALGEBRAIC[8] - ALGEBRAIC[19]<0.00000&&ALGEBRAIC[19] - ALGEBRAIC[5]<0.00000 ? 0.00000 : ALGEBRAIC[8] - ALGEBRAIC[19]<0.00000 ? - (ALGEBRAIC[19] - ALGEBRAIC[5])/CONSTANTS[3] : ALGEBRAIC[19] - ALGEBRAIC[5]<0.00000 ? (ALGEBRAIC[8] - ALGEBRAIC[19])/CONSTANTS[2] : (ALGEBRAIC[8] - ALGEBRAIC[19])/CONSTANTS[2] - (ALGEBRAIC[19] - ALGEBRAIC[5])/CONSTANTS[3]); RATES[2] = (ALGEBRAIC[19] - ALGEBRAIC[5]<0.00000&&ALGEBRAIC[5] - ALGEBRAIC[6]<0.00000 ? 0.00000 : ALGEBRAIC[19] - ALGEBRAIC[5]<0.00000 ? - (ALGEBRAIC[5] - ALGEBRAIC[6])/CONSTANTS[4] : ALGEBRAIC[5] - ALGEBRAIC[6]<0.00000 ? (ALGEBRAIC[19] - ALGEBRAIC[5])/CONSTANTS[3] : (ALGEBRAIC[19] - ALGEBRAIC[5])/CONSTANTS[3] - (ALGEBRAIC[5] - ALGEBRAIC[6])/CONSTANTS[4]); RATES[5] = (ALGEBRAIC[7] - ALGEBRAIC[8]<0.00000&&ALGEBRAIC[8] - ALGEBRAIC[19]<0.00000 ? 0.00000 : ALGEBRAIC[7] - ALGEBRAIC[8]<0.00000 ? - (ALGEBRAIC[8] - ALGEBRAIC[19])/CONSTANTS[2] : ALGEBRAIC[8] - ALGEBRAIC[19]<0.00000 ? (ALGEBRAIC[7] - ALGEBRAIC[8])/CONSTANTS[5] : (ALGEBRAIC[7] - ALGEBRAIC[8])/CONSTANTS[5] - (ALGEBRAIC[8] - ALGEBRAIC[19])/CONSTANTS[2]); } void computeVariables(double VOI, double* CONSTANTS, double* RATES, double* STATES, double* ALGEBRAIC) { ALGEBRAIC[2] = STATES[0]+STATES[1]; ALGEBRAIC[3] = CONSTANTS[13]*(exp( CONSTANTS[15]*(ALGEBRAIC[2] - CONSTANTS[14])) - 1.00000); ALGEBRAIC[4] = ALGEBRAIC[3]+CONSTANTS[8]; ALGEBRAIC[0] = (int)(VOI) % (int)(CONSTANTS[12]); ALGEBRAIC[1] = CONSTANTS[9]*exp( - CONSTANTS[10]*pow(( ALGEBRAIC[0]*CONSTANTS[6])/60.0000 - CONSTANTS[11], 2.00000)); ALGEBRAIC[18] = ALGEBRAIC[10]+ALGEBRAIC[4]; ALGEBRAIC[6] = CONSTANTS[34]*(STATES[3] - CONSTANTS[35]); ALGEBRAIC[7] = CONSTANTS[36]*(STATES[4] - CONSTANTS[37]); ALGEBRAIC[5] = CONSTANTS[32]*(STATES[2] - CONSTANTS[33]); ALGEBRAIC[8] = CONSTANTS[38]*(STATES[5] - CONSTANTS[39]); ALGEBRAIC[19] = ALGEBRAIC[15]+ALGEBRAIC[4]; ALGEBRAIC[20] = ALGEBRAIC[18] - ALGEBRAIC[19]; } void objfunc_0(double *p, double *hx, int m, int n, void *adata) { struct rootfind_info* rfi = (struct rootfind_info*)adata; #define VOI rfi->aVOI #define CONSTANTS rfi->aCONSTANTS #define RATES rfi->aRATES #define STATES rfi->aSTATES #define ALGEBRAIC rfi->aALGEBRAIC #define pret rfi->aPRET ALGEBRAIC[9] = p[0]; ALGEBRAIC[10] = p[1]; ALGEBRAIC[11] = p[2]; ALGEBRAIC[12] = p[3]; ALGEBRAIC[13] = p[4]; ALGEBRAIC[14] = p[5]; ALGEBRAIC[15] = p[6]; ALGEBRAIC[16] = p[7]; ALGEBRAIC[17] = p[8]; hx[0] = ALGEBRAIC[9] - (STATES[0] - ALGEBRAIC[11]); hx[1] = ALGEBRAIC[10] - ( ALGEBRAIC[1]*ALGEBRAIC[12]+ (1.00000 - ALGEBRAIC[1])*ALGEBRAIC[13]); hx[2] = ALGEBRAIC[12] - CONSTANTS[16]*(ALGEBRAIC[9] - CONSTANTS[17]); hx[3] = ALGEBRAIC[13] - CONSTANTS[18]*(exp( CONSTANTS[19]*(ALGEBRAIC[9] - CONSTANTS[20])) - 1.00000); hx[4] = ALGEBRAIC[14] - (STATES[1]+ALGEBRAIC[11]); hx[5] = ALGEBRAIC[15] - ( ALGEBRAIC[1]*ALGEBRAIC[16]+ (1.00000 - ALGEBRAIC[1])*ALGEBRAIC[17]); hx[6] = ALGEBRAIC[16] - CONSTANTS[21]*(ALGEBRAIC[14] - CONSTANTS[22]); hx[7] = ALGEBRAIC[17] - CONSTANTS[23]*(exp( CONSTANTS[24]*(ALGEBRAIC[14] - CONSTANTS[25])) - 1.00000); hx[8] = ALGEBRAIC[10] - (( ALGEBRAIC[1]*CONSTANTS[26]*(ALGEBRAIC[11] - CONSTANTS[27])+ (CONSTANTS[31] - ALGEBRAIC[1])*CONSTANTS[28]*(exp( CONSTANTS[29]*(ALGEBRAIC[11] - CONSTANTS[30])) - CONSTANTS[31])) - ALGEBRAIC[15]); #undef VOI #undef CONSTANTS #undef RATES #undef STATES #undef ALGEBRAIC #undef pret } void rootfind_0(double VOI, double* CONSTANTS, double* RATES, double* STATES, double* ALGEBRAIC, int* pret) { static double p[9] = {0.1,0.1,0.1,0.1,0.1,0.1,0.1,0.1,0.1}; double bp[9], work[LM_DIF_WORKSZ(9, 9)]; struct rootfind_info rfi; rfi.aVOI = VOI; rfi.aCONSTANTS = CONSTANTS; rfi.aRATES = RATES; rfi.aSTATES = STATES; rfi.aALGEBRAIC = ALGEBRAIC; rfi.aPRET = pret; do_levmar(objfunc_0, p, bp, work, pret, 9, &rfi); ALGEBRAIC[9] = p[0]; ALGEBRAIC[10] = p[1]; ALGEBRAIC[11] = p[2]; ALGEBRAIC[12] = p[3]; ALGEBRAIC[13] = p[4]; ALGEBRAIC[14] = p[5]; ALGEBRAIC[15] = p[6]; ALGEBRAIC[16] = p[7]; ALGEBRAIC[17] = p[8]; }