How to combine two graph from two function in another script?

조회 수: 1 (최근 30일)
cindyawati cindyawati
cindyawati cindyawati 2023년 6월 27일
댓글: Torsten 2023년 6월 28일
Hello, I have 2 function in another script. I want to combine the graph into 1 graph. Is it possible and how to combine it?

채택된 답변

Jon
Jon 2023년 6월 27일
Example
Run script 1
% Script 1
% define first curve
t1 = linspace(0,5);
y1 = t1 + 2*t1.^2;
Run script 2
% Script 2
% define second curve
t2 = linspace(0,5);
y2 = -t2 -3*t2.^2;
t1,y1,t2,y2 will be in the workspace. Now plot both curves either from the command line or in another script
plot(t1,y1,t2,y2)
  댓글 수: 5
cindyawati cindyawati
cindyawati cindyawati 2023년 6월 28일
I already try but I get error (error using plot, my vectors must be same length) @Torsten , my code like this
%Euler function
%input parameter
delta=50;
K1= 10^-4;
Ko=0.1;
n=3;
Oa=10;
Pa=100;
mu_1=10^-3;
mu_2=10^-3;
mu_3=10^-3;
mu_o=10^-4;
mu_p= 10^-5;
K2=5*10^-4;
K3=10^-3;
gamma=75;
%input initial condition
M11(1)=10;
M22(1)=0;
M33(1)=0;
O1(1)=0;
P1(1)=0;
%input for time
t(1)=0;
dt=0.01; %time interval
t=0:dt:100; %time span
%input empty array
T=zeros(length(t)+1,1); %empty array for t
M11=zeros(length(t)+1,1); %empty array for M1
M22=zeros(length(t)+1,1); %empty array for M2
M33=zeros(length(t)+1,1); %empty array for M3
O1=zeros(length(t)+1,1);
P1=zeros(length(t)+1,1);
for i = 1:length(t)
T(i+1)=T(i)+dt;
M11(i+1)=M11(i)+1./(1+exp(-T(i)));
M11(i+1) = M11(i)+(dt*(delta*M11(i+1)*(1-(M11(i+1)/gamma))-2*K1*M11(i+1)*M11(i+1)-M11(i+1)*(K2.*M22(i+1))-((Oa-n)*K3*M11(i+1)*M33(i+1))-((Pa-Oa)*Ko*M11(i+1)*O1(i+1))-(mu_1*M11(i+1))));
M22(i+1) = M22(i)+(dt*(K1*M11(i)*M11(i)-K2*M11(i)*M22(i))-(mu_2*M22(i+1)));
M33(i+1) = M33(i)+(dt*(K2*M11(i)*M22(i)-K3*M11(i)*M33(i)-mu_3*M33(i)));
O1(i+1) = O1(i)+(dt*(K3*M11(i)*M33(i)-Ko*M11(i)*O1(i)-mu_o*O1(i)));
P1(i+1) = P1(i)+(dt*(Ko*M11(i)*O1(i)-mu_p*P1(i)));
end
%RK function
tstart = 0;
tend = 100;
dt = 0.01;
T = (tstart:dt:tend).';
Y0 = [10 0 0 0 0];
f = @myode;
Y = fRK4(f,T,Y0);
M1 = Y(:,1);
M2 = Y(:,2);
M3 = Y(:,3);
O = Y(:,4);
P = Y(:,5);
plot (T,M11,'r-',T,M1,'r');
Error using plot
Vectors must be the same length.
legend('Euler', 'RK4');
xlabel('time (days)');
ylabel ('M1 (gr/ml)');
title('M1');
function Y = fRK4(f,T,Y0)
N = numel(T);
n = numel(Y0);
Y = zeros(N,n);
Y(1,:) = Y0;
for i = 2:N
t = T(i-1);
y = Y(i-1,:);
h = T(i) - T(i-1);
k0 = f(t,y);
k1 = f(t+0.5*h,y+k0*0.5*h);
k2 = f(t+0.5*h,y+k1*0.5*h);
k3 = f(t+h,y+k2*h);
Y(i,:) = y + h/6*(k0+2*k1+2*k2+k3);
end
end
function CM1 = myode (~,MM)
M1 = MM(1);
M2 = MM(2);
M3 = MM(3);
O = MM(4);
P = MM(5);
delta=50;
gamma=75;
K1= 10^-4;
K2=5*10^-4;
K3=10^-3;
Ko=0.1;
n=3;
Oa=10;
Pa=100;
mu_1=10^-3;
mu_2=10^-3;
mu_3=10^-3;
mu_o=10^-4;
mu_p= 10^-5;
CM1= zeros(1,5);
CM1(1) = (delta*M1*(1-(M1/gamma))-2*K1*M1*M1-M1*(K2*M2)-((Oa-n)*K3*M1*M3)-((Pa-Oa)*Ko*M1*O)-(mu_1*M1));
CM1(2) = (K1*M1*M1)-(K2*M1*M2)-(mu_2*M2);
CM1(3) = (K2*M1*M2)-(K3*M1*M3)-(mu_3*M3);
CM1(4) = (K3*M1*M3)-(Ko*M1*O)-(mu_o*O);
CM1(5) = (Ko*M1*O)-(mu_p*P);
end
Torsten
Torsten 2023년 6월 28일
tstart = 0;
tend = 100;
dt = 0.01;
T = (tstart:dt:tend).';
Y0 = [10 0 0 0 0];
f = @myode;
YRK4 = fRK4(f,T,Y0);
M1RK4 = YRK4(:,1);
M2RK4 = YRK4(:,2);
M3RK4 = YRK4(:,3);
ORK4 = YRK4(:,4);
PRK4 = YRK4(:,5);
YEuler = fEuler(f,T,Y0);
M1Euler = YEuler(:,1);
M2Euler = YEuler(:,2);
M3Euler = YEuler(:,3);
OEuler = YEuler(:,4);
PEuler = YEuler(:,5);
plot (T,M1Euler,'r',T,M1RK4,'b');
legend('Euler', 'RK4');
xlabel('time (days)');
ylabel ('M1 (gr/ml)');
title('M1');
function Y = fEuler(f,T,Y0)
N = numel(T);
n = numel(Y0);
Y = zeros(N,n);
Y(1,:) = Y0;
for i = 2:N
t = T(i-1);
y = Y(i-1,:);
h = T(i) - T(i-1);
Y(i,:) = y + h*f(t,y);
end
end
function Y = fRK4(f,T,Y0)
N = numel(T);
n = numel(Y0);
Y = zeros(N,n);
Y(1,:) = Y0;
for i = 2:N
t = T(i-1);
y = Y(i-1,:);
h = T(i) - T(i-1);
k0 = f(t,y);
k1 = f(t+0.5*h,y+k0*0.5*h);
k2 = f(t+0.5*h,y+k1*0.5*h);
k3 = f(t+h,y+k2*h);
Y(i,:) = y + h/6*(k0+2*k1+2*k2+k3);
end
end
function CM1 = myode (~,MM)
M1 = MM(1);
M2 = MM(2);
M3 = MM(3);
O = MM(4);
P = MM(5);
delta=50;
gamma=75;
K1= 10^-4;
K2=5*10^-4;
K3=10^-3;
Ko=0.1;
n=3;
Oa=10;
Pa=100;
mu_1=10^-3;
mu_2=10^-3;
mu_3=10^-3;
mu_o=10^-4;
mu_p= 10^-5;
CM1= zeros(1,5);
CM1(1) = (delta*M1*(1-(M1/gamma))-2*K1*M1*M1-M1*(K2*M2)-((Oa-n)*K3*M1*M3)-((Pa-Oa)*Ko*M1*O)-(mu_1*M1));
CM1(2) = (K1*M1*M1)-(K2*M1*M2)-(mu_2*M2);
CM1(3) = (K2*M1*M2)-(K3*M1*M3)-(mu_3*M3);
CM1(4) = (K3*M1*M3)-(Ko*M1*O)-(mu_o*O);
CM1(5) = (Ko*M1*O)-(mu_p*P);
end

댓글을 달려면 로그인하십시오.

추가 답변 (0개)

카테고리

Help CenterFile Exchange에서 Ordinary Differential Equations에 대해 자세히 알아보기

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by