数值分析常微分方程数值解实验题.doc

数值分析第四次上机练习实验报告常微分方程数值解实验题一、 问题的描述考虑一下问题：；初试条件为，。其精确解为；请分别用古典四级四阶显式Runge-Kutta方法和隐式二级四阶Runge-Kutta方法计算，计算区间取成，并与精确解比较。二、 方法描述Runge-Kutta方法是采用不同点上函数值的不同组合来提高方法的精度。又避免了函数f的偏导数计算。a) 古典四级四阶显式Runge-Kutta方法古典四级四阶显式Runge-Kutta方法相应的计算格式为b) 隐式二级四阶Runge-Kutta方法隐式二级四阶Runge-Kutta方法相应的计算格式为：三、 方案设计我们通过编写程序来进行运算，程序语言采用Matlab语言，运行环境为Matlab R2010b。在我们的程序文件中practice4.m文件为主程序文件， RK_ex44.m 为使用古典四级四阶显式Runge-Kutta方法计算的函数，RK_im24.m 为使用隐式二级四阶Runge-Kutta方法计算的函数。a) 古典四级四阶显式Runge-Kutta方法由function y=RK_ex44(f,a,b,y0,N,p)函数计算得到。该函数可计算函数f，计算区间上下限为a、b，初值为y0，N个点，p个未知数，使用古典四级四阶显式Runge-Kutta方法计算得到的值。显式计算方法非常便于实现，按照计算公式使用matlab语句编程实现即可。具体实现详见源程序。b) 隐式二级四阶Runge-Kutta方法由function y=RK_im24(A,c,a,b,y0,N) 函数计算得到。该函数可计算函数f，计算区间上下限为a、b，初值为y0，N个点使用隐式二级四阶Runge-Kutta方法计算得到的值。隐式计算方法不似显式那么便利，需要解非线性方程组。然后依然是迭代求解。具体实现详见源程序。四、 计算结果及误差分析 如上图所示为步长h=0.001时的仿真结果。图中蓝色线代表精确解，红色线代表古典显式四级四阶Runge-Kutta方法的解，绿色线（由于颜色叠加显示为黄色线）代表隐式二级四阶Runge-Kutta方法的解。由图可知，当步长h选择足够小时，两种方法均能较好的对实验结果进行模拟。下面对当选择不同的步长h时的实验结果进行讨论。首先讨论古典显式四级四阶Runge-Kutta方法。下述六幅图分别显示了当步长h=0.0015，h=0.0014，h=0.0013时，u(t)和v(t)的仿真结果。由图中可以看出，当步长不够小时，无法得到正确的仿真结果。而一旦选择足够小的步长，能够得到较好的仿真结果。当h=0.0015时当h=0.0014时当h=0.0013时当h=0.1时接下来讨论步长h对隐式二级四阶Runge-Kutta方法的仿真结果的影响。当h=0.01时可以看出，当h=0.1时，隐式二级四阶算法已能较好的精确v(t)，然而，因为u(t)值出现了阶跃，故并不能很好的对u(t)进行仿真。下表显示了h=0.1,h=0.08, h=0.06, h=0.04, h=0.02, 0.01时，u(t)的仿真图。并根据局部讨论的需要进行了局部的放大。由图中可以看出，随着步长的逐步缩减，在t=0处由于，u(t)的阶跃造成的误差在逐步缩减，直至可以忽略。五、 结论通过实验可以发现，显式RK方法虽然表达方式简明，步骤简单，然而需要较小的步长以及较长的计算时间才能达到较精确的解。相比之下，隐式RK方法需要的计算时间远小于显式方法，并且较大的步长就可以满足得到精确解的需求