（材料学专业论文）三维非线性有限元在热处理炉cad中的应用研究.pdf

上海交通大学硕士论文中文摘要 三维非线性有限元在热处理炉c a d 中的应用研究 摘要 现有的热处理炉c a d 系统虽然极大提高了绘图的效率，但普遍存在着温度场计算薄弱的缺 点，因此并未能明显提高炉子结构设计的科学性。正确求解热处理炉内温度分布应是热处理炉 c a d 系统的核心部分。为此本文将前人的研究成果和模型进行了综合和加工，建立了热处理炉通 用三维瞬态非线性温度场模型。该模型全面考虑了辐射、热物性参数和换热系数非线性变化对 温度场的影响，并使用伽辽金加权余量法推导了其有限元分析格式。 为适应形状复杂、存在多种异质材料边界的高度非线性温度场求解问题，本文在有限元分 析时采取了一些特殊处理方法。( 1 ) 网格划分。采用六面体单元划分网格既能节约计算资源， 满足计算精度，并且适合于材料物性参数相差很大而导致的异质单元计算；( 2 ) 时间步长控制。 采用向后差分的时间积分法，对于线性问题是无条件稳定的，而对于强非线性问题为了避免结 果发散，采用自适应步长控制功能；( 3 ) 迭代方法。经过大量的试算和分析，采用修正的全牛 顿一拉夫森迭代法在收敛速度与计算精度上同时取得可以接受的效果。 r 根据三维瞬态非线性温度场模型，对7 2k w 气体氮化炉温度场进行了模拟计算，计算机模 拟结果与实测结果吻合良好。根据模拟结果，从工艺指标和节能指标两个方面对数台热处理炉 进行了炉顶、炉壁、炉底及整体优化设计，可以得到以下结论： 1 、在炉盖热损失中，由保温箱侧面的耐热钢板所引起的散热在整个炉盖热损失中所占比例 最大，所以在设计时应该在满足强度要求的前提下，炉盖和料筐应尽量选用导热系数小 和薄的钢板，这样可以明显地降低炉顶的热损失。另外采用一些办法减少或挡掉炉膛对 炉盖的直接辐射，也可以明显地降低表面温升和减少炉顶的热损失。本文以7 2k w 井式 气体氮化炉为例，对比了改进前后的温度分布与热损失，实测数据表明可以减少热损失 达6 8 6 。 2 、由炉项填充材料所导出的热量只占整个炉顶热损失很小一部分。因此通过增加保温材料 厚度并不能有效减少炉顶的热损失，反而会增大蓄热量和生产成本。 3 、炉壁设计时，在耐火、隔热保温和强度能满足要求的前提下，应尽量选用热导率小、密 度轻、热容小的新型耐火材料与保温材料，可以大大减少热损失和缩短升温时间。例如 箱式电阻炉采用全硅酸钙材料时，对比常用轻质砖炉壁可以减少炉壁热损失达4 0 以上。 4 、本文对比分析两种炉子结构( 上承重和下承重结构) 的热损失，发现下承重结构由于重 质砖的存在，使炉子底部存在热短路，从而导致大量热损失。从本文的炉底优化实例可 上海交通大学硕士论文中丈摘要 以看山，常址的将炉础f i l l f u 炉底j i n 入油封中的人州升式渗碳炉f 承重结构。由r 热短路 会导致人姑的热越损火，刈1 t 能和炉湍均匀性很不利。冈此一般热处理炉设计时，上承 重炉底结构麻优先选川， 5 、旧炉的优化政造。 1 住小l 。t 芝炉1 二原纳 = f = 的甚础上，住原炉罐外加隔热屏( 内剃掣度较 小的材料铝筘) ，i 口以人k 减，b 炉啦的辐射热损火，降低炉外鼙温升，是一种简单易行、 低成本的高效1 y 能疗；：， 6 、要使炉内温度分布均匀泓i ? 小必彻考虑i u 加热体在炉膛内的分布位置。对丁人删升式 炉应实行分阿控湍 全文的创新点如i 、： 1 ) 建立了适合热处理炉的迎j ：维懈态雌埘生濡度场模艰该模玳全面考虑了辐射、热物胜参 数平换热系数1 r 线性变化对l & 温j 场| r , j 影响。 2 ) 通过研究解决丁王维瞬态雌e i t ! y t i ；m 无分析时单元选择的问题，在保证计算精度的前提了， 使计算姑人人降低：刘r i | t i l i 同利 | = ：i 物性参数棚筹根人而导致的异质单元问题世进行r 研 究井得剑解决，f 越j h j 刈计张中嫂嗽心 3 ) 将建立的热处理敖寻：獭，m 静川i u 处州炉他计优化，对热处理炉虚拟生产进行了探讨，亓 为解决热处理炉设计_ 热。i 踏、挖湍苫世汁难题捉供一些有毓的参节：) 之 疋，则a l 面净损失的辐射热为q l ( 一) ，a 2 面净得的辐射热为q 2 ( + ) ，其中： q 。( 一) = u - 一g 。h ，= 渺 占l 爿l 上海交通大学硕士论文第一章绪论 q 2 ( + ) - ( g 2 - j 2 蚺高萨 6 2 a 2 因为换热仅发生在a l 与爿2 之间，所以 q 卜2 = q l ( 一) = q 2 ( + ) 联立( 2 1 2 ) 、( 2 一1 3 ) 、( 2 - 1 5 ) 可以解出 q 。一2 = f 了e o t l - e o 百2 l j 一1 。一 s l 爿l4 if 1 2占2 4 2 对于两个无限大平行放置的灰体平面，因为a i = a 2 = a ，9 2 = q 2 l = 1 ，则有 。 对于两个同心圆柱彳l 和爿2 ( a l a 2 ) ，将式( 1 1 0 ) 代入，可得出 组一鱼亟：蔓k 墨! 监：+ 上 占l占1 ( 1 一1 2 ) ( 卜1 4 ) ( 1 _ 1 5 ) ( 卜1 6 ) 其中f 2 = 竺等。这是一个很重要的公式，在井式炉传热分析计算中，我们碰到的辐射换热基本 。讲 上都属于这种情况。 1 3 3 热处理炉三维瞬态非线性温度场的数值模拟 对于温度场的研究已有很多，有限差分法与有限元法都是十分成熟的，甚至边界元法也有 了一定的研究与应用。关于这些方法的基本原理和计算步骤及注意事项等等问题不再一一叙述。 这里仅仅叙述在几何形状，边界条件复杂的热处理炉温度场的研究现状及其难点。 热处理炉温度场模拟必须将三维问题，瞬态问题、热物性参数( 导热系数、密度、比热容) 与温度的关系、强迫对流和辐射边界条件等因素全盘考虑在内，而这方面的研究工作，以前由 于计算机硬件水平不高，计算方法还不完善，还很少，大部分研究工作只是考虑了这诸多因素 中的某些项，如文献 4 5 】使用有限差分法求解了井式渗碳炉的二维温度场模拟结果与实测结果 在形状比较复杂的炉子顶部和底部相差较大，这就是由于没有采用三维模型，不能正确反映热 1 9 坠。 坷一。一龟鱼q = 红4 上海交通大学硕士论文 第一章绪论 处理炉实际工作状态的一个例子。 温度场数值模拟中，非线性问题是一个难点，解决方法是否合适将会影响解的精度、收敛 速度和稳定性等。热处理炉温度场计算就是一个典型的非线性问题，其非线性的来源为以下几 个方面： 1 ) 表面综合换热系数中应包括与温度有关的辐射和对流两项，辐射和对流均是温度的高次 函数，在有限元法中我们称此类问题为高度非线性问题。 2 ) 筑炉材料的热物性参数( 导热系数、密度、比热容) 也与温度有关，随温度的变化而变 化。 通常非线性方程的解法有三类：迭代法、增量法及混合法。迭代法的主要优点是简单易行， 但缺点是不能保证收敛到精确解。增量法将载荷以微小增量形式逐渐施加，用分段线性化的办 法来近似非线性问题，其主要优点是收敛性好和适应性强，除个别情况外，它适应于各种类型 和各种程度的非线性问题。它主要缺点是很难事先估计载荷增量取多大，因此靠方法本身无法 判断其偏离精确解的近似程度，此外，它的估算量很大。混合法是增量形式的迭代法，综合了 迭代法和增量法的优点，计算精度更高，而且可以判断每一增量步终了时解的近似程度。综合 法成为当前非线性方程解法的主流，使用广泛的一些最新方法有牛顿- 拉斐逊法( f n r ) 、修正 牛顿拉斐逊法( m n r ) 以及准牛顿拉斐逊法( q nr ) 。 在热处理炉的温度场模拟计算中，由于物体温度变化很大，热物性参数也变化较大，而热 交换系数随表面温度的变化更大，这些都给非线性方程的解带来困难，如计算中采用简单迭代 法，则很容易出现发散现象，一般采用综合法中的最新的几种方法。随着大型商品化有限元软 件的出现，计算方法的成熟以及计算机软硬件水平的提高，为解决热处理炉的三维瞬态非线性 温度场成为可能。 1 4 本课题研究的主要内容 热处理炉是实现热处理工艺的重要设备。工业生产中要求热处理炉具有能保证产品质量、 保证生产效率、操作方便、安全可靠、节约能源等特性。因此如何利用发展的c a d 技术、热处 理数学模型和计算机模拟技术实现正确选择炉型、合理设计炉子结构，不断改进炉子热工性能 以实现热处理高效、节能，己成为现代热处理炉设计制造的重要研究课题之一。但该课题有待 于进一步深入研究，因为目前热处理c a d 大都是在机械c a d 的基础上开发的，虽然其绘图功 能强大，但模拟计算这一部分很薄弱。首先，大都局限于非常简单的二维形状，仍末涉足三维 问题，不能真正反映生产实际情况；其次，采用的模型所考虑的因素，有的不太科学，有的则 过分简化，需要进一步完善；最后，模拟计算所需材料性能参数，换热系数与表面温度的关系 等一系列数据严重缺乏与零乱，有待于进一步整理与研究。 本课题的主要任务是： 上海交通大学硕士论文 第一章绪论 l 、通过对热处理炉内的传热分析，建立热处理炉通用三维瞬态非线性温度场有限元分析模 型，该模型中综合考虑了辐射、强迫对流以及材料热物性参数等非线性因素对热处理炉内温度 场的影响： 2 、热处理炉三维瞬态非线性温度场的有限元求解是一个存在多种异质材料、高度非线性问 题，本文探讨了如何通过选择合适的单元、时间步长、迭代方法来满足其特有的有限元分析要 求。 3 、根据建立的热处理炉三维瞬态非线性温度场有限元分析的数学模型，在大型非线性有限 元分析软件m a r c 上对实际热处理炉温度场进行有限元模拟计算，并通过大量的实验测试，对数 值模拟的结果进行了验证、分析与讨论。 4 、根据有限元数值模拟的结果，进行热处理炉整体与局部传热分析。从满足工艺指标和节 能指标两个方面对现有的热处理炉提出优化措施，模拟与实测结果均证明优化措施的有效性， 为解决热处理炉设计中热短路、温度均匀区窄等设计难题提供一些有益的参考。 以上均是国家九五攻关项目“热处理c a d 系统及其应用”的一部分。本文工作的完成，说 明了利用计算机模拟技术进行热处理炉c a d ，不但可行，而且有效为热处理炉的虚拟生产与优 化改造工作打下了坚实的基础。 上海交通大学硕士论文第二章热处理炉传热分析及其温度场的有限元实现 第二章热处理炉传热分析及其温度场的有限元实现 2 1 热处理炉内热量传递的三种方式及其数学描述 研究热处理炉传热的基本任务是解决如何把电( 或其它发热体) 产生的热量有效地传递给 工件和如何减少炉子热损失的问题，因此有必要搞清楚在热处理炉内热量究竟以如何方式进行 传递。热量传递有三种基本的方式：传导、对流和辐射。热处理炉内尽管热传递很复杂，也不 外乎是两种或三种传热形式同时存在的综合。 1 ) 传导传热 热传导是在没有相对运动的介质中，由于温度梯度的存在，引起介质内部的能量传递。通 过某一截面的热流密度用著名的傅立叶定律，即： 口：一七旦l( 2 1 ) ah 式中：k 为介质导热系数，疗为截面外法线方向，t 为温度。 由傅立叶方程，可推出系统的能量守恒为： 詈- v ( k v r ) + q 沼：) 式中：p 为介质的密度，c 为比热，q 为单位容积的热产生率。 该方程就是导热的微分方程，它表达了温度丁和时间f 以及空间坐标等的关系。在三维直角 坐标下的表述方程为： 肛詈= 丢( k 篆 + 号 b 等 + 虱af ，k ：a 瑟t ) c z 在热处理炉中，炉衬各部分之间的热量传递主要靠导热方式来进行。 2 ) 对流换热 对流是指流体各部分发生相对位移冷热流体相互掺混所引起的热量传递方式。实际工作中 遇到的不是单纯的对流方式而是流体流过另一物体的表面时对流和导热联合作用的热量传递过 程，我们称之为对流换热。对流换热的热流密度可以用牛顿冷却定律计算，即： q=bt)(2-4) 式中：h 为对流换热系数，0 为流体温度，瓦为固体的边界温度。 在热处理炉温度场计算时，本文考虑了以下两种对流条件：其一为炉内的强迫对流，其二 上海交通大学硕士论文第二章热处理炉传热分析厦其温度场的有限元实现 为炉体与空气的自然对流。 3 ) 辐射传热 当物体温度高于绝对温度时，就会以电磁波的形式向外辐射能量。这种由于物体受热的作 用而发射的辐射能称之为热辐射，因热辐射而发生的热量传递称为辐射传热。在高温炉内，加 热元件通过辐射向工件传递的热量远远大于对流，而在真空炉内，所有的热交换都是通过辐射 进行的。辐射传热的计算相对传导和对流来说要更复杂。真实表面发射的热流密度用下式计算： p = 眦( 2 5 ) 式中：瓦为表面的绝对温度( k ) ，盯为斯蒂芬一玻耳兹曼常数b = 5 6 7 1 0 - 8 b 2 k 4 ， s 为表面发射率，其值在0 到1 之间。 在两个或多个表面之间的辐射换热过程中，辐射换热不仅取决于参加辐射换热的表面温度、 辐射性质( 如发射率、吸收率、反射率和透射率) 和辐射的方向和波长特性，而且还取决于表 面的几何形状和它们之间的空间相对位置( 空间角系数) ，甚至还牵涉到表面间介质的吸收、发 射和散射等性质。辐射换热的能量守恒方程是积分一微分方程，求解比较复杂。 为了使上述控制方程唯一地被确定下来，还需要相应的具体问题的初始条件和边界条件， 即必须附加边界与初始条件和式( 2 - 3 ) 联立求解。 温度场边界条件共分三类，图3 一l 是它们的示意图。 ( 1 ) 第一类边界条件是指物体某些边界上温度函数为己知，即： b = b 0 ，j ，= ，) ( 2 6 ) 式中，s 1 为物体第一类边界条件的表面面积，珏为已知壁面温度，它可以随位置与时间的变化 而变化。 ( 2 ) 第二类边界条件是指物体某些边界上s 2 热流密度g 己知，公式描述为： 七。譬l 品= 一绋( x ，y ，印) ( 2 - 7 ) d h 。 式中，热流密度方向是边界外法线n 的方向，一般认为热流量从物体内部流出，因此在计算时， 对第二类边界条件凡热流由物体向外流出者q s 取正，而流入者取负。第二类边界条件又称为自 然边界条件。 上海交通大学硕士论文第二章热处理炉传热分析及其温度场的有限元实现 节：英边晃条 q + 已知 界条件 以及 为巳 图2 - 1 热传导问题的三类边界条件 f i g 2 1t h r e ek i n do fb o u n d a r yc o n d i t i o n so fh e a tt r a n s f e rp r o b l e m ( 3 ) 第三类边界条件是指物体某些边界辐射与对流换热条件为已知，它其实是第二类边界 条件的特例。 对于对流： =亿一t)(2-8) 式中h 是对流系数，t 。是与固体产生对流流体参照面温度，例如气体流过固体边界的附面层温 度。 对于辐射： q s 3 = 咖。一巧) ( 2 - 9 ) 式中g = e l 2 ，是由两相互辐射物体黑度系数乘积，f 是与两辐射物体形状有关的平均角系数， 称为形状因子，盯。为斯蒂芬一玻耳兹曼常数( 盯o = 5 6 7 x 1 0 - 8 w m 二k 4 ) ，t f 辐射源的温度。 初始条件是指待求的非稳态传热问题在初始时刻整个区域内各点温度值是己知的，它即可 以是常量，也可以是空间坐标的函数。 7 1 ，。= t o ( x ，y ，z ) ( 2 1 0 ) 求解热处理炉温度场分布问题实际就是求解在边界条件( 2 - 7 ) 、( 2 - 8 ) 及( 2 - 9 ) 和初始 条件( 2 - 1 0 ) 下满足热传导方程( 2 - 3 ) 的场函数t ( x ，y ，z ，t ) 。 2 4 上海交通大学硕士论文 第二章热处理炉传热分析厦其温度场的有限元实现 2 2 三维瞬态温度场的有限元格式 从热传导问题的微分方程和边界条件导出其有限元求解方法有两类。一类是通过寻找热传 导方程和边界条件对应的能量泛函，使泛函求极值来推导有限元方程：另一类是通过使热传导 方程和边界条件取加权残差为零来近似导出。本文在温度场的有限元分析中，用伽辽金加权余 量法来推导三维瞬态温度场的有限元计算方程，并用伽辽金差分格式进行时间上的离散。 2 2 1 伽辽金加权余量法建立有限元方程 伽辽金加权余量法有限元分析非稳态温度场的基本步骤如下： l 、结构的离散化 将区域v 划分为e 个有限单元，每个单元有p 个结点。由于传导控制方程对定义域v 内每 一点都有效，因此对v 内任意子域和任意点都是有效的这样就可以考虑任意一个单元来确定 局部近似解。 2 、选择适当的插值函数 在每个单元内，对温度t 假发一个适当的插值形式，并把单元e 中的温度表示为： p t ex ，y ，z ，f ) = m b ，y ，：) 寸( 2 - 1 1 ) t = l n i 称为形函数，即温度的插值函数。显然单元内任意点的温度是用各结点的温度值加权而 成的，单元内任意点的坐标位置也可以用各个结点的坐标来表示，当温度插值函数与坐标位置 插值函数采用完全相同的形式时，这种有限元单元称这等参元，这是一种常用的离精度计算的 单元类型。 3 、单元矩阵和载荷向量的推导 利用伽辽金方法，取权函数等于形函数，有： 州如警 + 孙等 * 型a ：儿7 矿+ 腓- 蚜甜肌。 对上式的第一个积分式进行分步积分和格林变换可得 ( 2 一1 2 ) 上海交通大学硕士论文第二幸热处理炉传热分析厦其温度场的有限元实现 ，m 弘譬) d v = - v ，飘以和矿+ 抄( t 孙峦 协m ，f ，m 号( b 譬 d v = - v t ，警( 和y + 抄( b 势嬲 川 ，”孙譬 d p = - v t ，掣o zl 吒矿+ 卅：孙搬 于是方程( 2 - 1 2 ) 变为： 一心，警譬+ b 等警地警譬卜+ 外；譬”q 警矿t 譬斗诎州q i + q 2 - - 掣卜= 。 由于边界s 是由s l ，s 2 ，黾组成的，因为s i 上温度为已知值，对x ，y ，z 的导数应为零，所以表面 积分为零。在s ，s ，上应满足： h k 警”b 警矿t 譬斗批 一9 吨署吗一9 一z p 一乙b = 。 于是式( 2 - 1 6 ) 可以表达为下式： 一此，警警屿警警址警警卜 p 鼍啦一黔限乙) 吗+ 川蚴一譬卜= 。 ( 2 - 1 8 ) 上式是由伽辽金方法导出的热传导方程和边界条件的最后积分方程。在分析热传导时可以 根据具体情况分几类情况，分别由上式导出具体的控制方程。如线性稳态、非线性稳态、线性 稳态及非线性瞬态等。 把方程( 2 1 7 ) 表示为下述矩阵形式： 上海交通大学硕士论文第二章热处理炉传热分析及其温度场的有限元实现 式中 k ；弦( e ) ) + k ；弦。) ) + k ；弦。) 一 p o ) _ 0 蚓2 1 t i p 4 n 1 砂 p = 奴。， - 够， + 缸曲 陋】= 蚓 a n la n 2 a n p 叙缸叙 o n 】o n 2o n p 砂砂勿 a n l8 n 3矾p a za za z 【】= ( 。2 p ) 4 、集合单元方程得总的平衡方程组： ( 2 1 9 ) ( 2 2 0 ) ( 2 2 1 ) ( 2 2 2 ) ( 2 2 3 ) ( 2 2 4 ) ( 2 - 2 5 ) ( 2 - 2 6 ) ( 2 2 7 ) ( 2 2 8 ) ( 2 - 2 9 ) 砂 讼 p 防 九 口 k 肚妒 = i i k k p 咒 嬲 绞 盯 + 舻驴肌 饵 似 值 k o o 。l = -_ip 上海交通大学硕士论文第二章热处理炉传热分析及其温度场的有限元实现 【k 。弦 + k ：弦 + k 3 l i t _ p ( 2 3 0 ) 在热传导问题上，无论用变分法求泛函还是用加权余量法都可以推出上面的有限元总的平 衡方程组。 5 、引入s 上规定的边界条件与初始条件后，求解方程组( 2 - 3 0 ) 。 2 2 2 瞬态非线性温度场的时间积分 。 在求解瞬态温度场时，需要进行时间上的离散，一般用差分法将微分方程转化为线性方程 组。 口( 詈 ，小目旧。= 古亿以。) 差分方式有多种，按照上式中口( 0 s 目1 ) 来分，有 口= 1 2 时，得到c r a n k n i c o l s o n 格式 土2 f t , 曼o 三t ) ，+ ( 詈 。 = 古亿一i 一。) 口= 2 3 时，得到伽辽金( g a l e r k i n ) 格式 詈( 詈 ，+ ；( 乳。= 古亿也。) 此外，还有三点差分格式中的三点向后差格式 ( 乳= 击陬。褂u ( 2 3 1 ) ( 2 3 2 ) ( 2 3 3 ) ( 2 3 4 ) ( 2 3 5 ) ( 2 3 6 ) 式 一， 屯 一 黻 纯 黻 一 黼加黼 脚址粉知 脯 m 骺 一 删罚俐割，l 上海交通大学硕士论文第二幸热处理炉传热分析及其温度场的有限元实现 向前差分格式能够得到显式解，但是稳定性差，实用中要求选择的时间步长很小，目前已 很少采用。研究表明，向后差分格式、c r a n k n i c o l s o n 格式、g a l e r k i n 格式都是无条件稳定的， 而向前差分格式则是条件稳定的。在定步长计算中，当时间步长较大时，c r a n k n i c o l s o n 产生 的振荡强烈，三点差分格式也很强烈。g a l e r k i n 格式的振荡要轻得多，只有向后差分格式才能 无条件不振荡，只是精度稍差。 采用不同的时间差分格式，将得到不同的温度场有限元的集合方程。一般计算中均采用精 度较高而且无条件稳定的g a l e r k i n 格式。下面为不同差分格式对应的瞬态温度场的有限元基本 方程。 向后差分格式： ( k - 】+ k ：】+ 吉k ，】) 口) ，= 石1 k ，弦) ，。+ p ， c z s t ， c r a n k n i c o l s o n 格式 ( k - 】+ 区z 】+ 云医，】 f ) ，= ( 古略卜医- 】一医：】 丁 。，+ ( p ，+ p ，。 ) g a l e r k i n 格式 ( 2 - 3 8 ) ( z k 。卜z k ：】+ 丢医，】) f ) ，= ( 言k ，】一【足。卜医：】) 留) ，一。+ ( 2 p ，+ p 。，) ) 三点差分格式 ( 2 3 9 ) ( k - 】+ k z 】+ 互石3 【岛】) 留 。= 古k ，k f ) 。一p ) m 。，) + p ) ， c z a 。， 图2 2 为采用g a l e r k i n 格式进行瞬态温度场有限元分析的程序流程图。 2 3 有限元分析中几个特殊问题的处理 通过各国学者的共同努力，热处理的数学模型考虑的因素更完备，研究的范围更广泛，所 采用的数值方法也由过去的有限差分法逐步向有限元法，从二维平面向三维实体过渡，目前己 基本采用三维有限元法以适应生产实际情况、复杂的边界和高度非线性情况。但是其数值模拟 还有很多特殊的问题尚在摸索阶段或考虑的较少，以下分别进行讨论。 圭堡奎望查兰塑主丝查 苎三主垫竺墨生堡垫坌堑垦苎墨生望堕查! ! 垄壅墨 图2 - 2 g a l c r k i n 格式瞬态温度场计算的程序框图 f i g - 2 - 2p r o g r a mf r a m e o ft r a n s i e n t t e m p e r a t u r ef i e l dc a l c u l a t i o nu s i n gg a l e r k i nm e t h o d 上海交通大学项士论文 第二幸热处理炉传热分析厦其温度场的有限元实现 2 3 1 网格的划分 为了使模拟对象更接近实际以扩大其工业应用，需解决以下问题：1 ) 工件形状的复杂性： 2 ) 高度非线性。目前，随着大型商业有限元软件的普及，三维瞬态非线性温度场的计算机模拟 采用有限元法是必然发展趋势。网格划分是有限元计算的基础。由于热处理炉三维瞬态非线性 温度场的求解是一强非线性问题，再加上形体复杂( 三维) ，网格划分的优劣直接影响到计算结 果的精确度以及计算量。对于三维瞬态非线性温度场温度场有限元分析，一般用四面体单元或 六面体单元划分网格。 本文选取一正方体耐火材料( 0 1 m 0 1 m 0 i m ，k = o 1 6 5 + 0 o 0 0 1 9 t ，黑度取0 7 5 ， c 。= 0 9 1 9 6 + 0 0 0 0 2 5 t ) 进行单元选择试验。首先将此正方体均匀划分为6 8 5 个四面体单元( 图 2 - 3 ) 。施加固定温度、对流、辐射边界计算正方体的温度场，结果发现靠近棱边的单元中结点 温度与积分点温度相差较大。而结点温度反映了实体真实温度场的分布，若积分点温度与结点 温度相差过大，必将导致计算结果失真。分析发现导致四面体网格结点温度与积分点温度相差 过大的原因有：1 、四结点四面体单元只有一个积分点( 图2 - 5 ) ；2 、网格划分较大。由于本文 三维瞬态非线性温度场模拟中温度变化非常大，表面处四面体单元积分点不在表面上，导致积 分点温度与表面处结点温度差距很大。原则上网格划分越细，结点值与积分点值就越接近。但 一味采取过分细小的网格势必导致计算量急剧增大，受到计算资源的限制。 为了解决上述矛盾。用八结点六面体等参元来划分。八结点六面体等参元共有八个积分点， 和结点位置比较接近( 图2 6 ) 。同样用上述正方体来试验，划分为6 4 个八结点六面体单元( 图 2 4 ) ，边界条件施加也和前面一致，图2 7 与2 8 为计算对比结果。对比显示，不仅网格数较 四面体划分大大减少，计算时间相对缩短，结果却较为理想，结点温度与积分点温度相差无几， 满足计算要求。另外热处理炉所用耐火材料一般均是长方体形，用八结点六面体来划分单元可 以避免产生畸变单元。因此在三维瞬态非线性温度场的有限元分析中选择八结点六面体单元即 能节约计算资源，又能满足计算要求，适合于本文的有限元分析。 图2 - 3 四面体单元划分网格 f i g2 - 3m e s h i n gb yt e t r a h e d r a le l e m e n t 上海交通大学硕士论文第二章热处理炉传热分析厦其温度场的有限元实现 图2 - 4 六面体单元划分网格 f i g2 - 4m e s h i n gb yh e x a h e d r a le l e m e n t 图2 - 5 四面体单元示意图( + 为积分点) 2 f i g2 - 5s c h e m a t i cd i a g r a mo ft e t r a h e d r a le l e m e n t ( + i si n t e g r a t i o np o i n t ) 图2 - 6 六面体单元示意图( + 为积分点) f i g 2 - 6s c h e m a t i cd i a g r a mo fh e x a h e d r a le l e m e n t ( + i si n t e g r a t i o np o i n t ) 上海交通大学项士论又第二章热处理炉传热分析及其温度场的有限元实现 3 3 上海交通大学硕士论文第二章热处理炉传热分析度其温度场的有限元实现 2 3 2 时间步长选择 非线性分析中时间步长的选择关系到计算结果的收敛性和精度，以及分析所花费的代价( 主 要是指计算时间) 。在处理非线性问题时，可选择固定时间步长和自适应时间步长。当接近稳态 时，时间步长将增加；反之如果温度变化剧烈或边界条件高度非线性( 如强迫的对流、辐射) ， 时间步长将减少。 本文的温度场分析中采用向后差分的时间积分法，对于线性问题是无条件稳定的，对于强 非线性问题为了获得较得精确的结果，采用自适应步长控制功能。 2 3 3 非线性迭代方法的选择 按照几何、材料、边界非线性理论建立的最后有限元方程为非线性的，求解这些非线性方 程是有限元分析的重要部分。高精度、高效率和高度数值稳定性的非线性求解方法和快速的代 数方程组求解技术对于获得合理、准确的分析结果至关重要。与线性分析不同，在非线性分析 中很难找到一种适合各类非线性及各种非线性程度的解法。在非线性有限元的发展过程中出现 过种种解法，它们都有各自的适用范围，用得不好，解可能会发散。因此对这些解法应深刻了 解，并在使用中积累经验，以正确使用各种方法。适合于本文研究的三维瞬态非线性温度场的 迭代方法有以下两种： 全牛顿一拉夫森迭代法 需要在每次迭代中重新形成和分解切线刚度矩阵，具有两阶收敛性。收敛快，适用于高度 非线性问题，对几乎所有非线性问题都可得到良好的分析结果，但是计算时间需要较多。 修正的全牛棋拉夫森迭代法 每次迭代都采用增量步开始时的切线刚度矩阵，迭代过程中不需要重新形成刚度矩阵，一 个增量步只形成一次总刚度矩阵，因而每次迭代需要的计算时间少。但收敛速度比全牛顿一拉夫 森方法慢，这种迭代方法特别适用大规模中等非线性问题的求解。 本文经过大量的试算和分析认为：对于本文研究的三维瞬态非线性温度场的有限元分析， 全牛顿一拉夫森迭代法与修正的全牛顿拉夫森迭代法二者均可达到收敛，但全牛顿一拉夫森迭代 法的计算时间是修正的全牛顿一拉夫森迭代法的i 5 倍左右，而收敛速度与修正的全牛顿一拉夫 森迭代法相当，因此对于本文研究的三维瞬态非线性温度场的有限元分析我们采用修正的全牛 顿一拉夫森迭代法比较合适。 2 3 本章小结 热处理数学模型是开发高度知识密集型的热处理智能技术的关键，已日益为各国热处理界 所重视。目前，对热处理炉温度场数值模拟计算往往只限于二维模型，筑炉材料的热物性参数 和复杂的边界条件也大都简化处理，从而导致计算结果与实测数据误差较大。本节针对这一问 题，做了以下工作： 上海交通大学硕士论丈 第二幸热处理炉传热分析及其温度场的有限元实现 i 、住对热处理炉进行详尽f 热分忻的基础上，建立丁一个描述热处理炉的二维1 r 线性瞬 态传热数学模型，全面考虑了热物性参数和l 复杂边界条什( 强迫对流与辐射) 随温度的1 i 线性 变化对温度场的影响，升实现了，该1 r 线性方删的有限元解法。 2 、讨论了热处理炉= 维懈态1 r 线性有限元分析中的儿个特殊问题： 网格的选取问题。即对t 二维实体网格的划分，本文比较了选择四面体单元与八面体 单元对模拟计算结果的影响，发现本文所f i j i ：究的二维瞬态1 f 线性温度场的有限元分析中，川八 面体单元划分网格即能1 ，约i 1 + i f 7 ：铙源义能满足汁锋璎求： 时间步k 的选抒问题，本文的温度场分忻中采_ l i j 向斤彳莘分的时问积分法，对丁姐1 r 线 性问题为了获得较精确的结粜心采k _ i 适府步k 控制功能： 1 r 线性有限元法迭代疗法的选取题。本文经过分析利人彗的试算认为：对丁本文研 究的二维瞬态1 r 线性温度场的仃限元分析，全牛顿一拉大森迭代法与修止的全牛顿一拉大森迭代 法_ 二青均可达刨收敛，但企二| 二l i l n 一拉欠森迭代法的计算时问是修止的全牛顿一拉大森迭代法的lj 倍左右，而收敛述度与修止的仑牛蝴一拉大森迭代法相当，冈此对r 本文研究的二维瞬态1 r 线 性温度场的有限元分忻找们采川他j l 的全牛倾一拉大森迭代法比较合适。 ；激“糍嚣。+ | 叠i 、氮蕊藤描 上海交通大学硕士论文第三章7 2k w 井式气体氮化炉的温度场模拟 第三章7 2k w 井式气体氮化炉的温度场模拟 3 1 研究对象及实验方法 3 1 1 研究对象简介 本节以上海交通大学设计并与上海铁道大学工厂联合研制的7 2k w 井式气体氮化炉作为具 体算例，模拟计算其三维瞬态非线性温度场。该炉主要由炉壳、炉衬、加热器、马弗罐、导流 简、炉盖、炉盖升降机构、电气控制系统等组成。该炉的外形尺寸与材料分布图中给出了炉衬 各层的材料分布及尺寸，如图3 - 1 ( a 、b 、c 三点为实验测温点) ，表3 - 1 给出了各炉衬材料物 性参数。计算时认为炉内壁温度保持在工况温度t = 5 5 0 ；环境温度巧= 2 4 。炉壳外表 面与环境主要以自然对流与辐射进行热交换，炉内则为强迫对流与辐射传热，在模型中通过给 出表面对流系数、强迫对流系数和辐射函数对此进行了模拟。 3 1 2 实验方法 1 、测温点的位置确定 在实验中我们根据模拟计算结果在垂直方向选取三个具有代表性的不同高度( 如图4 - i 中 a 、b 、c 三点所示) ，在每个高度上我们由表及里测定不同r 位置的温度。 2 、测温方法 1 热处理炉表面温度用x m x 型表面温度计( 测量精度：0 5 0 1 ，测量范围：0 6 0 0 ) 进行测量： 2 测量炉衬内部温度时，在选定的高度用钻头在炉外壁上打一4 小孔，再将陶瓷套管一端 制成尖锐的锯齿状，由表及里逐渐将陶瓷套管打入炉衬内部，注意打入时使陶瓷套管尽量 与表面法线方向一致，打入过程中不断测量小孔的打入深度，并将钻下的耐火材料倒出。 当达到预定的测温位置时，将瓷管取出，然后在管内插入直径为巾3 的铠装k 型热电偶， 将炉子的半径减去热电偶伸入陶瓷套管的长度作为测温点径向位置。热电偶通过补尝导线 与电位计相连，用水银温度计测量补偿线冷端的温度，查表得到的对应的m v 值v l 。根据 u j 3 6 型电位差计所测得的m v 值v 2 加上v 1 查表即得到所测点的温度值。 3 实验注意事项： 测温过程中除了要注意热电偶的正确使用方法及一些常规注意事项外，还要注意我们模拟 的是炉子工况下的温度场分布，所以炉子应有定长的保温时间。本实验中在5 5 0 保温1 0 小 时后进行测量，炉内壁电热元件表面温度为5 8 0 c ，补偿线端的温度为2 4 c 。 第三章7 2k w 井式气体氯化炉的温度场模拟 图3 - 1 简化炉体与测温点( a 、b 、c ) 示意图 f i g 3 1s i m p l i f i e ds c h e m a o f7 2k w n i t r i d i n gf u r n a c e 表3 - i 筑炉材料的导热系数( k _ a + b t ) t a b l e3 - 1h e a tc o n d u c t i v i t yo ff u r n a c e c o n s t r u c t i n gm a t e r i a l s 材料名称标号a ( w ，m )b ( w m 2 ) 高铝砖+ 硅酸钙( 1 ：2 ) m 1o 5 9o0 0 1 8 6 硅酸钙m 2 0 0 5 60 0 0 0 1 保温材料 m 30 0 40 0 0 0 2 2 o n 0 6 砖 m 4o 1 6 50 0 0 0 1 9 耐火纤维 m 5o 0 20 0 0 0 1 2 钢 m 65 20 0 0 0 0 2 8 密封圈 m 7o 1 6 2 8o 3 7 上海交通大学硕士论文第三圭7 2k w 井式气体氮化炉的温度场模拟 表3 - 2 炉内对流换热系数表 t a b l e3 - 2c o n v e c t i v eh e a tt r a n s f e rc o e f f i c i e n t so fi n s i d es u r f a c ew i t hr e s p e c tt oi n t e r n a l g a sf l o w ( w ，m 2 ) 炉内气体流速炉内温度t ( ) ( m s )7 6 06 5 05 4 0 4 2 5 61 2 61 3 4 1 4 41 5 8 1 21 9 52 1 22 3 _ 32 5 9 1 82 6 32 8 83 2 03 6 1 2 43 3 13 7 14 0 74 2 3 表3 - 3 炉外表面的综合换热系数表 t a b l e3 - 3c o n v e c t i v eh e a tt r a n s f e rc o e f f i c i e n t so fe x t e r n a ls u r f a c e ( w m 2 ) 炉外壁温度 垂直壁面 水平壁面炉外壁温度 垂直壁面 水平壁面 ( ) 顶面底面 ( ) 顶面底面 2 59 01 07 66 01 2 21 4 09 9 3 51 0 21 1 68 48 01 3 41 5 2l o 8 4 51 0 8 1 2 38 81 0 01 4 71 6 71 1 9 3 2 数值模拟与分析 3 2 1 研究问题的简化 考虑到实际计算的可能性，作者在模拟计算时忽略诸如肋板、电机等对传热贡献不大的一 些次要因素，集中解决主要矛盾。对于炉内壁搁丝砖在硅酸钙上的炉衬部分，根据上海交通大 学凌均谷硕士论文，采用一种简化模型来处理，即认为硅酸钙层由两层均质材料组成：内层为 混合层，其导热系数取硅酸钙和高铝砖的加权平均导热系数( 根据体积百分比加权) ，外层为单 一的硅酸钙材料。 因为网格划分( 文中最小的网格密度取1 0m ) 的限制，设计中小于1 0m 的材料层，如炉 壳的钢板( 本为6 8n l l n ) ，均取为1 0m 。注意到虽然炉体直径由实际的1 7 4 0 改为1 7 6 0m ， 炉内壁的直径由1 1 0 0m 改为1 1 2 0m 来进行处理，但炉衬的总厚度没有改变，再对照炉子的 其它部分，基本与原来一致，所以作者认为在计算时作了相当的简化而使求解的问题有一定的 偏差，但偏差是很小的。根据炉体温度分布的实际情况，在炉盖与炉体结构过渡位置几个温度 急变处采用细分网格，以提高计算精度。图3 2 3 4 是温度场计算所使用的1 4 炉盖、炉底及 炉体中部有限元网格模型。 进一步考虑到计算机及计算软件的运算速度与内存的限制，在实际计算中，将整个炉分为 三部分分块进行计算。实验结果以及理论分析都表明，在大型井式炉的中间部分，等温线应该 是垂直分布的( 所谓井式炉的均温区) ，即只有半径方向而没有垂直方向的热量传递，所以本文 将炉项部分的下边界和炉底部分的上边界设为绝热边界，首先计算炉体顶部下边界温度值后传 递给炉体中部计算程序，最后计算炉体底部温度场分布。 3 8 上海交通大学硕士论文第三幸7 2k w 井式气体氮化炉的温度场模拟 注意：为了观看方便，对炉顶、炉底及炉体中部的三维有限元网格均调整了一定的角度。 图3 - 2炉盖1 4 有限元网格模型 f i g 3 - 2f i n i t ee l e m e n t m e s ho fo n eq u a r t e ro nf u r n a c et o p 图3 - 3 炉体中部1 4 有限元网格模型 f i g 3 3f i n i t ee l e m e n tm e s ho fo n e q u a r t e ro nf u r n a c em i d d l e 上海交通大学硕士论丈 第三章7 2k w 井式气体氛化炉的温度场模枞 囝3 - 4炉成1 7j 有限元网格模型 f i g3 4f i n i t ee l e m e n tm e s ho fo n eq u a r t e ro nf u r n a c eb o t t o m 3 2 2 数值模拟结果与分析 利用m a r c 7i 犬i 线。| + “一畦几坝”术蝌铺三蒂建立的三维 f 线陛瞬态传热数学模，m 全 部t 作在p i i5 0 0 ，内仃12 8 m 诎叭f 逍玎乐用1 1 9 汁算参数如表：j - 4 所不。 是j 一4炉项；盐度场模拟参数 t a b l e3 - 4p a r a m e t e r so fs i m u l a t i o no nt e m p e r a t u r eo ft o pf u r n a c e j 多野 n 圯毁i 点数单元类抖 炉体部分 炉顺部 0 17 01 2 5 6 4 八1 ，点八面体轴对称单元 炉中部 4 3 8 f 16 2 7 6 八1 ，点八面体轴对称单元 炉底部 0 5 7 3 4 5 2 八1 7 点八面体轴对称单元 表3 - 5 3 7 为轴向何。岢a 、叭( 的删1 试7 粜模拟汁算结果，模拟计算与实验测试数据的 比较幽可参址i 芏| 3 - 5 ：3 _ 7 。模拟炉体再。邗分j 维实体及一维剖面模拟温度场分布如图3 - 8 3 一】0 、 ( 为丁观石方便，对炉j | j ! 、炉- 0 ，乇炉休中 j 维有限元模拟结果均调整了一定n 勺角度。 表? 一j 铀向位置 模拟结果与实测数值 t a b l e3 - 5n u m e r i c a ls i m u l a t i o na n c le x p e r i m e n t a lr e s u l t sa tp o s i t i o nai nf i g1 狰同f ：= 7 ：置( m m ) 5 2 05 7 06 2 06 7 07 2 07 7 08 2 0 襄面 实洲结果( 。c ) 3 8 12 0 52 7 315 51 3 70 85 03 9 模拟结果( ) 4 133 ：52 9 91 7 41 4 51 l35 74 5 上海交通大学硕士论文第三章7 2k w 井式气体氮化炉的温度场模拟 表3 - 6 轴向位置c 模拟结果与实测数值 t a b l e3 - 6n u m e r i c a ls i m u l a t i o na n de x p e r i m e n t a lr e s u l t sa tp o s i t i o nci nf i g1 径向位置( m ) 6 3 06 7 07 2 07 7 08 2