当前位置:首页 » 《随便一记》 » 正文

Matlab之四阶龙格—库塔法方法:解常微分初值问题

0 人参与  2023年04月01日 11:21  分类 : 《随便一记》  评论

点击全文阅读


目录

1. 题目

2. 算法原理

3. 代码

4. 结果

4.1 运行结果

4.2 结果分析


【若觉文章质量良好且有用,请别忘了点赞收藏加关注,这将是我继续分享的动力,万分感谢!】

直接通过解题的方式进行学习,代入感更强

1. 题目

用经典四阶龙格库塔方法对初值问题\left\{\begin{array}{l} y^{\prime}=-20 y \\ y(0)=1 \end{array}\right.,步长分别取求解,观察稳定区间的作用。

2. 算法原理

某些常微分方程有解析解,但大多数都没有,因此需要进行数值解计算。

龙格—库塔法是利用f(x,y)在某些特殊点上的函数值的线性组合,来估算高阶单步法的平均斜率。

经典的龙格—库塔法是四阶的,也就是在(x_i,x_{i+1})中用四个点处的斜率来估计其平局斜率,构成四阶龙格—库塔公式

其准确解y(x)在一系列点xi处y(xi)的近似值yi的方法,yi称为数值解。经典的四阶龙格库塔法方程如下:

y_{i+1}=y_i+c_1 K_1+c_2 K_2+c_3 K_3+c_4 K_4

其中:

\left\{\begin{array}{l} K_1=h f\left(x_i, y_i\right) \\ K_2=h f\left(x_i+a_2 h, y_i+b_{21} K_1\right) \\ K_3=h f\left(x_i+a_3 h, y_i+b_{31} K_1+b_{32} K_2\right) \\ K_4=h f\left(x_i+a_4 h, y_i+b_{41} K_1+b_{42} K_2+b_{43} K_3\right) \end{array}\right.

其中的各个参数具体如下:

\begin{gathered} a_2=a_3=b_{21}=b_{32}=\frac{1}{2} \\ b_{31}=b_{41}=b_{42}=0 \\ a_4=b_{43}=1 \\ c_1=c_4=\frac{1}{6} \\ c_2=c_3=\frac{1}{3} \end{gathered}

其整合之后为:

\left\{\begin{array}{l} y_{i+1}=y_i+\frac{1}{6}\left(K_1+2 K_2+2 K_3+K_4\right) \\ K_1=h f\left(x_i, y_i\right) \\ K_2=h f\left(x_i+\frac{h}{2}, y_i+\frac{1}{2} K_1\right) \\ K_3=h f\left(x_i+\frac{h}{2}, y_i+\frac{1}{2} K_2\right) \\ K_4=h f\left(x_i+h, y_i+K_3\right) \end{array}\right.

其中h为步长。

3. 代码

clear;clc;for step = [0.1, 0.2]        x_0 = 0;    y_0 = 1;    num = floor(1/step);    n = 1;    X_output = [0];    Y_output = [1];    disp("y'= -20 * y")    while n <= num        x_1 = x_0 + step;        K_1 = step * fun(x_0,y_0);        K_2 = step * fun(x_0 + step/2, y_0 + K_1/2);        K_3 = step * fun(x_0 + step/2, y_0 + K_2/2);        K_4 = step * fun(x_0 + step, y_0 + K_3);        y_1 = y_0 + (K_1 + 2 * K_2 + 2 * K_3 + K_4) / 6 ;        X_output = [X_output x_1];        Y_output = [Y_output y_1];        x_0 = x_1;        y_0 = y_1;        n = n + 1;    end    figure()    plot(X_output,Y_output)    xlabel('x')    ylabel('y')    title(['Runge-Kutta4阶,步长为:', num2str(step)])    X_output    Y_output    clear X_output Y_outputend[x,y] = ode45('fun', [0:1], 1);figure()plot(x,y)xlabel('x')ylabel('y')title('自带函数求解结果')function dy = fun(x, y)dy = - 20*y;end

4. 结果

4.1 运行结果

Step = 0.1 时X_output =0    0.1000    0.2000    0.3000    0.4000    0.5000    0.6000    0.7000    0.8000    0.9000          1.0000Y_output =1.0000    0.3333    0.1111    0.0370    0.0123    0.0041    0.0014    0.0005    0.0002    0.0001    0.0000Step = 0.2时X_output =0    0.2000    0.4000    0.6000    0.8000    1.0000Y_output =1           5          25         125         625        3125

 

 

4.2 结果分析

用经典四阶龙格库塔方法求解,其求解结果与设置得步长有很大的相关性,步长设置合适时,其求解情况与真实值基本一致,趋于稳定。但步长加大时,其求解值与真实值相差太大。

​​​​​​​【若觉文章质量良好且有用,请别忘了点赞收藏加关注,这将是我继续分享的动力,万分感谢!】


点击全文阅读


本文链接:http://zhangshiyu.com/post/57757.html

<< 上一篇 下一篇 >>

  • 评论(0)
  • 赞助本站

◎欢迎参与讨论,请在这里发表您的看法、交流您的观点。

关于我们 | 我要投稿 | 免责申明

Copyright © 2020-2022 ZhangShiYu.com Rights Reserved.豫ICP备2022013469号-1