赵妍,令锋
(1.内蒙古工业大学理学院,内蒙古呼和浩特010051;2.肇庆学院计算机学院,广东肇庆526061)
熔化凝固问题的Keller盒式格式
赵妍1,2,令锋2
(1.内蒙古工业大学理学院,内蒙古呼和浩特010051;2.肇庆学院计算机学院,广东肇庆526061)
建立了热参数依赖于温度的熔化凝固问题的显热容模型的Keller盒式格式,运用冻结系数法证明了该格式的无条件稳定性,并给出了数值算例.
显热容法;盒式格式;稳定性
数值求解熔化凝固问题的显热容模型[1]的控制方程具有强非线性性,其热参数是依赖于温度的函数.在数值求解中,发生相变的温度区间的选取对计算结果有较大影响.为克服这一缺点,1985年,Hsiao对显热容模型做了改进[2].Keller盒式格式是Keller在1970年提出的一种差分格式[3].1984年,Meek和Norbury用Keller盒式格式求解了非线性移动边界问题[4];2008年,Hamzah等给出了Keller盒式格式的并行算法[5].本文中,笔者在前人工作的基础上,对熔化凝固问题的显热容模型做进一步改进,用Keller盒式格式离散其控制方程,并分析该格式的无条件稳定性,最后给出数值算例.
根据显热容法,熔化凝固问题的控制方程可表示为
初始条件为
边界条件为
式中:
C(T)是依赖于温度的有界间断函数,Tm是相变临界温度.
故Keller盒式格式的截断误差为R=O(τ2+h2).
文献[2]657-661的研究结果表明,直接应用式(4)会导致数值误差或解的不稳定性,为此,笔者依据文献[2]657-661的方法对热容量做进一步改进.由文献[2]657-661(见图1)知,在(i,j)点的热容量C(Ti,j)用(i-1,j),(i+1,j),(i,j-1),(i,j+1)4个点热容和的平均值来代替,即
依据同样的原理(见图2),在(j-1/2,n-1/2)点的热容量C(Tnj--11//22)用(j-1,n),(j,n),(j-1,n-1),(j,n-1)4个点热容和的平均值来代替,即
其中:Tn-1/2表示节点(j-1/2,n-1/2)处的温度,C(Tn-1/2,Tn)表示修正热容量.
图1 改进的显热容法的节点示意图
图2 Keller盒式格式的节点示意图
定理1差分格式(12)是无条件稳定的.
证用“冻结系数”法[6]讨论其稳定性.
使用记号(7),将式(12)展开得
其中:λ=τ/h2.
令Tnj=vneikjh,代入式(15)得
将式(16)整理后得
为验证本文算法的有效性,以下用Keller盒式格式求解存在解析解问题[7](18)~(21)的数值解:
其解析解为
由积分公式及式(18)与(19),式(21)中的第2式可以写为如下形式:
控制方程(25)定义在固定边界0<σ<1,0<τ<H上.
边界条件(19)变形为
初始条件(20)变形为
边界条件(21)中的第1式变形为
式(23)变形为
利用Keller盒式格式,控制方程(25)可以写为如下形式:
控制方程(30)离散为下列2个差分方程:
在σ=0处的边界条件可近似表达为
初始条件(27)离散为
边界条件(28)离散为
对两相界面的守恒条件(式(29))应用矩形公式,得
求解方程组(31)~(36)时,取Δσ=0.01,Δτ=0.05,N=100,H=5,J=100.
迭代过程的一般步骤如下:
1)先给出试探的s(τn+1);
2)用TDMA算法求解方程(31)~(35),求得vnj+1,unj+1(j=0,1,2,…,J-1);
3)根据求得的unj+1,由式(36)计算出新的s(τn+1);
4)用新的s(τn+1)重复步骤2),直到满足收敛条件
图3u(3,t)时数值解和精确解之差
图3 为点(3,t)的数值解和精确解之差.从图3中可以看出,本文得到的数值解非常接近于精确解,表明笔者提出的方法是求解熔化凝固问题的有效方法.
[1]BONACINA C,COMINI G,FASANO A,et al.Numerical solution of phase-change problems[J].Heat Mass Transfer,1973,16:1 825-1 832.
[2]HSIAO J S.An efficient algorithm for finite-difference analyses of heat transfer with melting and solidification[J].Numerical Heat Transfer,1985(8):653-666.
[3]KELLER H B.A new difference scheme for parabolic problems,in Numerical Solutions of Partial Differential Equations[M]. 2nd ed.New York:Academic Press,1971:327-350.
[4]MEEK P C,NORBURY J.Nonlinear moving boundary problems and a Keller box scheme[J].SIAM J Numer Anal,1984,21:883-893.
[5]HAMZAH N,ALIAS N,AMIN N S.The parallelization of the Keller box method on heterogeneous cluster of workstations[J]. Journal of Fundamental Sciences,2008(4):253-259.
[6]陆金甫,关治.偏微分方程数值解法[M].2版.北京:清华大学出版社,2004:57-59.
[7]FURZELAND R M.A comparative study of numerical methods for moving boundary problems[J].J Inst Maths Applics,1980, 26:411-429.
Keller Box Scheme for Melting and Solidification Problem
ZHAO Yan1,2,LING Feng2
(1.College of Sciences,Inter Mongolia University of Technology,Hohhot Inter Mongolia010051,China; 2.School of Computer Science,Zhaoqing University,Zhaoqing Guangdong526061,China)
The Keller box scheme for the mathematical model of melting and solidification problem with temperature-dependent heat parameters is established.The unconditional stability of the scheme is analyzed by using the freeze coefficient method.And a numerical example is presented.
apparent heat capacity method;Keller box scheme;stability
O241.82
A
1009-8445(2010)02-0001-05
(责任编辑:陈静)
2009-10-09
广东省自然科学基金资助项目(04011600)
赵妍(1983-),女,黑龙江哈尔滨人,内蒙古工业大学与肇庆学院联合培养硕士研究生.