% The function eval_deriv.m evaluates the right-hand side of a set
% of first-order differential equations which characterize the
% intact circulation with only linear elements at a desired time
% step.
%
% Function arguments:
% P - a 6x1 vector containing the current pressure values
% = [Pl; Pa; Pv; Pr; Ppa; Ppv]
% Qc - a 2x1 vector containing the current ventricular volume values
% - [Ql; Qr]
% Pbreathe - 4x1 vector containing current respiratory values
% [Qlu; Pth dPth Ppac (Palv)]
% th - current parameter values
% sumdeltaT - time surpassed in current cardiac cycle
%
% Function outputs:
% dP - a 6x1 vector containing the derivative of current
% pressure values
% = [dPl; dPa; dPv; dPr; dPpa; dPpv]
%
function dP = eval_deriv(P,Qc,Pbreathe,th,sumdeltaT);
% Computing flow rates from the current pressure values.
if (P(6) >= P(1))
qli = (P(6)-P(1))/th(20);
else
qli = 0;
end
if (P(1) >= P(2))
qlo = (P(1)-P(2))/th(15);
else
qlo = 0;
end
qa = (P(2)-P(3))/th(16);
if (P(3) >= P(4))
qri = (P(3)-P(4))/th(17);
else
qri = 0;
end
if (P(4) >= P(5))
qro = (P(4)-P(5))/th(18);
else
qro = 0;
end
qpa = (P(5)-P(6))/th(19);
% Computing ventricular elastances.
[El,dEl] = var_cap(th(1),th(2),th(24),sumdeltaT);
[Er,dEr] = var_cap(th(5),th(6),th(24),sumdeltaT);
% Computing dP from the above calculations and th parameter vector.
dP = [(El*(qli-qlo))+((dEl/El)*(P(1)-Pbreathe(2)))+Pbreathe(3);
((qlo-qa)/th(3))+(Pbreathe(3)/3);
(qa-qri)/th(4);
(Er*(qri-qro))+((dEr/Er)*(P(4)-Pbreathe(2)))+Pbreathe(3);
((qro-qpa)/th(7))+Pbreathe(3);
((qpa-qli)/th(8))+Pbreathe(3);
0;
0;];