定步长四阶经典公式 解决数值积分

2744557306

"定步长四阶经典公式"通常指的是数值积分中的四阶Runge-Kutta方法。这是一种常用的数值解常微分方程（ODE）的方法，其主要思想是通过逐步逼近来估计微分方程的解。
, O% m1 S/ T$ g定步长四阶经典公式是Runge-Kutta方法的一种，其中最常见的是经典的四阶Runge-Kutta方法。对于一个一阶常微分方程
[\frac{dy}{dt} = f(t, y)]
2 ~3 D; K7 G+ w& d1 g4 n5 A这个方法的迭代公式如下：
[k1 = h \cdot f(tn, yn)]
[k2 = h \cdot f(tn + \frac{h}{2}, yn + \frac{k1}{2})]
[k3 = h \cdot f(tn + \frac{h}{2}, yn + \frac{k2}{2})]
. M  l; @, x: ^; K# W5 l[k4 = h \cdot f(tn + h, yn + k_3)]
[y{n+1} = yn + \frac{1}{6}(k1 + 2k2 + 2k3 + k4)]
+ E# d' C! {+ W" o8 r& x# |其中，(tn) 是当前时间步，(yn) 是当前的解，(h) 是步长，(f(t, y)) 是微分方程右侧的函数。
这个方法的精度相对较高，因为它使用了函数 (f(t, y)) 在一个步长内的多个点上的信息。四阶Runge-Kutta方法在许多情况下被广泛应用，因为它相对简单且相对高效。

1. %四阶经典公式,微分方程为f.m- K4 f- G3 T  M3 y
   
2. 
   * w( L# E  o: a  Q: ?
3. if exist('f.m')==0                                           %在星号处输入文件名(把星号改为文件名)' G9 r3 O- a* Q; M  R# z
   
4.    disp('没有为方程创建名为f.m的函数文件,请参照下例建立它');7 ^. ~5 E1 R# N4 y/ N3 W: V* `. `  x  T
   
5.    disp('function z=f(x,y)');
   3 N9 o2 B; O0 G6 h0 c% G8 M6 I
6.    disp('z=y-2*x/y;');/ ^! G% m9 S5 E, \. S
   
7.    disp('并将该文件保存在work文件夹下');) X4 s. }/ Z8 R& o$ a
   
8. end 8 K9 a/ d\" D# u# q, z6 b6 y
   
9. 
   ( x! Y4 n1 `. E6 z
10. X1=input('请输入求解区间的左端点X1=');/ }% A% N6 r, y3 P
    
11. Y1=input('请输入微分方程的初始条件Y1=（X=X1时Y的值）');$ R( \; |  U\" \& k& _
    
12. Xn=input('请输入求解区间的右端点Xn=');
    4 h' L  d  [* a- Z7 Y* n
13. h=input('请输入求解步长h=');
    - a6 e; u; Y3 S+ B: F\" v
14. 
    ! J: V! }2 ~1 \, Y, |$ e
15. X=X1;! y6 {, P\" s/ Y
    
16. Y=Y1;                                                        %运算初始点
    \" h6 ~' S0 s: y
17. n=0;                                                         %节点序号变量置零
    , v5 Y& l9 {$ T# _
18. 
    , F\" p0 W0 n; v
19. while X<=Xn-h$ M! U+ q1 p$ ~0 s
    
20.     K1=f(X,Y);
    3 F5 g: _6 G! R: A
21.     K2=f(X+h/2,Y+K1*h/2);8 m% p4 M$ W3 ]' ~: V
    
22.     K3=f(X+h/2,Y+K2*h/2);
    $ H! g0 |' ~7 o! n$ R
23.     K4=f(X+h,Y+K3*h);- Y/ b: G$ Z5 G- a% J; J; b
    
24.     X=X+h;8 K4 ^3 Z: h  h
    
25.     Y=Y+h*(K1+2*K2+2*K3+K4)/6;                               %四阶标准的龙格-库塔公式* o# q+ J7 m9 a* F: k) q\" z
    
26.     n=n+1;                                                   %节点序号加1
    , n* {2 |% P% n1 Y. g
27. 
    7 m# m% i; d9 m& s' S( J5 j8 u
28.     fprintf('第%d个点的计算结果为X=%10.8f,Y=%10.8f\n',n,X,Y);4 y; [4 z( j\" W6 P. M
    
29.     plot(X,Y,'o')& G& q! d8 {2 |- J) D
    
30.     hold on
    \" D\" h8 m; m, n: F* B
31. end

复制代码

1. function z=f(x,y). D) s* g8 u5 E7 _\" T- z3 {$ v1 S0 Y
   
2. z=y-2*x/y;

复制代码

# W2 I, G# \; l6 n7 b* v! e6 Z- t