资源描述:
《四阶龙格库塔法的编程(赵)》由会员上传分享,免费在线阅读,更多相关内容在行业资料-天天文库。
1、例题一四阶龙格-库塔法的具体形式为:1.1程序:/*e.g:y'=t-y,t∈[0,1]/* y(0)=0/*使用经典四阶龙格-库塔算法进行高精度求解 /*y(i+1)=yi+(K1+2*K2+2*K3+K4)/6/*K1=h*f(ti,yi)/*K2=h*f(ti+h/2,yi+K1*h/2)/*K3=h*f(ti+h/2,yi+K2*h/2)/*K4=h*f(ti+h,yi+K3*h)*/#include#include#defineN20floatf(floatd,floatp)//要求解的微分方程的右部的函数e.g:t-y{flo
2、atderivative;derivative=d-p;return(derivative);}voidmain(){floatt0;//范围上限floatt;//范围下限floath;//步长floatnn;//计算出的点的个数,即迭代步数intn;//对nn取整floatk1,k2,k3,k4;floaty[N];//用于存放计算出的常微分方程数值解floatyy;//精确值floaterror;//误差inti=0,j;//以下为函数的接口printf("inputt0:");scanf("%f",&t0);printf("inputt:");scanf("%f",&
3、t);printf("inputy[0]:");scanf("%f",&y[0]);printf("inputh:");scanf("%f",&h);//以下为核心程序nn=(t-t0)/h;printf("nn=%f",nn);n=(int)nn;printf("n=%d",n);for(j=0;j<=n;j++){yy=t0-1+exp(-t0);//解析解表达式error=y[i]-yy;//误差计算printf("y[%f]=%fyy=%ferror=%f",t0,y[i],yy,error);//结果输出k1=h*f(t0,y[i]);//求K1k2=
4、h*f((t0+h/2),(y[i]+k1*h/2));//求K2k3=h*f((t0+h/2),(y[i]+k2*h/2));//求K3k4=h*f((t0+h),(y[i]+k3*h));//求K4y[i+1]=y[i]+((k1+2*k2+2*k3+k4)/6);//求y[i+1]t0+=float(h);i++;}}1.2结果:inputt0:0inputt:1inputy[0]:0inputh:0.2nn=5.000000n=5tiyiyyerror00000.20.0197360.0187310.0010050.40.0748130.0703200.004493
5、0.60.1583030.1488120.0094910.80.2646350.2493290.0153061.00.3893310.3678790.021452图一例题二(见《高等流体力学》P45页例1.6)给定速度场u=x+t,v=-y-t,求t=1时过(1,1)点的质点的迹线。解答:由迹线微分方程dx/dt=x+t,dy/dt=-y-t积分得x=Aet–t-1,y=Be-t–t+1t=1时过(1,1)点的质点有A=3/e,B=e最后得此质点的迹线方程为:x=3et-1–t-1,y=e-1-t–t+12.1程序#include#include