微分方程数值解法.docx
微分方程数值解法第四章 抛物型方程的有限差分法 l 模型问题 常系数线性抛物型方程初边值问题 2ì¶u¶u=a+f(x,t),0<t£T, (1.1)ï2¶xï¶tíu(x,0)=f(x),0<x<l, (1.3)ïu(0,t)=u(l,t)=0,0£t£T, (1.3)ïî12书中，f(x,t):=f(x)即与t无关。 作业讲解： P130： 将向前差分格式和向后差分格式作加权平均，得到下列格式： ujk+1-ujkt=ah2q(uj+1-2ujkk+1k+1+uj-1)+kkk+1(1-q)(uj+1-2uj+uj-1) 其中0£q£1，试计算其截断误差，并证明当q=24截断误差的阶最高。 12-112r时，解法： 本题的计算格式对应的微分方程是： ¶u¶t=a¶u¶x22,f=0 u(xj+1,tk)=u(xj,tk)+h+16h11203¶¶xu(xj,tk)+124h124h¶2¶22¶x44u(xj,tk)¶33¶xh5u(xj,tk)+55¶x6u(xj,tk)¶¶xu(xj,tk)+O(h)(1) u(xj-1,tk)=u(xj,tk)-h-16h11203¶¶xu(xj,tk)+124h124h¶2¶22¶x44u(xj,tk)¶33¶xh5u(xj,tk)+55¶x6u(xj,tk)¶¶xu(xj,tk)+O(h)(2) (1)+(2)得到： u(xj+1,tk)+u(xj-1,tk)=2u(xj,tk)+h+112h42¶22¶x6u(xj,tk)¶44¶xu(xj,tk)+O(h)即 1h2u(xj+1,tk)-2u(xj,tk)+u(xj-1,tk)=+112h2¶22¶xu(xj,tk)¶44¶xu(xj,tk)+O(h)4 (3) 同理： 1h2u(xj+1,tk+1)-2u(xj,tk+1)+u(xj-1,tk+1)=+112h2¶22¶xu(xj,tk+1)¶44¶xu(xj,tk+1)+O(h)4(4) u(xj,tk+1)=u(xj,tk)+t+12¶¶tu(xj,tk)3t2¶22¶tu(xj,tk)+O(t)¶22¶xu(xj,tk+1)=¶22¶xu(xj,tk)+t¶23¶x¶tu(xj,tk)+O(t) 2¶44¶xu(xj,tk+1)=¶44¶xu(xj,tk)+t¶45¶x¶tu(xj,tk)+O(t) 2由，得到： u(xj,tj+1)-u(xj,tj)t+12=¶¶tu(xj,tj)2t¶22¶tu(xj,tj)+O(t)将和代入，又得到： 1h2u(xj+1,tk+1)-2u(xj,tk+1)+u(xj-1,tk+1)=+¶11222¶xu(xj,tk)+t2¶3¶x¶t1122u(xj,tk)+ht2h¶444¶x2u(xj,tk)+¶45¶x¶tu(xj,tk)+O(t+h) 由于 æ¶öu(x,t)=u(x,t)çjkjk÷22¶x¶t¶xè¶tø4 æ¶2ö¶ça=u(xj,tk)÷=au(xj,tk)2ç24÷¶xè¶x¶xø¶3¶2¶2¶æ¶öu(xj,tk)=çu(xj,tk)÷2¶t¶tè¶tø23ö¶æ¶¶ça=u(xj,tk)÷=au(xj,tk)22ç÷¶tè¶x¶x¶tø ¶2=a2¶44¶xu(xj,tk)又由于要求网比r是常数，所以t=O(h)2，则有 1h2u(xj+1,tk+1)-2u(xj,tk+1)+u(xj-1,tk+1)=¶22¶xu(xj,tk)+(ta+24112h)2¶44¶xu(xj,tk)+O(t+h)u(xj,tk+1)-u(xj,tk)t+12=¶44¶¶tu(xj,tk)2ta2¶xu(xj,tk)+O(t)重写 1h2u(xj+1,tk)-2u(xj,tk)+u(xj-1,tk)=+112h4¶22¶xu(xj,tk)¶44¶xu(xj,tk)+O(h)4故我们最后求截断误差： R(xj,tj)=-aqu(xj,tk+1)-u(xj,tk)tu(xj+1,tk+1)-2u(xj,tk+1)+u(xj-1,tk+1)h+(1-q)2u(xj+1,tk)-2u(xj,tk)+u(xj-1,tk)h442=¶¶tu(xj,tj)+12ta2¶¶xu(xj,tj) -aq¶22¶xu(xj,tk)-aq(ta+¶22112h)1122¶44¶xhu(xj,tk) ¶44-a(1-q)2¶xu(xj,tk)-a(1-q)2¶xu(xj,tk) +O(t+h)=12at-qat-24224112qah-2112ah+2112qah2¶44¶xu(xj,tk)+O(t+h)=at=at221212-q-q-1h1212at¶x12r¶x¶44¶44u(xj,tk)+O(t+h)2424u(xj,tk)+O(t+h)所以当 at212-q-12-1112r¶x¶44u(xj,tk)=0 也即当q=12r时，截断误差阶可以达到最高阶O(t+h) 24 证明完毕 P136 1 1. 求证差分格式(4.1.17)当当0£q£1212£q£1时恒稳定 时稳定的充要条件是 12(1-2q)。 r£证明 首先将(4.1.17)写成如下形式： ujk+1-ujkt=ah2q(uj+1-2ujk+1k+1+uj-1)+k+1(1-q)(uj+1-2uj+uj-1)kkk-qruk+1j+1+(1+2qr)ujkj+1k+1-qruk+1j-1kjkj-1=(1-q)ru+(1-2(1-q)r)u+(1-q)ru这里网格比是：r=令 ath2æ0çç1S=ç0ççKèUk10O0k10OO1k2Lö÷0÷1÷÷÷0øuTkN-1=u(uL)则可以写成矩阵形式： (1+2qr)I-qrSUk+1=(1-2(1-q)r)I+(1-q)SUk从而 Uk+1=(1+2qr)I-qrS(1-2(1-q)r)I+(1-q)SU -1k其递增矩阵为 C=(1+2qr)I-qrS(1-2(1-q)r)I+(1-q)rS -1由于S的特征值为： S lj=2cosjph,j=1,2,.,N-1,h=1/N, 所以C的特征值为： l=Cj1-2(1-q)r+2(1-q)rcosjph1+2qr-2qrcosjph1-2(1-q)r(1-cosjph)1+2qr(1-cosjph)1-4(1-q)rsin2jph2=1+4qrsin1-cosA=2sin22jph2A2为使lj£1+Mt,就有 1-4(1-q)rsin-1£1+4qrsin22Cjph2£1jph2(&) 1)当 12£q£1时: 0£1-q£1/2 4qrsinÞ|1-4(1-q)rsin22jph2³4(1-q)rsin2jph2³0jph2|<1+4qrsin2jph2即恒成立,也就是当 定 12£q£1时，所讨论的格式恒稳2) 当 0£q£1/2时: 1/2£1-q£1,1-2q³0 由于 -1-4qrsin2jph2£1-4(1-q)rsin2jph2£1+4qrsin2jph2显然后一不等式恒成立 -2£4(2q-1)rsinÞ4(1-2q)rsinÞ(1-2q)r£1/2Þr£1/2(1-2q)22jph2jph2£2
