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

2744557306

"定步长四阶经典公式"通常指的是数值积分中的四阶Runge-Kutta方法。这是一种常用的数值解常微分方程（ODE）的方法，其主要思想是通过逐步逼近来估计微分方程的解。
- \* e! K% i2 Y定步长四阶经典公式是Runge-Kutta方法的一种，其中最常见的是经典的四阶Runge-Kutta方法。对于一个一阶常微分方程
[\frac{dy}{dt} = f(t, y)]
这个方法的迭代公式如下：
* \, S7 V: Y+ F' m! Z[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})]
" P0 X! I- j1 `* w6 c4 y[k4 = h \cdot f(tn + h, yn + k_3)]
* H6 L. {; u# z5 J" O- ~- }[y{n+1} = yn + \frac{1}{6}(k1 + 2k2 + 2k3 + k4)]
其中，(tn) 是当前时间步，(yn) 是当前的解，(h) 是步长，(f(t, y)) 是微分方程右侧的函数。
这个方法的精度相对较高，因为它使用了函数 (f(t, y)) 在一个步长内的多个点上的信息。四阶Runge-Kutta方法在许多情况下被广泛应用，因为它相对简单且相对高效。

1. %四阶经典公式,微分方程为f.m
   ( x, Q9 v3 G/ }' r& V6 ^4 j, d
2. 
   8 T% n$ j5 U, S8 Y7 c+ O+ R; a
3. if exist('f.m')==0                                           %在星号处输入文件名(把星号改为文件名)9 B/ o8 d; X8 X
   
4.    disp('没有为方程创建名为f.m的函数文件,请参照下例建立它');: J8 f% |5 X- V+ n
   
5.    disp('function z=f(x,y)');
   : {/ l/ X5 g1 @: {\" A; }6 v6 o
6.    disp('z=y-2*x/y;');1 L2 v+ K. |! G: v4 Y! Q1 X# l
   
7.    disp('并将该文件保存在work文件夹下');: h, P' N( b. U) @: @
   
8. end
   $ ?  H& W. K# q. E
9. 
   ' E! N0 X+ X) \# |8 h$ h, ]0 h
10. X1=input('请输入求解区间的左端点X1=');% j( p  z2 s7 t6 p
    
11. Y1=input('请输入微分方程的初始条件Y1=（X=X1时Y的值）');
    ) E3 ~; `; H) X0 ^) S) E
12. Xn=input('请输入求解区间的右端点Xn=');
    . n( S- p5 |6 n; Z\" m\" N3 i9 Q' B' _
13. h=input('请输入求解步长h=');
    9 g; P- E- k7 t
14. / k/ V6 \0 \, T1 |' k8 C
    
15. X=X1;
    0 W( ~+ Z2 ]0 [' O0 B
16. Y=Y1;                                                        %运算初始点\" ^5 ^* a4 n- `; L: f
    
17. n=0;                                                         %节点序号变量置零- I: {* g3 h3 g
    
18. 
    6 H) z  \5 K* T7 A& `
19. while X<=Xn-h
    7 W5 J% D. u  Q( Z
20.     K1=f(X,Y);. v5 T1 i' L3 R5 J
    
21.     K2=f(X+h/2,Y+K1*h/2);
    : e$ g& v# A9 P2 I; a' S. [\" f
22.     K3=f(X+h/2,Y+K2*h/2);\" v: o7 w- E6 H2 J1 D) V
    
23.     K4=f(X+h,Y+K3*h);
    ! L6 z6 ^2 d3 G! I' H0 |, [
24.     X=X+h;
    5 o- T9 @% R/ r8 j1 I- Q
25.     Y=Y+h*(K1+2*K2+2*K3+K4)/6;                               %四阶标准的龙格-库塔公式. w6 n* v- C2 U( ^5 I
    
26.     n=n+1;                                                   %节点序号加1% Q9 O/ N: w6 }8 y0 k* ]) j' y! m8 R
    
27. 
    # n/ p* z# e$ n- z- H, A
28.     fprintf('第%d个点的计算结果为X=%10.8f,Y=%10.8f\n',n,X,Y);- [3 y% ]\" W: {, S7 P- U
    
29.     plot(X,Y,'o')
    . K& {7 b' |0 P' n5 \
30.     hold on0 [! n, [+ Z. A$ S' U
    
31. end

复制代码

1. function z=f(x,y)4 u. L+ J/ U# x/ E! r
   
2. z=y-2*x/y;

复制代码