3d phase portrait

Asked by Aishah Malek

Aishah Malek (view profile)

on 28 Aug 2018
Latest activity Commented on by Aishah Malek

Aishah Malek (view profile)

on 2 Sep 2018
Hi, I want to plot a 3-d phase portrait for a set of 3 ODEs, i have used the following code and i get a basic plot, i was wondering how to add direction arrows and a mesh grid and why i only get single spirals for the trajectories. Here is the code i have used:
function f = cluster(t,y)
%BD
a=1;
b=1.2;
%equilibrium values
%c1 equilbrium value
f = zeros(size(y));
f(1) = -50*a*y(1)-b*y(1)+15*a*y(1)*y(2)+20*a*y(2)*y(3)+y(2)*b+9*a*y(2)^2+6*a*y(1)^2-60*a*y(2)-80*a*y(3)+24*a*y(2)*y(3)+16*a*y(3)^2;
f(2) = 10*a*y(1) - a*y(1)*y(2) -4*a*y(1)*y(3) -2*a*y(1)^2 -b*y(2) +3*a*y(2)^2 -10*a*y(2) +4*a*y(2)*y(3) +b*y(3);
f(3) = -2*a*y(1)*y(2) - 3*a*y(2)^2 -4*a*y(2)*y(3) +10*a*y(2)-b*y(3);
Plotting code:
[t,y] = ode45(@cluster,[0:0.01:1],[1 2 3]);
figure(1)
plot(t,y(:,3)); % plot of z(t) versus time
figure(2)
plot(t,y(:,1));
figure(3)
plot(y(:,1),y(:,3)); % plot of z versus x
figure(4)
plot3(y(:,1),y(:,2),y(:,3)); % 3D plot of trajectory
figure(5)
plot(y(:,1),y(:,2)); % plot of z versus x
figure(6)
plot(y(:,3),y(:,1)); I have computed the corresponding eigen values and vector points for specific equilbrium points in a separate file not sure if that will help?

Akshay Khadse (view profile)

on 31 Aug 2018
Can you elaborate on what is the "f" in your code above? According to me, for a phase portrait, "f" should be the gradients. However, you are plotting the solution of the differential equations, hence the single spirals.

Akshay Khadse (view profile)

on 31 Aug 2018
Edited by Akshay Khadse

Akshay Khadse (view profile)

on 31 Aug 2018

You can get 3D Phase Portraits by plotting the gradients against the co-ordinates using the " meshgrid ", and " quiver3 " functions.
" meshgrid " is used to generate the 2D or 3D grids and " quiver " or " quiver3 " is used to place arrows at these co-ordinates.
Creating meshgrid:
[x1,y1,z1] = meshgrid(-2:0.2:2,-2:0.2:2,-2:0.2:2);
u = zeros(size(x1));
v = zeros(size(y1));
w = zeros(size(z1));
t=0;
for i = 1:numel(x1)
Yprime = cluster(t,[x1(i); y1(i); z1(i)]);
u(i) = Yprime(1);
v(i) = Yprime(2);
w(i) = Yprime(3);
end
Plotting:
quiver3(x1,y1,z1,u,v,w); figure(gcf)
Please refer to the documentation for "quiver3" here for examples.

Aishah Malek

Aishah Malek (view profile)

on 2 Sep 2018
This makes sense , thankyou for you help