6.2 确定常微分方程数值解的Excel模板的制作
本节介绍用欧拉法、改进的欧拉法和龙格-库塔法确定初值问题数值解的方法，以Excel求解为例介绍模板的制作过程。
例6.1 求解初值问题：
初值问题（6.28）的正确解为
我们将区间[0，1]分成10等分，令h=0.1，xi=x0+0.1i（i=0，2，…，10）。
首先打开一个Excel工作簿，在这个工作簿中打开5个工作表，分别命名为欧拉法、改进的欧拉法、二阶龙格-库塔法、三阶龙格-库塔法和四阶龙格-库塔法。在每一个工作表内完成相应的操作。
1.欧拉法
由式（6.2）得到对应初值问题（6.28）的欧拉公式：
下面的所有制作过程以图6.1为参考。
图6.1 用欧拉法得到的数值解
第一步，输入步长和初始值。在单元格A2中输入步长“0.1”；在单元格A4、B4和C4中分别输入“0”“0”和“1”。
第二步，输入计算公式和误差。在单元格A5和B5中分别输入“A4+1”和“=B4+$A$2”；用式（6.30）在单元格C5中输入“=C4+$A$2*（C4-2*B4/C4）”；用式（6.29）在单元格D5中输入正确解的计算公式“=SQRT（1+2*B5）”；在单元格E5中输入误差的计算公式“=ABS（D5-C5）”。
第三步，填充。选中单元格区域A5：E5，把鼠标指针移动到单元格E5的右下角，待其变为填充柄（“
”符号）时，按住鼠标左键向下拖动到单元格E14，就得到区间[0，1]内10个点的数值解及其误差。
2.改进的欧拉法
由式（6.15）得到对应初值问题（6.28）的改进的欧拉公式：
下面的所有制作过程以图6.2为参考。
图6.2 用改进的欧拉法得到的数值解
第一步，输入步长和初始值。在单元格A2中输入步长“0.1”；在单元格A4、B4和D4中分别输入“0”“0”和“1”。
第二步，输入计算公式和误差。在单元格A5和B5中分别输入“A4+1”和“=B4+$A$2”；用式（6.31）在单元格C5和D5中分别输入计算y（0）1和y1的公式“=D4+$A$2*（D4-2*B4/D4）”和“=D4+$A$2/2*（（D4-2*B4/D4）+（C5-2*B5/C5））”；用式（6.29）在单元格E5中输入正确解的计算公式“=SQRT（1+2*B5）”；在单元格F5中输入误差的计算公式“=ABS（D5-C5）”。
第三步，填充。选中单元格区域A5：F5，把鼠标指针移动到单元格F5的右下角，待其变为填充柄（“
”符号）时，按住鼠标左键向下拖动到单元格F14，就得到区间[0，1]内10个点的数值解及其误差。
3.二阶龙格-库塔法
由式（6.23）得到对应初值问题（6.28）的二阶龙格-库塔公式：
(https://www.daowen.com)
下面的所有制作过程以图6.3为参考。
图6.3 用二阶龙格-库塔法得到的数值解
第一步，输入步长和初始值。在单元格A2中输入步长“0.1”；在单元格A4、B4和E4中分别输入“0”“0”和“1”。
第二步，输入计算公式和误差。在单元格A5和B5中分别输入“A4+1”和“=B4+$A$2”；用式（6.32）在单元格C5、D5和E5中分别输入计算K1，K12和y1的公式“=E4-2*B4/E4”“=（E4+$A$2/2*C5）-2*（B4+$A$2/2）/（E4+$A$2/2*C5）”和“=E4+$A$2*D5”；用式（6.29）在单元格F5中输入正确解的计算公式“=SQRT（1+2*B5）”；在单元格G5中输入误差的计算公式“=ABS（D5-C5）”。
第三步，填充。选中单元格区域A5：G5，把鼠标指针移动到G5的右下角，待其变为填充柄（“
”符号）时，按住鼠标左键向下拖动到单元格G14，得到区间[0，1]内10个点的数值解及其误差。
4.三阶龙格-库塔法
由式（6.25）得到对应初值问题（6.28）的三阶龙格-库塔公式：
下面的所有制作过程以图6.4为参考。
图6.4 用三阶龙格-库塔法得到的数值解
第一步，输入步长和初始值。在单元格A2中输入步长“0.1”；在单元格A4、B4和F4中分别输入“0”“0”和“1”。
第二步，输入计算公式和误差。在单元格A5和B5中分别输入“A4+1”和“=B4+$A$2”；用式（6.33）在单元格C5、D5、E5和F5中分别输入计算K1，K2，K3和y1的公式“=F4-2*B4/F4”“=（F4+$A$2/2*C5）-2*（B4+$A$2/2）/（F4+$A$2/2*C5）”“=（F4-$A$2*C5+2*$A$2*D5）-2*（B4+$A$2）/（F4-$A$2*C5+2*$A$2*D5）”和“=F4+$A$2/6*（C5+4*D5+E5）”；用式（6.29）在单元格G5中输入正确解的计算公式“=SQRT（1+2*B5）”；在单元格H5中输入误差的计算公式“=ABS（D5-C5）”。
第三步，填充。选中单元格区域A5：H5，把鼠标指针移动到单元格H5的右下角，待其变为填充柄（“
”符号）时，按住鼠标左键向下拖动到单元格H14，就得到区间[0，1]内10个点的数值解及其误差。
5.四阶龙格-库塔法
由式（6.27）得到对应初值问题（6.28）的三阶龙格-库塔公式：
下面的所有制作过程以图6.5为参考。
第一步，输入步长和初始值。在单元格A2中输入步长“0.1”；在单元格A4、B4和F4中分别输入“0”“0”和“1”。
第二步，输入计算公式和误差。在单元格A5和B5中分别输入“A4+1”和“=B4+$A$2”；用式（6.34）在单元格C5、D5、E5和F5中分别输入计算K1，K2，K3，K4和y1的公式“=F4-2*B4/F4”“=（F4+$A$2/2*C5）-2*（B4+$A$2/2）/（F4+$A$2/2*C5）”“=（F4-$A$2*C5+2*$A$2*D5）-2*（B4+$A$2）/（F4-$A$2*C5+2*$A$2*D5）”和“=F4+$A$2/6*（C5+4*D5+E5）”；用式（6.29）在单元格H5中输入正确解的计算公式“=SQRT（1+2*B5）”；在单元格I5中输入误差的计算公式“=ABS（D5-C5）”。
第三步，填充。选中单元格区域A5：I5，把鼠标指针移动到单元格I5的右下角，待变为填充柄（“
”符号）时，按住鼠标左键，向下拖动到单元格I14，就得到区间[0，1]内10个点的数值解及其误差。
图6.5 用四阶龙格-库塔法得到的数值解
6.误差比较图的绘制
在一个直角坐标系内绘制5种方法的误差比较图，横坐标为自变量x，纵坐标为误差，如图6.6所示，从图中可以观察到各方法的误差比价和优劣性。
图6.6 5种方法的误差比较图
