Matlab求微分方程的数值解

迷南。 2022-11-26 04:59 165阅读 0赞

注：首先计算微分方程的解析解，如果发现没有解析解，再用数值解

**一、Matlab中求微分方程的数值解函数**

[x,y]=solver('f',ts,x0,options)

1）x代表自变量

2）y代表函数值

3）solver代表求解函数，常用的为ode45->函数图像变化较为平稳，ode15s->函数图像存在突变（一个运行时间长时，换另一个）

![20200810163804170.png][]

![20200810163819291.png][]

4）ts=\[t0,t1\]，代表自变量函数的初值和终值，也可以指定以一个向量

5）x0表示函数的初始值

6）options（可选）：options=odeset('reltol','rt','abstol',at)，rt，at分别为设定的相对误差，默认为10e-3，10e-6

7）f表示由待解微分方程编写的m函数文件名

**二、一阶微分方程的数值解**

1、求解步骤：

1）化为标准型（等号左边为微分）如下图

![20200810164138427.png][]

2）带入函数求解

2、例子

例1（求方程）：

题目：![y-y\{^'\}=2x][y-y_2x],初始值为y(0)=3,x范围为\[0,2\]

解：

1）化为标准型：![y\{^'\}=y-2x][y_y-2x]

2）matlab求解如下：

2.1）新建一个m文件，用来编辑f函数

function dy = f(x,y) % dy是标准化后左边的式子，x是式子右边的自变量，y是因变量
    dy = y - 2*x;
    end

2.2）新建一个m文件用来计算结果

[x y] = ode45('f',[0,2],3); %x是自变量，y是因变量
    figure(1)
    plot(x,y,'*-')
    axis([0 2 -inf inf])
    grid on

3）运行结果

![watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3FxXzM2ODkwMzcw_size_16_color_FFFFFF_t_70][]

例2：（求方程组）

题目：

![\\left\\\{\\begin\{matrix\} y\{\_1^'\}=y\{\_2\}y\{\_3\}\\\\ y\{\_2^'\}=-y\{\_1\}y\{\_3\}\\\\ y\{\_3^'\}=-0.51y\{\_1\}y\{\_2\} \\end\{matrix\}\\right.][left_begin_matrix_ y_1_y_2_y_3_ y_2_-y_1_y_3_ y_3_-0.51y_1_y_2_ _end_matrix_right.]，自变量![x\\in\[0.4\\pi\]][x_in_0.4_pi]，初始值![y\{\_1\}(0)=0,y\{\_2\}(0)=y\{\_3\}(0)=1][y_1_0_0_y_2_0_y_3_0_1]

解：

1）化为标准型（题目已经为标准型）

2）matlab代码

2.1）新建一个m文件，用来编辑f函数

function dy = f(x,y) % dy是标准化后左边的式子，x是式子右边的自变量，y是因变量
        dy = zeros(3,1);
        dy(1) = y(2)*y(3);
        dy(2) = -y(1)*y(3);
        dy(3) = -0.51*y(1)*y(2);
    end

2.2） 新建一个m文件用来计算结果

[x, y] = ode45('f', [0 4*pi], [0 1 1]);  % 这里的y是一个有3列的矩阵哦！x是自变量，y是因变量
    plot(x, y(:,1), 'o', x, y(:,2), '*', x, y(:,3), '+') 
    legend('y1','y2','y3')  % 加上标注
    axis([0, 4*pi, -inf, +inf])  % 设置横坐标范围为0-4pi，纵坐标范围不需要设置，写成-inf到+inf

3）结果

![watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3FxXzM2ODkwMzcw_size_16_color_FFFFFF_t_70 1][]

**三、高阶微分方程的数值解**

1、步骤

1）标准化

高阶方程转化为一阶方程组合

2）按一阶求解

2、例子

![2020081017484445.png][]

1）标准化

![20200810174917230.png][]

2）matlab求解

![watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3FxXzM2ODkwMzcw_size_16_color_FFFFFF_t_70 2][]

[20200810163804170.png]: /images/20221124/7f7e5dd873a942afbb2ed6c026c6358a.png
[20200810163819291.png]: /images/20221124/26aee4cfac344686b94ce54b21488844.png
[20200810164138427.png]: /images/20221124/694cc8620b61498296003e6622b04548.png
[y-y_2x]: https://private.codecogs.com/gif.latex?y-y%7B%5E%27%7D%3D2x
[y_y-2x]: https://private.codecogs.com/gif.latex?y%7B%5E%27%7D%3Dy-2x
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3FxXzM2ODkwMzcw_size_16_color_FFFFFF_t_70]: /images/20221124/8f9ae1d07ba4486cbba47e61f100ae05.png
[left_begin_matrix_ y_1_y_2_y_3_ y_2_-y_1_y_3_ y_3_-0.51y_1_y_2_ _end_matrix_right.]: https://private.codecogs.com/gif.latex?%5Cleft%5C%7B%5Cbegin%7Bmatrix%7D%20y%7B_1%5E%27%7D%3Dy%7B_2%7Dy%7B_3%7D%5C%5C%20y%7B_2%5E%27%7D%3D-y%7B_1%7Dy%7B_3%7D%5C%5C%20y%7B_3%5E%27%7D%3D-0.51y%7B_1%7Dy%7B_2%7D%20%5Cend%7Bmatrix%7D%5Cright.
[x_in_0.4_pi]: https://private.codecogs.com/gif.latex?x%5Cin%5B0.4%5Cpi%5D
[y_1_0_0_y_2_0_y_3_0_1]: https://private.codecogs.com/gif.latex?y%7B_1%7D%280%29%3D0%2Cy%7B_2%7D%280%29%3Dy%7B_3%7D%280%29%3D1
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3FxXzM2ODkwMzcw_size_16_color_FFFFFF_t_70 1]: /images/20221124/98267ceb89854f05bf70e904a1cc50f1.png
[2020081017484445.png]: /images/20221124/f45f5d1c83fb4d14bc9a3ae236f9d290.png
[20200810174917230.png]: /images/20221124/1198c6d71e5344d18e7c69ebb9c67a85.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3FxXzM2ODkwMzcw_size_16_color_FFFFFF_t_70 2]: /images/20221124/980b4168ca6945829595cdbeda182439.png