《钢铁冶金过程计算机仿真》由会员分享，可在线阅读，更多相关《钢铁冶金过程计算机仿真（100页珍藏版）》请在金锄头文库上搜索。
1、钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真温度场数值模拟温度场数值模拟北京科技大学冶金学院刘青北京科技大学冶金学院刘青北京科技大学冶金学院刘青北京科技大学冶金学院刘青备瑟茂样里赤柜秩捞踩拨绽孔拍膝饲稗辗撩叶语雁矛道那吴丧劈冀迢人睹钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真1前言前言、钢铁冶金加工过程、钢铁冶金加工过程、钢铁冶金加工过程、钢铁冶金加工过程炼铁炼钢（铸造、锻造、焊接、挤、压、拉、拔、热处理）。专业设置。、数值模拟、数值模拟、数值模拟、数值模拟两个不同领域的现象，能用同一数学形式来描述，称这两个现象彼此是可模拟。模拟的方法是把一个领域内求解的问题过渡到另一个领域中去解决。随着计算
2、机和计算数学的发展，用计算机数值计算法求解问题的计算精度已经达到很接近解析解，此法称为计算机数值模拟。仿真的定义，与模拟的差异。、工艺优化、工艺优化、工艺优化、工艺优化 计算机数值模拟的最终目的是解决工艺优化设计问题。一旦全面实现了冶金加工过程的计算机数值模拟，材料的加工生产将会产生深刻的变革。CAD/CAM/CAE的发展，波音767飞机的设计和加工。多手珍抓瑞姆付压恋颐积吸代援戍仁释彭轮漓龄畴募为滤缕请潮余咙巷肋钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真2前言前言、计算机数值模拟的步骤、计算机数值模拟的步骤、计算机数值模拟的步骤、计算机数值模拟的步骤给定研究对象几何条件、物理条件、初始条件
3、、边界条件离散化处理将过程所涉及的区域在空间和时间上进行离散化处理。建立数值方程建立内节点和边界节点的数值方程。选择计算方法选择适当的计算方法求解线性代数方程组。编制计算机程序编制计算机程序，由计算机算出结果，并用数据、曲线、图形输出。优化工艺分析计算机模拟的结果，如果结果不理想，调整工艺参数，再进行计算机模拟，直到模拟结果为最佳结果。俯奇耍象笑罩育骂奉舜枚羊棵噪慈柏既粕蹋郎守怯已祸港消冻服广骄础闻钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真3前言前言、讲解的主要内容、讲解的主要内容通过一、二个例子，让大家亲自经历计算机数值模拟的一般步骤，从中掌握其基本原理和方法，为更好应用计算机数值模拟，优
4、化工艺参数打下基础。冶金过程中，比较常用、典型的数值模拟有：温度场数值模拟浓度场数值模拟组织与性能模拟猜雇溃窗拦鹏删膀场葡嗓连怪发双嫁巡眯喉贸嫩契山绝宋公梁微原桂雪眉钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真4典饰札瀑怂庄箭烦甸蒜勋斋冀琐呻帆寿迷营此榔起等荆善囚砾丹嫌汾歼瞩钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真51 温度场计算机数值模拟温度场计算机数值模拟1.1 1.1 传热的基本知识传热的基本知识传热的基本知识传热的基本知识1.1.1 1.1.1 传热的基本方式传热的基本方式传热的基本方式传热的基本方式导热导热属于接触传热，是连续介质就地传递热量，没有各部部分物质之间宏观的相对位移。
5、在不透明固体实体内部，由于各部分物质之间无法作宏观的相对位移，不透明无法传递辐射能，实体保证接触，所以只能依靠导热方式传递热量。导热的基本定律是傅立叶定律，即导热的比热流量q与温度梯度成正比，即： qgrad T=-grad T= - (W/m2) （1-1-1） n式中：q传热方向上单位时间、单位面积上所通过的热量，J/(s.m2)=W/m2材料的导热系数，W/(m.K)温度梯度，K/mn负号表示导热的热流量向温度低的方向传递，即与温度梯度的方向相反。比热流量是个向量，即它有大小和方向。挣菜研呛拿住糜馒父汲哟刚撂钱寝救岛昼辊床愈颧荤拍始蓉霞梨诵懊蕴笋钢铁冶金过程计算机仿真钢铁冶金过程计算机仿
6、真61.1.1 1.1.1 传热的基本方式传热的基本方式如果比热流量的分量和（X、Y、Z）坐标系相联系，那么X、Y、Z方向的热流量分量应是： qx= - , qy= - , qz= - （1-1-2） X Y Z比热流量 q=iqx+jqy+kqz （ 1-1-3）导热系数 物理意义：沿导热方向的单位长度上，温度减低,物质所容许通过的热流量。方向性：大多数液体和固体属于各向同性的物质。各向异性材料的导热系数具有方向性，如石墨。温度函数：值还随温度而变化。大多数金属的导热系数随温度的升高而降低。大多数液体（水和甘油除外）其导热系数随温度的升高而降低。 气体的导热系数随温度的升高而增加。 罢掖卿丝
7、淡鉴券横沾眉辫绅业叉隅喝虹插永嚣羌凛滋狄估零啄懈裕矫甥疗钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真71.1.1 1.1.1 传热的基本方式传热的基本方式对流对流是流体（气体和液体）中温度不同的各部分相互混合的宏观运动引起热量传递的现象。工程上最具有实际意义的是：相对运动着的流体与所接触的固体壁面之间的热量交换过程，一般称为对流换热。工程上在研究固体壁面和流体之间的对流换热时，除了高度稀薄的气体外，人们不去注意流体的单个质点，而把流体看成是连续介质。实际的流体总有粘性，流动时，受粘性和壁面摩擦的影响，在靠近壁面附近的流体将降低流速，在壁面上完全被滞止不动，即X=0时，=0，如图1-1所示。因此
8、，热量从壁面传给贴壁的那部分流体，将依靠导热。T(K)V(m/s)TwqwcTfVX(m)图1-1邻近壁面的流体速度分布和温度分布音蜀秤匈烙随件夷秤脑颊馏造旭灼惠靴烯鸦菱娘皇泰睁饵华迈辰灵占铁叁钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真81.1.1 1.1.1 传热的基本方式传热的基本方式流体和固体壁接触面上的“相”际热流密度为： qwc= - f x=+0 （1-1-4） X式中：qwc热流密度，W/m2 f流体的导热系数，W/(m.K) T流体的温度，K液体的导热系数较小，气体的导热系数更小，所以受热时，在贴壁处的流体温度势必沿着轴的反方向急剧升高。随着离壁面的距离的增加，流体的流动将带
9、走更多的热量，使轴方向的温度梯度连续下降，温度分布趋向平坦化。正是通过这种导热和对流的共同作用，使热量在流体内部得到传播，越临近壁面，导热越起主导作用。图1-1所示，假如厚度为的贴壁静止膜，膜内温度线性地从壁面温度TW降到远离壁面，尚未被加热的流体温度Tf，则 Tf-TW qwc= - f （1-1-5） 无界对流时壁面与流体的换热，钢锭与周围空气的对流换热属于这种情况。攻雍搁烹挫蓖湖讹姻獭党氟滞搔井子附皋高四料综亿彼醛理黑谣淫兼滴聋钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真91.1.1 1.1.1 传热的基本方式传热的基本方式流体在管和槽道内部的流动，称为有界对流，这时不存在远离壁面，尚未
10、被加热的流体温度，则采用截面平均温度作为流体的特征温度Tf，则qWC=aC(TW-Tf), W/m2 （1-1-6） 这就是所谓的“牛顿冷却定律”。式中：aC为放热系数，W/(m2.K)其实“牛顿冷却定律”并不是表达对流换热现象本质的物理定律。凡能影响流体流动的各种因素，包括流体的种类和状态、运动的速度、流道的形状和大小，以及固体壁面粗糙度等，都会对aC值产生影响。式1-1-6只不过形式地把放热过程的一切复杂性和计算上的困难，都转移到并且集中在放热系数这样一个物理量上罢了。aC代表流体和所接触的固体表面之间温度每相差，该流体与表面之间“相”际热流量的大小。运用式1-1-6可以进行对流换热的计算
11、。但由于对流换热的复杂性，该式中的放热系数aC需从相应的准则方程式求出。准则方程式是针对不同对流换热情况，在综合实验结果的基础之上，运用相似理论将表征某现象的物理量整理成一些相似准则，用因次分析法得到的各种类型的表达式。刑胡涎纱缚擒芬陈脸凉瘫稼涟瞎钧葫玻亚京施碌萨倍牙寅鉴辈馆掇嫩棠曳钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真10辐射只要温度高于绝对零度，任何物体都随时向周围空间辐射能力。辐射用斯蒂芬玻耳兹曼定律表达： E=，w/m2 1-1-7 式中：物体的辐射率或黑度,(0-1)；斯蒂芬玻耳兹曼常数或绝对黑度的辐射常数； 5.6710-8 W/m2 K4温度,K。实际上，辐射往往涉及二个物
12、体间辐射热交换。如果二个物体的表面温度不同，中间被空气所隔开，T1T2时，则相互辐射的结果，表面温度为T1的物体发射出去的辐射热超过了吸收来自表面温度为T2的物体辐射热，引起辐射换热的热流量则为： Q12= 12(1-2)F112 或q12= 12(1-2)12 1-1-8 式中：12 物体1与2综合黑度； F1物体1的表面积；12 物体的表面向外发射出去的辐射热量中，能投射到物体表面上的份额，称为角系数,(0-1)。1.1.1 1.1.1 传热的基本方式传热的基本方式曹瓢诉袄怕行颇拼继科士泽铣恒牵鲸韭移梳沙贮酷跨查辰脊质擎外肯兰窝钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真111.1.1 1
13、.1.1 传热的基本方式传热的基本方式壁面在气体中冷却，存在辐射换热和对流换热。考虑到壁面与气体之间存在着辐射换热，其热流密度为Qr=arF(Tw-Tf) , w 或qwr=ar(Tw-Tf) , w/m2 1-1-9式中：qwr单位面积的辐射量，w/m2 ar辐射放热系数，w/(m2.k)Tw辐射物体表面温度，k Tf透明的气体介质的特征温度，k考虑到壁面与气体之间还存在着对流换热，其热流密度为qwc=ac(Tw-Tf) , w/m2 1-1-10 由壁面传走的总热流密度qw应是qwr和qwc二者之和qw=ac(Tw-Tf)+ar(Tw-Tf) 1-1-11 令a0=ac+ar 则qw=a0
14、(Tw-Tf) ,w/m2 1-1-12用辐射放热系数ar，可以形式地把辐射换热折合成对流换热，用总放热系 奏荡谋栖写锣盏盲低苫太戴维墅化限阁诉锐乎梧灵昨醒仿均释向鹊僚呢竭钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真121.1.1 1.1.1 传热的基本方式传热的基本方式 数a0兼顾辐射换热的影响，从而有利于简化复杂传热的分析和计算。如远离表面的外界表面温度趋于环境温度Tf，并且12=1时，由式1-1-8得 qwr= 12(w-f) 1-1-13由式1-9得：qwr=ar(w-f) 1-1-14 由式1-1-13和式1-1-14得： 12(w-f)= ar(w-f) 1-1-15 因w-f=(
15、w-f)(w3+w2f +wf2 +f3) 设m=1/2(w+f)，(w3+w2f +wf2 +f3) 4m3 ar 4 12 m3,w/(m2.k) 1-1-16从此可见，(w-f)降为零时，ar并非零值，而是以4 12 m3。随着温差的扩大和平均温度m的升高，ar值将迅速增加。由于ac随温差的变化较小，在高温范围和大温差情况下，ar有可能成为a0的主要组成部分。ar与气体的运动状况无关，而ac与气体速度的降低而减小，在气体自然对流的情况下，ac 5,w/(m2.k)，即使m只有300K，ar就可能占a0的一半左右。 岭搂娘篮稿扭肆择集矽巷足筹财荔堪介匠洛琢诈往甄饭回擦信迂袖邹危依钢铁冶金过
16、程计算机仿真钢铁冶金过程计算机仿真131.1.2 1.1.2 导热的偏微分方程式导热的偏微分方程式钢锭一般属于各向同性的物质，其加热或冷却过程数学模拟计算依据的基本数学模型，是不稳定导热偏微分方程。下面讨论各向同性材料导热方程式的建立。直角坐标系直角坐标系直角坐标系直角坐标系导热导热导热导热偏微分方程偏微分方程偏微分方程偏微分方程导热偏微分方程的建立，是通过考察处于导热过程中的物质的微元体积（xyz)的能量平衡来建立。如图1-2所示。在时间内，通过六个面的导热所获得的能量，加上微元体内产生的内热源热量，要等于微元体积内物质积蓄热量的改变，即温度升高或降低。xzydQx+xdQxdQzdQydQ
17、z+zdQy+y图1-2直角坐标系导热方程式的微元体矫缩板镭镀扳搐峡呻搭啸禁隅垮图泅古茸布褐吠啡颈瞄辙饼预幕板肮洁羊钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真141.1.2 1.1.2 导热的偏微分方程式导热的偏微分方程式图1-2中，微元控制体尺寸x、y、z，按照傅立叶导热定律，在X方向流入微元体左表面的热流可表示为： TdQx=-(yz)，W 1-1-17 X 式中： 导热系数，W/(m.k)yzX方向微元体表面积，m2T X方向的温度梯度，k/mX在X方向流出微元体右表面的热流，可以应用泰勒级数展开，保留级数的第一项和第二项而得到： TdQx+x= dQx + (dQx)x=dQx +
18、-(yz) x X X X T = dQx - ( )xyz1-1-18 X X 在X方向导热的净热流为 T dQx - dQx+ = ( )xyz1-1-19 X X 蝎落链郝庞闪吟卿痕些冠撞依瓶弧剐冉裙惕戒库屡勘廓祷奏诣炔椿雪声惨钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真151.1.2 1.1.2 导热的偏微分方程式导热的偏微分方程式用同样的方法，可以得出Y,Z方向与式-19相类似的导热的净热流方程式，即 T dQy - dQy+y = ( )xyz1-1-20 Y Y T dQz - dQz+z = ( )xyz1-1-21 Z Z 在三个坐标方向净热流的总和为： T T T ( )+
19、 ( ) + ( ) xyz1-1-22 X X Y Y Z Z如果单位时间、单位空间所产生的热量为Q(x,y,z,)，那么微元体的发热量为：Q(x,y,z,)xyz1-1-23由于导热传进微元体的净热流式1-22和微元体内产生的热量式1-23一起用于增大微元体的内能。微元体的内能的增加表现在微元体能量存储随时间的变化率，即 T Cp xyz 1-1-24 辑陶伎靶哗淘顿秀酝贰耙沈松瞪愿彩渗掣埋项何陵波毖卯繁围年民域纷牵钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真161.1.2 1.1.2 导热的偏微分方程式导热的偏微分方程式式中：密度，kg/m3 Cp比热，J/(kg.k) 时间，s 对微元
20、体进行能量平衡，使能量存储的时间变化率与由导热引起的流入微元体的净热流和微元体内产生的热量之和相等，得下式： T T T T Cp = ( )+ ( ) + ( )+Q 1-1-25 X X Y Y Z Z 圆柱坐标系圆柱坐标系圆柱坐标系圆柱坐标系导热导热导热导热偏微分方程偏微分方程偏微分方程偏微分方程实际问题常常涉及柱面对称问题，且边界条件给定在一个表面上，因此表面具有一个坐标保持不变的性质。在这种情况下，采用圆柱坐标系是合适的。图1-3圆柱坐标系导热方程式的微元体肮远脑骇寺陀绸匈儡臻蝉手澜钨主赦忧坞峨佐荧们近溯焙略深肢橇写限崇钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真171.1.2 1.
21、1.2 导热的偏微分方程式导热的偏微分方程式图1-3所示圆柱坐标系，直接按内外两个圆弧面和其它四个平面组成的微元体积，在导热过程中热量必须按收支平衡的原则导出，此时，微元体积为：dv=(r d)dzdr 1-1-26 沿内圆弧面流入微元体积的热流：TdQr=-(r ddz) 1-1-27 r沿外圆弧面流出微元体积的热流： TdQr+dr= dQr + (dQr)dr = dQr + -(r ddz)dr r r r T 1 T = dQr + (-r) ddzdr = dQr- (r) dv 1-1-28 r r r r r 沿平面流入微元体积的热流：TdQ=-(dr dz) 1-1-29r沿
22、+d平面流出微元体积的热流： TdQ+d= dQ + (dQ)rd = dQ + - (drdz) rd r r r踏枢尝邻散欲拦故痊柒陋拨罗疙祖刀版探肌牛羊在政炭茨青雌悔仲故熊绅钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真181.1.2 1.1.2 导热的偏微分方程式导热的偏微分方程式 1 = dQ - () dv 1-1-30 r2 沿z面流入微元体积的热流： TdQz=-(r ddr) 1-1-31 z沿z+dz面流出微元体积的热流： TdQz+dz= dQz + (dQz)dz = dQz + -(r ddr)dz z z z T = dQz - () dv 1-1-32 z z 根据
23、直角坐标系导热微分方程推导的思路，可得到： 1 T 1 T T (r) + () + () + Q = Cp 1-1-33 r r r r2 z z 逸者器蹿埔邑肥刺煮速吠帝俩藉溪赋恭图菠纬糊匙赔酪瞄螟眨势木柑瞧淹钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真191.1.2 1.1.2 导热的偏微分方程式导热的偏微分方程式球坐标系球坐标系球坐标系球坐标系导热导热导热导热偏微分方程偏微分方程偏微分方程偏微分方程对于球坐标系，如图1-4所示。图1-4球坐标系导热方程式的微元体由内、外两个球面、两个圆弧面和两个平面所组成的微元体积为： dv=(r d)（rSind)dr 1-1-34监粳添亏聪酱空存薯
24、耿然捍懒乡薪闰钝厚浇朽局惺怂胞竿祥俄遥蜡堂陛挫钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真201.1.2 1.1.2 导热的偏微分方程式导热的偏微分方程式沿半径方向流入微元体积的热流： TdQr=-(r Sind.rd) 1-1-35 r沿半径方向流出微元体积的热流： TdQr+dr= dQr + (dQr)dr = dQr + -(rSind.rd)dr r r r T 1 T = dQr+(-r2)Sind.ddr=dQr- (r2)dv 1-1-36 r r r2 r r 沿方向流入微元体积的热流：TdQ=-(rSind.dr) 1-1-37 r沿方向流出微元体积的热流： TdQ+d=
25、dQ + (dQ)rd = dQ + -(rSinddr)rd r r r T T= dQ- (Sin)rddr.rd=dQ- (Sin)dv r2 Sinr2 1-1-38 游辩桓烂滓丧马柠柬爬娘塘畏敛翘射隆瞥亢墅成腺俄瀑涉钡撮互镶抓械我钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真211.1.2 1.1.2 导热的偏微分方程式导热的偏微分方程式沿方向流入微元体积的热流：TdQ=-(dr.rd) 1-1-39 (rSin)沿方向流出微元体积的热流： dQ+d= dQ + (dQ).(rSin)d (rSin) T = dQ + (-dr.rd)(rSin)d (rSin) (rSin) T =
26、 dQ - ()dv 1-1-40 (r2Sin2)同理可整理得： T T T(r2)+(Sin)+()+Q = Cp r2r r r2Sin r2Sin 1-1-41 憎翘箍匹底憨豹氟祟岂郎肉萄菱悔雕使朔供窘贰惨拇总模君述次俭五域我钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真221.2 1.2 导热方程的有限差分解法导热方程的有限差分解法求解不稳定导热偏微分方程的数值解法，主要有：有限差分解法、有限元法、边界元法三类。边界元法正在研究和完善之中，目前常用的是有限差分解法和有限元法。我们专门讨论有限差分解法的数学基础，数值方程的建立，差分方程的稳定性和收敛性等问题。有限差分解法是用差分方程近似
27、地代替微分方程，建立差分方程有直接法和能量平衡法两种。1.2.11.2.1直接代换法直接代换法直接代换法直接代换法直接代换法是从微分形式出发，用差商直接代换微商的办法建立差分方程。1.2.1.1数学基础微商和差商的定义若T(x)是连续函数，则其导数为：dT(x+x)-T(x)=lim=lim1-2-1dxx0xx0x疗陵砾特赂是鲸沸蔫饭拱哥党章桔谰褂传洗棺缄那广刀烁府讼还粪傻奥吼钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真231.2.1.11.2.1.1数学基础数学基础1-2-1式右边/x是有限的差商，x和都不为零。式左边d/dx是T/x当x0时的差商，称之为微商。在x没有达到零之前，T/x只
28、是dT/dx的近似。实际应用x0。如果把T/x趋于dT/dx的过程认为是近似向精确过渡，那么，用T/x代替dT/dx，则两者的差值T/x-dT/dx表示差商代替微商的偏差。误差多大？需要做误差分析，才能大胆地应用。误差分析假设函数T(x)在x时的值为T(x)，在x+x时的值为T(x+x)，如果函数T(x)在x处的各阶导数存在，则按照泰勒级数展开，T(x)与T(x+x)的关系如下式所示：dT(x)2d2T(x)ndnTT(x+x)=T(x)+x+1-2-2dx2!dx2n!dxn整理后得：TT(x+x)-T(x)dTxd2T(x)n-1dnT=+1-2-3xxdx2!dx2n!dxn从上式可知，
29、T(x)在x处的差商T/x等于函数T(x)在x处的各阶导数的线性组合，只能是近似地等于差商。两者之间也必然有偏差。崖俄拽探馁仗孰氖层施篮设猴元疲曲倪企啪观啮嫂彭惹冻导睡林推飘钝屯钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真241.2.1.11.2.1.1数学基础数学基础图1-2-1表示了一阶差商与一阶微商之间的关系。用不同方法得到的差商去代替微商，它们带来的误差是不同的。即向前差商：dT/dxT(x+x)-T(x)/x1-2-4向后差商：dT/dxT(x)-T(x-x)/x1-2-5中心差商：dT/dx1/2T(x+x)+T(x)/x+T(x)-T(x-x)/x=T(x+x)-T(x-x)/2
30、x1-2-6图1-2-1差商与微商剥正杠予楼菩稍甘感彰踊好埂诞为迈动苍役芬寄演区詹番赞嗣清妥未苗懊钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真251.2.1.11.2.1.1数学基础数学基础按照泰勒级数展开，T(x)与T(x+x)的关系如下式所示：dT(x)2d2T(x)ndnTT(x+x)=T(x)+x+1-2-7dx2!dx2n!dxn整理后得：T(x+x)-T(x)dTxd2T-=+=O(x)1-2-8xdx2!dx2即向前差商的偏差是截去了泰勒级数展开式中的高阶项而引起的，常称“截断误差”，其截断误差为与x同级的小量O(x)。同理dT(x)2d2T(x)3d3TT(x-x)=T(x)-
31、x+-+1-2-9dx2!dx23!dx3整理后得：T(x)-T(x-x)dTxd2T-=-+=O(x)1-2-10xdx2!dx2即向后差商的截断误差为与x同级的小量O(x)。青轿逻瀑记蚕论辖雏周冯胳评薪引闭咎宾鞘桃瓶胶缉邯纤囱夏缝乱原押屠钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真261.2.1.11.2.1.1数学基础数学基础由式1-2-7减式1-2-9，将2xdT/dx移至等式左边，两边再除以2x，得：T(x+x)-T(x-x)dT(x)2d3T-=+=O(x)21-2-112xdx3!dx3即中心差商的截断误差为与(x)2同级的小量O(x)2。当x固定时，上述一阶差商一般仍是x的函数
32、，对它们还可以求差商。这种一阶差商的差商称为二阶差商，它是二阶微商的近似，常用向前和向后差商来二阶微商，即d2TT(x+x)-T(x)T(x)-T(x-x)-/xdx2xxT(x+x)-2T(x)+T(x-x)=1-2-12(x)2由式1-2-7和式1-2-9相加，经简化后得：T(x+x)-2T(x)+T(x-x)d2T-=O(x)21-2-13(x)2dx2即二阶差商的截断误差为与(x)2同级的小量O(x)2。拘叠球梯荷绸龙耐插那已牺悍娄抛嘶洽屿趾镍挝好焕拘彼搏尉埋彤汤垒脆钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真271.2.1.2 1.2.1.2 建立内节点差分方程建立内节点差分方程一维
33、系统一维系统假定有一高宽无限(即高宽方向上无热流)，厚度为L的平板，T=f(x,)即温度是x方向位置和时间的函数，所谓一维系统是指几何空间为一维。初始时刻=0，T=T0，为了简化，考虑无内热源，、Cp均为常数。选取网格点间距x和时间步长将研究对象离散化。显式差分格式 T T 一维不稳定导热方程为 Cp = ( ) X X 该方程在区域0，0xL内全部点都成立，如图1-2-2所示。将方程应用于内节点i可写成： T p 2T p Cp() = ( ) 1-2-13 i X2 i上式等号左端的微分式用温度对时间的一阶向前差商来近似，即：T p Tp+1i -Tpi () = + O() 1-2-14
34、 i 上式等号右端的微分式用温度对空间的二阶差商来近似，即：T p Tpi+1-2Tpi+Tpi-1 () = +O(x)2) 1-2-15 x2 i (x)2碗荷测碾钢誓浩所欲报闹暮钮差票弄矩尚谗量栈猿施强陵潦趣栋箕印舍梗钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真281.2.1.2 1.2.1.2 建立内节点差分方程建立内节点差分方程图1-2-2一维显式差分将式1-2-14和式1-2-15代入式1-2-13得： Tp+1i Tpi Tpi+1-2Tpi+Tpi-1 Cp+O()= + O(x)2) 1-2-16 (x)2若在上式中去掉O()和O(x)2),整理得：Tp+1i= Tpi+F0
35、(Tpi+1-2Tpi+Tpi-1) 1-2-17 式中F0=/Cp(x)2=导热速率/蓄热速率，称F0为傅立叶数。闭斑喜狱拟刀示人伶块桐美捞抚笔斯斌丢厨乒酝增御柠臼启继基侧咕朱支钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真291.2.1.2 1.2.1.2 建立内节点差分方程建立内节点差分方程F0的数值小意味着加热或冷却此物体所需要的时间长，反之，所需要的时间短，F0是一个无因次数字。当初始条件和边界条件已知时，用式1-2-17就可模拟区域内各节点随时间增长的温度值Tpi(i=2,3,，N-1；p=1,2,3,）隐式差分格式一维隐式差分如图1-2-3所示，将一维不稳定导热微分方程应用于内节点
36、i，则：Tp2Tp+1Cp()=()1-2-18iX2i式1-2-18和1-2-13相比，式的左边完全一样，温度对时间的一阶偏微商，仍用一阶向前差商来近似，而式1-2-18和1-2-13右边有所不同。式1-2-13中温度对距离的二阶偏微商是对应时刻p的，在用差商近似微商时，用p时刻的T值；而式1-2-18中，温度对距离的二阶偏微商是对应时刻p+1的，用差商近似微商时，用p+1时刻的T值。即式1-2-18相应的差分方程为：Tp+1iTpiTp+1i+1-2Tp+1i+Tp+1i-1Cp+O()=+O(x)2)1-2-19(x)2若在上式中去掉O()和O(x)2),整理得：Tp+1i=Tpi+F0
37、(Tp+1i+1-2Tp+1i+Tp+1i-1)1-2-20绑逢笛狙番夜蹲廖止素禹肿嵌比朔华酿裕武旭播滩扦蜗硅俱疵烦再拂铭喧钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真301.2.1.2 1.2.1.2 建立内节点差分方程建立内节点差分方程式中F0=/Cp(x)2。比较式1-2-17和1-2-20，可以看出显式差分格式的突出优点是每个节点方程都可以独立求解，整个计算过程十分简便。但是，的选取要受到限制，有时为了满足差分格式稳定性条件，可能选得很小，使计算工作量加大。隐式差分格式的最大优点是，它对任意F0值都是稳定的。这种稳定是绝对的，即不受边界条件、x、热物参数的影响，于是可以选的较大，计算速
38、度加快。但是，对于节点i，只从式1-2-20不能独立求解，必须涉及Tp+1i+1、Tp+1i、Tp+1i-1的联立线性代数方程组才能求解，也就是说，它含有三个未知数。时间步长每前进一步，从坐标0xL，网格点1iN整个区域的每个点，上述方程都要列出一次（见图1-2-3)。因此，向一个新的时间步长每移动一步就必须解一个方程组。当按顺序列出这些方程时，除了要的第一点和最后一点的方程只有两个未知数外，其余每一个方程都含有三个未知数，于是方程是三对角的。对于这种情况，已有适用的求解程序，后面将讨论。关于差分方程稳定性将在后面讨论。虚涂是谭琉娘撮滓掠辐迂庄输疮聂粳圭绢凰蓄鲍汉崖钞容满误卜献竖贾陶钢铁冶金过
39、程计算机仿真钢铁冶金过程计算机仿真311.2.1.2 1.2.1.2 建立内节点差分方程建立内节点差分方程图1-2-3一维隐式差分二维系统二维系统显式和隐式差分格式建立的方法和两种差分格式的特点在前面讨论过，对二维系统同样适用。为简略起见，在此只讨论二维系统显式差分方程在此只讨论二维系统显式差分方程的建立的建立。为使问题简化，仍然假定热物性值为常数，考虑无内热源。首先把所讨论的区域离散化，如图1-2-4所示。门嚷蜡骸堂娱间冀生朋廊玖爽便象奶裤咸湾屯客骨谆驭郝悸别郸惫嘶丛葛钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真321.2.1.2 1.2.1.2 建立内节点差分方程建立内节点差分方程网线：用
40、平行于X、Y轴的直线（网线），进行空间离散化节点：网线与网线的交点步长：节点间的距离图1-2-4二维均质物体的分割x与y可以是不变的常量，即等步长，也可以是变量（即在区域内的不同处是不同的），即变步长。如果区域内各点处的温度梯度相差很大，则在温度变化剧烈处，网格布得密些，在温度变化不剧烈处，网格布得疏些。至于网格多少，步长取多少为宜，要根据计算精度与计算工作量等因素而定。对于0xL1和0y1/2，由上式可见，当p时刻i节点温度Tpi越大，因为(1-2F0)为0时，全部内节点的温度Tp+1i的值总是处于Tpi+1，Tpi，Tpi-1 三个值的最大值与最小值之间的某个中间数值。这一事实显然是符合热
41、力学第二定律的。为什么隐式差分方程无条件稳定？将隐式差分方程稍加整理，可得：僻徘衅趾疽杯性加锋篮暮寒拓猴肆垮皆艇坟紫釉漂巍纷赎蹦脐燃衅舜惮恰钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真541.2.51.2.5差分方程的稳定性和收敛性差分方程的稳定性和收敛性Tpi= Tp+1i+F0(2Tp+1i-Tp+1i+1-Tp+1i-1) 1-2-75假定p+1时刻，在区域内某一节点i处取得区域内最大温度值，即Tp+1i为这一时刻区域内的最大温度,Tp+1i Tp+1i+1 ；Tp+1iTp+1i-1,式1-2-75等号右边括号内三项的代数和必然大于零，从而必然有TpiTp+1i。也就是说，在在p p时
42、时刻，区域内的最大温度必然大于刻，区域内的最大温度必然大于p+1p+1时刻区域内的最大温度时刻区域内的最大温度。依此类推，必然将最大温度值推到初始条件。假定p+1时刻，在区域内某一节点i处取得区域内最小温度值，即Tp+1i为这一时刻区域内的最小温度,Tp+1i Tp+1i+1 ；Tp+1iTp+1i-1,式1-2-75等号右边括号内三项的代数和必然小于零，从而必然有TpiTp+1i。按照上面的分析方法可知，整个过程的温度最小值，必然出现在初始条件中。总之式1-2-75不管F0为何值，它的运算逻辑都是符合热力学原理的。即隐即隐式差分格式是无条件稳定的。式差分格式是无条件稳定的。由上作出结论：隐式
43、差分格式是无条件稳定的，显式差分格式是有条件稳隐式差分格式是无条件稳定的，显式差分格式是有条件稳定的。定的。这一结论，原则上对二维、三维系统也是适用的原则上对二维、三维系统也是适用的。显式差分格式的稳定性条件是要求Tpi项的系数0。由式1-2-17（一维）、式1-2-25（二维）、式1-2-28（三维）和式1-2-52（二维对流边界）等差分方程可导出具体的稳定性条件如下：内部节点：到谢蟹举昔钧太句类鄂肌酚账宏蚤块衷巴孵紊朗蛇儡轻丙刑轴水岛矫肿她钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真551.2.5 1.2.5 差分方程的稳定性和收敛性差分方程的稳定性和收敛性一维直角坐标：F01/2二维直角
44、坐标：F01/4三维直角坐标：F01/6边界节点（设为对流边界）11一维直角坐标：F0二维直角坐标：F02(1+acx/)2(2+acx/)1三维直角坐标：F02(3+acx/)分析ac=0得到内部节点判据；当其它条件不变时，ac增大，收敛条件要求更小，即更严格。根据这一原则，其它边界条件下边界节点差分方程的稳定性也可以相应导出。对于圆柱坐标系、球坐标系的显式差分方程，也要按上述原则分别导出具体的稳定性条件。总之，用显式差分格式数值求解导热问题，离散化处理时，时间步长总之，用显式差分格式数值求解导热问题，离散化处理时，时间步长和和空间步长空间步长x的选取是受到稳定性条件制约的。的选取是受到稳定
45、性条件制约的。一般的习惯是先定空间步长，再计算时间步长。馁鞍歹难枣皿狙到易矽主泄因进魄吗洼俄惹黔良频疙乌担门宴白乡级瀑窒钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真561.2.6 1.2.6 二维交替隐式格式二维交替隐式格式显式差分格式的优点：每个节点方程能独立求解，整个计算过程十分简单。缺点：为满足稳定性条件，空间步长和时间步长的选取必须遵守下列关系：(x)2Cp一维系统：(F01/2)2(x)2Cp二维系统：(F01/4)4(x)2Cp三维系统：(F01/6)6例如，一维系统，x=1mm，=7840kg/m3，Cp=610J/(kg.)，=30J/(s.m.)，求得：0.079,s。即1h
46、约5万次。这就表明：为了满足精度要求，当网格划分较细时，时间步长大大减小，当维数愈高，在同样的空间步长条件下，时间步长也愈小。因此，为了为了得到稳定而又比较精确的答案，就要进行大量的计算得到稳定而又比较精确的答案，就要进行大量的计算。为了消除稳定性时间步长的限制，可以采用隐式格式求解为了消除稳定性时间步长的限制，可以采用隐式格式求解。但是，在每惠垫揪颂棕泥举菠冲嫌鸳褒戎斌绰丘溅罩摔诡剑看毫副斯淀腥谴犬措豁拒钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真571.2.61.2.6二维交替隐式格式二维交替隐式格式一个时间层上，对二维系统就要解一个五对角的N*M个未知量的代数方程组(N为x方向的网格数，
47、M为y方向的网格数)，即移动一个时间步长，就必须解整个方程组，每一个方程都有五个未知数和一个对角占优元素。这样一来，就没有解三对角方程组那样方便。因此，需要建立新的差分格式，希望这一差分格式能满足如下要求：无条件稳定；有合理精度的解；所产生的代数方程组易于求解。交替方向隐式方法正是满足这些要求的一种差分格式交替方向隐式方法正是满足这些要求的一种差分格式。（简称IAD）1.2.6.11.2.6.1差分方程的建立差分方程的建立差分方程的建立差分方程的建立交替方向隐式方法需要对给定时间步长推导出二组有限差分方程。这些方程是显式和隐式项的混合，对于前半个时间步长，x方向各项是隐式形式，y方向各项是显式
48、形式；而在下半个时间步长时则相反。若把二维导热微分方程式应用于节点(i,j)，可以写成如下形式：第一个1/2 T p 2T p+1/2 2T p Cp() = () +() 1-2-76 i,j x2 i,j y2 i,j Tp+1/2i,j Tpi,j Tp+1/2i+1,j-2Tp+1/2i,j+Tp+1/2i-1,j Cp()=+ /2(x)2 Tpi,j+1-2Tpi,j+Tpi,j-1 1-2-77 (y)2边骄饰睛铭孵熟突掘佛歌摈缅沼运埔疏敦磅牡渊今灌彭里燎尚篮扇样彭担钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真581.2.6 1.2.6 二维交替隐式格式二维交替隐式格式第二个1/
49、2 T p+1/2 2T p+1/2 2T p+1 Cp() = () + () 1-2-78 i,j x2 i,j y2 i,j Tp+1i,j Tp+1/2i,j Tp+1/2i+1,j-2Tp+1/2i,j+Tp+1/2i-1,j Cp()=+ /2(x)2 Tp+1i,j+1-2Tp+1i,j+Tp+1i,j-1 1-2-79 (y)2从上可见：第p+1/2层是一个过渡时间层，通过式1-2-77计算过渡时间层的值，然后用式1-2-79计算第p+1时刻的值。由于现在有一种迅速解三对角线矩阵的算法，所以交替方向隐式法是有效的。1.2.71.2.7三对角线矩阵的解法三对角线矩阵的解法三对角线
50、矩阵的解法三对角线矩阵的解法前面介绍了显式、隐式和交替方向隐式法等差分格式，从数学方面来看，它们都是一个线性代数方程组。线性代数方程组的求解一般有两类方法线性代数方程组的求解一般有两类方法：一类是直接方法直接方法，如三对角线矩阵的解法、高斯消去法、矩阵分解法；另一类是迭代方法迭代方法，如高斯塞德尔迭代、超松弛迭代、强隐式迭代。一般说来，直接方法只要进行有限此运算就可以得到方程的解。它的缺点是必须有较大的计算机存储量和花费较多的计算时间，而且在使用时，肋滚僧凰钧窄帘遂翻穿爱影锦肿舱拴战轮论泛讫侧苹烟停推摔肮围奇封侩钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真591.2.7 1.2.7 三对角线矩
51、阵的解法三对角线矩阵的解法还必须注意抑制舍入误差的影响，可以采用增加计算数字的位数，或利用主元消去等手段。迭代法没有舍入误差积累的问题，一般需要较少的计算机存储量，计算时间也比较节省。但是，如何能确保迭代收敛和收敛得比较快，这在使用中必须注意解决。当离散化方程是三节点关系式时，高斯消去法就转换成人们所常用的三对角线矩阵法，这种方法也称为托马斯(Thomas)算法或TDMA、追赶法。它之所以称为三对角线矩阵算法，是由于用矩阵形式写出这些三节点关系式时，所有系数都是沿矩阵的三对角线排列的。设离散化方程是如下形式的三节点关系式：-aiTi-1+biTi-ciTi+1=di1-2-80i=1,2,N，
52、其中i=1和i=N分别为边界节点。由于T0和TN+1没有定义，故取a1=0和CN=0。因此，当i=1时，式1-2-80变成T1和T2之间的关系式，即可以用T2表示T1，T1=f1(T2)。（一个方程二个未知数）当i=2时，式1-2-80变成T1、T2和T3之间的关系式，可以用T2表示T3，T2=f2(T3)。（以上二个方程三个未知数）获番态局竹牛车千蚕窑谢况怀吓彩绸裙夸充唾吮溜宁对孪蒜蚜座禄岩禁难钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真601.2.7 1.2.7 三对角线矩阵的解法三对角线矩阵的解法这种替换关系，可以一直进行到TN用TN+1来表示，即TN=fN(TN+1)，由于CN=0(T
53、N+1的系数)，于是可以求得TN的值。（以上增加一个方程，没有增加未知数，所以可求出TN）而TN-1可从TN求得，TN-2可从TN-1求得，T可从T3求得，T1可从T2求得，这就是三对角线矩阵解法的实质。具体公式推导如下：设待求的递推替换关系式为：Ti=PiTi+1+Qi1-2-81则Ti-1=Pi-1Ti+Qi-11-2-82将式1-2-82代入式1-2-80，得-ai(Pi-1Ti+Qi-1)+biTi-ciTi+1=di，经整理得：Cidi+aiQi-1Ti=Ti+1+1-2-83bi-aiPi-1bi-aiPi-1比较式1-2-83和式1-2-81得Cidi+aiQi-1Pi=，Qi=
54、1-2-84bi-aiPi-1bi-aiPi-1式1-2-84是个递推关系式，因为它们是用Pi-1，和Qi-1来表示Pi，和Qi的。梧李畜至凸拒颜兰案哉闻肮磺喉拥轰束差橇特爬兆卓饿佑魔冰闻拾湃砂乱钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真611.2.71.2.7三对角线矩阵的解法三对角线矩阵的解法根据以上分析，其求解步骤是：取i=1，因为a1=0，所以P1=C1/b1，Q1=d1/b1；对于i=2,3,N-1，用式1-2-84分别求出Pi和Qi；当i=N时，因为CN=0，所以PN=0，由式1-2-81得TN=QN，即求出点的温度；对于i=N-1,N-2,2,1，用式1-2-81分别求出TN-
55、1，T2，T1。1.2.7.11.2.7.1三对角线矩阵解法程序三对角线矩阵解法程序三对角线矩阵解法程序三对角线矩阵解法程序TurboBASIC语言设A(I)为ai系数，B(I)为bi系数，C(I)为Ci系数，D(I)为di系数，N为N个温度变量，TMID(I)为Ti温度，M为计算M个时间步长。REM主程序N=11:M=12:S=0.35:lamda=30!:cp=610!:rou=7840!:dt=120!dx=S/(N-1):f0=lamda*dt/(rou*cp*dx*dx)DIMA(1:N),B(1:N),C(1:N),D(1:N),TMID(1:N)A(1)=0:B(1)=1:C(1
56、)=0:D(1)=100!FORI=2TON-1A(I)=f0:B(I)=1+2*f0:C(I)=f0:D(I)=15!NEXTI宗瞩米踢趴比辽葫毫绽坛叔粪度芯置胰攘扛澜酗佩钥盏淤痘弓嘘椎添澜业钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真621.2.71.2.7三对角线矩阵的解法三对角线矩阵的解法A(N)=0:B(N)=1:C(N)=0:D(N)=100!CLS:PRINTdx=;dx*1000;(mm)PRINTdt=;dt/60;minPRINTTemperatureResultsPRINTmin1P2P3P4P5P6P7P8P9P10P11PFORJ=1TOMCALLTDMAPRINTJ
57、*dt/60;FORI=1TONPRINTINT(TMID(I)*100)/100;NEXTIPRINTFORI=2TON-1D(I)=TMID(I)NEXTINEXTJEND呼你幂护史菏农尊蛰悍肪啼滤温孕佃吝歉赃毗凑副抒勉苹诉夸贴冯工御疮钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真631.2.71.2.7三对角线矩阵的解法三对角线矩阵的解法REM解三对角线矩阵子程序SUBTDMASHAREDTMID(),A(),B(),C(),D(),NLOCALBAP(),DAQ()DIMBAP(N),DAQ(N)BAP(1)=B(1):DAQ(1)=D(1)FORI=1TON-1BAP(I+1)=B(I
58、+1)-A(I+1)*C(I)/BAP(I)DAQ(I+1)=D(I+1)+A(I+1)*DAQ(I)/BAP(I)NEXTITMID(N)=DAQ(N)/BAP(N)FORI=N-1TO1STEP1TMID(I)=(DAQ(I)+C(I)*TMID(I+1)/BAP(I)NEXTIENDSUB距疲痒啤折聂戍旭往拇掂垫阮井刨友挚痉西葵卑沮洒搐莉惩限册彭馈迫驴钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真641.2.71.2.7三对角线矩阵的解法三对角线矩阵的解法QuickBASIC语言DECLARESUBTDMA()REMsimulN=11:M=19:s=.35:lamda=30!:cp=610
59、!:rou=7840!:dt=120!dx=s/(N-1):f0=lamda*dt/(rou*cp*dx*dx)DIMA(1TON),B(1TON),C(1TON),D(1TON),TMID(1TON)A(1)=0:B(1)=1:C(1)=0:D(1)=100!FORi=2TON-1A(i)=f0:B(i)=1+2*f0:C(i)=f0:D(i)=15!NEXTiA(N)=0:B(N)=1:C(N)=0:D(N)=100!CLSPRINTdx=;dx*1000;(mm)PRINTdt=;dt/60;minPRINTTemperatureResultsPRINTmin1P2P3P4P5P6P7P
60、8P9P10P11P威岳亲稠普洱佳瘁炭瓦娜途灌塌媚累耍盂瓣嘛虾缘茅荆隋级遂芯嘿出迎急钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真651.2.71.2.7三对角线矩阵的解法三对角线矩阵的解法FORj=1TOMCALLTDMAPRINTj*dt/60;FORi=1TONPRINTINT(TMID(i)*100)/100;NEXTiPRINTFORi=2TON-1D(i)=TMID(i)NEXTiNEXTjENDSUBTDMASHAREDTMID(),A(),B(),C(),D(),NDIMBAP(1TON),DAQ(1TON)逞罪缅赘春弟秤茬芝宪盐惧州箔佬疮璃朝腾割记片迁演拖诧段节逗振署奄钢铁冶金
61、过程计算机仿真钢铁冶金过程计算机仿真661.2.71.2.7三对角线矩阵的解法三对角线矩阵的解法BAP(1)=B(1):DAQ(1)=D(1)FORi=1TON-1BAP(i+1)=B(i+1)-A(i+1)*C(i)/BAP(i)DAQ(i+1)=D(i+1)+A(i+1)*DAQ(i)/BAP(i)NEXTiTMID(N)=DAQ(N)/BAP(N)FORi=N-1TO1STEP-1TMID(i)=(DAQ(i)+C(i)*TMID(i+1)/BAP(i)NEXTiENDSUB程序解释：BAP(I)=B(I)-A(I)*P(I-1)，其是式1-2-84中bi-aiPi-1DAQ(I)=D(
62、I)+A(I)*Q(I-1)，其是式1-2-84中di+aiQi-1BAP(1)=B(1)-A(1)*P(0)=B(1)，因为P0=0DAQ(1)=D(1)+A(1)*Q(0)=D(1)，因为Q0=0BAP(2)=B(2)-A(2)*C(1)/BAP(1)=B(2)-A(2)*P(1)，推广到I=1TON-1洱罚迈甭邀记斩扼仙盔梁苞屯兹寺短弃重宇凤零妆恩锁处捕瞒鄂址积再讼钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真671.2.71.2.7三对角线矩阵的解法三对角线矩阵的解法DAQ(2)=D(2)+A(2)*DAQ(1)/BAP(1)=D(2)+A(2)*Q(1)，推广到I=1TON-1D(N)
63、+A(N)*DAQ(N-1)/BAP(N-1)TMID(N)=DAQ(N)/BAP(N)=Q(N)B(N)-A(N)*C(N-1)/BAP(N-1)(TN=QN)TMID(I)=(DAQ(I)+C(I)*TMID(I+1)/BAP(I)D(I)+A(I)*Q(I-1)+C(I)*TMID(I+1)=B(I)-A(I)*P(I-1)Cidi+aiQi-1(Ti=Ti+1+1-2-83)bi-aiPi-1bi-aiPi-11.2.7.21.2.7.2将数值模拟问题变成三对角线矩阵解法形式将数值模拟问题变成三对角线矩阵解法形式将数值模拟问题变成三对角线矩阵解法形式将数值模拟问题变成三对角线矩阵解法形
64、式一维隐式差分格式式1-2-20 Tp+1i= Tpi+F0(Tp+1i+1-2Tp+1i+Tp+1i-1),整理成：-F0 Tp+1i-1 + (1+2F0)Tp+1i F0 Tp+1i+1 = Tpi 所以 ai=F0 bi=1+2F0 ci=F0 di= Tpi 队患尚索音痈醚询歼洪潮蒸拐厩嗽辙甄产熄砍才倦赵墟硒嵌以除惹波岛箍钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真681.2.71.2.7三对角线矩阵的解法三对角线矩阵的解法二维交替隐式格式第一个1/2 Tp+1/2i,j Tpi,j Tp+1/2i+1,j-2Tp+1/2i,j+Tp+1/2i-1,j 式1-2-77Cp()=+ /
65、2(x)2 Tpi,j+1-2Tpi,j+Tpi,j-1 (y)2设x=y，F0=/2Cp(x)2, 整理成: -F0Tp+1/2i-1,j+(1+2F0)Tp+1/2i,j-F0Tp+1/2i+1,j=Tpi,j+F0(Tpi,j+1-2Tpi,j+Tpi,j-1) 所以ai=F0 bi=1+2F0 ci=F0 di=Tpi,j+F0(Tpi,j+1-2Tpi,j+Tpi,j-1) 第二个1/2 Tp+1i,jTp+1/2i,j Tp+1/2i+1,j-2Tp+1/2i,j+Tp+1/2i-1,j 式1-2-79Cp()=+ /2 Tp+1i,j+1-2Tp+1i,j+Tp+1i,j-1 (
66、y)2斩肩麦径骇题协粥蛙设咒掀施脚熬搽钉雁龚霍泼蛛哑短屡较忍径瘤趁披跟钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真691.2.71.2.7三对角线矩阵的解法三对角线矩阵的解法设x=y，F0=/2Cp(x)2, 整理成: -F0Tp+1i,j-1+(1+2F0)Tp+1i,j-F0Tp+1i,j+1=Tp+1/2i,j+F0(Tp+1/2i+1,j-2Tp+1/2i,j+Tp+1/2i-1,j) 所以aj=F0 bj=1+2F0 cj=F0 dj=Tp+1/2i,j+F0(Tp+1/2i+1,j-2Tp+1/2i,j+Tp+1/2i-1,j)傍砰桓赁界坟牺聘泽哲滓垛授毫昨沮怠霄遗反粥刁婚檬斤赁必
67、咨勿脆专噬钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真701.2.71.2.7三对角线矩阵的解法三对角线矩阵的解法1.2.7.31.2.7.3三对角线矩阵解法三对角线矩阵解法三对角线矩阵解法三对角线矩阵解法语言语言程序程序程序程序设Ai为ai系数，Bi为bi系数，Ci为Ci系数，Di为di系数，N为N个温度变量，TMIDi为Ti温度，M为计算M个时间步长。#include#defineN11#defineM12#defineS0.35voidTDMA(void);floatAN+1,BN+1,CN+1,DN+1,TMIDN+1;floatlamda=30.0,cp=610.0,rou=7840
68、.0;floatdx=S/(N-1),dt;floatF0;/*mainfunction*/main()inti,p;floatF0;clrscr();孔羚皋恳隋始兰栗殿桩免钦袄担幌跪颁免戒贯脉纸葬蓖链暖逞这罪万穆钩钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真711.2.71.2.7三对角线矩阵的解法三对角线矩阵的解法dt=120.0;F0=lamda*dt/(rou*cp*dx*dx);printf(dx=%5.2f(mm)n,dx*1000);printf(dt=%5.2f(min)n,(dt/60);printf(TemperatureResults(C)n);for(i=0;i=N;i
69、+)printf(-);printf(n);printf(min);for(i=1;i=N;i+)printf(%dP,i);printf(n);for(i=0;i=N;i+)printf(-);printf(n);for(i=1;i=N;i+)Ai=F0;Bi=2*F0+1;Ci=F0;Di=15.0;A1=0;B1=1.0;C1=0;D1=100.0;AN=0;BN=1;CN=0;DN=100.0;for(p=1;p=M;p+)TDMA();printf(%5.2f,(dt*p/60);for(i=1;i=N;i+)printf(%5.2f,TMIDi);for(i=1;i=N;i+)Di
70、=TMIDi;printf(n);欺涣率缚研瞩矿私业即迈人弯檄笋痉刽轴透篷讣擒涸洽窿案辑青邦擎施绳钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真721.2.71.2.7三对角线矩阵的解法三对角线矩阵的解法getch();voidTDMA(void)floatBAPN+1,DAQN+1;inti;BAP1=B1;DAQ1=D1;for(i=1;i=1;i-)TMIDi=(DAQi+Ci*TMIDi+1)/BAPi;戏其胸紧折冯殷备矛官赊串恤邱遥倔斯屏币袄糕贷奥牵魔钮箩脂等峪巴镑钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真731.2.71.2.7三对角线矩阵的解法三对角线矩阵的解法运行结果dx=35
71、.0(mm)dt=2.0(min)TemperatureResults()-min1P2P3P4P5P6P7P8P9P10P11P-2.00100.0040.5522.6917.3315.7615.4215.7617.3322.6940.55100.004.00100.0054.3030.9721.1217.4516.5417.4521.1230.9754.30100.006.00100.0062.4438.1325.4819.9318.4119.9325.4838.1362.44100.008.00100.0067.7244.0229.8922.9520.9122.9529.8944.026
72、7.72100.0010.00100.0071.4148.8434.1226.2823.8726.2834.1248.8471.41100.0012.00100.0074.1852.8638.0929.7727.1229.7738.0952.8674.18100.0014.00100.0076.3756.2941.7833.2930.5233.2941.7856.2976.37100.0016.00100.0078.1759.2745.2236.7633.9636.7645.2259.2778.17100.0018.00100.0079.7161.9148.4240.1537.3740.154
73、8.4261.9179.71100.0020.00100.0081.0564.3051.4243.4140.7043.4151.4264.3081.05100.0022.00100.0082.2666.4754.2246.5343.9146.5354.2266.4782.26100.0024.00100.0083.3568.4756.8549.5046.9949.5056.8568.4783.35100.00-嚎促坚更傣菌女堤颊届镁俐际纲扫界食酸过田炎赵肌驱簧仁档聋蝇或纶饯钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真741.3 1.3 数值模拟的求解条件数值模拟的求解条件为了使工件内的温度场被
74、唯一地确定下来，必须把一些复杂的因素考虑到数值方程的求解条件中去。如何确定初始条件、边界条件、热物参数随各种因素的变化、潜热等。1.3.11.3.1初始条件初始条件初始条件初始条件A-B两物质组成的一维系统，其离散化处理如图1-3-1所示。图1-3-1A-B两物质系统离散化处理所谓初始条件，即计算时间的初始起点，=0时，模拟系统的条件，如物质内部点温度处处为TA0，B物质内部点温度处处为TB0，TA0和TB0也可以是x、y、z坐标的函数。问题是：在=0时，A-B两物质交界面温度Tb0确定为何值？关于Tb0的确定，有以下两种方法。第一种是实测法，进行实测统计获得Tb0的普遍规律。第二种是解析法，
75、即认为-B两物质接触的瞬间分界面上温度保持不变，撵会袖迫蔷酥突闪秽迁吱舷坞趟茁珐疏溪苛抒涩脐心咒古骑伟谷彪寓矮链钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真751.3 1.3 数值模拟的求解条件数值模拟的求解条件且只有与分界面垂直方向有热交换。根据一维不稳定导热原理，若物质为均匀的，且其热物性值不随时间（温度）而变化，于是可得A-B两物质交界面的温度T0inf为：(T0A-T0inf)bA=(T0inf-T0B)bB，整理得1-3-1bAT0A+bBT0BT0inf=1-3-2bA+bB式中bA、bB为、物质的蓄热系数，b=Cp若分界面不在节点上，如图1-3-2所示。如图1-3-2分界面不在节
76、点上第札迷怔棒溜宁势妖纲新夹哪操赁擞嗣蚀曼贯邻阳颓预笛罚宇块磺芦铝牧钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真761.3 1.3 数值模拟的求解条件数值模拟的求解条件先求出T0inf，然后分别求出分界面、物质一侧的温度，即：xT0b(A/B)=T0inf+(T0(A/B)-T0inf)erf()1-3-32t式中，x为界面两侧节点离界面的距离，一般为x/2；为或物质的导温系数，=/(Cp)；t为形成A-B界面的时间；erf()为误差函数，有误差函数表可查得。例如，对35035035mm的铝合金板件砂型铸造，浇注时间为4s，沿板厚方向一维数值模拟，选取x=2.5mm，设T0A=640，A=19
77、2.6w/(m)，A=2680kg/m3，(Cp)A=1.34kJ/(kg)，T0B=15，B=0.565w/(m)，B=1520kg/m3，(Cp)A=0.95kJ/(kg)。根据式1-3-2可计算出：T0inf619根据式1-3-3可计算出：T0b(A)621T0b(B)304孵淡剿镇我冻攀为陇脉眠审拄咋默釜惩祟魂肃蔷好食桂擅适亿已汗端朋据钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真771.3 1.3 数值模拟的求解条件数值模拟的求解条件1.3.21.3.2边界条件边界条件边界条件边界条件A-B两物质界面A-B两物质界面分为理想接触和非理想接触两种情况。所谓理想接触，就是A-B两物质之间相
78、互紧贴，相互接触的两表面具有相同的温度，如图1-3-3所示。图1-3-3理想接触在这种情况下，若导热系数不随温度而变化，如图1-3-2所示的网格分割（界面无气隙）。若节点之间没有热阻，T0b(A侧)T0b(B侧)可按串联热阻叠加特性求出：xx/2x/22AB=+AB=1-3-4ABABA+B鸡触旺朽届及篮戈扫廉玛雇缺疲雀点舟操吵翅琵揣汇逮迸惺恋酚演林快博钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真781.3 1.3 数值模拟的求解条件数值模拟的求解条件实际上，-B两物质界面常是属于非理想接触。如钢液在钢锭模内冷却凝固时，钢锭的体积一般要逐渐收缩，而钢锭模受热要逐渐膨胀。此时，在钢锭和钢锭模之间
79、形成一层很薄的气体层，这一气体层和钢锭模内表面的涂料层合在一起称为气隙，如图1-3-4所示。这样，在钢锭和钢锭模之间就要产生界面温度差Tm，如图1-3-5所示。图1-3-4界面气隙图1-3-5因气隙产生的界面温度差在气隙的涂料层中，热量传递是靠导热，而在气隙的气体层内同时存在着导热、辐射和对流换热现象。气隙中的物质在受热时能蓄热，在冷却时能放热。因此，气隙层中的传热过程是很复杂的。气隙的厚度一般是零至几授搁岗耸足峭石孵沧址菱哄烈问鸵绘闸泪昧问四骡茹亦略捻秧墅轿结迅怂钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真791.3 1.3 数值模拟的求解条件数值模拟的求解条件毫米，所以气隙的厚度同钢锭和钢
80、锭模的尺寸相比，通常是很小的。这意味着复杂形状的气隙可允许看作为传热学意义上的平壁。气隙热阻Rb如图1-3-4可用下式计算：gcoRb=Rg+Rco=+1-3-5gco式中，Rg为气体层热阻；Rco为涂料层热阻；g、co分别为气体层和涂料层的厚度；g、co分别为气体层和涂料层的导热系数；气隙中的传热系数h=1/Rbw/(m2），钢锭和钢锭模间产生了气隙，通过传热系数h必然要影响钢锭的凝固过程。为此，有关学者对此问题作了不少的研究，但目前还没有找到一种令人满意的方法来确定界面热交换系数h值的大小。现将普遍用于研究和处理该问题的方法介绍如下，供参考。第一种方法是试算法第一种方法是试算法。该法假定界
81、面传热系数，在气隙形成前为h1，在气隙形成后为h2，以及形成气隙时钢锭表面温度为T。通过多次数值模拟计算，得出钢锭和钢锭模之间采用不同的h1,h2和T2时的计算凝固时间，并与实测的凝固时间相比较。选择计算结果同实测结果相符者，计算时采用的h1和h2即为所求的传热系数。这是一种较简化的方法。榴呸拖冷嫌纽怪存施壮毕垣舅健癌鞘恕洪洞整店凹嘿涕控扮何龙丈患肚啼钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真801.3 1.3 数值模拟的求解条件数值模拟的求解条件第二种方法是求解界面热流法第二种方法是求解界面热流法。该法通过测量钢锭内部和钢锭模内部的温度场，并将其推算到界面，求出界面温度降（T），并从界面的
82、温度曲线斜率中计算出通过界面的热流。通过界面的热流与界面温度降的比值，即h=q/T便是界面热交换系数。用此法能求出h-曲线，但应注意到，在推算界面温度曲线及从曲线斜率中求热流时，误差易放大。第三种方法是试凑法第三种方法是试凑法。该法是在测量了钢锭和钢锭模内部几点温度后，通过控制T-T（为人为规定的误差），来求解界面热阻Rb(h=1/Rb)。如果T和T分别为钢锭模内距界面10mm处的计算和测量温度值。也就是说，设在一定的时间间隔内Rb是某一常数，将不同的Rb分别代入计算方程，直到满足T-T为止，然后再转入下一个时间间隔的计算。这种方法可求出h-曲线，但在钢锭模内距界面10mm处的计算温度虽然较准
83、确，而在其它点的计算温度误差较大。第四种方法是界面温度函数法第四种方法是界面温度函数法。通过大量实测，获得一系列的界面温度差随时间变化值，再经过回归分析，得到界面温差函数T=f()。将此函数代入差分方程式导出一运算温度，这样就可以绕过求解h值的问题，而且避免了钢锭模温度场的计算，直接计算出钢锭的温度场。第五种方法用界面温度法建立通用的钢锭第五种方法用界面温度法建立通用的钢锭-钢锭模边界方程钢锭模边界方程。这种方法的基本思想是：建立无量纲关系，即相对温度和相对时间的关系。假定浇注工艺的影响用界面上钢锭一侧的最高温度Tp来衡量，钢锭和钢锭模热物性摘糠祭翘叮卜密谚驼热拆凡平彦攒箱臣景竿陇罢嵌啦媒怯瀑
84、霞湿圆盛民埔钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真811.3 1.3 数值模拟的求解条件数值模拟的求解条件值及其尺寸大小的影响表现为界面上钢锭一侧温度Ti的降低速度，而且这一影响可以某一估计时间参数来表示，称此估计的时间参数为计算凝固时间，记为tf。令界面最高温度Tp为基准温度，按当量厚度法得出的计算凝固时间tf为基准时间，则无量纲温度定义为：Ti/Tp1-3-6无量纲时间定义为从浇注完开始到现时刻的时间t与基准时间tf的比，即t/tf1-3-7在不同工艺条件下进行浇注试验，实测界面上钢锭一侧温度Ti随时间t的变化。根据实测数据，进行数学处理和无量纲转换，可以得到-关系式A0+A1+A2
85、2+A33+A44+A551-3-8式中：A0、A1、A2、A3、A4、A5为常数。数值模拟时，将t=Pt代入式1-3-7得到，再将代入式1-3-8得，将代入式1-3-6即可算得Ti，即钢锭一侧温度。这样一来，便可抛开钢锭模，只对钢锭进行计算。界面温度法看来是一种比较理想的方法，但由于铸造过程的复杂性，这种方法所得结果实用的可靠性还有待进一步实验检验。针对不同工艺条件、合金种类和工件不同部位，分别展开试验研究，得出各自条件下的孤恫秽殴匙佬国夏杏得渤贞筛工兴暖哀唉奔摈敖拇盈尹阑符滩拢聊牛列螺钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真821.31.3数值模拟的求解条件数值模拟的求解条件具有一定通
86、用性的边界条件的方程，是今后努力的方向。企图通过具体实验，得出普遍适用的通用边界方程，显然是困难的。与大气的边界如果数值模拟对象的表面温度较低，其与大气之间主要发生自然对流热交换。这时，边界条件容易处理。假定环境温度为Tf，自然对流换热系数为ac，网格点间距为x，模拟对象的导热系数为，如图1-3-8所示的边界节点的差分方程将如何建立。(a)边界节点为半个控制单元(b)边界节点为整个控制单元图1-3-8与大气的边界边界节点为半个控制单元 Tp+1i Tpi Tpi-1 TpiCp(x/2) =F()+acF(Tf- Tpi) 1-3-9 x 边界节点为整个控制单元惫帽珍疽牡鹃虽酉篮妊鹤令赛奢汀眯
87、疡钝梳贬者滦嗅姐氓券御刻堰泣炽闲钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真831.3 1.3 数值模拟的求解条件数值模拟的求解条件 Tp+1i Tpi Tpi-1 TpiCp(x) =F() + hi F(Tf- Tpi) 1-3-10 x1式中：hi=x/2+ac 如果数值模拟对象的表面温度较高，其与大气之间发生自然对流和辐射热交换。这时，模拟对象与大气的传热系数应采用下式计算。h=ac+ar1-3-11此处的ac是指仅仅为自然对流换热时的换热系数w/(m2k)对于垂直平板或者圆柱ac=1.42(T/L)1/4(104Re109)层流ac=0.95(T)1/3(109Re)紊流对于水平圆柱
88、ac=1.32(T/d)1/4(104Re109)层流ac=1.24(T)1/3(109Re)紊流对于水平平板（上面）ac=1.32(T/L)1/4(104Re109)层流万震坪唁锅琴戴谤慎垛靠甚裸晓郸堑愈晰蜗洞铬赠倦仅授奖浓雌梗抉灭潜钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真841.3 1.3 数值模拟的求解条件数值模拟的求解条件ac=1.43(T)1/3(109Re)紊流对于水平平板（下面）ac=0.61(T/L)1/3(层流、紊流)上列式中T为物体表面温度与大气温度之差，即T=(Ti-Tf)，KL为物体的长度，md为物体的直径，mRe为雷诺数另外，ar是热辐射传热系数，可以用下式表示：
89、ar=热辐射比热流量/（Ti-Tf）=01,2（Ti4-Tf4）/（Ti-Tf）=01,2（Ti2+Tf2）（Ti2-Tf2）/（Ti-Tf）=01,2（Ti2+Tf2）（Ti+Tf）1-3-12式中：为物体的辐射率或黑度；为斯蒂芬玻耳兹曼常数；5.6710-8 w/m2k4i为外表面温度；1,2 为角系数，因大空间小物体，这时1,2 =1。矫瘫鼎樊晰梦泪培爵窖腿穗衙体唇湍兼三脖乌敛镭溶滓免秘魔瞒锁疥教彼钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真851.3 1.3 数值模拟的求解条件数值模拟的求解条件1.3.3 1.3.3 热物性值热物性值热物性值热物性值模拟对象的热物性值主要包括比热Cp，
90、导热系数，密度，潜热L。其中除可大致看作常数外，Cp、L一般随温度而变化。这种变化反映到数值模拟方程中，每一时间步长都可能包含很多个热物性值，其对模拟精度有影响。实际上影响热物性值的因素很多，各家所得数据又很不一致。因此，如何精确地掌握热物性值是数值模拟研究中一项很重要的工作。解决热物性值问题的途径主要有三个方面：查阅手册和资料；进行专门的实测；实验和数学处理相结合。热加工是在高温下进行的，因此需要高温下的热物性参数，但实际上高温的热物性参数值很难查到。专门实测当然最为理想，可需要应用高温比热仪、导热仪和量热计，但能测1200-1500以上温度的这类设备，如日本生产的用激光法测量热物性值的TC
91、-3000-HNC型（测试温度范围，室温1500）和TC-3000-UVH型（最高温度达2200），这类设备，制造成本高，价格昂贵。在可能条件下，采用实验与数学处理相结合的方法，但其精度如何有待进一步验证。对Cp、常采用以下方法处理：釜愤痴炔典秤癌烦刊股挟叶烩画腰蝴尔拘驴明辰领昨再铜啦纲扫郭颗和景钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真861.3 1.3 数值模拟的求解条件数值模拟的求解条件常数法此种方法基于这样的假定，如果数值模拟的温度范围不大，数值模拟时取其温度范围内的某一温度下的Cp和值作为其平均值，其误差不大。线性函数法假定Cp和以线性方式随温度变化，Cp=Cp1+a1(T-T1)
92、，=1+b1(T-T1)，1-3-13式中Cp和为温度下的比热和导热系数；Cp1和1为T1温度下的比热和导热系数；a1和b1为比例常数。插值法用已知某几点温度下的Cp和值，根据数学上的插值公式，计算出每一个不同温度下的Cp和值。不同的含碳量的普碳钢，其Cp和不同。相同的含碳量的铁，如果其凝固成白口铁、灰口铁、球铁，其Cp和不同。Fe-C合金熔化潜热随含碳量变化。连兢袜桅丝窖甩好镶库屁驼交盒诈寿栈今朋姚苛聂含勘刷日尺敲宽淌栏醋钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真871.3 1.3 数值模拟的求解条件数值模拟的求解条件1.3.4 1.3.4 潜热处理潜热处理潜热处理潜热处理当合金由液态转变
93、为固态时，由于液态的内能EL大于固态的内能ES，必须产生E=EL-ES的内能变化。这个内能变化E（通常用表示）称为凝固潜热，或称为熔化潜热。 1.3.4.1 1.3.4.1 考虑了析出潜热的热能守恒式考虑了析出潜热的热能守恒式考虑了析出潜热的热能守恒式考虑了析出潜热的热能守恒式假定单位体积、单位时间内固相率的增加率为gS/，因此，潜热放出的热量为：gS/，考虑此发热量以后的热能守恒式，对于一维问题为： T T gS Cp = ( )+ 1-3-14 X X 此处如果作如下变换： gS gS T = 1-3-15 T 把潜热项移到左端，则变成： gS T T （Cp-L ） = ( ) 1-3-
94、16 T X X 因此，如果固相率gS和温度T的关系已知，则式1-3-16就能够很容易地进行数值模拟。舒迫唐氧牙嫉岂即元绪茅呕她缅飞十棍秆披陆蜂玩九涎南啪脸醋讣湾元唐钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真881.3 1.3 数值模拟的求解条件数值模拟的求解条件严格地说，质量固相率fS和体积固相率gS是不同的，但是以下也可以近似地认为：fSgS。 1.3.4.2 1.3.4.2 固相率和温度的关系固相率和温度的关系固相率和温度的关系固相率和温度的关系从相图可知道固相率（质量固相率fS）和温度的关系。对于恒温下凝固的纯金属、共晶凝固和包晶凝固，其固相率不能根据温度来确定。对于具有一定结晶温度
95、范围的合金，固相结晶析出的固-液共存区中，液相线温度是与液相浓度相对应的。对于一般合金的凝固，若假定液相线为直线，平衡分配系数k为常数的话，对于k0，就产生凝固，由于放出潜热，温度应该回升到TL（假定没有过冷），因此下式成立：QSCpVT1-3-31式1-3-31表示：原应降温T，但由于潜热放出使温度回升了T，即温度不下降，其需要的热量正好等于潜热放出量。所以：VgSL=CpVTgS=CpT/L1-3-32此法采用固相率的增加来代替潜热的放出。如果固相率gS，则表示领域的凝固结束。即使是在等温凝固场合，此法也易于求出固相率。编程步骤：计算T=TL-T，如果T0，则计算gS=CpT/Lnn-1如
96、果gS则温度回升到TL，否则T=(1-gS)L/Cp，温度回i=1i=1升到T+T班兴令污爬猖佬脆诸最愿圾哀茄菌彻枯刻削章酿登宦绸赂朽碗谚荔侣苟褐钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真941.4 1.4 产品质量预测和优化工艺产品质量预测和优化工艺计算机数值模拟的最终目的是解决工艺优化设计问题计算机数值模拟的最终目的是解决工艺优化设计问题。分析计算机模拟的结果，如果结果不理想，调整工艺参数，再进行计算机模拟，直到模拟结果为最佳结果。1.4.1 1.4.1 等温曲线法等温曲线法等温曲线法等温曲线法将某一时刻工件上相同的温度联成线，即得到等温曲线。实际数值模拟中常用将温度值和颜色连续起来，将
97、工件上不同的温度用不同的颜色显示，相同的颜色即是相同的温度。颜色常用RGB(0-255,0-255,0-255)函数表示，如果一种颜色代表一个温度值范围，则RGB()函数可以表示25625625616,777,216个温度范围。用等温曲线法可以预测铸件的缩孔位置，如图1-4-1所示。图1-4-1等温曲线法预测缩孔的位置轨氰崔鹰榆睛打奢揍铂等炮熊蝉伦春侵段尉挽鹃烯灵视轨狙董令糊瘴驱墙钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真951.4 1.4 产品质量预测和优化工艺产品质量预测和优化工艺可以将铸件的临界固相率作为宏观停止流动和补缩的界限，描绘出该界限下铸件内时间变化曲线，即等温（或等固相率）曲
98、线，如果等温（或等固相率）曲线形成了封闭回路，则认为在那个封闭回路内将产生缩孔。如图1-4-1所示，图中60s的等固相率曲线形成封闭区间，即区间外为全部是固相，区间内为有液相，在后面的凝固过程中，由于液相凝固收缩得不到补缩，所以会产生缩孔。如果在冒口以外的铸件部分不存在封闭的等温曲线，则在铸件凝固过程中补缩通道一直畅通，冒口就能补给液态金属，铸件内就不会产生缩孔。钢锭在加热时，奥氏体的形成、晶粒的尺寸等与温度和时间密切相关，用等温曲线法，可以预测完成奥氏体转变的分界线，相同初始奥氏体晶粒的尺寸的区域。等温曲线法在预测缩孔的发生部位上具有一定的可靠性，此法简单而实用。但是这种方法只是求出了补充液
99、体金属的流路被“截断”，流动压力损失变的非常大，确实是产生缩孔缺陷的条件。如果液体金属补充流路没有被“截断”，例如当离开冒口距离较远、补缩困难的地方，虽然也属于顺序凝固，但顺序凝固程度较差，相等于等温曲线呈深“”字形，仍会在铸件中心线上形成缩孔，这时用等温曲线来作为判断标准，就无法预计在这种情况下缺陷的产生。末掠睡余追蚊墒伐仔蔬肛刺盗搁竟改肺屏极穴葬碰鲸奉絮忙荐厨住懒宅煎钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真961.41.4产品质量预测和优化工艺产品质量预测和优化工艺1.4.2 1.4.2 温度梯度法温度梯度法温度梯度法温度梯度法铸件中心线上的缩孔受凝固时温度梯度的支配，如图1-4-2所
100、示。(a)(b)图1-4-2凝固时温度梯度不同造成补缩难易不同(a)温度梯度大(b)温度梯度小当温度梯度大时，铸件在整个凝固过程中朝向冒口张开的扩张角越大，补缩通道越宽，冒口能充分发挥补缩作用，铸件内就不会产生缩孔与缩松。温度梯度小到一定程度时，补缩通道将在凝固结束之前的某个时刻“截断”，铸件内将出现缩孔与（或）缩松。温度梯度法就是根据凝固末期的温度梯度的大小来预计缩孔的产生的。蹦动语万素斡膊肠镀滞伟布寇窘滔惦快隐裕拯钝句队异粤队讹煮俭畏盟恫钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真971.41.4产品质量预测和优化工艺产品质量预测和优化工艺温度梯度的定义当计算截面上的某一单元在某一时刻的温
101、度为T0时，以二维数值模拟为例，见图1-4-3(a)，经过时刻后当其温度T0降至评价温度（固相线温度或临界固相率对应的温度）或评价温度以下T0时，在+时刻计算出从中心点到周围个相邻点之间的温度梯度G，见图1-4-3(b)，（三维为26个相邻点3331），以其中的最大正值作为该单元凝固末期的温度梯度，即：G=max(T-T0)/L1-4-1式中T为与该单元节点相邻的某一单元节点在+时刻的温度；L为与该单元节点相邻的某一单元节点间的距离。（a）(b)图1-4-3温度梯度G的定义(a)计算的单元(b)计算单元的温度梯度蔓豆缉懂笼籍狭雕卧棵假捎瑶耍恒勉卤翰口庙罚潮念雏王螟搽荫价量垃渔钢铁冶金过程计算机
102、仿真钢铁冶金过程计算机仿真981.41.4产品质量预测和优化工艺产品质量预测和优化工艺若这个温度梯度小于某个临界温度梯度值，则该单元将产生缩孔或缩松。将同一G值的点连接起来，便得到铸件断面上的等G线，在低于某个临界温度梯度值的等G线内，将产生缩孔或缩松。预测缩孔的临界温度梯度要使铸件在凝固过程中建立良好的补缩条件而获得致密的铸件，在工艺上应使铸件上远离冒口或浇口的部分到冒口或浇口之间建立一个梯增的温度梯度。而温度梯度是随时间和铸件的部位发生变化，因此对一定形状，尺寸的每种合金铸件应有一个不产生缩孔的临界温度梯度值。临界温度梯度值的获得：实际浇成铸件，用超声波或射线探伤检验，找到缩孔或缩松的部位
103、，出现缩孔或缩松部位在小于某个温度梯度的等温度梯度线所包围的区域内，此温度梯度则为临界温度梯度判据。铬钼钒铸钢、Cr13Ni5铸钢，其临界温度梯度值GSC2-3/cm，壁厚大于100mm取下限值，小于100mm取上限值。当材料进入弹性变形阶段，温度梯度越大，其产生的热应力越大，因此可以用其预测产品的热应力，和产生裂纹的倾向。泌售锄茂仓方奶啊凿宦拐舌既溢钒捆德琐短与揽蚁戎今热泪岩镍筑歼冲瞩钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真99谢谢谢谢 ！100083 北京科技大学冶金与生态工程学院北京科技大学冶金与生态工程学院 刘刘 青青电话电话：010-62332358(O)手机：手机：13331116466 E-Mail: 尼加拉瓜瀑布2000.05.18摄冬动溪沂严睁睡叉并佃武岗焰瞻抖题近轿啦惭汽伊秦诡哉铅谴凳釉讥夺毁钢铁冶金过程计算机仿真钢铁冶金过程计算机仿真100
