四阶龙格-库塔法.docVIP

2013-2014（1）专业课程实践论文 题目：四阶龙格—库塔法 一、算法理论 由定义可知，一种数值方法的精度与局部截断误差 有关，用一阶泰勒展开式近似函数得到欧拉方法，其局部截断误差为一阶泰勒余项，故是一阶方法，完全类似地若用p阶泰勒展开式 进行离散化，所得计算公式必为p阶方法，式中 由此，我们能够想到，通过提高泰勒展开式的阶数，可以得到高精度的数值方法，从理论上讲，只要微分方程的解充分光滑，泰勒展开方法可以构造任意的有限阶的计算公式，但事实上，具体构造这种公式往往相当困难，因为符合函数的高阶导数常常是很烦琐的，因此，泰勒展开方法一般不直接使用，但是我们可以间接使用泰勒展开方法，求得高精度的计算方法。 首先，我们对欧拉公式和改进欧拉公式的形式作进一步的分析。 如果将欧拉公式和改进的欧拉公式改写成如下的形式： 欧拉公式 改进的欧拉公式 ， ， 。 这两组公式都是用函数在某些点上的值的线性组合来计算的近似值，欧拉公式每前进一步，就计算一次的值。另一方面它是在处的一阶泰勒展开式，因而是一阶方法。改进的欧拉公式每前进一步，需要计算两次的值。另一方面它在处的泰勒展开式与在处的泰勒展开式的前三项完全相同，因而是二阶方法。这启发我们考虑用函数在若干点上的函数值的线性组合来构造计算公式。构造时，要求计算公式在处的泰勒展开式，与微分方程的解在处的泰勒展开式的前面若干项相同，从而使计算公式打到较高的精度。这样，既避免了计算函数的偏导数的困难，又提高了计算方法的精度，这就是龙格——库塔方法的基本思想。 二、算法框图 = 三、算法程序 程序代码： #include stdio.h #include conio.h float func(float x,float y) { return(2*x*y); } float runge_kutta(float x0,float xn,float y0,int n) { float x,y,y1,y2,h,xh; float d1,d2,d3,d4; int i; x=x0; y=y0; h=(xn-x0)/n; for(i=1;i=n;i++) { xh=x+h/2; d1=func(x,y); d2=func(xh,y+h*d1/2.0); d3=func(xh,y+h*d2/2.0); d4=func(xh,y+h*d3); y=y+h*(d1+2*d2+2*d3+d4)/6.0; x=x0+i*h; } return(y); } int main() { float x0,xn,y0,e; int n; printf(\ninput n:\n); scanf(%d,n); printf(input x0,xn:\n); scanf(%f%f,x0,xn); printf(input y0:\n); scanf(%f,y0); e=runge_kutta(x0,xn,y0,n); printf(y(%f)=%6.6f,y0,e); } 四、算法实现 四阶经典龙格——库塔方法算法如下： 输入区间等分数，初值。 输出在的个点处得近似值。 置。 计算 ， ， 置，输出。 若，置，转（4）；否则，停机。 例1．利用四阶龙格——库塔方程公式计算 的数值，取步长。 解：运行结果 例2．利用4阶龙格——库塔方程公式计算 取步长。 解：运行结果 始 输入 输出k,x,y k=n 终