How to write a phase plane in matlab?
조회 수: 48 (최근 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
채택된 답변
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.
추가 답변 (0개)
참고 항목
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!