一维热传导方程的前向 、紧差分格式

中南林业科技大学本科课程论文学院：理学院专业年级：09信息与计算科学一班课程：偏微分方程数值解法论文题目：一维热传导方程的前向Euler和紧差分格式指导教师：陈红斌2012年7月一维热传导方程的前向、紧差分格式学生姓名：唐黎学号：20093936分工：程序编写，数值例子学生姓名：何雄飞学号：20093925分工：格式建立，资料收集学生姓名：汪霄学号：20093938分工：文档编辑，资料整理学生姓名：毛博伟学号：20093931分工：公式编辑，查找资料学生姓名：倪新东学号：20093932分工：数据分析，查找资料学生姓名：何凯明学号：20093924分工：数据分析，查找资料一维热传导方程的前向、紧差分格式目录1引言..............................................12物理背景.............................................13网格剖分..........................................24．1.1向前Euler格式建立.............................24.1.2差分格式的求解...............................44.1.3收敛性与稳定性...................................44.1.4数值例子......................................74.2.1紧差分格式建立................................104.2.2差分格式求解..................................124.2.3数值例子......................................13总结..................................................17参考文献...........................................18附录...............................................19一维热传导方程的前向、紧差分格式第1页共25页1引言本文考虑的一维非齐次热传导方程的定解问题：22(,),0,0,uuafxtxltTtx(,0)(),0,uxxxl(0,)(),(1,)(),0.uttutttT其中a为正常数，(,),(),(),()fxtxtt为已知函数，(0)(0),(1)(0).目前常用的求解热传导方程的差分格式有前向Euler差分格式、向后Euler差分格式、Crank-Nicolson格式、Richardson格式[1,2,3]．本文将给出前向Euler格式和紧差分格式，并给出其截断误差和数值例子．2物理背景热传导是由于物体内部温度分布不均匀,热量要从物体内温度较高的点流向温度较低的点处.以函数,,,uxyzt表示物体在t时刻,,MMxy处的温度,并假设,,uxyz关于,,xyz具有二阶连续偏导数,关于t具有一阶连续偏导数.,,kkxyz是物体在,,Mxyz处的热传导系数,取正值.设物体的比热容为,,ccxyz,密度为,,xyz.根据Fourier热传导定律,热量守恒定律以及Gauss公式得,uuuuckxkktxxyyzz如果物体是均匀的,此时,kc以及均为常数.令2kac,上式方程化为22222222,tuuuuaauxyz若考虑物体内有热源,其热源密度函数为,,FFxyz,则有热源的热传导方程为2,,,,tuaufxyzt其中Ffc.一维热传导方程的前向、紧差分格式第2页共25页3网格剖分取空间步长Nlh/和时间步长MT/，其中MN,都是正整数．用两族平行直线),1,0(Njjhxj和),,1,0(Mkktk将矩形域}0,0{TtlxG分割成矩形网格，网格节点为),(kjtx．记),(kjkjtxuu．以hG表示网格内点集合，即位于开矩形G的网点集合；hG表示所有位于闭矩形的网点集合；hhhGG是网格界点集合.引进如下记号:0maxiim，21()mijh，2111()miijhh，2211，分别称为无穷范式（一直范式）2范数（2L范数，平均范数），差商的2范数（差商的2L范式）和1H范式.4.1.1向前Euler格式建立定义h上的网格函数{0,0},kiUUimkn其中(,),kiikUuxt0,im01.kn在结点处考虑微分方程(3.1-1),有22(,)(,)(,),11,01.ikikikuuxtaxtfxtimkntx(3.2)将224112241(,)[(,)2(,)(,)](,)12ikikikikikkuhuxtuxtuxtuxttxhx一维热传导方程的前向、紧差分格式第3页共25页242114(,),12kxiikkiikihuUtxxx和2121(,)[(,)(,)](,)2ikikikiikuuxtuxtuxtxtt212(,),2ktiiikkikkuDUxttt代入(3.2),得到224224(,)(,)(,),212kktixiikiikikkuahuDUaUfxtxttx11,01.imkn(3.3-1)注意到初边值条件(3.1-2)和(3.1-3),有0(,0)(),0,iiiUuxxim(3.3-2)0(),kkUt(),1.kmkUtkn(3.3-3)在(3.3-1)~(3.3-3)中略去小量项224(1)24(,)(,),212ikiikikkuahuRxttx(3.4)并用kiu代替kiU，得到如下差分格式2(,),kktixiikDuaufxt11,01.imkn(3.5-1)0(),iiux0,im(3.5-2)0(),kkut1.kn(3.5-3)称(1)kR为差分格式的局部截断误差。记241240101001max{max,max},212xxtTtTuauctx(3.6-1)则有(1)21(),ikRch11,01.imkn(3.6-2)一维热传导方程的前向、紧差分格式第4页共25页4.1.2差分格式的求解记2=/rah，称r为步长比。差分格式(3.5-1)可写为111(12)()(,),kkkkiiiiikururuufxt01,01.imkn上式表明第k+1层上的值由第k层上的值显示表示出来。若已知第k层的值{|0}kiuim，则由上式就可直接得到第k+1层上的值+1{|0}kiuim。有时也称(3.5-1)为古典显格式。可把古典显格式写成矩阵形式1111012221222111112(,)12(,)12(,)12(,)kkkkkkkkkmmmkkkkmmmkmrruufxtrurrruufxtrrruufxtrruufxtru4.1.3收敛性与稳定性收敛性设{(,)|0,0}ikuxtimkn为定解问题（3.1-1）~（3.1-3）的解，{|0,0}kiuimkn为差分格式（）（）的解，则当1/2r时，有210max(,)(),0kikiimuxtucThkn，其中1c由（3.6-1）定义证明记(,),0,0.kkiikieuxtuimkn将（3.3-1）~（3.6-1）分别于（3.5-1）~（3.5-3）相减，得到误差方程2(1),11,01,kktixiikDeaeRimkn00,0,ieim00,0,1.kkmeekn当1/2r时，有1(1)2211110max()(),1.kkikimieRckhcThkn一维热传导方程的前向、紧差分格式第5页共25页证明完毕.稳定性如果在应用差分格式(3.5-1)~(3.5-3)时，计算右端函数(x,t)ikf有误差gki，计算初值(x)i有误差i，则实际得到的是如下差分方程的解。2(,),kkktixiikiDvavfxtg11,im01,kn0(),iiivx0,im(3.8)0(),kkvt(),kmkvt1.kn令,kkkiiivu0,im0,kn将(3.5-1)~(3.5-3)与(3.8)相减，可得摄动方程组2,kkktixiiDag11,im01,kn0,ii0,im(3.9)00,k0,km1.kn当1/2r时，有10,kkllg1.kn(3.10)上式说明当和10niig很小时，误差1maxkkn也很小。摄动方程组(3.9)和差分方程(3.5-1)~(3.5-3)的形式完全一样。上述结果可叙述如下。当1/2r时，差分格式(3.5-1)~(3.5-3)关于初值和右端项在下述意义下是稳定的：设{|0,0}kiuimkn为差分方程组2,kkktixiiDuauf11,01,imkn0,iiu0,im00,0,kkmuu1kn的解，则有一维热传导方程的前向、紧差分格式第6页共25页100,kklluuf1.kn下面来考虑1/2r的情况。此时必存在0h当0hh时，2(1)1sin.22mhr于是2(1)14sin1.2mhr设0,kig11,01,imknsin(1),iimx0,im则（3.9）为20,kktixiDa11,01,imkn0sin(1),iimx0,im00,k0,km1.kn可以验证其解为20(1)14sin,2kkiimhr0,im0.kn由此易知20(1)14sin,2kkmhr20(1)14sin.2kkmhr由于当k时，2(1)14sin,2kmhr所以不论初始误差多么小，均会使解有较大的误差。我们得到如下结论。定理3.1.4当1/2r时，差分格式(3.5-1)~(3.5-3)是不稳定的。由定理3.1.3和定理3.1.4知，当步长比1/2r时，差分格式(3.5-1)~(3.5-3)是不稳定的；当步长比1/2r时，差分格式(3.5-1)~(3.5-3)是不稳定的。我们把这种稳一维热传导方程的前向、紧差分格式第7页共25页定性称为条件稳定性。实际计算时选取步长比必须满足1/2r，即2/1/2ah4.1.4数值例子应用向前Euler格式(3.5-1)-(3.5-3)计算定解问题220,01,01uuxttx(,0),01xuxex1(0,),(1,),01ttuteutet上述定解问题的精确解为(,)xtuxte.部分节点处数值解、精确解和误差的绝对值(1/10,1/200