求解微分方程的时候，如果不能将求出结果的表达式，则可以对利用数值积分对微分方程求解，获取数值解。龙格库塔方法是的一种常用的数值解法。今天主要介绍二阶和四阶的龙格库塔方法实例求解微分方程。
二阶的龙格库塔公式为：
四阶的龙格库塔公式为：
实例
求解微分方程dy/dt=-y+t+1,y(0)=1,t的取值为0到2，步长h=0.1,用欧拉法、二阶和四阶的龙格库塔方法求解微分方程并将结果与y(t)=exp(-t)+t比较。
主程序
clc; clear all; close all; h = 0.1;%步长 y0 = 1;%初值 t = 0:h:1;%x范围 y = ex我爱线报网p(-t)+t;%真解 n = length(t); numy = zeros(1,n); f2 = -y0; numy(1) = y0; %二阶龙格库塔参数初始化 runge2_y = zeros(1,n); runge2_y0 = y0; runge2_y(1) = y0; k1 = -y0+1; k2 = -(y0+h*k1); %四阶龙格库塔参数初始化 range4_y = zeros(1,n); range4_y(1) = y0; range4_y0 = y0; fori =我爱线报网2:n %欧拉法计算 numy(i) = euler1(y0,h,f2); y0 = numy(i); f2 = f1(t(i),y0); %二阶龙格库塔计算 runge2_y(i) = runge2(k1,k2,h,runge2_y0); runge2_y0 = runge2_y(i); k1 = ff(t(i),runge2_y0); k2 = ff(t(i),runge2_y(i)+h*k1); %四阶龙格库塔计算 runge4_y(i)=runge4(t(i),range4_y0,h); r我爱线报网ange4_y0 = runge4_y(i); end %绘图 figure; plot(t,y,r-,linewidth,1); hold on; plot(t,numy,b-,linewidth,1); plot(t,runge2_y,c-,linewidth,1); plot(t,runge4_y,g-,linewidth,1); xlabel(t); grid on; title(Euler法,二阶和四阶龙格库塔法求系统的输出响应); ylabel(输出响应y(t)); legend(真我爱线报网解,Euler法,二阶龙格库塔法,四阶龙格库塔法,location,southeast); text(0.08,0.5,四阶龙格库塔法); text(0.5,1,二阶龙格库塔法); text(0.1,1.1,真解y(t)); text(0.9,1.3,Euler法); wucha_euler = (numy-y).^2; wucha_2 = (runge2_y-y).^2; wucha_4 = (runge4_y-y).^2; disp(Euler法误差平方:); wucha_euler disp(二阶龙格库塔法误差平方:); wuch我爱线报网a_2disp(四阶龙格库塔法误差平方:); wucha_4自定义函数1
function y = euler1(y0,h,f2) %%输入参数 y0表示 t=0时 y的取值 即初值 %h表示步长 %f 表示函数值 %输出y表示方程的响应y y=y0+h*f2; end自定义函数2
function f= f1(x,y) f = -y+x+1;%微分方程 dy/dt=-y+t+1 初值y(0)=1 微分方程右边的剩余部分构成的函数 end自定义函数3
function f= ff(x,y) f = -y+x+1;%微分方程 dy/dt=-y+t+1 初值y(0)=1微分方程右边的剩余部分构我爱线报网成的函数end自定义函数4
function yk_1 = runge2(k1,k2,h,yk) yk_1 = yk+h*(k1+k2)/2; end自定义函数5
function yk_1 = runge4(t,y0,h) k1 = h*ff(t,y0); k2 = h*ff(t,y0+k1/2); k3 = h*ff(t,y0+k2/2); k4 = h*ff(t,y0+k3); yk_1 = y0+(k1+2*k2+2*k3+k4)/6; end结果
结果 Euler法误差平方: wucha_euler = 0 0.0110 0.0097 0.0086 0.0076 0.0067 0.0058 0.0051 0.0044 0.0039 0.0034 二阶龙格库塔法误我爱线报网差平方: wucha_2 = 0 0.0030 0.0030 0.0029 0.0029 0.0029 0.0029 0.0028 0.0028 0.0028 0.0028 四阶龙格库塔法误差平方: wucha_4 = 1.0000 0.0000 0.0001 0.0002 0.0003 0.0004 0.0005修改时间t为0到10的结果为
作 者 | 郭志龙
编 辑 | 郭志龙
校 对 | 郭志龙本文内容来源于网络，仅供参考学习，如内容、图片有任何版权问题，请联系处理，24小时内删除。
给力项目：线报网会员可免费下载 加入会员友情提醒： 请尽量登录购买，防止付款了不发货！
QQ交流群：226333560 站长微信：qgzmt2
温馨提示：本站提供的一切软件、教程和内容信息都来自网络收集整理，仅限用于学习和研究目的；不得将上述内容用于商业或者非法用途，否则，一切后果请用户自负，版权争议与本站无关。用户必须在下载后的24个小时之内，从您的电脑或手机中彻底删除上述内容。如果您喜欢该程序和内容，请支持正版，购买注册，得到更好的正版服务。我们非常重视版权问题，如有侵权请邮件与我们联系处理。敬请谅解！
