在 Matlab 中求解 ODE

Resolving ODE in Matlab

我正在尝试编写一个程序来解决以下系统 ode23:

2y’ + z’–y + 2z = 0
y’ + 3z’ –3y +z = 0

初始值:

y(0) = 1
z(0) = 0

及解析解:

y= cos(x)
z= sin(x)

但是当我更改变量 4 时:

function dy = eqdif2(t,y)
%2y’ + z’ –y + 2z = 0;
%y’ + 3z’ –3y +z = 0
% y(0) = 1, Z(0) = 0
% y=y(1), z=y(2), y'=y(3), z'=y(4)
dy = [-2*y(3)+y(1)-2*y(2);3*y(1)-y(2)-3*y(4)];

我遇到了 ode23 的问题,仅报告了 2 个解决方案:

clc,clear;
yp = [1 0];  %initial values
options = odeset('RelTol', 1e-4);
[t,y]= ode23('eqdif2',[0 20],yp,options); 
ya=cos(x);
za=sin(x);
figure;
plot(t,y(:,1),'-');
figure;
plot(t,ya,t,za);

我们开始...

ODE 是线性的,我们知道,有一个代数解等等。要使用 ode23,您需要一个 explicit 集,格式为 -- ps。注意因变量 x 和自变量 t 的向量形式用法:

x'=f(t,x)

在这种情况下,它被简化为琐碎但无用的设置:

x(1)' = x(2) x(2)' = -x(1)

求解为:

f=@(t,x)([-x(2);x(1)]);
y0=[1 0];
% f(0,y0) %Check IC = y'
[t, x] = ode23(f,[0 pi*2],y0);
plot(t,x)

但是,如果我们需要或想要,或者我们被迫直接解决系统问题 - 即隐式-,我们应该将集合定义为 - 同样,g 可以是方程的向量:

g(t,x,x')=0

通过使用 ode15i - i 代表 "implicit":

g=@(t,x,dx)([ 2*dx(1)+dx(2)-x(1)+2*x(2); ...
    dx(1)+3*dx(2)-3*x(1)+x(2)]);
y0=[1;0];
dy0=[0;1]; % kind of magic to find...
% g(0,y0,dy0) %Check IC = 0
[t, x] = ode15i(g, [0 pi*2],y0,dy0);
plot(t,x)

两个选项的解决方案相同...