偏微分方程数值解例题答案

二、改进的Euler方法梯形方法的迭代公式(1.10)比Euler方法精度高,但其计算较复杂,在应用公式(1.10)进行计算时,每迭代一次,都要重新计算函数),(yxf的值,且还要判断何时可以终止或转下一步计算.为了控制计算量和简化计算法,通常只迭代一次就转入下一步计算.具体地说,我们先用Euler公式求得一个初步的近似值1ny,称之为预测值,然后用公式(1.10)作一次迭代得1ny,即将1ny校正一次.这样建立的预测－校正方法称为改进的Euler方法:预测:),,(1nnnnyxhfyy校正:)].,(),([2111nnnnnnyxfyxfhyy(1.15)这个计算公式也可以表示为11(,),(,),1().2pnnncnnpnpcyyhfxyyyhfxyyyy例1取步长0.1h，分别用Euler方法及改进的Euler方法求解初值问题d(1),01,d(0)1.yyxyxxy解这个初值问题的准确解为()1(21)xyxex.根据题设知).1(),(xyyyxf(1)Euler方法的计算式为)],1([1.01nnnnnyxyyy由1)0(0yy,得,9.0)]101(1[1.011y,8019.0)]9.01.01(9.0[1.09.02y这样继续计算下去，其结果列于表9.1.(2)改进的Euler方法的计算式为110.1[(1)],0.1[(1)],1(),2pnnnncnpnpnpcyyyxyyyyxyyyy由1)0(0yy，得110.1[1(101)]0.9,10.1[0.9(10.10.9)]0.9019,1(0.90.9019)0.900952pcyyy20.900950.1[0.90095(10.10.90095)]0.80274,0.900950.1[0.80274(10.20.80274)]0.80779,1(0.802740.80779)0.805262pcyyy这样继续计算下去，其结果列于表9.1.表9.1Euler方法改进的Euler方法准确值nxnyny)(nxy0.10.90000000.90095000.90062350.20.80190000.80526320.80463110.30.70884910.71532790.71442980.40.62289020.63256510.63145290.50.54508150.55761530.55634600.60.47571770.49055100.48918000.70.41456750.43106810.42964450.80.36108010.37863970.37720450.90.31454180.33262780.33121291.00.27418330.29235930.2909884从表9.1可以看出,Euler方法的计算结果只有2位有效数字,而改进的Euler方法确有3位有效数字,这表明改进的Euler方法的精度比Euler方法高.例2试用Euler方法、改进的Euler方法及四阶经典R-K方法在不同步长下计算初值问题1)0(,10),1(ddyxxyyxy在0.2、0.4、0.8、1.0处的近似值，并比较它们的数值结果.解对上述三种方法,每执行一步所需计算)1(),(xyyyxf的次数分别为1、2、4。为了公正起见，上述三种方法的步长之此应为4:2:1。因此，在用Euler方法、改进的Euler方法及四阶经典R-K方法计算0。2、0。4、0。8、1。0处的近似值时，它们的步长应分别取为0。05、0。1、0。2，以使三种方法的计算量大致相等。Euler方法的计算格式为)].1([05.01nnnnnyxyyy改进的Eluer方法的计算格式为).(21)],1([1.0)],1([1.011cpnpnpncnnnnpyyyyxyyyyxyyy四阶经典R-K方法的计算格式为)]2.0)(2.0(1)[2.0()],22.0)(22.0(1)[22.0()],22.0)(22.0(1)[22.0(),1(),22(62.0334223112143211kyxkykkyxkykkyxkykyxykkkkkyynnnnnnnnnnnnnn初始值均为1)0(0yy，将计算结果列于表9.2.表9.2Euler方法（步长h=0.05）改进的Euler方法（步长h=0.1）四阶经典R-K方法（步长h=0.2）准确解nxnynyny)(nxy0.20.80318660.80526320.80463630.80463110.40.62717770.63256510.63146530.63145290.60.48255860.49055100.48919790.48918000.80.36930360.37863970.37722490.37720451.00.28274820.29235930.29100860.2909884从表9.2可以看出，在计算量大致相等的情况下，Euler方法计算的结果只有2位有效数字，改进的Euler方法计算的结果有3位有效数字，而四阶经典R-K方法计算的结果却有5位有效数字，这与理论分析是一致的。例1和例2的计算结果说明，在解决实际问题时，选择恰当的算法是非常必要的。需要指出的是Runge-Kutta方法的基于Taylor展开法，因而要求解具有足够的光滑性。如果解的光滑性差，使用四阶Runge-Kutta方法求得数值解的精度，可能不如改进的Euler方法精度高。因此，在实际计算时，要根据具体问题的特性，选择合适的算法。一、应用向前欧拉法和改进欧拉法求由如下积分20xtyedt所确定的函数y在点x=0.5,1.0,1.5的近似值。解：该积分问题等价于常微分方程初值问题2'(0)0xyey其中h=0.5。其向前欧拉格式为2()100ihiiyyhey改进欧拉格式为22()2(1)10()20ihihiihyyeey将两种计算格式所得结果列于下表iix向前欧拉法iy改进欧拉法iy000010.50.50.4447021.00.889400.7313731.51.073340.84969二、应用4阶4步阿达姆斯显格式求解初值问题'1(0)1yxyy00.6x取步长h=0.1.解：4步显式法必须有4个起步值，0y已知，其他3个123,,yyy用4阶龙格库塔方法求出。本题的信息有：步长h=0.1；结点0.1(0,1,,6)ixihii；0(,)1,(0)1fxyxyyy经典的4阶龙格库塔公式为11234(22)6iihyykkkk1(,)1iiiikfxyxy121(,)0.051.0522iiiihkhkfxyxyk232(,)0.051.0522iiiihkhkfxyxyk433(,)0.11.1iiiikfxhyhkxyk算得11.0048375y,21.0187309y,31.0408184y4阶4步阿达姆斯显格式1123(5559379)24iiiiiihyyffff1231(18.55.93.70.90.243.24)24iiiiiyyyyyi由此算出4561.0703231,1.1065356,1.1488186yyy三、用Euler方法求'1,0101xyeyxxy问步长h应该如何选取，才能保证算法的稳定性？解：本题,1xfxyeyx,0,01xyfxyex本题的绝对稳定域为111xhhe得02xhe，故步长应满足02,00.736heh求梯形方法111[(,)(,)]2kkkkkkhyyfxyfxy的绝对稳定域。证明：将Euler公式用于试验方程'yy，得到11[]2kkkkhyyyy整理11(1)22kkhhyy设计算ky时有舍入误差,0,1,2,kk，则有11(1)22kkhh据稳定性定义，要想1kk，只须1122hh因此方法绝对稳定域为复平面h的整个左半平面（？），是A-稳定的。五、对初值问题'(0)1yyy01x证明：用梯形公式111[(,)(,)]2nnnnnnhyyfxyfxy求得的数值解为22nnhyh并证明当步长0h时，ny收敛于该初值问题的精确解xnye证明：由梯形公式，有1111[(,)(,)][]22nnnnnnnnnhhyyfxyfxyyyy整理，得122nnhyyh由此递推公式和初值条件，有02222nnnhhyyhh[0,1]x，则有在区间0,0,1x上有nxxnh，步长xhn，由前面结果有02222022limlimlim1222lim12xnhnnnhxhhhxhhhyhhheh由x的任意性，得所证。六、对于微分方程'(,)yfxy，已知在等距结点0123,,,xxxx处的y的值为0123,,,yyyy，h为步长。试建立求4y的线性多步显格式与与隐格式。解：取积分区间24[,]xx，对'(,)yfxy两端积分：442242(,)xxxxyxyxdyfxydx对右端(,)fxy作123,,xxx的二次插值并积分4242021112222233(,)[()(,)()(,)()(,)]xxxxfxydxlxfxylxfxylxfxydx112233123((,)(,)(,))337hfxyfxyfxy得到线性4步显格式42113123()337yyhfff若对右端在34,xx两点上作线性插值并积分，有424201331144(,)[()(,)()(,)]xxxxfxydxlxfxylxfxydx442(,)hfxy由此产生隐格式42442,yyhfxy七、证明线性多步法111(3)()2nnnyhffnn-1n-2（y-y）-y存在的一个值，使方法是4阶的。解：由本题的公式，有111(3)()2nnnyhffnn-1n-2（y-y）+y11()nnnTyxhy234(4)5[()'()''()'''()()()]2!3!4!nnnnnhhhyxhyxyxyxyxOh1[(()())(2)(3)(''())]2nnnnnyxyxhyxhhyyxh234(4)5[()'()''()'''()()()]2!3!4!nnnnnhhhyxhyxyxyxyxOh234(4)5()(()'()''()'''()()())2!3!4!nnnnnnhhhyxyxhyxyxyxyxOh234(4)5(2)(2)(2)(()2'()''()'''()()())2!3!4!nnnnnhhhyxhyxyxyxyxOh23(4)51(3)('()'()''()'''()()())22!3!nnnnnhhhyxyxhyxyxyxOh