实验5 常微分方程的数值解法.docVIP

课程 名称 计算方法 实验项目 名 称 常微分方程的数值解法 实验项目类型 验证 演示 综合 设计 其他 指导 教师 胡小兵 成 绩 √ 1．实验目的： （1）学会四阶龙格-库塔方法的使用 （2）设计出相应的算法，编制相应的函数子程序 （3）会用这些函数解决实际问题 2．实验内容 （1）分别取h=0.05，N=10；h=0.025，N=20；h=0.01，N=50，用四阶龙格-库塔方法求解微分方程初值 问题：y’=-50y，y(0)=10 （2）某跳伞者在t=0时刻从飞机上跳出，假设初始时刻的垂直速度为0，且跳伞者垂直下落。已知空气阻力为F=cv2，其中c为常数，v为垂直速度，向下方方向为正。写出此跳伞者的速度满足的微分方程；若此跳伞者的质量为M=70kg，且已知c=0.27kg/m，利用四阶龙格-库塔公式计算t=20s的速度（取h=0.1s） 3.实验原理 (1) format long; F=50*y; a=0; b=0.5; h=0.01; n=(b-a)/h; X=a:h:b; Y=zeros(1,n+1); Y(1)=10; for i=1:n x=X(i); y=Y(i); K1=h*eval(F); x=x+h/2; y=y+K1/2; K2=h*eval(F); x=x; y=Y(i)+K2/2; K3=h*eval(F); x=X(i)+h; y=Y(i)+K3; K4=h*eval(F); Y(i+1)=Y(i)+(K1+2*K2+2*K3+K4)/6; end %准确解 temp=[]; f=dsolve(Dy=50*y,y(0)=10,x); df=zeros(1,n+1); for i=1:n+1 temp=subs(f,x,X(i)); df(i)=double(vpa(temp)); end disp( x RK 准确解); disp([X,Y,df]); %pic figure; plot(X,df,k*,X,Y,--r); grid on; title(四阶经典R--K法解常微分方程); legend(准确解,四阶经典R--K法); 实验结果及分析 1.h=0.05 ； 2.h=0.025； 3.h=0.01时, (2)dv/dt=Mg-cv^2 x R--K 准确解 0 0 0 0.10000000000000 0.99987144716508 0.99987144840507 0.20000000000000 1.99897206046570 1.99897206295078 2.40000000000000 22.36755004004259 22.36755008131588 2.50000000000000 23.16756935550245 23.16756939922576 4.30000000000000 35.04076235589000 35.04076244476733 4.40000000000000 35.56007546772376 35.56007555873195 5.70000000000000 41.10978904603041 41.10978915610126 5.80000000000000 41.45244872371784 41.45244883451885 6.80000000000000 44.32804778430844 44.32804789677157 6.90000000000000 44.56602894632287 44.56602905843579 7.90000000000000 46.54050880231677 46.54050890691594 8.00000000000000 46.70212238522891 46.70212248874960 9.30000000000000 48.34520236041540 48.34520244684039 9.4000000000000