微分方程组的四阶龙格库塔公式求解matlab版
微分方程组的四阶龙格库塔公式,求解matlab版的代码。下面是我的代码,运行不出来。
function varargout=rungekutta(varargin)
clc,clear
x0=0;xn=1.0;y0=1;h=0.05;%h为步长
[y,x]=rungekutta4(x0,xn,y0,h);
Function z=f(x,y);
z=y-2*x/y;
function [y,x]=rungekutta4(x0,xn,y0,h)
x=x0:h:xn;
n=(xn-x0)/h;
y1=x;
y1(1)=y0;
for i=1:n %龙格库塔法的算法
K1=f(x(i),y1(i));
K2=f(x(i)+h/2,y1(i)+h/2*K1);
K3= f(x(i)+h/2,y1(i)+h/2*K2);
K4= f(x(i)+h,y1(i)+h*K3);
y1(i+1)=y1(i)+h/6*(K1+2*K2+2*K3+K4);
yy1(i)=(1+2*x(i+1)) ^0.5;
error(i)=y1(i+1)-yy1(i);
yy(i+1)=yy1(i);
e(i+1)=error(i);
end
y=y1;
n=(xn-x0)/h;
fprintf(‘i x(i) yi数值解 y(i)真值 误差\n’);
fprintf(‘-------------------------------------------------------------\n’);
for i=1:n
fprintf(‘%2d %12.4f %14.8f %14.8f %12.8f\n’,i, x(i+1), y(i+1), yy(i+1), error(i));
end
plot(x,e,‘r.’);%画出误差曲线
xlable(‘x轴’),ylable(‘误差’);
title(‘步长为0.05时的误差曲线’); 返回小木虫查看更多
一、几个自定义函数位置错误。应为这样放置
1、主程序 ,
function rungekutta()
clc,clear
x0=0;xn=1.0;y0=1;h=0.05;%h为步长
[y,x]=rungekutta4(x0,xn,y0,h);
end
2、 四阶龙格库塔函数程序
function [y,x]=rungekutta4(x0,xn,y0,h)
。。。
end
3、z函数程序
function z=f(x,y);
z=y-2*x/y;
end
二、xlable(‘x轴’),ylable(‘误差’);这句代码中函数书写错误。应为
xlabel('x轴'),ylabel('误差');
三、修改后运行得到如下结果
什么结果?大神。
运行结果:
i x(i) yi数值解 y(i)真值 误差
-------------------------------------------------------------
1 0.0500 1.04880886 1.04880885 0.00000001
2 0.1000 1.09544514 1.09544512 0.00000002
3 0.1500 1.14017546 1.14017543 0.00000004
4 0.2000 1.18321600 1.18321596 0.00000005
5 0.2500 1.22474493 1.22474487 0.00000006
6 0.3000 1.26491113 1.26491106 0.00000007
7 0.3500 1.30384056 1.30384048 0.00000008
8 0.4000 1.34164088 1.34164079 0.00000009
9 0.4500 1.37840498 1.37840488 0.00000011
10 0.5000 1.41421368 1.41421356 0.00000012
11 0.5500 1.44913781 1.44913767 0.00000014
12 0.6000 1.48323985 1.48323970 0.00000015
13 0.6500 1.51657526 1.51657509 0.00000017
14 0.7000 1.54919353 1.54919334 0.00000019
15 0.7500 1.58113904 1.58113883 0.00000021
16 0.8000 1.61245178 1.61245155 0.00000023
17 0.8500 1.64316793 1.64316767 0.00000026
18 0.9000 1.67332034 1.67332005 0.00000028
19 0.9500 1.70293895 1.70293864 0.00000031
20 1.0000 1.73205115 1.73205081 0.00000034
图像上传不上来
,
function rungekutta()
clc,clear
x0=0;xn=1.0;y0=1;h=0.05;%h为步长
[y,x]=rungekutta4(x0,xn,y0,h);
End
function [y,x]=rungekutta4(x0,xn,y0,h)
x=x0:h:xn;
n=(xn-x0)/h;
y1=x;
y1(1)=y0;
for i=1:n %龙格库塔法的算法
K1=f(x(i),y1(i));
K2=f(x(i)+h/2,y1(i)+h/2*K1);
K3= f(x(i)+h/2,y1(i)+h/2*K2);
K4= f(x(i)+h,y1(i)+h*K3);
y1(i+1)=y1(i)+h/6*(K1+2*K2+2*K3+K4);
yy1(i)=(1+2*x(i+1))^0.5;
error(i)=y1(i+1)-yy1(i);
yy(i+1)=yy1(i);
e(i+1)=error(i);
end
function z=f(x,y);
z=y-2*x/y;
end
y=y1;
n=(xn-x0)/h;
fprintf('i x(i) yi数值解 y(i)真值 误差\n');
fprintf('-------------------------------------------------------------\n');
for i=1:n
fprintf(‘%2d %12.4f %14.8f %14.8f %12.8f\n’,i, x(i+1), y(i+1), yy(i+1), error(i));
end
plot(x,e, 'r. ');%画出误差曲线
xlabel('x轴'),ylabel('误差');
title('步长为0.05时的误差曲线');
大神,我修改完,还是运行不出来啊。你能不能直接把你的完整代码发给我,我自己找找错误。
怎么发给你?