求助:用龙格库塔解方程

wangh

从来没接触过龙格库塔啊,今天导师居然让我们用这个方法解这个方程,不会啊,谁可以帮助我解下,以前的帖子搜了下,看不懂^^!
y''''=4+x+x*x*x     在x=0时 : y=1 , y'=1 , y''=1 ,  y'''=1
谁能给我个完整的程序吗?最好是用C编程.谢谢大虾们咯.

心灯

网上搜索徐世良的 C语言算法程序集？ （本论坛不提供讨论下载共享的问题，还请注意）

另外，没有接触过，但是你懂了原理，自己实现也很快的，前提是你有C语言的基础。

风花雪月

1. /*
   
2. /*编程者：刘艮平
   
3. /*完成日期：2004.5.29
   
4. /*e.g:y'=y-2x/y,x∈[0,0.6]
   
5. /*　　y(0)=1
   
6. /*使用经典四阶龙格-库塔算法进行高精度求解　
   
7. /* y(i＋1)＝yi＋( K1+ 2*K2 +2*K3+ K4)/6
   
8. /* K1=h*f(xi,yi)
   
9. /* K2=h*f(xi+h/2,yi+K1/2)
   
10. /* K3=h*f(xi+h/2,yi+K2/2)
    
11. /* K4=h*f(xi+h,yi+K3)
    
12. */
    
13. #include <stdio.h>
    
14. #include <math.h>
    
15. 
    
16. float f(float den,float p0)  //要求解的微分方程的右部的函数 e.g: y-2x/y
    
17. {
    
18.     float rus;
    
19. //   den=w/(W0+sl);
    
20. //   rus=k*A*f/Wp*sqrt(RTp)-(k-1)*.....
    
21.     rus=p0-2*den/p0;
    
22.     return(rus);
    
23. }
    
24. //float fden()
    
25. //{
    
26. //}
    
27. 
    
28. void main()
    
29. {
    
30.     float x0; //范围上限
    
31.     int x1;   //范围下限
    
32.     float h;  //步长
    
33.     int n;    //计算出的点的个数
    
34.      float k1,k2,k3,k4;
    
35.     float y[3]; //用于存放计算出的常微分方程数值解
    
36.     int i=0;
    
37.     int j;
    
38. 
    
39. //以下为函数的接口
    
40.     printf("intput x0:");
    
41.     scanf("%f",&x0);
    
42. 
    
43.     printf("input x1:");
    
44.     scanf("%f",&x1);
    
45. 
    
46.     printf("input y[0]:");
    
47.     scanf("%f",&y[0]);
    
48. 
    
49.     printf("input h:");
    
50.     scanf("%f",&h);
    
51. 
    
52. // 以下为核心程序
    
53.     n=((x1-x0)/h);
    
54.     n=3;
    
55.    
    
56.     for(j=0;j<n;j++)
    
57.     {   
    
58.       k1=h*f(x0,y[i]); //求K1
    
59.       k2=h*f((x0+h/2),(y[i]+k1/2)); //求K2
    
60.       k3=h*f((x0+h/2),(y[i]+k2/2)); //求K3
    
61.       k4=h*f((x0+h),(y[i]+k3)); //求K4
    
62. 
    
63.       y[i+1]=y[i]+((k1+2*k2+2*k3+k4)/6); //求yi+1
    
64.       x0+=float(0.2);
    
65.       printf("y[%f]=%f\n",x0,y[i+1]);
    
66.       ++i;     
    
67.     }
    
68. }

复制代码

这种东西网上遍地都是，建议今后自己好好找找