1、计算机在材料科学与工程中的应用计算机在材料科学与工程中的应用第四章第四章 材料科学与工过程中典型材料科学与工过程中典型 物理场物理场的的数值模拟数值模拟引言引言引言材料制备和加工主要内容：4.1 温度场的数学模型及求解 4.2 浓度场的数学模型及求解4.3 应力场的数学模型简介微分方程描述物理、化学、力学过程求解（定解条件）4.1 温度场的数学模型及求解4.1 温度场的数学模型及求解4.1 温度场的数学模型及求解加热设备加热设备冷却过程冷却过程加热方式加热方式电阻炉、燃料炉、浴炉、流化床炉、真空炉等电阻炉、燃料炉、浴炉、流化床炉、真空炉等 感应加热、电子束加热、激光表面处理、离子轰击感应加热、
2、电子束加热、激光表面处理、离子轰击加热等各种冷却介质的冷却性能和各种冷却方式高效加热等各种冷却介质的冷却性能和各种冷却方式高效率、节能的新加热方法。率、节能的新加热方法。各种冷却介质冷却，各种冷却方式。各种冷却介质冷却，各种冷却方式。材料学中传热学重要课题材料学中传热学重要课题2 24.1 温度场的数学模型及求解4.1.1 温度场的基本知识传热方式传热方式导热（热传导）导热（热传导）1 1）定义）定义 指温度不同的物体各部分或温度不同的两物体间直接接触时，依靠分子、原子及自由电子等微观粒子热运动而进行的热量传递现象。导热可以在固体、液体、气体中发生。2 2）导热的特点：）导热的特点： 必须有温
3、差；物体直接接触；依靠分子、原子及自由电子等微观粒子热运动而传递热量；不发生宏观的相对位移。3 3）导热机理）导热机理 气体气体：气体分子分子不规则热运动时 相互碰撞的结果。 导电固体导电固体：自由电子自由电子运动(类似气体分子)。 非导电固体非导电固体：晶格结构晶格结构的振动(弹性波)。 液体液体：很复杂。类似于气体；类似于非导电固体4.1 温度场的数学模型及求解4.1.1 温度场的基本知识导热基本定律导热基本定律Fourier(傅里叶傅里叶)定律定律单位时间内通过给定截面所传递的热量，正比例于垂直于该截面法线方向上的温度变化率，而热量传递的方向与温度升高的方向相反。qx热流密度 ；材料热传
4、导系数（热导率）W/(mK)；负号表示传热方向与温度梯度方向相反（P84）。l热导率是指具有单位温度差（1K）的单位厚度的物体(1m)，在它的单位面积上(1m2)、每单位时间(1s)的导热量(J)。热导率表示材料导热能力大小；物性参数。l不同材料的导热系数值不同，即使同一种材料导热系数值与温度等因素有关。金属材料最高，良导电体，也是良导热体，液体次之，气体最小。对于铁、碳钢和低合金钢，其值随温度的增加而下降；对于高合金钢（不锈钢、耐热钢），其值随温度的增加而增加。 x方向上的温度梯度4.1 温度场的数学模型及求解4.1.1 温度场的基本知识定义定义：流体中（气体或液体）温度不同的各部分之间，由
5、于发生相对的宏观运动宏观运动而把热量由一处传递到另一处的现象。伴随着热传导：伴随着热传导：热对流必然同时伴随着热传导伴随着热传导，自然界不存在单一的热对流。对流换热对流换热：流体与温度不同的固体壁间接触时的热量交换过程。C)(mW2 分类分类：有相变/无相变对流换热；自然对流/强制对流/沸腾换热及凝结换热传热方式传热方式热对流热对流1 1）定义与特征）定义与特征4.1 温度场的数学模型及求解4.1.1 温度场的基本知识)(fwtthq)(wftthq2 2）对流换热基本定律）对流换热基本定律牛顿冷却公式C)(mW2 h h 对流换热系数对流换热系数加热时：加热时：冷却时：冷却时：v h 的物理
6、意义： 对流换热系数，比例系数(表面传热系数)。 表面传热系数的大小与传热过程中的许多因素有关。它不仅取决于物体的物性、换热表面的形状、大小相对位置，而且与流体的流速有关。一般地，就介质而言：水的对流换热比空气强烈；就换热方式而言：有相变的强于无相变的；强制对流强于自然对流。对流换热研究的基本任务：用理论分析或实验的方法推出各种场合下表面换热导数的关系式。4.1 温度场的数学模型及求解4.1.1 温度场的基本知识1 1）定义）定义物体通过电磁波电磁波来传递能量的方式称为辐射辐射；因热的原因发出热辐射能热辐射能的现象称为热辐射热辐射，高温物体失去热量而低温物体得到热量，这种传热方式叫辐辐射传热射
7、传热。辐射热流量用斯蒂芬斯蒂芬- -波尔兹曼定律波尔兹曼定律表示。2 2）条件：条件：温差、发射电磁波。3 3）取决于：）取决于：两物体空间位置（角度系数）、表面特性（黑度）。 物体的温度越高、辐射能力越强；若物体的种类不同、表面状况不同，其辐射能力不同传热方式传热方式热辐射热辐射4.1 温度场的数学模型及求解4.1.1 温度场的基本知识4）特点）特点I.辐射换热是一个动态过程，处于热平衡时，辐射换热量为0.II.不需要冷热物体直接接触：即：不需要中间介质的存在，在真空中就可以传递能量，在真空中辐射能的传递最有效。III. 在辐射换热过程中伴随着能量形式的转换：物体热力学能电磁波能物体热力学能
8、。IV. 是一种双向热流同时存在的换热过程。无论温度高低，物体都在不停地相互发射电磁波能、相互辐射能量；高温物体辐射给低温物体的能量大于低温物体辐射给高温物体的能量;总的结果是热由高温传到低温。V.物体的辐射能力与其温度性质有关。区别于热传导、热对流的基本特点。 4.1 温度场的数学模型及求解4.1.1 温度场的基本知识24mW TEbb斯蒂芬-玻尔兹曼定律（Stefan-Boltzmann law） 黑体黑体：能全部吸收投射到其表面辐射能的物体。 或称绝对黑体。 （Black body）吸收率等于1的物体 黑体的辐射能力与吸收能力在同温度的物体中最强。1）黑体向外发射的辐射能： 黑体表面的绝
9、对温度（热力学温度） bET b)K(mW105.6742-8 K 绝对黑体辐射力 斯蒂芬-玻尔兹曼常数，2）实际物体辐射能力（低于同温度黑体）：24mW TEb 实际物体表面的发射率（黑度），0-1；与物体的种类、表面状况和温度有关。 4.1 温度场的数学模型及求解4.1.1 温度场的基本知识温度场温度场4.1 温度场的数学模型及求解4.1.2 温度场的数学模型导热微分方程的建立导热微分方程的建立建立导热微分方程，可以揭示连续温度场随空间坐标和时间变化的内在联系。建立导热微分方程，可以揭示连续温度场随空间坐标和时间变化的内在联系。理论基础：理论基础：傅里叶定律傅里叶定律 + + 能量守恒方程
10、能量守恒方程 温度场：温度场：Temperature fieldTemperature field。物质系统内各个点上温度的物质系统内各个点上温度的集合称为温度场，它是时间集合称为温度场，它是时间和空间坐标的函数和空间坐标的函数 。),(zyxft （导入微元体的总热流量）（导入微元体的总热流量）+ +（微元体内热源生成热）（微元体内热源生成热） = =（导出微元体的总热流量）（导出微元体的总热流量）+ +（微元体内能增加）（微元体内能增加）dxdyztQdxdzytQdydzxtQzyxdxdydzzttzQdxdzdyyttyQdydzdxxttxQdzzdyydxx)()()(dxdyd
11、ztcp微元体内能的增量dxdydzqQ.微元体内热源生成热微元体内热源生成热: 其中：其中：、c 、Q、微元体的密度、比热容、单位时间内单位体积内热源的生成微元体的密度、比热容、单位时间内单位体积内热源的生成热及时间。热及时间。 导入微元体的总热流量导入微元体的总热流量Q Q1 1=Q=Qx x+Q+Qy y+Q+Qz z ； 导出微元体的总热流量导出微元体的总热流量Q Q2 2=Q=Qx+dxx+dx+ Q+ Qy+dyy+dy+ Q+ Qz+dzz+dzcttttqxxyyzz笛卡尔坐标系中三维非稳态导热微分方程的一般表达式笛卡尔坐标系中三维非稳态导热微分方程的一般表达式4.1 温度场的
12、数学模型及求解4.1.2 温度场的数学模型4.1 温度场的数学模型及求解4.1.1 温度场的基本知识v笛卡尔坐标系中三维非稳态导热微分方程物理意义：反映了物体的温度随时间和空间的变化关系。能量方程是目前温度场数值模拟中普遍使用的描述方程，它不仅适用于固体，也适用于流体。其中，为材料的密度（kg/m3）；cp为材料的比热容（J/(kgK)）；为时间（s）；分别为材料沿x，y，z方向的热导率（W/(m*K)）；q为材料内部的热源密度（W/kg）。上式中，第一项为体元升温需要的热量；右侧第一、二和三项是由x，y和z方向流入体元的热量；最后一项体元内热源产生的热量。v微分方程的物理意义：微分方程的物理
13、意义：体元升温所需的热量应该等于流入体体元升温所需的热量应该等于流入体元的热量与体元内产生的热量的总和。元的热量与体元内产生的热量的总和。4.1 温度场的数学模型及求解4.1.2 温度场的数学模型v若若边界条件边界条件和和内部热源密度内部热源密度Q Q不随时间变化不随时间变化则经过一定时间后物则经过一定时间后物体内部各点的温度将达到平衡，即有稳态热传导方程：体内部各点的温度将达到平衡，即有稳态热传导方程：根据系统有无内热源、是否导热过程为稳态导热，以及一维、二维和三根据系统有无内热源、是否导热过程为稳态导热，以及一维、二维和三维的情况，可进行相应的简化。维的情况，可进行相应的简化。三维稳态热传
14、导方程为：三维稳态热传导方程为：4.1 温度场的数学模型及求解4.1.2 温度场的数学模型0tttqxxyyzz0tttqxxyyzzxttcqxx二维稳态热传导方程为：二维稳态热传导方程为：一维非稳态热传导方程为：一维非稳态热传导方程为：一维稳态热传导方程为：一维稳态热传导方程为：二维非稳态热传导方程为：二维非稳态热传导方程为：4.1 温度场的数学模型及求解4.1.2 温度场的数学模型0ttqxxyytttcqxxyy0 xtqxxPpCQztytxtct.222222)(讨论：讨论：1 1）直角坐标下有内热源的非稳态导热微分方程）直角坐标下有内热源的非稳态导热微分方程为常数时为常数时: :
15、 其中其中 称扩散系数（热扩散率称扩散系数（热扩散率）。）。 pca)(222222ztytxtctp0222222ztytxtQ2 2）直角坐标下）直角坐标下无内热源无内热源的稳态导热微分方程，的稳态导热微分方程，=0 =0 ，且，且 时时 3 3）常物性、稳态、无内热源，若）常物性、稳态、无内热源，若为常数时，且属稳态，即：为常数时，且属稳态，即：0t4.1 温度场的数学模型及求解4.1.2 温度场的数学模型4.1 温度场的数学模型及求解4.1.2 温度场的数学模型zzryrxsincos圆柱坐标下有内热源的非稳态导热微分方程：圆柱坐标下有内热源的非稳态导热微分方程：圆柱坐标：圆柱坐标：把
16、直角坐标的把直角坐标的xyxy平面变换为极坐标，而平面变换为极坐标，而z z轴不变轴不变。圆柱坐标与直角坐标之间的变换关系圆柱坐标与直角坐标之间的变换关系：4.1 温度场的数学模型及求解4.1.2 温度场的数学模型222222211ptttttqcrrrrzcPpCQtrtrrrtrct.2222222)sin1)(sinsin1)(1(cossinsincossinrzryrx球坐标下有内热源的非稳态导热微分方程：球坐标下有内热源的非稳态导热微分方程：球面坐标与直角坐标之间的变换关系球面坐标与直角坐标之间的变换关系：4.1 温度场的数学模型及求解4.1.2 温度场的数学模型初始条件和边界条件
17、初始条件和边界条件 1 1）初始条件）初始条件初始条件是指求解问题的初始条件是指求解问题的初始温度场初始温度场，也就是在零时刻温度场的，也就是在零时刻温度场的分布。它可以是均匀的，此时有分布。它可以是均匀的，此时有也可以是不均匀的，各点的温度值已知或者遵从某一函数关系。也可以是不均匀的，各点的温度值已知或者遵从某一函数关系。00TTt00( , , )tTT x y z 4.1 温度场的数学模型及求解4.1.2 温度场的数学模型2）边界条件边界条件是指物体表面或者边界与周围环境的热交换情况边界条件是指物体表面或者边界与周围环境的热交换情况，通常有三类重要的边界条件。（1）第一类边界条件 第一类
18、边界条件是指物体边界上的温度分布函数已知，表示为 或wsTTtzyxTTws,（2）第二类边界条件 第二类边界条件是指边界上的热流密度已知，表示为： 或 n为物体边界的外法线方向，并规定热流密度的方向与边界的外法线方向相同。（3）第三类边界条件 又称为对流边界条件，是指物体与其周围环境介质间的对流传热系数k和介质的温度已知，表示为：wssqnTqtzyxqnTqwss,fTTknT4.1 温度场的数学模型及求解4.1.2 温度场的数学模型Tf,h第一类边界qwTw第二类边界第三类边界边界上的温度分布函数已知边界上的温度分布函数已知边界上的热流密度边界上的热流密度对流传热系数和介质温度已知对流传
19、热系数和介质温度已知4.1 温度场的数学模型及求解4.1.2 温度场的数学模型4.1 温度场的数学模型及求解4.1.3 温度场的有限差分求解l导热问题的求解是对导热微分方程在已知边界条件和初始条件下积分求解，即解析解。但目前由于数学上的困难，在工程实际中许多问题还不能采用分析解法进行求解，如物体的几何形状比较复杂及边界形状不规则，材料的物性常数随温度变化等。近年来，随着计算机技术和计算技术的迅速发展，数值方法已经得到广泛应用并成为有力的辅助求解工具，已发展了许多用于工程问题求解的数值计算方法。数值积分法、有限差分法和有限元法离散数学为基础计算机作为工具4.1 温度场的数学模型及求解4.1.3
20、温度场的有限差分求解有限差分法 1）把原来物体内随时间和空间连续分布的温度问题转化为在时间和空间领域内有限个离散点的温度值问题，再用这些离散点上的温度值去逼近连续的温度分布。 2）有限差分是以差分代替微分，差商代替微商，建立以节点温度为未知量的线性代数方程组，然后求解得到各节点温度的近似值。yxdydx(i,j)(i+1,j)(i-1,j)(i,j+1)(i,j-1)4.1 温度场的数学模型及求解4.1.3 温度场的有限差分求解网格划分网格划分绝绝热热TwTf，hqwi,ji,j-1i-1,ji+1,jxy0L2L1i,j+1xxxii1yyyii102222yTxT网格划分网格划分二维稳态无
21、内热源导热二维稳态无内热源导热 教材教材P92P921 14.1 温度场的数学模型及求解4.1.3 温度场的有限差分求解4.1 温度场的数学模型及求解4.1.3 温度场的有限差分求解4.1 温度场的数学模型及求解4.1.3 温度场的有限差分求解4.1 温度场的数学模型及求解4.1.3 温度场的有限差分求解4.1 温度场的数学模型及求解4.1.3 温度场的有限差分求解 nininnnnnncccTTTaaaaaaaaa11212222111211CTA，其中：二维稳态导热二维稳态导热1差分方程建立差分方程建立 CTA构成矩阵形式构成矩阵形式4.1 温度场的数学模型及求解4.1.3 温度场的有限差
22、分求解二维稳态问题求解二维稳态问题求解1123A100E860B100C100D100H200G300F400 Eg1.Eg1. 边界点温度已知，且区域内无内热源。利用有限边界点温度已知，且区域内无内热源。利用有限差分方法来计算节点差分方法来计算节点1 1、2 2、3 3的温度。的温度。)(411,1, 1, 1,jijijijijiTTTTT0404043232312121TTTTTTTTTTTTTTTFEDGCHBA）节点）节点）节点 Eg.Eg. 边界点温度已知，且区域内无内热源。利用有限边界点温度已知，且区域内无内热源。利用有限差分方法来计算节点差分方法来计算节点1 1、2 2、3 3
23、的温度。的温度。1360400440043232121TTTTTTTCTCTCT400240160321解方程得解方程得4.1 温度场的数学模型及求解4.1.3 温度场的有限差分求解4.1 温度场的数学模型及求解4.1.3 温度场的有限差分求解实际工作中遇到的导热问题通常为非稳态导热，其特点是温度不仅随空间坐标的变化而变化，而且还随时间的变化而变化。因此温度场的分布与时间和位置两个因素有关。非稳态问题的求解原理、离散化方法和主要求解步骤与稳态问题的求解类似，但由于非稳态导热中增加了时问变量，因此，在差分格式、解的特性以及求解方法上都要复杂一些。如在区域离散化中，不仅包括空间区域的离散化，还有时
24、间区域的离散化。 一块无限大平板（如图一块无限大平板（如图4-11所示），其一半厚度为所示），其一半厚度为L=0.1m，初始温度，初始温度t0=1000，突然将其插入温度，突然将其插入温度t=20的流体介质中。平板的导热系的流体介质中。平板的导热系数数=34.89W/m，密度，密度=7800kg/m3，比热，比热c=0.712 J/kg，平板与，平板与介质的对流换热系数为介质的对流换热系数为h=233W/m2，求平板内各点的温度分布。，求平板内各点的温度分布。310无限大平板非稳态导热无限大平板非稳态导热 非稳态问题的有限差分格式非稳态问题的有限差分格式24.1 温度场的数学模型及求解4.1.
25、3 温度场的有限差分求解1）数学描述，由于平板换热关于中心线是对称的，仅对平板一半区域进行计算即可。坐标x的原点选在平板中心线上，因而一半区域的非稳态导热的数学描述为：2200,0,0,ttaxtttxxtxLh ttx该数学模型的解析解为:02coscossinsin210FnnnnnnneLxtttt20LaFniBctg/hLBi其中,为方程的根，。(4.34)(4.35)(4.36)(4.37)(4.38)4.1 温度场的数学模型及求解4.1.3 温度场的有限差分求解表4-1给出了在平板表面(x=L)处由上式计算得到的不同时刻的温度值。表4-1平板表面各不同时刻温度值时时 间（间（S）
26、12345678910温温 度度()981.84974.47968.88964.20960.11956.14953.08949.97947.07944.342）微分方程的离散 一维非稳态导热指的是空间坐标是一维的。若考虑时间坐标，则所谓的一维非稳态导热实际上是二维问题(见下图)，即：有时间坐标和空间坐标x两个变量。但要注意，时间坐标是单向的，就是说，前一时刻的状态会对后一时刻的状态有影响，但后一时刻的状态却影响不到前一时刻，如图示出了以x和为坐标的计算区域的离散，时间从=0开始，经过一个个时层增加到K时层和K+1时层。平面区域的时间和空间离散4.1 温度场的数学模型及求解4.1.3 温度场的有
27、限差分求解平面区域的时间和空间离散对于i节点，在K和K+1时刻可将微分方程写成下面式子：221122KKiiKKiittaxttax将上式的左端温度对时间的偏导数进行差分离散为：111KKKiiiKKKiiittttttKit1Kit这两个式子的右端差分式完全相同，但在两个式子中却有不同含义。对1式，右端项相对i点在K时刻的导数是向前差分向前差分。而在2式中，右端项的向后差分向后差分。 是i点在K+1时刻的导数(4.39)(4.40)(4.41)(4.42)4.1 温度场的数学模型及求解4.1.3 温度场的有限差分求解将式(4.41)和(4.42)分别代入式(4.39)和(4.40)，并将式(
28、4.39)和(4.40)右端关于x的二阶导数用相应的差分（二阶中心差分二阶中心差分）代替，则可得到下列显式和隐式显式和隐式两种不同的差分格式： 显式:KiKiKiKifTTffTT111)21 (K=0,1,2, , i=2,3,N-1)全隐式：KiKiKiKiTfTfTfT11111211(K=0,1,2, , i=2,3,N-1)以上两式中的2xaf 从式(4.43)可见，其右端只涉及K时刻的温度，当从K=0（即=0时刻）开始计算时，在，在K=0K=0时等号右端都是已知值，因而直接可计算出时等号右端都是已知值，因而直接可计算出K=1K=1时刻各点时刻各点的温度（显示）的温度（显示）。由K=
29、1时刻的各点的温度值，又可以直接利用式(4.43)计算K=2时刻的各点的温度，这样一个时层一个时层的往下推，各时层的温度都能用式(4.43)直接计算出来，不要求解代数方程组。而对于式(4.44) 等号右端包含了与等号左端同一时刻但不同节点的温度，因而必须通过求解包含了与等号左端同一时刻但不同节点的温度，因而必须通过求解代数方程组才能求得这些节点的温度值（隐式）。代数方程组才能求得这些节点的温度值（隐式）。(4.43)(4.44)4.1 温度场的数学模型及求解4.1.3 温度场的有限差分求解显式：对两时间层格式，知道n时刻各空间层上函数值而推的下一时刻的值。而隐式格式中，知道n时刻某空间层上函数
30、值，而不能得出待求的下一时刻各空间层上之值，必须联立求解一个与网格点数目相同的方程组。显式格式可以直接求解，缺点是计算稳定性差。隐式格式的优点是稳定性好，只是求解麻烦。4.1 温度场的数学模型及求解4.1.3 温度场的有限差分求解3)边界条件的离散 对于式(4.36)和(4.37)所给出的边界条件，可以直接用差分代替微分，也可以用元体平衡法给出相应的边界条件，亦有显式和隐式之分。通常，当当内部节点采用显式时，边界节点也用显式离散内部节点采用显式时，边界节点也用显式离散；内部节点用隐式时，边界节内部节点用隐式时，边界节点亦用隐式点亦用隐式。边界节点的差分格式是显示还是隐式，取决于如何与内部节点的
31、差分方程组合。用K+1时刻相应节点的差分，代替式(4.36)和(4.37)中的微分，可得到边界节点的差分方程：txhtxhtttKNKNKK111121111(4.45)最终的离散格式显式： 初始值： 0ttii=1,2,3,N) KiKikiKitfftftt21111(i=2,3,N-1) (4.46)(4.47)4.1 温度场的数学模型及求解4.1.3 温度场的有限差分求解txhtxhtttKNKNKK11112 1111其中K=0,1,2,(4.48)隐式：初始值： 00tti1211KKttKiKiKiKitftftft11111211(i=2,3,N-1) 11111KKNNh x
32、ttth x其中K=0,1,2, (4.49)(4.50)(4.51)(4.52)在用隐式差分计算时，每个时层都需要迭代求解代数方程组（4.49）-式（4.52）。在每个时层计算时，都要先假定一个温度场（一般取上一时层的温度场为本时层的初始场），然后迭代计算直至收敛。显式差分格式：每个节点方程可以独立求解，但需要考虑稳定性。)0 , 0( ,122Lxtxtninitxt1221()nnniiitttO得到：221122)()(2xOxtttxtnininini(4.53)2()(1121nininininitttxtt220)()(xcxFpnininininittxttt12111)(24
33、.1 温度场的数学模型及求解4.1.3 温度场的有限差分求解OxyvL2L1O)(3exp(2rrhQQm二维不稳态焊接热传导问题求解二维不稳态焊接热传导问题求解2 Eg2.Eg2.薄板焊接中移动热源薄板焊接中移动热源二维焊接离散化模型二维焊接离散化模型( (Qm-Qm-最大热源密度；最大热源密度；h h- -板厚；板厚；r-离开热源中心距离离开热源中心距离)0222222)0,()1(1)(TyxTKQyTxTTvyrr初始条件：简单传热学计算机分析简单传热学计算机分析环境温度初始温度；换热系数；换热绝热aTTaaaTTyTkLxLyTTyTkLxyTTxTkLyLxxTkLyx012221
34、1)(,0 ,)(,0 , 0)(,0 ,)(,0 , 0二维不稳态焊接热传导问题求解二维不稳态焊接热传导问题求解2 Eg2.Eg2.薄板焊接中移动热源薄板焊接中移动热源二维焊接离散化模型二维焊接离散化模型边界条件边界条件OxyvL2L1O简单传热学计算机分析简单传热学计算机分析二维不稳态焊接热传导问题求解二维不稳态焊接热传导问题求解2根据二维不稳定导热方程根据二维不稳定导热方程(1)以及初始条件、边界条件以及初始条件、边界条件,可建立差分可建立差分方程方程(2)kQyTTTxTTTTTnjinjinjinjinjinjinjinjinji,21,1,2,1,1,1,2221）（根据差分方程以
35、及各边界节点差分方程根据差分方程以及各边界节点差分方程,在满足稳定条件下就可求在满足稳定条件下就可求出不同时刻的温度分布。出不同时刻的温度分布。简单传热学计算机分析简单传热学计算机分析二维不稳态焊接热传导问题求解二维不稳态焊接热传导问题求解2时间步长时间步长DS=0.01s；热源作用时间热源作用时间S1=0sSe=5s=0.0008cal/cm2s；Qm=2500cal/cm；板厚板厚H=1cm；焊速焊速v=0.05cm/s节点数节点数N1=N2=10,T0=Ta=20,比热比热C=0.16cal/g;密度密度=7.82g/cm2;导热系数导热系数k=0.1cal/cms,节点间距节点间距DX
36、=0.2cm输入参数输入参数简单传热学计算机分析简单传热学计算机分析二维不稳态焊接热传导问题求解二维不稳态焊接热传导问题求解2计算框图计算框图START读参数判断差分稳定？按差分方程计算各时刻温度分布按时间为0.5sn的温度存入文件SS2？S=S+dSEND重设步长YNYN简单传热学计算机分析简单传热学计算机分析二维不稳态焊接热传导问题求解二维不稳态焊接热传导问题求解20.5s0.5s温度场温度场简单传热学计算机分析简单传热学计算机分析二维不稳态焊接热传导问题求解二维不稳态焊接热传导问题求解23.0s3.0s温度场温度场简单传热学计算机分析简单传热学计算机分析MatlabPDE二维不稳态焊接热
37、传导求解二维不稳态焊接热传导求解3区域设置区域设置单击单击工具，在窗口拉出一个矩形，双击矩形区域，工具，在窗口拉出一个矩形，双击矩形区域，在在ObjectDialog对话框输入对话框输入Left为为0，Bottom为为0，Width为为2,Height为为2。与默认的坐标相比，图形小的看不见，所以要调整与默认的坐标相比，图形小的看不见，所以要调整坐标显示比例。方法：选择坐标显示比例。方法：选择Options-AxesLimits,把把X，Y轴的自动选项打开。轴的自动选项打开。设置设置Options-Application为为HeatTransfer（设置程序应（设置程序应用热传输模型）用热传输
38、模型）简单传热学计算机分析简单传热学计算机分析MatlabPDE二维不稳态焊接热传导求解二维不稳态焊接热传导求解3边界条件设置边界条件设置单击单击使边界变红色，然后分别双击每段边界，打开使边界变红色，然后分别双击每段边界，打开BoundaryConditions对话框，设置边界条件（根据边界对话框，设置边界条件（根据边界条件）。所有的边界都为条件）。所有的边界都为Neumann条件。条件。简单传热学计算机分析简单传热学计算机分析MatlabPDE二维不稳态焊接热传导求解二维不稳态焊接热传导求解3方程类型设置方程类型设置单击单击，在，在PDESpecification对话框中设置对话框中设置方程
39、类型为方程类型为Parabolic（抛物型），（抛物型），rho(密度）为密度）为7.82，C（比热）为（比热）为0.16，k（导热系数）为（导热系数）为0.1，Q（热源）为（热源）为4000*exp(-3*(x.2+(y-0.4*t).2)/0.49)简单传热学计算机分析简单传热学计算机分析MatlabPDE二维不稳态焊接热传导求解二维不稳态焊接热传导求解3网格设置网格设置初值和误差设置初值和误差设置单击Solve菜单中Parameters选项，打开SolveParameters对话框，输入Time为0:0.5:5，u(t0)为20,其他不变。解方程解方程单击 ，或者加密网格，单击 。单击
40、。简单传热学计算机分析简单传热学计算机分析7)7)整理数据整理数据 单击单击Mesh-Export MeshMesh-Export Mesh输出输出 p e t p e t 的数值，单击的数值，单击Solve-ESolve-Export Solutionxport Solution输出输出u.u.回到回到MatlabMatlab主窗口执行下面两条命令：主窗口执行下面两条命令：u1=p,u(:,7) %u1=p,u(:,7) %将节点坐标和其在将节点坐标和其在3s3s时的温度组成新矩阵时的温度组成新矩阵u2=sortrows(u1,3) %u2=sortrows(u1,3) %将将u1u1按温度
41、值大小升序排列。按温度值大小升序排列。u1=p,u(:,4) %u1=p,u(:,4) %将节点坐标和其在将节点坐标和其在1.5s1.5s时的温度组成新矩阵时的温度组成新矩阵u2=sortrows(u1,3) %u2=sortrows(u1,3) %将将u1u1按温度值大小升序排列。按温度值大小升序排列。MatlabPDE二维不稳态焊接热传导求解二维不稳态焊接热传导求解31.5s时时温度场温度场简单传热学计算机分析简单传热学计算机分析MatlabPDE二维不稳态焊接热传导求解二维不稳态焊接热传导求解33s时时温度场温度场简单传热学计算机分析简单传热学计算机分析4.2 浓度场的数学模型及求解4.
42、2 浓度场的数学模型及求解v 任何不均质的材料，在定的热力学条件下，都将趋向于均匀化。 譬如，通过扩散退火可以改善因凝固带来的成分不均匀性这是在合金中分布不均匀的溶质原子从高浓度区域向低浓度区域运动(扩散)的结果。所以固态中的扩散本质是在扩散力(浓度、电场、应力场等的梯度)作用下，原于定向、宏观的迁移。这种迁移运动的结果是使系统的化学自由能下降。 材料的扩散现象在工程中广泛存在，如压力加工时的动态恢复再结晶，双全属板的生产、焊接过程，热处理中的相变，化学热处理以及粉末冶金的烧结等。 扩散理论的研究主要方面之是宏观规律的研究，它重点讨论扩散物质的浓度分布与时间的关系，根据不同条件建立一系列的扩散
43、方程，并技其边界条件求解。用计算机数值计算方法代替传统的、复杂的数学物理方程对浓度场问题进行研究已成为发展的趋势。4.2 浓度场的数学模型及求解上坡扩散上坡扩散l 上坡扩散是沿着浓度升高的方向进行扩散，即由低浓度向高浓度方向扩散，使浓度发生两级分化。l 例如，奥氏体向珠光体转变的过程中，碳原子由浓度较低的奥氏体向浓度较高的渗碳体扩散，就是上坡扩散。l 例如，将碳含量相近的碳钢（Wc=0.441%）与硅钢（Wc=0.478%）对焊在一起，在1050进行长时间(13d)加热后，焊接面处两侧的碳浓度分布如图。硅的存在推高了碳的活度和化学位，从而驱使碳从含硅的一侧向另一侧上坡扩散。 l 发生浓度低的向
44、浓度高的方向扩散，产生成分的偏聚而不是成分的均匀化。硅钢侧硅钢侧碳钢侧碳钢侧4.2 浓度场的数学模型及求解v 菲克第一定律菲克第一定律 稳定浓度场模型如图4-21和4-22所示 ，x轴上两单位面积1（A）和2（B），间距dx，面上原子浓度为CA、CB， 则平面1到平面2上原子数n1=C1dx，平面2到平面1上原子数n2=C2dx，若原子平均跳动频率f, d时间内跳离平面1的原子数为n1fd，跳离平面2的原子数为n2fd。4.2 浓度场的数学模型及求解v 菲克第一定律(一定时间内，浓度不随时间变化浓度不随时间变化dC/d=0)，单位时间内通过垂直于扩散方向的单位截面积的扩散物质流量（扩散通量）与
45、该面积处的浓度梯度成正比。定义：组分i每单位时间通过单位面积的质量传输量正比于浓度梯度。定义式： v D扩散系数。 负号表示质量传输的方向与浓度梯度的方向相反。J扩散通量，g/cm2*s 式中负号表明扩散通量的方向与浓度梯度方向相反。可见，只要存在浓度梯度，就会引起原子的扩散，物体的扩散系数单位：m2/s，物理意义：单位传质量相当于单位浓度梯度下的扩散传质通量。影响因素：物体的种类，物体的结构，温度、压力等。 CJDx 4.2 浓度场的数学模型及求解菲克第二定律菲克第二定律v 解决溶质浓度随时间变化的情况解决溶质浓度随时间变化的情况，即dc/dt0。v 如图：单元模型v 则单元体积中溶质积累速
46、率为：v 两个相距dx垂直x轴的平面组成的微体积，J1、J2为进入、流出两平面间的扩散通量，扩散中浓度变化为：C12CdxJJxxcDJ)(1 由由Fick第一定律：第一定律：4.2 浓度场的数学模型及求解v 即第二个面的扩散通量为第一个面注入的溶质与在这一段距离即第二个面的扩散通量为第一个面注入的溶质与在这一段距离内溶质浓度变化引起的扩散通量之和。若内溶质浓度变化引起的扩散通量之和。若D不随浓度变化，则不随浓度变化，则2122()CCCdxJJDDdxxxx 22()CCDx菲克第二定律在三维直角坐标系下的形式：菲克第二定律在三维直角坐标系下的形式：222222()CCCCDxzz4.2 浓
47、度场的数学模型及求解v 渗碳浓度场的有限差求解渗碳浓度场的有限差求解 渗碳是最常用的一种化学热处理工艺，工件渗碳后表面具有高硬度和耐磨性，心部具有良好的强韧性。工件渗碳后的碳浓度分布决定了工件的性能，求解菲克第二定律偏微分方程，能预测工件渗碳后的碳浓度分布。实际工件的形状是多种多样的，其表面形状可归为平面、凸柱面、凹柱面和球面四类。许多工作表明，工件的表面形状和曲率半径影响渗碳后的碳浓度分布，对不同形状的工件或表面，应采用不同的坐标求解。数学描述：工件内任一时刻任何位置的碳浓度变化规律由菲克第二定律偏微分方程描述： 222222()CCCCDxyz4.2 浓度场的数学模型及求解v 当工件表面为
48、平面，可将碳在工件内的扩散看成只在垂直于工件表面的一维方向上进行，忽略扩散系数D随碳浓度的变化，考虑到气体渗碳过程主要由碳向工件表面的传递工件表面的传递和碳在工件内部的扩散内部的扩散两部分组成，则工件内的碳浓度分布可归结为求下列定解问题，直角坐标系下的方程：2200 xgsxCCDxCDCCxCx4.2 浓度场的数学模型及求解v 半径为R，长径比较大的圆柱形工件，可认为碳在工件内只沿半径方向进行扩散，利用坐标变换，归结为求下列定解问题，极坐标系下的方程：2210mr Rgsr rCcCDrrrCDCCrCr4.2 浓度场的数学模型及求解v 半径为R的球形工件内的碳浓度分布可归结为求下列定解问题
49、，球坐标系下的方程：2220mr Rgsr rCCCDrrrCDCCrCr4.2 浓度场的数学模型及求解v 在一维扩散条件下，边界条件为：在一维扩散条件下，边界条件为：0SsCCCxCDJxCDCCx外层内层外层外层利用有限差分法求解时，一般分两步进行：利用有限差分法求解时，一般分两步进行：(1)将连续函数将连续函数C=f(x，)离散化离散化, 将将x-平面划分为如图平面划分为如图4-23所示的网格，所示的网格，图中图中x、分别代表距离步长和时间步长。两条平行线的交点称为节点分别代表距离步长和时间步长。两条平行线的交点称为节点，并以有限个节点上的函数值，并以有限个节点上的函数值C(xi，tn)
50、代替连续函数代替连续函数Cf(x，)。为简便起。为简便起见，将见，将C(xi，n)写为写为Ci,n,即表示在即表示在n时刻时刻xi处的浓度值；同理可用处的浓度值；同理可用Ci,n+1表示表示在在n+1(n+)时刻时刻xi+1处的浓度值处的浓度值 。(2)用差分代替微分用差分代替微分，对每个节点用差分代替微分对每个节点用差分代替微分，此时在，此时在xi处可作下列处可作下列代换：代换：1nniiCCC（4.78）4.2 浓度场的数学模型及求解v 其浓度对空间的二阶偏导数同样也用二阶中心差分替代：其浓度对空间的二阶偏导数同样也用二阶中心差分替代：111211112222212nnnnnniiiiii
