Paste code

Era = (Emaxra - Eminra) * yra + Eminra  
Vrar = (1 - yra) * (Vrard - Vrars) + Vrars


Erv = (Emaxrv - Eminrv) * yrv + Eminrv
Vrvr = (1 - yrv) * (Vrvrd - Vrvrs) + Vrvrs


Ela = (Emaxla - Eminla) * yla + Eminla  
Vlar = (1 - yla) * (Vlard - Vlars) + Vlars


Elv = (Emaxlv - Eminlv) * ylv + Eminlv 
Vlvr = (1 - ylv) * (Vlvrd - Vlvrs) + Vlvrs


dVla_dt = Fpv - Fla


Fla = (Plac-Plvc)/Rla  


Pla = (Vla - Vlar) * Ela - Px2 * (1 / ((Vla / Vx8) - 1)) 
Plac = Pla + Ppcdc


dVlv_dt = Fla - Flv 


Flv = (Plvc-Paopc)/Rav


Plv = afs_con2 * (Vlv - Vlvr) * Elv - Px2 * (1 / ((Vlv / Vx8) - 1))
Plvc = Plv + Ppcdc


dVra_dt = Fvc - Fra + Fcorev  


Fra = (Prac - Prvc) / Rra  


Pra = (Vra - Vrar) * Era - Px2 * (1 / ((Vra / Vx8) - 1))  
Prac = Pra + Ppcdc


dVrv_dt = Fra - Frv 


Frv = (Prvc-Ppapc)/Rpuv


Prv = afs_con2*(Vrv-Vrvr)*Erv - Px2*(1/((Vrv/Vx8)-1))
Prvc = Prv + Ppcdc

COutput = (Flv - COutput)/COtau
SV = COutput/HR


Rsa = Rsao + (Kr*(4*F_vaso)) + (Kr*(Vsa_max/Vsa)**2)
Rvc = (KR*(Vmax_vc/Vvc)**2) + Ro

Paopc_t = (Flv-Vaop_t-Faop-Fcorao)*((1/Ccorao)+((Px2/(1))*exp(Vcorao/Vx1)/(exp(Vcorao/Vx1)-1)**2)) + PplcFOL_t
Paop = Paopc - PplcFOL

MAP_t = (Psa-MAP)/COtau

Paodc = (Rtaod*Rcrb*Faop - Rtaod*Rcrb*Faod + (Vaod*Rcrb/Caod) - Rcrb*Px2/((Vaod/Vx8)-1) + Pvcc*Rtaod)/(Rcrb+Rtaod)
Paod = Paodc - Pbs

PaodFOL_t = (Paod-PaodFOL)

Psa_a = Kc*( ((Vsa-Vsa_o)/Do) 
Psa_p = Kp1*(tau_p*(Vsa-Vsa_o)) + Kp2*(Vsa-Vsa_o)
Psa = F_vaso*Psa_a + (1-F_vaso)*Psa_p

Psap = Vsap/Csap - Px2*(1/((Vsap/Vx8)-1))
Psc  = Vsc/Csc - Px2*(1/((Vsc/Vx8)-1))
Psv  = -Kv*((Vmax_sv/Vsv)-0.99)


Pvc = D2+K2*exp(Vo/Vmin_vc)+K1*(Vvc-Vo) - Px2/(exp(Vvc/Vx8)-1) + D2+K2*exp(Vvc/Vmin_vc) - Px2/(exp(Vvc/Vx8)-1)


Pvcc = Pvc + Pplc

Fcrb = (Paodc-Pvcc)/Rcrb
Fsap = (Psap-Psa)/Rsap
Fsa = (Psa-Psc)/Rsa
Fsc = (Psc-Psv)/Rsc
Fsv = (Psv-Pvcc)/Rsv
Fvc = (Pvcc-Prac)/Rvc

Vaop_t = (Paopc - (Vaop/Caop) + Px2*(1/((Vaop/Vx8)-1)) - PplcFOL)/Rtaop
Vaod_t = Faop - Faod - Fcrb
Vsap_t = Faod - Fsap
Vsa_t = Fsap - Fsa
Vsc_t = Fsa - Fsc
Vsv_t = Fsc - Fsv
Vvc_t = Fsv + Fcrb - Fvc

Faop_t = (Paopc - Faop*Raop - Paodc) / Laop
Faod_t = (Paodc - Faod*Raod - Psap) / Laod

Vtot = Vra+Vrv+Vpap+Vpad+Vpa+Vpc+Vpv+Vla+Vlv+Vaop+Vaod + Vsa+Vsap+Vsc+Vsv+Vvc+perifl+Vcorcirc
VBcirc = Vtot-perifl

SysArtVol = Vaop + Vaod + Vsa + Vsap
SysVenVol = Vsv + Vvc
PulArtVol = Vpap + Vpad + Vpa
PulVenVol = Vpv


Ppapc1 = (Rpuv*(Vpap/Ctpap) - Rpuv*(Px2*(1/((Vpap/Vx8)-1))) + Prvc*Rtpap - Rpuv*Rtpap*Fpap + Pplc*Rpuv)/(Rtpap+Rpuv)
Ppapc2 = (Vpap/Ctpap + Pplc - Rtpap*Fpap - Px2*(1/((Vpap/Vx8)-1)))


Ppapc =  Ppapc1 + Ppapc2

Ppap = Ppapc-Pplc

Ppadc = Vpad_t*Rtpad + Pplc + Vpad/Ctpad - Px2*(1/(Vpad/Vx8)-1))
Ppad = Ppadc-Pplc

Ppa = Vpa/Cpa - Px2*(1/((Vpa/Vx8)-1))
Ppac = Ppa + Pplc

Ppc = Vpc/Cpc - Px2*(1/((Vpc/Vx8)-1))
Ppcc = Ppc + Pplc

Ppv = Vpv/Cpv - Px2*(1/((Vpv/Vx8)-1))
Ppvc = Ppv + Pplc

Fps = (Ppac-Ppvc)/Rps
Fpa = (Ppac-Ppcc)/Rpa
Fpc = (Ppcc-Ppvc)/Rpc
Fpv = (Ppvc-Plac)/Rpv

Vpap_t = Frv - Fpap
Vpad_t = Fpap - Fpad
Vpa_t = Fpad - Fps - Fpa
Vpc_t = Fpa - Fpc
Vpv_t = Fpc+Fps-Fpv

Fpap_t = (Ppapc - Ppadc - Fpap*Rpap)/Lpa
Fpad_t = (Ppadc - Ppac - Fpad*Rpad)/Lpad


Nbr = (K*(Paod) + (a1*K*PaodFOL_t)) - ((a2+a)*Nbr_t) - (a2*a*(Nbr_t_t))


N_hrv_t = (-N_hrv + (K_hrv*Nbr(t-L_hrv)))/T_hrv

F_hrv = a_hrv + (b_hrv / ((tau_hrv*(N_hrv-No_hrv)) + 1.0))
b_hrv = 1 - a_hrv


N_hrs_t = (-N_hrs + (K_hrs*Nbr(t-L_hrs)))/T_hrs

F_hrs = a_hrs + (b_hrs / ((tau_hrs*(N_hrs-No_hrs)) + 1.0))
b_hrs = 1 - a_hrs


N_con_t = (-N_con + (K_con*Nbr(t-L_con)))/T_con

F_con = a_con + (b_con / ((tau_con*(N_con-No_con)) + 1.0))


N_vaso_t = (-N_vaso + (K_vaso*Nbr(t-L_vaso)))/T_vaso

b_vaso = 1 - a_vaso
F_vaso = a_vaso + (b_vaso / ((tau_vaso*(N_vaso-No_vaso)) + 1.0))

HRcont = (h1 + (h2*F_hrs) - (h3*F_hrs**2) - (h4*F_hrv) + (h5*F_hrv**2) - (h6*F_hrv*F_hrs))
afs_con = amin + (Ka*F_con)
bfs_con = bmin + (Kb*F_con)


Rc = Kc_air_t + Rco + Kc_air * (Vcmax / Vc)**2 + Rco


Rs = As * (Ks * (VA - RV) / (Vstar - RV)) + Bs
Ru = Au + Ku * (dV_t)

dV_t = (Vcw - dV) / tau  

Pl = Al * (Kl * VA) + Bl
Pcw = A * (2 *  * RespR * (t - tresp)) - B


Pc = Ac - Bc * (Vc / Vcmax - 0.7)**2 + Ac-Bc*(0.5 - 0.7)^2-Bcp*ln(Vcmax/Vc - 0.999)


Ppl = Pcw
Pve = Vve / Cve
PplmmHg = Ppl
Pmouth = Pbs
Pplc = Ppl
PplcFOL_t = (Pplc - PplcFOL) / 0.001


RespRcont = C/(2*PI)

Pcc = Pplc + Pc
Ps = Pcc - PA  

Vcw_t = VA_t + Vc_t

Qdotco = (Pmouth - Pcc) / (Ru + Rc)
Qdotup = Qdotco
Qdotve = Pve / Rve
Qdotsm = Ps / Rs

Qdotsm - Alvflux - VA_t = 0
VA_t - Vve_t - Qdotve = 0
Vc_t = Qdotco - Qdotsm

Alvflux = (Pstp * Tbody / ((PA + 760) * Tstp)) * (O2flux + CO2flux)  

PA = Ppl + Pve + Pl
PD = Qdotup * Ru
PDc = Qdotco * Rc + Pcc
Pup = Qdotup * (Ru + Rc)


SHbO2_ao = (PO2_ao / P50_O2)**nH / (((PO2_ao / P50_O2)**nH) + 1) 
SHbO2_pa = (PO2_pa / P50_O2)**nH / (((PO2_pa / P50_O2)**nH) + 1)  
SHbO2_aoISR = SHbO2_ao * 100

dSHbO2dPO2_ao = nH * (PO2_ao**-1) * (PO2_ao / P50_O2)**nH / (((PO2_ao / P50_O2)**nH) + 1)**2
dSHbO2dPO2_pa = nH * (PO2_pa**-1) * (PO2_pa / P50_O2)**nH / (((PO2_pa / P50_O2)**nH) + 1)**2

CtO2_ao = alphaO2 * PO2_ao + CHb * SHbO2_ao * Hcrit  
CtO2_pa = alphaO2 * PO2_pa + CHb * SHbO2_pa * Hcrit
CtO2_isf = alphaO2 * PO2_isf


DL_O2 = ((0.397)+((0.0085)*PO2_ao)-((0.00013 )*(PO2_ao**2)))+((5.1*exp(-7))*(PO2_ao^3))+((Vpc/Vpcmax)**(1/2))*((0.397)+((0.0085)*PO2_ao)-((0.00013 )*(PO2_ao^2)))+((5.1*exp(-7))*(PO2_ao^3))


V_O2 = 0 + PO2_isf + (Vcytox * tweight * alphaO2 * PO2_isf) / (Kcytox + (alphaO2 * PO2_isf))  
O2flux = DL_O2 * (PA_O2 - PO2_ao) 


PO2_ao_t = (COutput * (CtO2_pa - CtO2_ao) + DL_O2 * (PA_O2 - PO2_ao)) / (Vpc * (alphaO2 + Hcrit * CHb * dSHbO2dPO2_ao))
PO2_pa_t = (COutput * (CtO2_ao - CtO2_pa) + PS * tweight * alphaO2 * (PO2_isf - PO2_pa)) / ((Vsc + Vcorcap) * (alphaO2 + Hcrit * CHb * dSHbO2dPO2_pa))
PO2_isf_t = PS * tweight * (PO2_pa - PO2_isf) / Visf - V_O2 / (alphaO2 * Visf)

SHbCO2_ao = PCO2_ao / P50_CO2 / (PCO2_ao / P50_CO2 + 1)
SHbCO2_pa = PCO2_pa / P50_CO2 / (PCO2_pa / P50_CO2 + 1)

dSHbCO2dPCO2_ao = (PCO2_ao**-1) * (PCO2_ao / P50_CO2) / (PCO2_ao / P50_CO2 + 1)**2
dSHbCO2dPCO2_pa = (PCO2_pa**-1) * (PCO2_pa / P50_CO2) / (PCO2_pa / P50_CO2 + 1)**2

CtCO2_ao = alphaCO2 * PCO2_ao + CHb * SHbCO2_ao * Hcrit
CtCO2_pa = alphaCO2 * PCO2_pa + CHb * SHbCO2_pa * Hcrit
CtCO2_isf = alphaCO2 * PCO2_isf


DL_CO2 = 16.67 * (Vpc / Vpcmax) ** (1/2) 


V_CO2 = RQ * V_O2


CO2flux = DL_CO2 * (PA_CO2 - PCO2_ao)*22400

PCO2_ao_t = (COutput * (CtCO2_pa - CtCO2_ao) + DL_CO2 * (PA_CO2 - PCO2_ao)) / (Vpc * (alphaCO2 + Hcrit * CHb * dSHbCO2dPCO2_ao))
PCO2_pa_t = (COutput * (CtCO2_ao - CtCO2_pa) + PS * tweight * alphaCO2 * (PCO2_isf - PCO2_pa)) / ((Vsc + Vcorcap) * (alphaCO2 + Hcrit * CHb * dSHbCO2dPCO2_pa))
PCO2_isf_t = PS * tweight * (PCO2_pa - PCO2_isf) / Visf + V_CO2 / (alphaCO2 * Visf)

Pao_O2 = (Pao - PH2O) * r_Pao_O2
Pao_CO2 = (Pao - PH2O) * r_Pao_CO2
Pao_N2 = (Pao - PH2O) - Pao_O2 - Pao_CO2


PD_O2_t = (Qdotup * Pao_O2 - Qdotco * PD_O2) / VD + Qdotup + (Qdotup * PD_O2 - Qdotco * PC_O2) / VD


PC_O2_t = (((Qdotco*PD_O2 - Qdotsm*PC_O2 - PC_O2*Vc_t)+(Qdotco*PD_O2 - Qdotsm*PA_O2 - PC_O2*Vc_t)+(Qdotco*PC_O2 - Qdotsm*PC_O2 - PC_O2*Vc_t)+(Qdotco*PC_O2 - Qdotsm*PA_O2 - PC_O2*Vc_t))/Vc)+Qdotup


PD_CO2_t = (Qdotup * Pao_CO2 - Qdotco * PD_CO2) / VD

Pmouth_CO2 = PD_CO2


PC_CO2_t = ((Qdotco*PD_CO2 - Qdotsm*PC_CO2 - PC_CO2*Vc_t)+(Qdotco*PD_CO2 - Qdotsm*PA_CO2 - PC_CO2*Vc_t)+(Qdotco*PC_CO2 - Qdotsm*PC_CO2 - PC_CO2*Vc_t)+(Qdotco*PC_CO2 - Qdotsm*PA_CO2 - PC_CO2*Vc_t))/Vc+Qdotup
PD_N2 = (PDc + Pao) - PD_O2 - PD_CO2
PC_N2 = (Pcc + Pao) - PC_O2 - PC_CO2
PA_N2 = (PA + Pao) - PA_O2 - PA_CO2


DL_N2 = ((Vpc/Vpcmax)) ** (1/2)*0.25


GO2isf = ((Vcytox / (Kcytox + CtO2_isf))  + (Vcytox / (Kcytox + CtO2_isf))

PBC_pc_t = (Fpc / Vpc) * (PBC_sc - PBC_pc) + kp5 * CtCO2_ao * (Cheme - PBC_pc) * (SHbO2_ao / (1 + H_ao / K3bgh) + (1 - SHbO2_ao) / (1 + H_ao / K2bgh)) - kp5 * PBC_pc * H_ao * (SHbO2_ao / K6bgh + (1 - SHbO2_ao) / K5bgh)
PBC_sc_t = (COutput / (Vsc + Vcorcap)) * (PBC_pc - PBC_sc) + kp5 * CtCO2_pa * (Cheme - PBC_sc) * (SHbO2_pa / (1 + H_pa / K3bgh) + (1 - SHbO2_pa) / (1 + H_pa / K2bgh)) - kp5 * PBC_sc * H_pa * (SHbO2_pa / K6bgh + (1 - SHbO2_pa) / K5bgh)

HCO3_ao_t = (Fpc / Vpc) * (HCO3_pa - HCO3_ao) + CFb * (kp1 * CtCO2_ao - km1 * HCO3_ao * H_ao / K1bgh)
HCO3_pa_t = (COutput / (Vsc + Vcorcap)) * (HCO3_ao - HCO3_pa) + CFb * (kp1 * CtCO2_pa - km1 * HCO3_pa * H_pa / K1bgh)
HCO3_isf_t = (PS_HCO3 * tweight / Visf) * (HCO3_ao - HCO3_isf) + CFt * (kp1 * CtCO2_isf - km1 * HCO3_isf * H_isf / K1bgh)

DVRO11 = 1 + dSHbO2dPO2_ao * O2cap * Hbconc / alphaO2
DVRO21 = 1 + dSHbO2dPO2_pa * O2cap * Hbconc / alphaO2

Hout_t =  H_ao + (H_ao + tauHout
H_ao_t = (Fpc / Vpc) * (H_pa - H_ao) + (2.303 / Betab) * H_ao * (CFb * (kp1 * CtCO2_ao - km1 * H_ao * HCO3_ao / K1bgh) + 1.5 * (PBC_pc_t) - 0.6 * alphaO2 * DVRO11 * (PO2_ao_t)) - Hout_t
H_pa_t = (COutput / (Vsc + Vcorcap)) * (H_ao - H_pa) + (2.303 / Betab) * H_pa * (CFb * (kp1 * CtCO2_pa - km1 * H_pa * HCO3_pa / K1bgh) + 1.5 * (PBC_sc_t) - 0.6 * alphaO2 * DVRO21 * (PO2_pa_t))

H_isf_t = (PS_H * tweight / Visf) * (H_ao - H_isf) + (2.303 / Betat) * H_isf * CFt * (kp1 * CtCO2_isf - km1 * H_isf * HCO3_isf / K1bgh) + (facid * GO2isf * tweight * CtO2_isf / Visf)

pH_ao = -(H_ao / HpNorm)
pH_pa = -H_pa / HpNorm) 
pH_isf = -(H_isf / HpNorm)


Ppcd = K_pcd * ((Vpcd - Vpcd_o) / phi_pcd) - Px2 * (1 / ((Vpcd / Vx75) - 1))
Ppcdc = Ppcd + Pplc
Vpcd = Vrv + Vra + Vlv + Vla + perifl + Vmyo + Vcorcirc


Pcorisfc = ((Plvc - Ppcdc) / 2)

Fcorao = (Pcoraoc - Pcoreac) / Rcorao
Fcorea = (Pcoreac - Pcorlac) / Rcorea
Fcorla = (Pcorlac - Pcorsac) / Rcorla
Fcorsa = (Pcorsac - Pcorcapc) / Rcorsa
Fcorcap = (Pcorcapc - Pcorsvc) / Rcorcap
Fcorsv = (Pcorsvc - Pcorlvc) / Rcorsv
Fcorlv = (Pcorlvc - Pcorevc) / Rcorlv
Fcorev = (Pcorevc - Prac) / Rcorev

Pcorao = Paop
Pcorea = Vcorea / Ccorea - Px1 * (1 / ((Vcorea / Vx1) - 1))
Pcorla = Vcorla / Ccorla - Px1 * (1 / ((Vcorla / Vx1) - 1))
Pcorsa = Vcorsa / Ccorsa - Px1 * (1 / ((Vcorsa / Vx1) - 1))
Pcorcap = Vcorcap / Ccorcap - Px1 * (1 / ((Vcorcap / Vx1) - 1))
Pcorsv = Vcorsv / Ccorsv - Px1 * (1 / ((Vcorsv / Vx1) - 1))
Pcorlv = Vcorlv / Ccorlv - Px1 * (1 / ((Vcorlv / Vx1) - 1))
Pcorev = Vcorev / Ccorev - Px2 * (1 / ((Vcorev / Vx8) - 1))

Pcoraoc = Paopc
Pcoreac = Pcorea + Ppcdc
Pcorlac = Pcorla + Pcorisfc
Pcorsac = Pcorsa + Pcorisfc
Pcorcapc = Pcorcap + Pcorisfc
Pcorsvc = Pcorsv + Ppcdc
Pcorlvc = Pcorlv + Ppcdc
Pcorevc = Pcorev + Ppcdc

Vcorao_t = Flv - Vaop_t - Faop - Fcorao
Vcorea_t = Fcorao - Fcorea
Vcorla_t = Fcorea - Fcorla
Vcorsa_t = Fcorla - Fcorsa
Vcorcap_t = Fcorsa - Fcorcap
Vcorsv_t = Fcorcap - Fcorsv
Vcorlv_t = Fcorsv - Fcorlv
Vcorev_t = Fcorlv - Fcorev

Vcorcirc = Vcorao + Vcorea + Vcorla + Vcorsa + Vcorcap + Vcorsv + Vcorlv + Vcorev