Trying to graph two ode45 functions on a single plot
1 visualización (últimos 30 días)
Mostrar comentarios más antiguos
hans hans
el 27 de Abr. de 2019
Comentada: darova
el 27 de Abr. de 2019
Hello I am working with this windkessel model with different parameters
%% Defining Physiological Parameters
clc,clear;
% Defining parameters for 2 element model
R = 1.08; %systemic peripheral resistance (mmHg/cm^3/sec)
C = 1.1; %systemic arterial compliance in (cm^3/mmHg)
% Heart rate and cycle times
HR = 60; %beats per minute
Tc = 60/HR; %Length of one cardiac cycle (in seconds)
Ts = (2/5)*Tc; %Length of systole is 40% of the cardiac cycle
InCond = 80; %pressure during diastole (80mmHg)
Timespan=[0 Tc];
twoelemmodel = @(t,P) (1/C).*(QF(t,HR)-(P./R));
[t,P] = ode45(twoelemmodel,Timespan,InCond);
plot(t,P,'r')
title('2 Element Windkessel Model (Physiological Parameters)');
xlabel('time (in seconds)');
ylabel('Pressure (in mmHg)');
%% %% Defining Theoretical Parameters
% Defining parameters for 2 element model
R = 1.3; %systemic peripheral resistance (mmHg/cm^3/sec)
C = 1.1; %systemic arterial compliance in (cm^3/mmHg)
% Heart rate and cycle times
HR = 120; %beats per minute
Tc = 60/HR; %Length of one cardiac cycle (in seconds)
Ts = (2/5)*Tc; %Length of systole is 40% of the cardiac cycle
InCond = 80; %pressure during diastole (80mmHg)
Timespan=[0 Tc];
twoelemmodel = @(t,P) (1/C).*(QF(t,HR)-(P./R));
[t,P] = ode45(twoelemmodel,Timespan,InCond);
plot(t,P,'b')
title('2 Element Windkessel Model (Theoretical Parameters)');
xlabel('time (in seconds)');
ylabel('Pressure (in mmHg)');
Functions I am using:
% Functions
function Qout = QF(t,HR)
HR = 60;%Heart Rate in beats per minute
Tc = 60/HR; %Length of one cardiac cycle (in seconds)
Ts =(2/5)*Tc; %Length of systole is 40% of the cardiac cycle
if HR == 60
Qo = 260;
else
Qo = 520;
end
if t > Ts
Qo = 0;
end
Qout = Qo.*sin((2.*pi.*t/Tc));
function dqdt = dQf(t,HR)
Tc = 60/HR; %Length of one cardiac cycle (in seconds)
Ts = (2/5)*Tc; %Length of systole is 40% of the cardiac cycle
if HR == 60
Qo = 260;
else
Qo = 520;
end
if t > Ts
Qo = 0;
end
dqdt = (2.*pi/Tc).*Qo.*cos((2.*pi.*t)/Tc);
I can't seem to get the two grpahs to plot together.
Hold on functions seem to work, but the plots seem identical.
Any idea?
Respuesta aceptada
Star Strider
el 27 de Abr. de 2019
For the first integration, save the outputs as:
[t{1},P{1}] = ode45(twoelemmodel,Timespan,InCond);
and the second as:
[t{2},P{2}] = ode45(twoelemmodel,Timespan,InCond);
and delete the existing plot calls, then plot the integration results as:
figure
plot(t{1}, P{1})
hold on
plot(t{2}, P{2})
hold off
grid
legend('Physiological Parameters', 'Theoretical Parameters')
title('2 Element Windkessel Model')
xlabel('time (in seconds)');
ylabel('Pressure (in mmHg)');
1 comentario
Walter Roberson
el 27 de Abr. de 2019
Or just
plot(t{1}, P{1}, 'r', t{2}, P{2}, 'b')
Más respuestas (0)
Ver también
Productos
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!