基于MATLAB实现固定步长的四阶龙格库塔方法求解常微分方程初值问题

本文介绍使用四阶龙格-库塔(RK)方法求解特定微分方程的MATLAB代码实现,并通过实例展示如何计算一暂态电路问题的数值解。同时,给出了该问题的精确解析解,并比较了数值解与解析解的差异。

摘要生成于 C知道 ,由 DeepSeek-R1 满血版支持, 前往体验 >

代码

构造函数

function [yn,xn] = RKorder4(f,a,b,h,y0)
%a,b为上下界;h为步长;y0为初值。
%返回yn、xn向量,并绘制曲线
yn=zeros(1,((b-a)/h)+1);
yn(1)=y0;
xn=a:h:b;
for i=1:((b-a)/h)
    k1=h*f(xn(i),yn(i));
    k2=h*f(xn(i)+h/2,yn(i)+k1/2);
    k3=h*f(xn(i)+h/2,yn(i)+k2/2);
    k4=h*f(xn(i)+h,yn(i)+k3);
    yn(i+1)=yn(i)+(1/6)*(k1+2*k2+2*k3+k4);
end
end

实例

对于某暂态电路的初值问题:
{didt=622sin314t−20ii(0)=0 \left\{\begin{matrix}\frac{\mathrm{d} i}{\mathrm{d} t}=622sin314 t-20i \\ i(0)=0 \end{matrix}\right. {dtdi=622sin314t20ii(0)=0
 f(tn,in)=622sin314t−20i\ f(t_n,i_n)=622sin314 t-20i f(tn,in)=622sin314t20i
取h=0.001,在[0,0.01]的时间区间内求解微分方程,得出电流值。
在命令窗口输入

a=0;b=0.01;h=0.001;
y0=0;
f=@(t,i) 622*sin(314*t)-20*i;
[yn,xn] = RKorder4(f,a,b,h,y0);
plot(xn,yn);grid on

在这里插入图片描述
又由暂态电路理论分析得,此电路的解析解为:
i(t)=4882724749(e−20t+10157sin314t−cos314t) i(t)=\frac{48827}{24749}\left( e^{-20t} +\frac{10}{157}sin314t-cos314t\right) i(t)=2474948827(e20t+15710sin314tcos314t)
输入命令,将解析解和数值解的曲线绘制在同一图内:

i=@(t) (48827/24749).*(exp(-20.*t)+(10/157)*sin(314.*t)-cos(314.*t));
plot(linspace(0,0.01,1000),i(linspace(0,0.01,1000)),'r',xn,yn,'b');grid on

在这里插入图片描述

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值