定步长四阶经典公式解决数值积分 - 数学建模社区-数学中国

%四阶经典公式,微分方程为f.m; U' Q/ r, b3 d- k- g9 k! K" q
9 n3 k9 T( U+ Z9 L3 u
if exist('f.m')==0 %在星号处输入文件名(把星号改为文件名)
0 _; h5 W/ W6 y. z0 r. f* @5 M
disp('没有为方程创建名为f.m的函数文件,请参照下例建立它');, W+ i6 p; P i( W* F# c
disp('function z=f(x,y)');; X X! J7 a) x, ^: Y' V8 c
disp('z=y-2*x/y;');8 r5 L5 G: [2 J1 d% ^# F
disp('并将该文件保存在work文件夹下');
3 P ^3 b3 a* h/ B4 T- W& s' N
end / P, ~4 x& r" ]( P. p8 t
6 i M7 n l: m, S
X1=input('请输入求解区间的左端点X1=');$ n$ A" j- ^8 h& P
Y1=input('请输入微分方程的初始条件Y1=（X=X1时Y的值）');
" D. p, ?2 B* Q
Xn=input('请输入求解区间的右端点Xn=');
% E7 W, I* t- q
h=input('请输入求解步长h=');) t9 k. c U- K: l1 a; e
5 C% W) k. i$ S) x- U7 e
X=X1;( o" B2 a9 ~! \* H J# k, a
Y=Y1; %运算初始点* y I$ @: t7 t3 m/ [
n=0; %节点序号变量置零4 D& C5 ]5 E/ [6 e& m9 g0 u5 m. H- b# T
& D* I0 Z" m3 Y, s: n) O3 c" p
while X<=Xn-h I- V( Z% G7 Z+ ~
K1=f(X,Y);
% g9 n/ k, q s
K2=f(X+h/2,Y+K1*h/2);, a& A9 n. d) @7 y! g( ?8 Y# |
K3=f(X+h/2,Y+K2*h/2);0 y! d- k4 e5 L, |: O
K4=f(X+h,Y+K3*h);
1 `- _( r- }) X ~8 c0 r
X=X+h;
A# L0 v, `/ l9 J3 f( Q
Y=Y+h*(K1+2*K2+2*K3+K4)/6; %四阶标准的龙格-库塔公式" F* D# J; Q) s& [5 [
n=n+1; %节点序号加1
' ~8 A8 g$ T/ p6 T3 _
- Z$ ~! m1 [9 @: s8 T' N: f* n7 [
fprintf('第%d个点的计算结果为X=%10.8f,Y=%10.8f\n',n,X,Y);/ [& d3 B+ G* g( }
plot(X,Y,'o')! E @+ m2 }% y
hold on. v. m+ d& f1 j" t% B
end

复制代码