3 dimensional plotting by varying the initial conditions

9 views (last 30 days)
I have a set of four differential equations(DE) and I want to see how one of them varies with the other two. I can solve the equation and use plot3 to see the trajectory, but i want to do more and get more information out of it by some kind of 3d plots. I am stuck and hence don't have a lot to show my attempt but here is what i did. The four variables in the system are N, T, A and I. I want to vary the initial conditions of only two of them in some way and solve the DE by fixing the third one. Is it possible i can get a 3d plot of say, how N changes with variation in T and I.
I can change the IC completely and use a for loop and get a set of lines in 3D but i was looking for some sort of 3D surface plot.
My attempt: (Honestly this is not much of an attempt but rather what I have done)
timerange= 0:0.5:200;
IC= [1,0.1,0,0];%initial conditions
[t,y] =ode45(@(t,y) fn(t,y),timerange, IC);
% For 3d plotting defining a meshgrid to specify Initial condition IC
T0=linspace(0,1,100);
I0=linspace(0,1,100);
[T0,I0]= meshgrid(T0,I0);
%Can't figure out how to procceed from here.
plot3(y(1,:),y(3,:),y(2,:));
xlabel('N')
ylabel('I')
zlabel('T')
grid on
xlabel('Time')
ylabel('Fraction of Population')
title('Evolution of cells')
grid on
function rk1 =fn(~,y)
r = 0.60516;
K = 8.97523;
A0 = 0.4;
gammaA = 0.04;
eps = 0.00379;
rho = 0.02733;
alpha1 = 2.15877;
c1 = 0.02718;
I0 = 0.3;
gammaI = 0.0208;
m = 0.08006;
q = 0.61760;gammaT=4;
T0 = 0.4;
alpha = (gammaA*m)/(A0);
beta = (gammaT*q)/(T0);
n= y(1);
A= y(2);
T= y(3);
I= y(4);
rk1(1)= r*n*(1- n/K)+ alpha*A*n -eps*n*I -beta*n*T;
rk1(2) = A0 - gammaA*A;
rk1(3) = T0- gammaT*T;
rk1(4) = I0 + (rho*I*n)/(alpha1+n) -c1*I*n - gammaI*I;
rk1=rk1(:);
end
Any sort of help or an hint would be very helpful. Thank you.

Answers (1)

darova
darova on 28 Apr 2020
Try quiver3 and streamline
clc,clear
cla
r = 0.60516;
K = 8.97523;
A0 = 0.4;
gammaA = 0.04;
eps = 0.00379;
rho = 0.02733;
alpha1 = 2.15877;
c1 = 0.02718;
I0 = 0.3;
gammaI = 0.0208;
m = 0.08006;
q = 0.61760;gammaT=4;
T0 = 0.4;
alpha = (gammaA*m)/(A0);
beta = (gammaT*q)/(T0);
A = 1;
[n,T,I] = meshgrid(0:.2:1);
dn = r*n.*(1- n/K)+ alpha*A.*n -eps*n.*I -beta*n.*T;
% rk1(2) = A0 - gammaA*A;
dT = T0 - gammaT*T;
dI = I0 + (rho*I.*n)./(alpha1+n) -c1*I.*n - gammaI*I;
[X,Z] = meshgrid(0:.2:1);
quiver3(n,T,I,dn,dT,dI,2)
streamline(n,T,I,dn,dT,dI,X,X*0+1,Z)
axis vis3d
xlabel('n axis')
ylabel('T axis')
zlabel('I axis')
  4 Comments

Sign in to comment.

Products


Release

R2019b

Community Treasure Hunt

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

Start Hunting!