急!!!求matlab 用四阶龙格-库塔法求解常微分方程

发布网友发布时间：2022-05-05 06:43

我来回答

共1个回答

热心网友时间：2023-10-13 13:50

建立.m文件

---------------------------------------------

function theta=dan(t,X)

x=X(1);

dx=X(2);

ddx=-sin(x);

theta=[dx;ddx];

----------------------------------------------

命令窗口输入

>> [t,Y]=ode45(@dan,[0 6],[pi/3 -1/2]);

>> plot(t,Y(:,1),'ro-',t,Y(:,2),'bv-');

>> legend('\theta-t','d\theta-t')

自编龙格库塔

------------------------------------------------

function [y,z]=Runge_kutta(a,b,y0,z0,h)

x=a:h:b;

y(1)=y0;

z(1)=z0;

n=(b-a)/h+1;

for i=2:n

K(1,1)=f1(x(i-1),y(i-1),z(i-1));

K(2,1)=f2(x(i-1),y(i-1),z(i-1));

K(1,2)=f1(x(i-1)+h/2,y(i-1)+K(1,1)*h/2,z(i-1)+K(2,1)*h/2);

K(2,2)=f2(x(i-1)+h/2,y(i-1)+K(1,1)*h/2,z(i-1)+K(2,1)*h/2);

K(1,3)=f1(x(i-1)+h/2,y(i-1)+K(1,2)*h/2,z(i-1)+K(2,2)*h/2);

K(2,3)=f2(x(i-1)+h/2,y(i-1)+K(1,2)*h/2,z(i-1)+K(2,2)*h/2);

K(1,4)=f1(x(i-1)+h,y(i-1)+K(1,3)*h,z(i-1)+K(2,3)*h);

K(2,4)=f2(x(i-1)+h,y(i-1)+K(1,3)*h,z(i-1)+K(2,3)*h);

y(i)=y(i-1)+h/6*(K(1,1)+2*K(1,2)+2*K(1,3)+K(1,4));

z(i)=z(i-1)+h/6*(K(2,1)+2*K(2,2)+2*K(2,3)+K(2,4));

end

y(2)

plot(y,'r') %θ-t图

hold on

plot(f1(x,y,z),'g') % dθ-t图

hold on

plot(f2(x,y,z),'b-')%d2θ-t图

%f1.m

function f1=f1(x,y,z)

f1=z;

%f2.m

function f2=f2(x,y,z)

f2=-sin(y);

-----------------------------------------

Runge_kutta(0,6,pi/3,-1/2,0.02)

热心网友时间：2023-10-13 13:50

建立.m文件

---------------------------------------------

function theta=dan(t,X)

x=X(1);

dx=X(2);

ddx=-sin(x);

theta=[dx;ddx];

----------------------------------------------

命令窗口输入

>> [t,Y]=ode45(@dan,[0 6],[pi/3 -1/2]);

>> plot(t,Y(:,1),'ro-',t,Y(:,2),'bv-');

>> legend('\theta-t','d\theta-t')

自编龙格库塔

------------------------------------------------

function [y,z]=Runge_kutta(a,b,y0,z0,h)

x=a:h:b;

y(1)=y0;

z(1)=z0;

n=(b-a)/h+1;

for i=2:n

K(1,1)=f1(x(i-1),y(i-1),z(i-1));

K(2,1)=f2(x(i-1),y(i-1),z(i-1));

K(1,2)=f1(x(i-1)+h/2,y(i-1)+K(1,1)*h/2,z(i-1)+K(2,1)*h/2);

K(2,2)=f2(x(i-1)+h/2,y(i-1)+K(1,1)*h/2,z(i-1)+K(2,1)*h/2);

K(1,3)=f1(x(i-1)+h/2,y(i-1)+K(1,2)*h/2,z(i-1)+K(2,2)*h/2);

K(2,3)=f2(x(i-1)+h/2,y(i-1)+K(1,2)*h/2,z(i-1)+K(2,2)*h/2);

K(1,4)=f1(x(i-1)+h,y(i-1)+K(1,3)*h,z(i-1)+K(2,3)*h);

K(2,4)=f2(x(i-1)+h,y(i-1)+K(1,3)*h,z(i-1)+K(2,3)*h);

y(i)=y(i-1)+h/6*(K(1,1)+2*K(1,2)+2*K(1,3)+K(1,4));

z(i)=z(i-1)+h/6*(K(2,1)+2*K(2,2)+2*K(2,3)+K(2,4));

end

y(2)

plot(y,'r') %θ-t图

hold on

plot(f1(x,y,z),'g') % dθ-t图

hold on

plot(f2(x,y,z),'b-')%d2θ-t图

%f1.m

function f1=f1(x,y,z)

f1=z;

%f2.m

function f2=f2(x,y,z)

f2=-sin(y);

-----------------------------------------

Runge_kutta(0,6,pi/3,-1/2,0.02)

急!!!求matlab 用四阶龙格-库塔法求解常微分方程

K(1,4)=f1(x(i-1)+h,y(i-1)+K(1,3)*h,z(i-1)+K(2,3)*h); K(2,4)=f2(x(i-1)+h,y(i-1)+K(1,3)*h,z(i-1)+K(2,3)*h); y(i)=y(i-1)+h/6*(K(1,1)+2*K(1,2)+2*K(1,3)+K(1,4)); z(i)=z(i-1)+h/6*(K(2,1)+2*K(2,2)+2*K(2,3)+...

随机（正弦）振动

正弦振动多用于找出产品设计或包装设计的脆弱点。看在哪一个具体频率点响应最大（共振点）；正弦振动在任一瞬间只包含一种频率的振动，而随机振动在任一瞬间包含频谱范围内的各种频率的振动。由于随机振动包含频谱内所有的频率，所以样品上的共...

Matlab四阶龙格库塔法求解常微分方程

用Matlab四阶龙格库塔法求常微分方程可以按照以下方法去实现。1、首先建立自定义微分方程函数 function f = ode_fun(x,y)f=y+2*x/y^2;end 2、然后用四阶龙格库塔法求其数值解 figure(2)y0=[1]; %初值y（0）=1 h=0.1;a=0;b=5;[x,y] = runge_kutta(@(x,y)ode_fun(x,y),y0...

求编程达人帮忙用matlab编程用龙格库塔方法解微分方程

功能：用四阶Runge-Kutta 法求解常微分方程 --- function R=Rungkuta4(f, a, b, n, ya)f:微分方程右端函数句柄 a,b:自变量取值区间的两个端点 n:区间等分的个数 ya:函数初值y(a)R=[x',y']:自变量X 和解Y 所组成的矩阵 h=(b-a)/n;x=zeros(,n+1);y=zeros(1,n+1);x=a:h...

四阶常微分龙格—库塔法求解常微分方程的初值问题。 Matlab程序算例:用...

clear, clc%清除内存中的变量 %数值解y=inline('x*exp(x)+2*x-1');y(1) % 四阶龙格库塔法y0=[-1 3 2];[x1,y1] = ode45(@fun,[0,1],y0); y1(end,1) % y(1)的值function dy=fun(x,y)dy=zeros(3,1);dy(1)=y(2);dy(2)=y(3);dy(3)=y(3)+y(2)-y(1)+...

怎么用MATLAB写龙格库塔四阶方程

四阶龙格-库塔是求解常微分方程（常微分方程组）精度最高的一种数值方法。四阶龙格-库塔迭代公式为根据四阶龙格-库塔迭代公式，可用matlab软件编写得到如下计算代码：调用方法：clc;close all;y0=1;h=0.1;a=0;b=1;[x,y] = runge_kutta(@func,y0,h,a,b);disp(' x y')[x',y']...

matlab求解常微分方程/偏微分方程

在MATLAB中，利用RK方法（四阶龙格-库塔方法）求解非刚性常微分方程，首先定义函数，然后使用ode45函数求解。例如，对于方程f = -2y + 2x^2 + 2x，通过以下MATLAB代码实现求解：// doty.m function f = doty(x, y)f = -2*y + 2*x^2 + 2*x;end // main.m [x, y] = ode45('doty...

MATLAB编写程序用四阶龙格库塔法求解常微分方程组,自己写了算的出错...

function [t,x]=rk4(funname,t0,t1,x0,dt)t=[];x=[];while t0<t1 t=[t;t0];x=[x;x0];k1=funname(t0,x0);t0=t0+dt/2;k2=funname(t0,x0+dt*k1/2);k3=funname(t0,x0+dt*k2/2);t0=t0+dt/2;k4=funname(t0,x0+dt*k3);x0=x0+dt/6*(k1+2*k2+2*k3+k4);t...

急求MATLAB编程源代码用四阶龙格库塔法解如下微分方程 y'=y-2x/y...

以下另存为文件 myrk4.mfunction [x,y]=myrk4(ufunc,y0,h,a,b)%参数: 函数名称，初始值向量，步长，时间起点，时间终点n=floor((b-a)/h);%求步数 x(1)=a;%时间起点 y(:,1)=y0;%赋初值%按龙格库塔方法进行求解 for ii=1:n x(ii+1)=x(ii)+h; k1=ufunc(x(ii),y...

急急急!用MATLAB按四阶龙格库塔法求解微分方程组,急用于毕业设计!

在Matlab下输入：edit，然后将下面两行百分号之间的内容，复制进去，保存 function dxdt=ode_Miss_ghost(t,x)分别用x(1),x(2),x(3),x(4)代替N1,P1,N2,P2 N1=x(1);P1=x(2);N2=x(3);P2=x(4);K=2;tau_c=3e-9;tan_p=6e-12;beta =5e-5;delta=0.692;eta =0.0001;fm ...

matlab龙格库塔法解微分方程

function [Y] = RK45(t,X,f,h)K1=f(t,X);K2=f(t+h/2,X+h/2*K1);K3=f(t+h/2,X+h/2*K2);K4=f(t+h,X+h*K3);Y=X+h/6*(K1+2*K2+2*K3+K4);end 以上是4阶龙格库塔法的代码:自己写函数，存为f.m function dxdt = f (t,x)dxdt(1)=exp(x(1)*sin(t))+...