필터 지우기
필터 지우기

How to indentify the variables when plotting after solving a coupled system of ODEs?

조회 수: 1 (최근 30일)
Hello, i had a system of coupled ODEs and i think i have solved it, the problem is that the variable Y which gives the output, i dont know how to tell what column is the response to the x(t), y(t), z(t) or their first derivatives. Help would be much appreciated. Here is the code
syms x(t) y(t) z(t) p(t) Y
p=sin((0.5+3.5*t)*t);
D1y = diff(y,t);
D2y = diff(y,t,2);
D1x = diff(x,t);
D2x = diff(x,t,2);
D1z = diff(z,t);
D2z = diff(z,t,2);
eqns=[D2x+2*x-1*y==p,D2y+(0.5*10*exp(-10*t))*D1y-x+2*y-z==0,D2z-y+2*z==0];
yode = odeToVectorField(eqns);
Yodefcn = matlabFunction(yode, 'Vars',[t Y]);
tspan = (0:0.001:10);
Y0 = [0 0,0 0,0 0];
[T,Y] = ode45(Yodefcn, tspan, Y0);
plot(T, Y)

채택된 답변

Star Strider
Star Strider 2017년 3월 18일
The easiest way is to ask for a second output from your odeToVectorField call. It tells you the substitutions the function made to define each element in the returned column:
[yode,Sbs] = odeToVectorField(eqns);
here, the substitutions are:
Sbs =
y
Dy
x
Dx
z
Dz
The columns of ‘Yodefcn’ exactly reproduces those returned by odeToVectorField, and ode45 returns them in that order.
That can be used with some other assignments to automatically produce a legend that labels every line in the plot.
The Code
syms x(t) y(t) z(t) p(t) Y
p=sin((0.5+3.5*t)*t);
D1y = diff(y,t);
D2y = diff(y,t,2);
D1x = diff(x,t);
D2x = diff(x,t,2);
D1z = diff(z,t);
D2z = diff(z,t,2);
eqns=[D2x+2*x-1*y==p,D2y+(0.5*10*exp(-10*t))*D1y-x+2*y-z==0,D2z-y+2*z==0];
[yode,Sbs] = odeToVectorField(eqns);
Yodefcn = matlabFunction(yode, 'Vars',[t Y]);
tspan = (0:0.001:10);
Y0 = [0 0,0 0,0 0];
[T,Y] = ode45(Yodefcn, tspan, Y0);
plot(T, Y)
Lgndc = sym2cell(Sbs);
Lgnds = regexp(sprintf('%s\n', Lgndc{:}), '\n', 'split');
legend(Lgnds(1:end-1), 'Location','N')
The Plot

추가 답변 (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