필터 지우기
필터 지우기

How to write a phase plane in matlab?

조회 수: 48 (최근 30일)
Moira
Moira 2020년 4월 24일
댓글: Moira 2020년 4월 30일
I'm trying to write a phase plane for ODE.
The equation goes like: da/dt = -2*k1*a^2 + k2*b; db/dt = 2*k1*a^2 - k2*b, where k1 & k2 are const.
I want it to have arrows and lines and stuff like this:
I'm still very new to matlab so it's quite painful to read long long unfamiliar codes. If you can put detailed annotations with codes I'd appreciate it very much!
Thank you!
  댓글 수: 3
BALAJI KARTHEEK
BALAJI KARTHEEK 2020년 4월 24일
use quiver for arrows
Moira
Moira 2020년 4월 24일
The initial condition: a = 1, b = 0

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

채택된 답변

Raunak Gupta
Raunak Gupta 2020년 4월 30일
편집: Raunak Gupta 2020년 4월 30일
Hi,
I have compiled some code which essentially plot the same figure you have. For plotting the straight lines, you must choose the starting points in either a or b. I am choosing those in variable y10 and y20 as seen from the graph. Finally, for getting the orange line you need to check where both da/dt and db/dt are zero. Since here both equations are the same, I have solved one and plotted the result at the end. You may put the title and other legend text as per required.
k1 = 3;
k2 = 1;
% Y is essentially [da/dt,db/dt]
f = @(t,Y) [-2*k1*(Y(1))^2 + k2*Y(2); 2*k1*(Y(1))^2 - k2*Y(2)];
% for plotting the vector field using quiver
y1 = linspace(0,1,21);
y2 = linspace(0,1,21);
[x,y] = meshgrid(y1,y2);
u = zeros(size(x));
v = zeros(size(x));
t=0;
for i = 1:numel(x)
% Calculating gradient value for a and b
Yprime = f(t,[x(i); y(i)]);
u(i) = Yprime(1);
v(i) = Yprime(2);
end
quiver(x,y,u,v,'r');
hold on
% Plotting the integrated value
for y10 = [0.3 0.5 0.7 0.8 1]
[ts,ys] = ode45(f,[0,50],[y10;0]);
plot(ys(:,1),ys(:,2),'k-')
end
for y20 = [0.35 0.65 0.8]
[ts,ys] = ode45(f,[0,50],[0;y20]);
plot(ys(:,1),ys(:,2),'k-')
end
grid on
% solve the equation for which gradient is zero and plotting it
syms a b
eqn = -2*k1*a^2 + k2*b == 0;
fimplicit(eqn, [0 1 0 1]);
hold off
Hope it is understandable.
  댓글 수: 1
Moira
Moira 2020년 4월 30일
It's so helpful!!! Thank you so much!!!

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

추가 답변 (0개)

카테고리

Help CenterFile Exchange에서 2-D and 3-D Plots에 대해 자세히 알아보기

Community Treasure Hunt

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

Start Hunting!

Translated by