ANSYS Workbench有限元分析实例详解（热学和优化）

版权信息

书名：ANSYS Workbench有限元分析实例详解（热学和优化）

ISBN：978-7-115-64355-1

本书由人民邮电出版社发行数字版。版权所有，侵权必究。

您购买的人民邮电出版社电子书仅供您个人使用，未经授权，不得以任何方式复制和传播本书内容。

我们愿意相信读者具有这样的良知和觉悟，与我们共同保护知识产权。

如果购买者有侵权行为，我们可能对该用户实施包括但不限于关闭该帐号等维权措施，并可能追究法律责任。

版　　权

著    周 炬 苏金英

责任编辑 胡俊英

人民邮电出版社出版发行　　北京市丰台区成寿寺路11号

邮编　100164 　电子邮件　315@ptpress.com.cn

网址　http://www.ptpress.com.cn

读者服务热线：(010)81055410

反盗版热线：(010)81055315

内容提要

本书以联系和对比的方式系统且全面地讲解了ANSYS Workbench热学分析及优化设计过程中的各种问题，从工程实例出发，侧重解决热学计算问题和优化设计流程与工程问题。

本书内容分为4章：第1章介绍热学有限元分析的基本概念；第2章详细介绍ANSYS Workbench软件中Mechanical模块的热分析，包含热传导、热对流、热辐射、瞬态和非线性热分析、扩散、热固耦合分析；第3章详细介绍ANSYS Workbench软件中Fluent模块的热分析，包含热传导、移动热源、对流和辐射、多孔介质的蒸发、凝华结霜、热分析标量方程、流固耦合传热分析；第4章详细介绍ANSYS Workbench软件中的优化设计过程，重点讲解拓扑优化、增材制造和试验优化设计的基本原理及分析过程。

本书内容丰富新颖，融入了现代工程教育思想，旨在帮助读者提升运用工程方法解决实际问题的能力。本书主要面向ANSYS Workbench软件的初级和中级用户，可供机械、材料、土木、能源、汽车交通、航空航天、水利水电等专业的本科生、研究生、教师、工程技术人员和CAE爱好者学习参考。

前  言

CAE（Computer Aided Engineering，计算机辅助工程）是利用计算机辅助求解复杂工程和产品结构各项性能和优化设计等问题的一种近似数值分析方法，存在于工程的整个生命周期。ANSYS软件是最经典的CAE软件之一，应用广泛。近些年，ANSYS公司收购了多款顶级流体、电磁类软件，并重点发展ANSYS Workbench平台。与ANSYS经典软件界面相比，ANSYS Workbench 软件具有一目了然的分析流程图，整个分析就像在做一道填空题。在 ANSYS 12.0 版本之后，越来越多的用户转向使用ANSYS Workbench，同时有关ANSYS Workbench软件的参考书也越来越多。

本书以先进性、科学性、实用性、服务性为原则，在表达风格上力求通俗、简洁、直观。本书采用对比的方式详细讲解ANSYS热学分析过程中的各种问题，例如Mechanical热分析模块和CFD（Computer Fluid Dynamics，计算流体力学）热分析模块的区别和联系，稳态热平衡与瞬态热平衡的后处理差异，并以工程实例的方式帮读者分析问题；也采用外延拓展的方式详细介绍ANSYS热学分析和优化设计相关的计算，例如以热学计算原理衍生扩散分析，以模态分析原理演绎热模态分析，以拓扑优化同步扩展增材制造分析，并对同类问题进行适当扩展，以处理实际工程问题。书中不仅详细介绍了操作流程，还清晰阐述了“为什么要这样操作？同类的问题该有怎样的分析思路？”正所谓“入木三分方能见微知著”。书中内容结合相关理论知识，从实际应用出发，文字通俗易懂、深入浅出，引领读者轻松掌握 ANSYS Workbench软件的热学分析和优化设计方法。

全书共4章，第1章介绍了热学有限元分析的基本概念，并对比了FEM（Finite Element Method，有限元法）和FVM（Finite Volume Method，有限体积法）；第2章详细说明了ANSYS Workbench软件中Mechanical模块的热分析，包含热传导、热对流、热辐射、瞬态和非线性热分析、扩散、热固耦合分析，特别讲解了稳态和瞬态热平衡的后处理；第3章详细说明了ANSYS Workbench软件中Fluent模块的热分析，包含热传导、移动热源、对流和辐射、多孔介质的蒸发、凝华结霜、热分析标量方程、流固耦合传热分析，重点对Fluent中UDF的应用进行了介绍；第4章详细介绍了ANSYS Workbench软件中的优化设计过程，重点讲解了拓扑优化、增材制造和试验优化设计的基本原理及分析过程，包含相关优化设计方法、增材制造前处理、后处理及优化和反演设计。

本书内容新颖，紧密联系实际工程问题。例如，利用冰箱对不同初始温度的水进行冷却，通过查看结冰时间以验证Mpemba效应；对包含多层多孔介质的蒸锅进行蒸发计算，以验证蒸锅哪一层温度更高等。此外，本书还利用Excel软件完成相关任务，例如基本热学有限元计算和优化设计，接触热阻的参数优化设计，Fluid116、Combine37单元热计算，Material Designer模块的元胞法计算，Radiosity和AUX12辐射算法比较，热接触和Joint连接计算，热屈曲和热模态分析，Fluent凝华分析，基于Fluent的流固耦合传热分析，流体拓扑优化分析，增材制造前处理、工艺模拟、材料分析及过程模拟，以直接优化实现网格参数化得到网格划分基本规则“厚三圆十”，反演设计，以及可用于数字孪生的ROM（Reduced Order Model，降阶模型）分析。

本书主要面向ANSYS Workbench软件的初级和中级用户，可作为机械、材料、土木、能源、汽车交通、航空航天、水利水电等专业的本科生、研究生和专业教师的自学和教学用书，也可供相关领域从事产品设计、仿真和优化设计等工作的工程技术人员及广大CAE工程师学习参考。

本书配套提供全书案例的模型文件，读者直接在ANSYS Workbench 2000及以上版本的软件中打开或导入即可使用。本书配套资源可在异步社区网站或QQ群“CAE基础与提高交流”（601859149）的群公告中查看并下载。

本书由周炬、苏金英撰写。在写作过程中得到了丁德馨教授、雷泽勇教授、邱长军教授、李必文教授、唐德文教授、朱红梅教授的悉心指导，在此深表感谢！本书还得到了湖南省学位与研究生教育改革课题“湖南省研究生高水平教材（93YSM001）”的资助。

由于时间仓促，加之本书内容新、专业性强且作者水平有限，书中难免存在不足之处，恳请广大读者批评指正。

作者  

2024年5月

资源与支持

资源获取

本书提供如下资源：

● 案例素材；

● 彩图文件；

● 本书思维导图；

● 异步社区7天VIP会员。

要获得以上资源，您可以扫描下方二维码，根据指引领取。

提交错误信息

作者和编辑尽最大努力来确保书中内容的准确性，但难免会存在疏漏。欢迎您将发现的问题反馈给我们，帮助我们提升图书的质量。

当您发现错误时，请登录异步社区（https://www.epubit.com），按书名搜索，进入本书页面，单击“发表勘误”，输入错误信息，单击“提交勘误”按钮即可（见下图）。本书的作者和编辑会对您提交的错误信息进行审核，确认并接受后，您将获赠异步社区的100积分。积分可用于在异步社区兑换优惠券、样书或奖品。

与我们联系

我们的联系邮箱是contact@epubit.com.cn。

如果您对本书有任何疑问或建议，请您发邮件给我们，并请在邮件标题中注明本书书名，以便我们更高效地做出反馈。

如果您有兴趣出版图书、录制教学视频，或者参与图书翻译、技术审校等工作，可以发邮件给我们。

如果您所在的学校、培训机构或企业，想批量购买本书或异步社区出版的其他图书，也可以发邮件给我们。

如果您在网上发现有针对异步社区出品图书的各种形式的盗版行为，包括对图书全部或部分内容的非授权传播，请您将怀疑有侵权行为的链接发邮件给我们。您的这一举动是对作者权益的保护，也是我们持续为您提供有价值的内容的动力之源。

关于异步社区和异步图书

“异步社区”是由人民邮电出版社创办的IT专业图书社区，于2015年8月上线运营，致力于优质内容的出版和分享，为读者提供高品质的学习内容，为作译者提供专业的出版服务，实现作者与读者在线交流互动，以及传统出版与数字出版的融合发展。

“异步图书”是异步社区策划出版的精品IT图书的品牌，依托于人民邮电出版社在计算机图书领域的发展与积淀。异步图书面向IT行业以及各行业使用IT的用户。

第1章　关于热分析的基本解析

热以能量形式从高温区域流向低温区域，这种能量转移与温差（温度梯度）和高低温区域相对应。其中涉及的主要物理量为热量（能量，单位：J）、热流率（功率，单位：W）、温度（基本量纲，单位：开氏度K、摄氏度℃、华氏度℉、兰氏度°R等）。

热量传递有3种途径：传导、对流和辐射。传导的机理是在微观尺度下，分子受到输入热量作用而被激发，导致振动振幅增加，与相邻分子相互碰撞并传递了部分能量，即热量由高温部位传递给低温部位，传递介质为固体和流体。对流的机理是在流体（液体或气体）被加热时，靠近热源的流体体积增大，温度升高，密度降低而上升，而温度较低/密度较高的流体又填补到热源区域，周而复始，热量即通过热流体运动实现传递，传递介质为流体。辐射的机理是温度高于绝对零度的所有物质都以电磁波的形式散发热量，且传播时不需要介质，热量即通过电磁波进行传递。

由传递介质可知，热分析时不仅涉及固体，还涉及流体。因此在ANSYS Workbench软件中既存在Mechanical模块（各类Thermal分析），又存在CFD模块（Fluent、CFX等）。其计算流程虽不尽相同，但都可以完成传导、对流和辐射分析。

1.1　Mechanical模块与CFD模块的热分析计算

下面以一个简单例子说明各模块进行热分析的流程。已知一个电熨斗的截面如图 1-1-1所示，热传导系数为60.5 W·m⁻¹·℃⁻¹，板厚为0.005 m，置于环境温度为25℃的空间内，电功率为1 200 W，考虑辐射在内的表面换热系数为80 W·m⁻²·℃⁻¹，求稳态条件下电熨斗上下平面的温度。

图1-1-1　电熨斗的截面

解：在笛卡儿坐标系内，三维非稳态传热通用微分方程为：

式中，T为温度，t为时间，ρ为密度，c为比热容，k为热传导系数，x、y、z为位置坐标，Q为内热源（未标注单位，均为国际单位制单位）。

对于一维（仅x向）、稳态和无内热源模型，方程简化为：

其通解为：

第一类边界条件为（给定壁面温度）：

式中，为输入热通量，(W·m⁻²)。

第二类边界条件为（给定壁面热通量）：

式中，为环境温度，h为表面换热系数，为厚度，。

将C₁和C₂代入通解方程，则：

当 m时，上表面温度为629.91℃；当 m时，下表面温度为625.94℃。

下面用ANSYS Workbench的Steady-State Thermal、CFX和Fluent模块分别进行计算，以验证其计算结果。

1．建立分析流程

图1-1-2所示为热分析流程。其中包括A框架结构的Spaceclaim建模模块、B框架结构的Steady-State Thermal稳态热分析、C框架结构的划分网格、D框架结构的CFX热分析、E框架结构的划分网格（复制C框架结构的划分网格）、F框架结构的Fluent热分析。

注意

CFX和Fluent流体分析流程有多种形式，本例采用C3（Mesh）选项与CFX的D2（Setup）选项建立关联；E框架由C框架复制而来，具体操作是右击C框架，选择Duplicate；E3（Mesh）选项与Fluent的F2（Setup）选项建立关联。

2．在Steady-State Thermal模块定义材料

在B2处双击鼠标左键，进入Engineering Data界面，新建材料为al，定义热传导系数（即Isotropic Thermal Conductivity）为60.5 W·m⁻¹·℃⁻¹，如图1-1-3所示。

图1-1-2　热分析流程

图1-1-3　材料定义

注意

对于稳态热分析，至少需要输入热传导系数这项材料参数。

3．稳态热分析之前处理

在B4处双击鼠标左键，进入Mechanical前处理。

首先单击Geometry→Geom\实体，在下方的Details窗口中展开Material选项，将Assignment设置为al，即选择定义的al材料。然后右击Mesh→Insert→Sizing，选择整个实体，在下方的Details窗口中的Type项中选择Element Size，将Element Size设置为2 mm，如图1-1-4所示。此处定义单元尺寸是为了保证在同一网格尺度下对比计算结果。

图1-1-4　网格划分

4．稳态热分析之边界条件

Analysis Settings各选项全部采用默认设置。边界条件设置如图1-1-5所示。选择电熨斗的上平面，对其加载Heat Flow（热流），数值为1 200 W；选择电熨斗的下平面，对其加载Convection（对流），其中Film Coefficient（换热系数）设置为80 W·m⁻²·℃⁻¹，Ambient Temperature（环境温度）设置为25℃。

图1-1-5　边界条件设置

5．稳态热分析之后处理

单击Steady-State Thermal→Solution，分别查看Temperature和Total Heat Flux项，其中Temperature的后处理结果如图 1-1-6 所示。计算结果显示电熨斗上平面（定义热流面）的温度为629.92℃，电熨斗下平面（定义对流面）的温度为625.947℃。查看Total Heat Flux结果为48 075.74 W·m⁻²。

对比之前的计算结果，数值之间的差异原因仅在于面积取值精度，如果将软件计算的热流密度代入重算，可认定两种计算结果完全一致。

图1-1-6　稳态热分析之Temperature的后处理结果

6．CFX热分析之划分网格

在C3处双击鼠标左键，进入Mesh划分网格。在下方的Details窗口中将Physics Preference设置为CFD，Solver Preference设置为CFX（以CFX求解器进行CFD计算），Element Order设置为Linear，如图1-1-7所示。再插入Sizing，选择整个实体，在下方的Details窗口中的Type项中选择Element Size，将Element Size设置为2 mm，与图1-1-4定义一致。

图1-1-7　CFX网格划分

注意

在Physics Preference项可选择Mechanical、Electromagnetics、CFD、Explicit、Nonlinear Mechanical、Hydrodynamics和Custom 7种物理场形式。表1-1-1列出本系列书籍相关物理场网格划分的默认设置参数。其中Element Order（单元的阶）表示选择线性（一次）单元或二次单元；Straight Sided Elements（直边单元）表示选择二次单元时，中节点是否与首尾两节点呈直线，一般为保证有限元模型与实体模型的匹配精度，都定义为No；Sizing Options（尺寸调整）分别基于网格大小调整（Use Adaptive Sizing）或曲率角度调整（Capture Curvature）；Transition（过渡）/Growth Rate（增长率）表示相邻网格尺寸的变化程度，例如1.2的增长率即为后续网格比前面网格的尺寸大20%；Span Angle Center（跨度角中心）/Curvature Normal Angle（曲率法向角）表示一个单元跨越的模型角度，越小则有限元模型与实体模型的匹配精度越高；Smoothing（平滑）表示网格节点自动匹配实体模型轮廓的程度，平滑度越高则有限元模型与实体模型的匹配精度越高；Collision Avoidance（避免冲突），需要在Inflation-View Advanced Options选项中选择Yes才能查看，因为曲面膨胀很可能导致相邻处冲突，产生棱锥形网格，进而造成网格质量下降，所以默认使用Stair Stepping（梯步）法，可以减少此类网格，如果采用Layer Compression（层压缩）法则效果更好；Transition Ratio（膨胀-过渡比）表示膨胀层体积变化比，但是特别注意除了CFX的过渡比为0.77，其余均是0.272，这是因为Fluent等流体软件以单元中心方式计算膨胀比，单元控制体积与单元一致，而CFX则以节点中心方式计算膨胀比，单元控制体积是以点为基准的各个单元形心连接构成区域，所以数值上存在差异，但处理结果相似；Rigid Body Behavior（刚体行为）表示对刚体进行全体划分网格处理（Full Mesh）或者进行接触区域部分面划分网格处理（Dimensionally Reduced）。

表1-1-1　默认网格参数设置

划分完网格后，建议再右击Mesh，选择Update，会出现“The mesh translation to CFX/Fluent was successful”的提示，这样可以方便后续操作。

如图1-1-8所示，采用Named Selections分别给电熨斗模型面定义不同名称，以方便后续分析，其中上面定义名称为“heatflow”，下面定义名称为“convection”，周边 5 个面集合定义名称为“wall”。Named Selections的具体操作参见《ANSYS Workbench有限元分析实例详解（静力学）》^[1]，其功能不仅对于流体分析极其重要，而且特别适用于优化设计等批量处理功能。

[1]　《ANSYS Workbenoh有限元分析实例详解（静力学）》，周炬、苏金英著，ISBN：978-7-115-44631-2。本书后边提到的与此同名的书均是指的同一本。

图1-1-8　定义名称

7．CFX热分析设置

在D2 Setup处双击鼠标左键，进入CFX分析模块。由于Mesh项已经与CFX的Setup项建立关联，所以Outline窗口下Mesh处出现B4（体模型默认名称）和之前定义的名称；在1区Materials下双击Aluminium项，弹出Details of Aluminium对话框，Basic Settings选项卡中全部采用默认设置，Material Properties选项卡中仅修改Transport Properties区域中的Thermal Conductivity为60.5 W·m⁻¹·K⁻¹（与题设参数匹配）；在2区双击Analysis Type项，定义Option为Steady State（稳态，与题设条件一致），先右击Default Domain，将其删除，再右击Flow Analysis1，插入Domain（创建域1），在Domain1的Basic Settings选项卡中定义Location为B4（体模型），在Domain Type处选择Solid Domain（固体域），其余项采用默认设置，在Domain1的Solid Models选项中，Heat Transfer-Option项选择Thermal Energy，其余项采用默认设置，如图1-1-9所示。

图1-1-9　CFX分析设置（1）

相关选项说明

在Details of Aluminium对话框中的Basic Settings选项卡中，Option下拉列表框中的Pure Substance选项用于对纯物质创建物化参数（如黏度、密度或摩尔质量等）；Fixed Composition Mixture选项基于固定质量分数创建固定成分混合物，且在空间或时间模拟过程中不允许改变；Variable Composition Mixture选项用于创建在空间和时间模拟过程中质量分数允许变化的混合物；Homogeneous Binary Mixture选项基于饱和度创建平衡相变混合物；Reacting Mixture选项用于创建化学反应混合物；Hydrocarbon Fuel选项用于定义碳氢燃料。

在Details of Domain 1 in Flow Analysis 1对话框中的Solid Models选项卡中，Heat Transfer区域的Option下拉列表框中的None选项表示不考虑传热的分析；Isothermal选项表示基于特定温度条件下，不考虑传热的分析；Thermal Energy选项表示不考虑流体动能所产生的热量变化，均为低速状态的传热分析，例如低于0.3马赫的不可压缩流体传热；Total Energy选项表示考虑流体动能所产生的热量变化。

如图 1-1-10 所示，右击左侧Domain 1项，依次插入边界条件Boundary 1、Boundary 2和Boundary 3，在3区中双击Boundary 1，在Details of Boundary 1 in Domain 1 in Flow Analysis 1对话框中的Basic Settings选项卡中，Location下拉列表框中选择“convection”，Boundary Details选项卡的Option下拉列表框中选择Heat Transfer Coefficient，Heat Trans.Coeff.定义为80 W·m⁻²·K⁻¹（与题设参数匹配），Outside Temperature定义为25℃（与题设参数一致）。同样的方法设置Boundary 2和Boundary 3。Boundary 2的Location下拉列表框选择“heatflow”，Option下拉列表框选择Heat Flux，Heat Flux in设置为48 075.7 W·m⁻²（为保证精度，采用稳态热分析计算结果数值）；Boundary 3的Location下拉列表框选择“wall”，Option下拉列表框选择Heat Flux，Heat Flux in定义为0 W·m⁻²。在4区双击Solver Control，在Details of Solver Control in Flow Analysis 1对话框的Basic Settings选项卡中将Max. Iteration设置为200，其余全部采用默认设置。

图1-1-10　CFX分析设置（2）

相关选项说明

因为仅对固体模型进行热分析，所以Boundary Type均为Wall，在Option下拉列表框中Adiabatic选项为绝热边界，即热通量为0，本例第三个边界条件也可以采用此项；Fixed Temperature、Heat Flux和Heat Transfer Coefficient选项均与Mechanical热分析边界条件一致。

使用耦合求解器时，CFX对稳态分析求解采用伪瞬态法，一般在Details of Solver Control in Flow Analysis 1对话框的 Basic Settings选项卡的Convergence Control区域，将Max. Iterations设置为100～200个迭代步即可保证收敛。如果在200步内还没有收敛，则需要考虑调整Timescale（时间尺度），而不是增加迭代步数。

Timescale Control计算总时间和每步计算时间达到实现计算收敛的目的，CFX分为3种设置方法：（1）Auto Timescale（自动），即软件依据模型自动计算网格长度和平均速度之比，如果同时将Timescale Factor设置为0.3则更好；（2）Local Timescale Factor（局部因子），即针对速度差异较大且网格尺寸一致的分析，软件以局部因子乘上Timescale以控制收敛；（3）Physical Timescale，一种又好又快的收敛设置，对于通用模型，该值一般为网格长度和平均速度之比的三分之一；对于旋转机械模型，该值为0.1/ω～1/ω（ω为转速，单位为rad/s）；对于自然对流模型，该值为（l为垂直于温度梯度的单元长度，单位为m；为热膨胀系数，单位为℃⁻¹；g为重力加速度，单位为m·s⁻²；ΔT为温差，单位为℃）。本例只涉及固体传热分析，可通过Solid Timescale Control进行调整，因为固体的传热能量方程只依据热传导系数、密度、比热容等参数，计算非常稳定，所以固体的Timescale参数一般为流体的Timescale参数的100倍。

Convergence Criteria用于判定是否收敛且计算结束，其优先级高于Max. Iterations设置。Residual Type可设置为MAX（最大）或RMS（均方根）。Elapsed Wall Clock Time Control用于控制计算时间，如果计算花费超过设置计算时间，则不管是否收敛均停止运算。Interrupt Control基于CEL表达式进行计算终止判断。其他设置默认即可。

在D3 Solution处双击鼠标左键，进入CFX求解模块，如图1-1-11所示。对1区采用默认设置（说明：本例不需要设置Double Precision，双精度一般用于域尺寸很大、网格纵横比较大、边界条件差异较大等计算，以保证细节参数足够精确。本例选中Double Precision复选框是为了方便进行计算精度对比），2区对Initial Values采用默认设置（说明：该处与分析设置定义的初始值匹配，对于稳态分析，并不需要必须设定初始值，但是设定合理的初始值对于计算收敛有帮助），Initialization Option选项中的4个选项依次为Automatic（尽量不要使用）、Cached Solution Data、Current Solution Data、Initial Conditions（建议优先使用），最后单击Start Run，即进行计算。

图中显示为T-Energy的均方根残差收敛曲线，在图表中右击，在出现的3区菜单中选择Monitor Properties，在弹出的对话框中可插入各种残差收敛曲线；单击4区New Monitor图标，Type选择Plot Monitor，还可以插入其他检测图形，因为热分析中最为重要的是热能不平衡值（Imbalance），所以必须输出。由5区数据可知，边界条件1为热量输出（对流散热），其值为1 192.3 J；边界条件2为热量输入（加载热流），其值为1 200 J（1 200 W × 1 s，稳态单位时间），因为热分析必须遵守能量守恒原则，所以本次计算差值为 7.645 1 J，不平衡值为7.645 1/1 200 = 0.637 1%（能量差值除以能量输入总量）。

图1-1-11　CFX求解设置

CFX收敛说明

CFX以迭代残差达到某一设定值作为收敛判定的依据。迭代残差是指相邻两次迭代计算之间同一物理量的差值，相关概念可参考《ANSYS Workbench有限元分析实例详解（静力学）》。因为固体热分析只以能量残差进行收敛判定，所以计算到低于图1-1-10中定义的Residual Target为1E−4时即完成计算，并不完全遵守设定的最大迭代步数。

在D4 Results处双击鼠标左键，进入CFX后处理模块，如图1-1-12所示。首先为了温度读数方便，在1区依次单击Variables→Temperature，在Units下拉列表框中选择C，即采用℃单位；在2区右击User Locations and Plots，使用快捷菜单中的命令插入Contour 1，如图设置后，可得上平面温度为627.69℃，下平面温度为623.56℃。

图1-1-12　CFX后处理模块

8．Fluent热分析之划分网格

在E3处双击鼠标左键，进入Mesh划分网格。由于网格和命名选择设置已经完成，仅需要将Solver Preference修改为Fluent（以Fluent求解器进行CFD计算），如图1-1-13所示。

图1-1-13　Fluent网格划分

9．Fluent热分析设置

在F2 Setup处双击鼠标左键，保持默认设置，单击Start进入Fluent分析模块，如图1-1-14所示。建议此刻回到图1-1-2所示的热分析流程界面进行存盘设置，如此可以保证临时存取文件不占用系统盘空间。在1区单击Units按钮，在弹出的Set Units对话框中修改温度单位为℃；在2区右击Energy，在弹出的快捷菜单中选择On（开启能量模型，热分析必须设置）；为设置方便，不必新建材料，直接修改软件默认材料的参数，在3区双击Materials→Solid→aluminum，修改Thermal Conductivity为60.5 W·m⁻¹·K⁻¹（与题设参数匹配）；在 4 区展开 Cell Zone Conditions，右击默认模型，在弹出的快捷菜单中选择Type→Solid（将域设置为固体，且采用软件默认的固体材料Aluminum参数），其余均为默认设置。

图1-1-14　Fluent分析设置（1）

如图1-1-15所示，在Boundary Conditions→wall项中依次双击convection、heatflow和wall（依据命名选择自动创建边界条件项），在5区中Thermal选项卡的Convection处将Heat Transfer Coefficient设置为80 W·m⁻²·K⁻¹（与题设参数匹配），Free Stream Temperature设置为25℃（与题设参数一致）；在6区中将Thermal选项卡中的Heat Flux设置为48 075.7 W·m⁻²（为保证精度，采用稳态热分析计算结果数值）；在7区中将Thermal选项卡中的Heat Flux设置为0 W·m⁻²。其余全部采用默认设置。

图1-1-15　Fluent分析设置（2）

相关选项说明

Fluent的两种求解类型分别为Pressure-Based和Density-Based，其中Pressure-Based（压力基）主要适用于计算0～3马赫范围的大部分工程问题；Density-Based（密度基）主要用于大于3马赫的工程问题，或者冲击波相互作用等特定工况条件。

Fluent热边界条件共有7种，其中Heat Flux（热通量）、Temperature（温度）、Convection（对流）与Mechanical热分析和CFX边界条件一致；Radiation（辐射）类型在第3章详细描述；Mixed（混合）类型合并了对流和辐射条件；Via System Coupling（通过系统耦合）用于在ANSYS Workbench中Fluent与其他模块进行系统耦合分析，即传输不在Fluent中定义的热分析参数；Via Mapped Interface（通过映射面）用于部件接触存在穿透或间隙时，自动创建映射面以传递热。

在边界条件中有Wall Thickness和Shell Conduction选项，其中Wall Thickness只需要定义厚度，不需要生成网格，热量只沿模型厚度的法向传递；Shell Conduction不仅定义厚度，还根据层数（layer）生成虚拟网格，热量可向任何方向传递，如图1-1-16所示。

（a）Wall Thickness

（b）Shell Conduction

图1-1-16　Wall Thickness和Shell Conduction选项的区别

10．Fluent热分析求解设置

如图1-1-17所示，双击1区Solution→Methods，设置求解方法，在Scheme处采用默认的SIMPLE法；双击2区Monitors→Residual，设置残差极限，设置energy的Absolute Criteria为1E−6；双击3区Initialization，进行初始化，采用默认的Hybrid Initialization法；双击4区Run Calculation，设置Number of Iterations为500。其余均保持默认设置。

图1-1-17　Fluent求解设置

相关选项说明

采用Pressure-Based求解器后，可以选择的算法有SIMPLE、SIMPLEC、PISO和Coupled，其中SIMPLE为默认算法，适用大多数不可压缩问题；SIMPLEC对于四边形和六面体网格模型更加有优势；PISO通常用于瞬态分析；Coupled的研究对象通常为可压缩流体或考虑浮力与结构运动耦合的不可压缩流体。

Fluent通过指定初始值来提高迭代计算稳定性并加速收敛。其中Hybrid初始化方法为Fluent默认的初始化方法，表现为所有单元初始值非均匀化；Standard初始化方法为最简单的方法，表现为所有单初始值均匀化。

与CFX类似，当能量残差低于定义的Energy的Absolute Criteria值（即1E−6）时即完成计算。

11．Fluent热分析后处理

如图1-1-18所示，右击1区Contours创建“contour-1”，在Contours of区域的下拉列表框中分别选择Temperature和Static Temperature，在Surfaces下方选择convection、heatflow和wall，再单击Compute按钮，即可得到最高温度为630.12℃，最低温度为625.79℃^[2]；双击2区Reports→Fluxes，在Flux Reports对话框中选中Total Heat Transfer Rate单选按钮，同样选择convection、heatflow和wall，再单击Compute按钮，即可得到输入和输出热功率，以及两者的差1.91 W。

[2]　此处的“最高温度”和“最低温度”相比于图1-1-18所示的参数进行了四舍五入处理，保留至小数点后两位。

12．计算误差分析

以上3个模块的计算采用了同样尺度的网格，基本均保持默认设置，3个模块的计算结果对比如表1-1-2所示。

图1-1-18　Fluent后处理设置

表1-1-2　3个模块的计算结果对比

由表1-1-2可知，3个模块的计算误差均小于1%，导致误差的原因是计算时输入和输出的热量不平衡。其中Mechanical模块的精度最高，如果去除输入数值误差，其与理论计算值一致，Fluent模块的精度其次，CFX模块的精度最低，但是不能就此简单地评价3个模块的计算精度。原因有二：其一，Mechanical和CFD（含Fluent和CFX）模块的比较应该考虑算法和机时对精度的影响；其二，CFD模块默认设置不同导致精度差异。

下面以不修改CFD迭代步数为条件，仅修改CFD相关默认设置，继续对比。

注意

对于所有已求解完成的流体分析模型，在重新导入外部设置前，都建议右击Solution，在出现的快捷菜单中选择Clear Generated Data，否则容易出现冗余设置。

在D2 Setup处双击，进入CFX设置项，在Solver Control项中修改Residual Target参数1E−4为1E−6，这样才与Fluent计算中定义Energy的Absolute Criteria为1E−6相匹配。重新计算，如图1-1-19所示，在图表中右击，在弹出的快捷菜单中选择Monitor Properties，在弹出的对话框中选中IMBALANCE→Domain 1→T-Energy Imbalance（%）In Domain 1，还可监控不平衡值曲线，其最终输入和输出能量不平衡值为0%；查看后处理结果可知上平面温度为629.92℃，下平面温度为625.95℃，该结果与Mechanical模块的计算结果一致。

注意

因为设定的能量收敛残差为1E−6，所以整个计算过程残差均大于设定值，虽然没有实现规定的残差收敛，但是仍按图1-1-10定义的Max. Iteration为200完成计算。对于此类热分析，以不平衡值作为收敛的判定标准，即要求输入和输出能量平衡。

图1-1-19　CFX相关设置

在F2 Setup处双击鼠标左键进入Fluent设置项，双击Solution→Methods，在弹出的对话框中的Scheme处将SIMPLE修改为Coupled，软件会自动选中Pseudo Transient，这样才与CFX的耦合求解器相匹配。重新计算6步即收敛，如图1-1-20所示，查看后处理结果可知上平面温度为629.89℃，下平面温度为625.43℃，输入和输出热功率相差1.74 W，计算精度高于默认的SIMPLE算法。

图1-1-20　Fluent相关设置

对比两者结果可知，CFX计算误差已经可以忽略不计，Fluent仍然存在很微小的误差，但是即便CFX在70余步（以不平衡值判定）达到收敛，计算所需机时仍长于Fluent。

如何提高CFX的计算效率？

（1）修改Timescale参数

如图1-1-21所示，双击Solver→Solver Control，在弹出的对话框的Solid Timescale下拉列表框中选择Physical Timescale，将Solid Timescale设置为10 000 s，注意Residual Target依然采用默认的1E−4，重新计算3步即收敛，输入和输出能量不平衡值为0.031 4%；查看后处理结果可知上平面温度约为630.00℃，下平面温度约为 625.87℃，与理论值相比误差约为0.012%。

图1-1-21　CFX相关设置（1）

注意

Solid Timescale默认为100 s，对于固体热分析该参数都适用，该值较小，使计算趋势更容易收敛，但所需机时更长；如果该值较大，使计算趋势更容易发散，但所需机时更短，对于固体热分析可以设置500～1 000 s以加快收敛。本例设置为10 000 s仅仅为了对比清晰，实际分析时不建议采用。

（2）修改Conservation Target参数

如图1-1-22所示，双击Solver→Solver Control，在弹出的对话框中选中Conservation Target复选框，Value设置为0.000 1，注意Solid Timescale默认为Auto Timescale，Residual Target默认为1E−4，重新计算20步即收敛，输入和输出能量不平衡值为0.006 8%（0.006 8%<0. 01%，以不平衡值为收敛目标）；查看后处理结果可知上平面温度为629.90℃，下平面温度为625.92℃，与理论值相比误差约为0.001 5%。

图1-1-22　CFX相关设置（2）

注意

该选项以不平衡值为收敛目标，精度和速度都可以很好地满足计算要求，建议选用。此外，Timescale项可以与Conservation Target项同时设置，但Timescale项优先级别更高。

如何提高Fluent的计算精度？

（1）修改Solution Controls的Energy参数

如图1-1-23所示，双击Solution→Methods，在弹出的对话框中将Scheme选择为Coupled并选中Pseudo Transient项，然后双击Solution→Controls，在弹出的对话框中将Energy设置为1，重新计算4步即收敛，输入和输出能量不平衡值为0.031 9；查看后处理结果可知上平面温度为630.01℃，下平面温度为625.85℃，与理论值相比误差约为0.015%。

图1-1-23　Fluent相关设置（1）

注意

此处表示采用Pressure-Based求解器时，通过调整亚松弛因子控制收敛，该参数意义为相邻两次迭代计算时某物理量变化的大小。通过对比图1-1-20可得，迭代次数较之前少了两步即可收敛，但计算精度相差无几。这是因为最终的收敛解与亚松弛因子无关，该参数只调整收敛所需的迭代次数。对于稳态固体传热，可以设为1以提高收敛效率。

（2）修改Solid Time Scale的Pseudo Time Step参数

如图1-1-24所示，保留Solution Methods对话框的设置，双击Solution→Controls，在弹出的对话框中将Energy设置为0.75，修改Run Calculation对话框中Solid Time Scale区域的Time Step Method为User-Specified，Pseudo Time Step设置为100 000，重新计算6步即收敛，输入和输出能量不平衡值为1.92；查看后处理结果可知上平面温度为629.79℃，下平面温度为625.33℃，与理论值相比误差约为0.019%。

注意

该项设置与CFX模块的Timescale参数设置几乎一致，区别为CFX模块中其默认值为100 s，Fluent模块中其默认值为1 000 s，按照图1-1-21设置统一放大100倍后，Fluent的计算精度并没有显著提升。

同理修改Solid Time Scale的Time Scale Factor参数。如图1-1-25所示，保留之前设置，在Run Calculation对话框中设置Solid Time Scale的Time Scale Factor为100，重新计算6步即收敛，输入和输出能量不平衡值为1.74；查看后处理结果可知上平面温度为629.89℃，下平面温度为625.43℃，与理论值相比误差约为0.009%。

图1-1-24　Fluent相关设置（2）

图1-1-25　Fluent相关设置（3）

注意

通过调整该值可以略微提高精度，一般设置为5倍进行调试，数值定义较大不能再提高精度，而且会使计算趋势发散。本例定义为100仅仅为了清晰对比，实际分析时不建议采用。

综上所述，对于CFD固体稳态热分析，CFX模块可以通过修改收敛残差、修改Solid Time Scale值、以不平衡值为收敛目标等方法提高计算精度和效率；Fluent模块与CFX模块基本类似，但是由于两者默认初始值不同，计算精度和效率略有差异。此外，CFX模块的很多选项需要使用者调整，而Fluent模块由于在历次升级中更新频繁，其默认选项基本不需要过多修改。当然，提高精度最直接的方法就是在控制收敛残差尽量小的情况下增加迭代步数。

1.2　Mechanical模块与CFD模块在热分析中的区别

通过前例可知，Mechanical模块与CFD模块均可以进行热分析，两者的区别如下。

1．计算原理不同

Mechanical模块采用有限单元法（FEM），Fluent模块采用有限体积法（FVM）的单元中心（cell-centered）方式，CFX模块采用有限体积法的节点中心（node-centered）方式。

注意

在FEM中将划分网格形成的块称为Element，在FVM中称为Cell。

下面以通量平衡基础方程说明FEM与FVM的区别：

式中，u为任意守恒物理量，为该物理量的通量，即单位面积、单位时间内通过域控制面的流动量。

采用有限元形式求解此类偏微分方程，并不拘泥强形式微分方程解，但至少满足积分点上的条件，即弱形式微分方程解。因此对模型域进行积分建立方程：

式中，为有限元计算加权试函数，V为域体积。

为了控制通量与模型内部向量的联系，对采用散度定律，则：

式中，为域的边界，n为边界上的法向量，S为域内封闭曲面的面积。

同时

则

代入模型域的积分方程，得：

该式即为FEM对偏微分求解的基本方程。

当试函数与形函数一致时，即为伽辽金法。形函数只与几何、物理量相关，一般设为多项式，也可以定义为常数，当时，则上式为：

该式即为FVM对偏微分求解的基本方程。

由此可知，FVM只是FEM求解中的一个特例，不同之处在于，FEM选取离散后有限数量的形函数，FVM选取离散后有限数量的体积域，分别使其满足控制域内的基本方程。

图 1-2-1 表示：对于内部节点，只有周围 6 个网格剖面线的单元（Element）与域的域积分有关联；对于边界上（）的节点，只有周围 3 个单元（Element）与域的域积分有关联，只有周围2 个斜线剖面线的单元（Element）与的域积分有关联。

图1-2-1　FEM离散原理

图1-2-2表示：每个单元（Cell）都被视为一个域，不管是内部还是边界上的单元（Cell），都与域的域积分有关联。

图1-2-3表示：在内部以每条边的中点加上每个单元（Cell）的中心构成一个域；在边界处除了边的中点和单元（Cell）的中心，再以边界进行封闭构成一个域（如果存在边界层，且为三角形网格，对于直角和钝角三角形单元，则以三角形最长边的中点加边的中点构成域）。不管是内部还是边界上的域，都与域的域积分有关联。

图1-2-2　FVM单元中心方式离散原理

图1-2-3　FVM节点中心方式离散原理

对比FEM和FVM离散原理可知，FEM可以保证全局通量守恒，但不能保证局部通量守恒。换言之，在域边界处通量保证平衡，但是难以控制内部局部通量。以热分析为例，纯传导分析易于计算，FEM可以采用较粗的二次单元完成计算且保证内部热通量平衡；但是针对对流、辐射为主导的热分析，虽然可以在边界处得到精度较高的结果，但是内部热量流动过程较难处理。

FVM可以保证局部通量守恒，但是在边界上的解难以保证精度，即便在边界处划分更多的网格也不行，除了降低效率，还会影响边界上整体与全域内局部的匹配。因此在热分析中，FVM在边界处通常都存在一些热不平衡量误差，但是针对对流、辐射为主导的热分析，可以很容易得到内部热量流动过程。

2．分析的侧重点不同

Mechanical模块利用单位等效类比方式将结构有限元分析扩展到热学有限元分析，即将结构分析中的位移自由度变为热学分析中的温度自由度，同理将结构中的其他物理量变为热学中的物理量，如表1-2-1所示。

表1-2-1　结构与热学物理量对应关系及单位

热传导计算公式为：

式中，q为热通量，Q为热量，A为截面积，k为热传导率，d为温差面之间的距离，T为温度。

对流计算公式为：

式中，q为热通量，h为对流换热系数（W·m⁻²·K⁻¹），T为温度。

辐射计算公式为：

式中，Q为热量，为辐射率（0～1），为斯特藩-玻尔兹曼常数（约为5.67 × 10⁻⁸ W·m⁻²·K⁻⁴），A为截面积，F为辐射角系数，T为温度。

对于传导过程，采用FEM可以快捷、准确地完成计算，完全不需要CFD计算；对于对流过程，如果对流换热系数为已知条件，FEM比CFD更方便计算，但是对流换热系数往往比较难获得一个较为准确的平均值，即便定义对流换热系数相对于时间或温度的函数，也存在一定的误差，所以对于强制对流形式，对流换热系数相对容易获得，可以采用CFD为主FEM为辅的计算形式，对于自然对流形式，对流换热系数较难获得，可以采用CFD计算；对于辐射过程，关键难度在于获得任意表面的辐射角系数，FEM将其过程大大简化，仅定义辐射率和Hemicube分辨率控制辐射所有参数，因此如果需要考虑更多的辐射因素，例如光谱依赖性、半透明固体的吸收、辐射热通量的平滑分布、太阳辐射等，采用CFD更加可靠。

下面以自然对流形式简要说明CFD分析流程。已知一个地暖区域，其二维截面积如图 1-2-4所示，下边温度为315 K，周围环境温度为290 K，求稳态条件下该区域的对流换热系数。

图1-2-4　对流计算模型

对流换热由高温模型边界外流体运动产生，将较热的流体从边界层表面输送出去，该区域空间由较冷的流体补充，以至循环。这种微观运移简化为宏观描述，即为对流。本例中没有其他外界因素促使该区域内空气流动，为自然对流。自然对流中流体循环的动力仅来源于热流和冷流的密度差异，重力起到了非常重要的作用。

在自然对流计算过程中，采用瑞利数作为流动特征评估参数，不采用雷诺数。瑞利数计算公式为：

（，层流；，湍流）

式中，（格拉晓夫数），用于定义浮力与黏性力之比，其中为热膨胀系数，单位为℃⁻¹，g为重力加速度，单位为m·s⁻²，L为垂直于温度梯度的长度，单位为m，ΔT为温差，单位为℃，为流体运动黏度，单位为m²·s⁻¹；（普朗特数），用于定义动量扩散与温度扩散之比，其中为流体热扩散系数，单位为m²·s⁻¹，在计算的这个公式中k为流体热传导系数，单位为W·m⁻¹·K⁻¹，为流体密度，单位为kg·m⁻³，C为比热容，单位为J·kg⁻¹·K⁻¹。

本例的物理参数如表1-2-2所示。

表1-2-2　本例的物理参数（均采用国际单位）

解：

将数值代入，则，流动形式为层流。

Nu（努塞尔特数）计算公式为：，式中参数如表1-2-3所示。

表1-2-3　努塞尔特数计算公式中的参数（均为系数，无单位）

本例中取C = 0.27，n = 1/4，则Nu = 34.93；对流换热系数计算公式为，其中L对于竖板模型为模型上垂直于温度梯度的长度，对于平板模型为面积与周长之比，本例中L = 0.167 m，则h = 5.062 W·m⁻²·K⁻¹。

下面用ANSYS Workbench的Fluent模块进行计算，以验证其计算结果。

（1）建立分析流程

如图1-2-5所示，建立热分析流程。其中包括A框架结构的划分网格、B框架结构的Fluent热分析。

图1-2-5　热分析流程

注意

本例采用A3（Mesh）选项与Fluent的B2（Setup）选项建立关联的分析流程，也可以采用Geometry组合Fluent（with Fluent Meshing）模块的分析流程。Fluent（with Fluent Meshing）模块可以在开始时定义划分网格和求解计算的硬件资源，在单独使用Fluent进行大型工程分析时建议采用。

（2）前处理

建模过程省略。

如图1-2-6所示划分网格，其中Details窗口中的Physics Preference处选择CFD，Solver Preference处选择Fluent；插入Sizing项对整个面定义网格尺寸，在Geometry处选择整个面，Type处选择Element Size，Element Size设置为0.005 m；插入Inflation定义边界层，在Geometry处选择整个面，Boundary处选择模型的4条边，Maximum Layers设置为10。

划分网格后，在Details窗口的Quality→Mesh Metric处查看网格质量。在Mesh Metric处选择Aspect Ratio（网格纵横比）。

注意

可以在Smoothing处选择High项进一步提高网格质量。

图1-2-6　划分网格

注意

对于自然对流分析，必须定义边界层，一般定义5～10层，网格质量重点为纵横比，必须小于40。

如图1-2-7所示定义命名选择，以便Fluent定义边界条件。其中上边线命名为“top”，下边线命名为“bottom”，左边线命名为“left”，右边线命名为“right”，整个面命名为“air-face”。

图1-2-7　命名选择

（3）Fluent热分析流程

在B2 Setup处双击鼠标左键，保留默认设置（2D模型），单击Start进入Fluent分析模块，如图1-2-8所示，在1区选中Gravity复选框，在Y方向定义重力加速度为−9.81（自然对流必须设置）；在2区设置Energy为On，Viscous为Laminar（层流）；在3区双击Materials→Fluid→air，在弹出的对话框中修改Density为boussinesq，设置数值为1.140 5 kg·m⁻³，设置Thermal Expansion Coefficient的数值为0.003 448 27 K⁻¹（切记不要遗漏），其余参数均为默认设置，修改完成单击Change/Create按钮保存；在4区Cell Zone Conditions处选择Fluid（定义为流体域，软件可能默认为固体域）；在5区Boundary Conditions处双击，在弹出的对话框中单击Operating Conditions…按钮，在出现的对话框中设置Operating Pressure为0 Pa，Operating Temperature为290 K。

图1-2-8　Fluent分析设置（1）

如图1-2-9所示，在Boundary Conditions→wall项中依次双击bottom、left、right、top（依据命名选择自动创建边界条件项），对“bottom”边在Thermal选项卡中设置Temperature为315 K；对“left”和“right”边，在Thermal选项卡中设置Heat Flux为0；对“top”边，在Thermal选项卡中设置Temperature为290 K。其余全部采用默认设置。

图1-2-9　Fluent分析设置（2）

注意

对于自然对流分析，必须定义重力。在定义流体密度参数时，一般采用Incompressible ideal gas（设定密度，推荐）或Boussinesq（单击Operating Conditions…按钮，在弹出的对话框的Operating Temperature处设置温度，温差在20%以内）形式，就计算精度而言，Boussinesq计算结果会导致羽流宽度变窄，存在一定误差。本例仅展示计算流程，所以采用Boussinesq。因为本计算模型为开放区域，所以Operating Pressure设置为0 Pa。

（4）Fluent热分析求解设置

如图1-2-10所示，在1区双击Solution→Methods，在弹出的对话框中，Scheme下拉列表框选择Coupled，Gradient下拉列表框选择Least Squares Cell Based，Pressure下拉列表框选择PRESTO！，Momentum下拉列表框选择First Order Upwind，Energy下拉列表框选择First Order Upwind，选中Pseudo Transient复选框；在2区双击Solution→Controls，在弹出的对话框中所有选项均保留默认设置；在3区双击Monitors→Residual，在弹出的对话框中，将Continuity的Absolute Criteria项设置为1E−5，其余项均为默认设置。

图1-2-10　Fluent分析设置（3）

注意

对于自然对流分析，因为流场和能量为强耦合，所以Scheme下拉列表框选择Coupled，如果采用SIMPLE也可以计算收敛，但计算效率和精度均低于Coupled。同样对于自然对流，Pressure下拉列表框只能选择PRESTO!或者Body-Force Weighted，其中前者应用范围更广一些，主要用于压力梯度较大的计算中，且计算效率高于后者。First Order Upwind与Second Order Upwind的主要区别在于收敛速度，前者容易收敛，但精度略逊于后者。选中Pseudo Transient（伪瞬态）复选框主要为了提高较大纵横比模型的收敛性能。

Solution Controls用于定义求解稳定性的亚松弛因子，该值与最终收敛结果无关，Fluent默认值可以用于大多数场合，如有必要可以调节（调小，提高收敛稳定性；调大，提高收敛速度）。

Fluent稳态分析收敛判断基于3点：①计算残差小于设定的标准残差；②计算域内的质量、动量、能量和标量达到平衡；③某个需要监控的物理量达到稳定状态。对于本例而言，如果仅以小于标准残差或系统平衡判定收敛，则容易出现很大的计算误差。在自然对流中，“left”边和“right”边的温度保持稳定才是该稳态分析的关键。如图1-2-11所示，先定义监控项，右击Report Definitions项，在出现的快捷菜单中选择New→Surface Report→Area Weighted Average…，在弹出的对话框中进行设置，Name文本框采用默认命名，在Field Variable下方的下拉列表中分别选择Temperature…和Static Temperature，其中“report-def-0”对应的Surfaces选择为“left”，“report-def-1”对应Surfaces选择为“right”。

图1-2-11　Fluent收敛设置（1）

如图1-2-12所示定义收敛监控项，右击Monitors→Report plots项，在出现的快捷菜单中选择New，在弹出的对话框中进行设置，Name文本框采用默认命名，分别单击Add按钮添加“report-def-0”和“report-def-1”，其余项均为默认设置。

图1-2-12　Fluent收敛设置（2）

单击Run Calculation，将Number of Iteration设置为500，其余项采用默认设置，计算完成后的监控图如图1-2-13所示。可以看到，在400步左右，监控左右两边线温度达到稳定。

图1-2-13　Fluent计算收敛

（5）Fluent热分析后处理

自然对流计算完成，先查看速度与温度分布云图，如图1-2-14所示。右击Contours项，在出现的快捷菜单中选择New，在弹出的对话框中进行设置，Contour Name文本框采用默认命名，在Contours of…下方的下拉列表框中选择“Velocity…”“Velocity Magnitude”和“Temperature…”“Static Temperature”，分别创建“contour-1”和“contour-2”，在Surfaces下方的列表框内选择“surface_body”，单击Compute按钮，最后单击Save/Display按钮即可。采用类似操作，生成速度向量图，如图1-2-15所示。

图1-2-14　速度和温度分布云图

图1-2-15　速度向量图

由图可见，从“bottom”边界到“top”边界的自然对流速度和热边界层的运移趋势，在“right”边界旁边形成浮力羽流，在“left”边界旁边形成浮射流。浮力羽流形成是因为流体初始动量很小（图中显示速度较慢，0.1 m·s⁻¹左右），但是主要受浮力作用（密度差，图中显示温差跨度约为模型的整个高度）继续向上流动和扩散，其扩散面形如羽毛。浮射流形成是因为流体初始动量较大（图中显示速度较快，0.25 m·s⁻¹左右），同时也受浮力作用（图中显示温差跨度约为模型高度的一半），在动量占优势的情况下，向下流动和扩散，最终动能耗散，又形成羽流，构成一个循环。

如图1-2-16所示，双击Results→Reports→Surface Integrals，在弹出的Surface Integrals对话框中，在Report Type下拉列表框中选择Area-Weighted Average，Field Variable下方的下拉列表框分别选择“Wall Fluxes…”和“Surface Heat Transfer Coef.”，在Surfaces下方的列表框内选择bottom、left、right、top，单击Compute按钮，可得对流换热系数为5.038 56 W·m⁻²·K^−1[3]，与之前的计算结果5.062 W·m⁻²·K⁻¹相差约0.46%。

[3]　因为Fluent的计算结果为对象的输出值，所以选取图中对应数值的绝对值。

图1-2-16　对流换热系数计算

综上所述，对流换热系数的获得虽然可以采用公式计算（或查表），但是对于复杂形貌、不能简单套用公式计算的流场区域，CFD模块可以很方便地计算自然对流状态的热传递过程并求出对流换热系数。但是该系数表现为区域内均值，不能反映区域内细节，即便FEM热计算时可将其作为输入已知边界条件，仍然对区域内的热分布细节估计不足，存在一定的系统误差。

1.3　热分析有限元计算的原理

1.1节和1.2节通过解析法对相关热分析进行了计算，但是解析法毕竟只适用一些特定模型。为解决复杂模型的热分析，通常采用有限元法进行计算，其优势主要表现为：复杂模型的适应性、热固耦合的连续性、各种非线性问题（非线性材料、非线性边界条件等）的可操作性。下面通过3个例子分别讲述有限元法进行热分析的计算原理。

1.3.1　复合棒材稳态热传导的有限元计算

已知一个直径为20 mm复材圆棒，如图1-3-1所示。其中长度为100 mm的棒材热传导系数为100 W·m⁻¹·℃⁻¹，长度为200 mm的棒材热传导系数为20 W·m⁻¹·℃⁻¹，长度为400 mm的棒材热传导系数为80 W·m⁻¹·℃⁻¹，圆棒周向绝热，左侧输入热流4 000 W·m⁻²，右侧保持温度为80℃，求稳态条件下棒材的温度分布。

图1-3-1　复合棒材模型

注意

为保证3个零件模型共节点，在Spaceclaim目录树上单击SYS，在其下的Properties中展开Analysis，将Share Topology设置为Merge。

1．建立单元热传导矩阵

对于一维（仅x向）、稳态、截面积一致和无内热源的模型，其第一类边界条件为（给定壁面温度）：

其中，，即热流出口温度与入口温度之差；，即单元长度；A为单元截面积；k为热传导系数。

对于热量入口，上式变换为：

根据能量守恒定律，热量出口为：

列成单元矩阵为：

2．建立整体单元热传导矩阵

简化计算，将每种材料的棒材仅划分为一个单元，从左至右节点编号依次为1、2、3、4，则每个单元的热传导矩阵为：

；；

整体热传导矩阵由各单元热传导矩阵叠加而成：

对于本例，模型内部没有热源和能量损失，则，上面的矩阵简化为：

计算可得T₁ = 144℃、T₂ = 140℃、T₃ = 100℃、q₄ = −4 000 W·m⁻²。

提示

矩阵计算用Excel也可以轻易完成，不需要另装其他软件。以本计算为例，矩阵形式为：

简化为：

如图 1-3-2 所示采用 Excel 完成矩阵计算，先计算逆矩阵，框选 G91:I93 格，在插入函数中输入 = Minverse(B91:D93)，按Ctrl + Shift + Enter组合键即可得到逆矩阵；再框选M91:M93格，在插入函数中输入 = MMULT(G91:I93,K91:K93)，按Ctrl + Shift + Enter组合键即可得到T₁、T₂、T₃的结果。同理计算MMULT(B91:E94,M91:M94)，即可得到q₄的结果。

图1-3-2　Excel矩阵计算

注意

本计算是基于每种材料的棒材仅划分为一个单元，如果划分为两个单元，同理计算可得T₁ = 144℃、T₂ = 142℃（位于原1单元中间）、T₃ = 140℃、T₄ = 120℃（位于原2单元中间）、T₅ = 100℃、T₆ = 90℃（位于原3单元中间），此结果正确；但是q₂ = −4.5E−11、q₃ = −2.3E−11、q₄ = 0、q₅ = 0、q₆ = 1.13E−11，这与模型内部没有热源和能量损失有关。

3．Steady-State Thermal模块计算

如图1-3-3所示，分别创建a、b、c三种材料，其Isotropic Thermal Conductivity分别定义为100、20、80。

图1-3-3　创建材料并定义材料参数

单击Geometry→Geom\SYS\a、b、c，在下方的Details窗口中展开Material选项，将Assignment分别对应设置为a、b、c，即选择定义的a、b、c材料。

Mesh设置保持默认。

边界条件设置如图1-3-4所示。其中选择a模型的左端面定义加载Heat Flux，数值为4 000 W·m⁻²；选择c模型的右端面定义加载Temperature，数值为80℃。

图1-3-4　边界条件定义

计算完成后，分别插入Temperature（Geometry设置为All Bodies）、Temperature2（Geometry设置为a、b模型交界面）、Temperature3（Geometry设置为b、c模型交界面）和Total Heat Flux四项，后处理的计算结果如图1-3-5所示。

图1-3-5　后处理的结果

由图可知，系统最高温度为144℃，位于模型最左侧；a、b模型交界面温度为140℃；b、c模型交界面温度为100℃，总热通量为4 000.1 W·m⁻²。与整体热传导矩阵计算结果相比，仅在总热通量上有0.002 5%的误差。

4．伽辽金法求解稳态热传导分析

对于物理问题一般都可以采用微分方程进行表达，但是并不是所有方程都可以通过泛函的极值方法进行求解。当采用加权余量法时，如果试函数与形函数一致，即为伽辽金法。伽辽金法有广泛的适用性，有限元软件一般都基于该方法进行计算。下面采用二阶单元对本模型进行计算。

由于二阶单元为3个节点，基于等参化建立插值函数（不是单元的形函数）为：

利用插值条件，可得：

求出a_i、b_i、c_i，代入插值函数，则：

令，则（取）、（取）、（取），插值函数为：

根据伽辽金法，热传导矩阵中的任意一项为：

　 i,m = 1, 3

又因为，对于任一单元，其中单元长度，，则：

，

热传导矩阵中任一项为：

　　i,m = 1,2,3，即

每种材料的棒材仅划分为一个单元，则a材料的热传导矩阵为：

同理，b、c材料的热传导矩阵分别为：

整体热传导矩阵为（每个材料单元的最后一个节点与下一个材料单元的第一个节点共点）：

则：

其中已知：Q₁/A = 4 000，Q_2-6 = 0，T₇ = 80。

计算可得：T₁ = 144℃、T₂ = 142℃、T₃ = 140℃、T₄ = 120℃、T₅ = 100℃、T₆ = 90℃（均无误）、q₂ = 0 W·m⁻²、q₃ = −1.7E−11 W·m⁻²、q₄ = −1.1E−11 W·m^-2、q₅ = 4.2E−12 W·m⁻²、q₆ = 0 W·m⁻²、q₇ = −4 000 W·m⁻²。对比线性单元，每个材料棒材划分两个单元的结果：q₂ = −4.5E−11、q₃ = −2.3E−11、q₄ = 0、q₅ = 0、q₆ = 1.13E−11，不论量级和数值，二次单元精度均高于线性单元。Excel计算如图1-3-6所示。其中逆矩阵为对传导矩阵前6 × 6矩阵求逆，热通量第五项、第六项计算基于高斯消元法（O79 = −G79*H81，O80 = −G80*H81），温度为逆矩阵与热通量矩阵之积。

图1-3-6　Excel矩阵计算

1.3.2　复合棒材对流状态的有限元计算

实际工程中，纯粹绝热状态是不存在的，往往伴随对流、辐射状态。在1.3.1小节的基础上，本小节再对模型表面增加空气对流，如图1-3-7所示，其对流换热系数为1 W·m⁻²·℃⁻¹，左侧输入热流4 000 W·m⁻²，右侧保持温度为80℃，求稳态条件下棒材的温度分布。

图1-3-7　复合棒材对流模型

1．伽辽金法求解稳态热传导分析

依据1.3.1小节的内容，采用伽辽金法求解本模型。为使计算简便，每种材料的棒材仅划分为一个单元，并且采用线性单元计算，同时假设任意截面的温度一致，忽略对流对模型截面的温度影响。

由于线性单元为两个节点，建立插值函数为：

利用插值条件，可得：

求出a_i、b_i，代入插值函数，其中单元长度，则：

依据伽辽金法，传热矩阵中的任意一项为：

　 i,j = 1, 2；P为截面周长

则单元传热矩阵为：

式中前部为传导矩阵项，后部为对流矩阵项。

传热平衡方程为：

对右式第一项进行分步积分，则：

式中第一项为内热源项、第二项为温度梯度项（对应热通量）、第三项为环境对流温度项。为环境温度。

依据单元传热矩阵，则a材料传热矩阵为：

同理，b、c材料热传导矩阵分别为：

整体传热矩阵为（每个材料单元的最后一个节点与下一个材料单元的第一个节点共点）：

本例在传热平衡方程中，由于无内热源项，所以只计算第二项和第三项，本模型仅在a模型入口和c模型出口存在热通量边界条件，则第二项为，第三项对流边界条件为，即

计算可得T₁ = 120.3℃、T₂ = 116.7℃、T₃ = 89.9℃、q_out = 1 910 W·m⁻² 。

2．Steady-State Thermal模块计算

采用与1.3.1小节一样的材料定义和网格划分，边界条件设置如图1-3-8所示。其中选择a模型的左端面定义加载Heat Flux，数值为4 000 W·m⁻²；选择c模型的右端面定义加载Temperature，数值为80℃；选择a、b、c三个模型的三段外圆面定义加载Convention，其中Film Coefficient数值定义为1 W·m⁻²·℃⁻¹，Ambient Temperature数值定义为80℃（环境温度）。

图1-3-8　边界条件定义

计算完成后，分别插入Temperature（Geometry设置为All Bodies）、Temperature2（Geometry设置为a、b模型交界面）、Temperature3（Geometry设置为b、c模型交界面）和Total Heat Flux四项，后处理的计算结果如图1-3-9所示。

图1-3-9　后处理的结果

由图1-3-9可知，系统最高温度为120.62℃，位于模型最左侧，与传热矩阵计算温度结果120.3℃存在 0.25%的误差；a、b模型交界面温度为 117.01℃，与传热矩阵计算温度结果116.7℃存在 0.34%的误差；b、c模型交界面温度为 90.223℃，与传热矩阵计算温度结果 89.9℃存在 0.33%的误差；输出热通量为 1 914.1 W·m⁻²，与整体传热矩阵计算热通量结果1 910 W·m⁻²存在0.21%的误差。与理论计算结果几乎一致。

1.3.3　棒材瞬态热传导的有限元计算

已知一个直径为20 mm，长度为200 mm，密度为2 000 kg·m⁻³的圆棒，如图1-3-10所示。其热传导系数为100 W·m⁻¹·℃⁻¹，比热容为500 J·kg⁻¹·℃⁻¹，左侧迅速升温到100℃，并保持该温度，右侧保持在环境温度为20℃，求瞬态条件（100 s）下棒材的温度分布。

图1-3-10　棒材瞬态计算模型

1．伽辽金法求解瞬态热传导分析

采用伽辽金法求解本模型。为使计算简便，棒材划分为均匀分布的5个单元，相邻两节点相距40 mm，并且采用线性单元计算。本例重点讲解瞬态迭代计算过程，传热矩阵推导不再复述。

对于只有热传导的瞬态分析，传热矩阵包含热传导矩阵和热容矩阵。其中单元传导矩阵为：

整体传导矩阵为：

单元热容矩阵为：

整体热容矩阵为：

对于瞬态热分析，其系统平衡方程为，式中为温度对时间的导数。对于本例分析，模型的首尾温度均为定值，且系统内没有内热源，所以

，，

基于有限差分，，式中为时间步长。则系统平衡方程为：

依据上式列出递推公式：

两边同时乘以热容矩阵的逆矩阵，上式即为，该式即为瞬态分析温度迭代计算公式。

将相应矩阵代入整体平衡方程，则：

因为模型首尾温度恒定，所以对中间节点温度进行求解，则矩阵简化为：

因为上式存在8个未知数，所以调用迭代计算公式，当时，。，则：

迭代计算为（前20步取，后8步取）：

……

　　……

4个节点的瞬态温度分析结果如图1-3-11所示。

图1-3-11　棒材中间4个节点瞬态温度分析结果

由图可知，4 个节点瞬态分析温度趋势均表现为：先逐渐升温，然后保持稳定，这与实际工况表现一致。但是在90～100 s，温度曲线出现拐点，这是因为步长设置过大；同时在0～10 s，3 个节点温度低于 20℃，且部分出现了温度下降趋势，这明显与实际不符，这是因为理论计算中数值迭代采用了恒定热容矩阵。

2．Transient Thermal模块计算

如图 1-3-12 所示，创建 a 材料，其 Density 定义为 2 000 kg·m⁻³，Isotropic Thermal Conductivity定义为100 W·m⁻¹·℃⁻¹，Specific Heat定义为500 J·kg⁻¹·℃⁻¹。

图1-3-12　创建材料并定义材料参数

单击Geometry→Geom\SYS\a，在下方的Details窗口中展开Material选项，将Assignment定义为a，即选择定义的a材料。

在后处理模块中定义截面，以方便与解析解进行对比。新建4个坐标系，如图1-3-13所示。其中Origin项中Define By采用默认的Global Coordinates，Origin X设置为0 m，Origin Y依次设置为0.04 m、0.08 m、0.12 m、0.16 m，Origin Z设置为0 m；Orientation About Principal Axis项中Axis设置为Z，Define By设置为Global Y Axis。即X轴采用默认方向，Z轴采用全局坐标系的Y轴方向。这是因为后续截面定义必须基于局部坐标系的XY平面。

图1-3-13　新建坐标系

创建截面如图1-3-14所示。在Construction Geometry处分别新建4个截面（Surface、Surface 2、Surface 3和Surface 4），在每个截面的Coordinate System项中分别选择之前定义的Coordinate System、Coordinate System 2、Coordinate System 3、Coordinate System 4。具体操作参见《ANSYS Workbench有限元分析实例详解（静力学）》。

图1-3-14　创建截面

Mesh设置保持默认。

边界条件及分析设置如图1-3-15所示。其中选择模型的左端面定义加载Temperature，数值为100℃；选择模型的右端面定义加载Temperature 2，数值为20℃。单击Analysis Settings，在下方的Details窗口中将Step End Time设置为100 s，Auto Time Stepping设置为Program Controlled（默认），Time Integration设置为On（默认）。具体操作参见《ANSYS Workbench有限元分析实例详解（动力学）》^[4]。

[4]　 《ANSYS Workbench有限元分析实例详解（动力学）》，周炬、苏金英著，ISBN：978-7-115-51065-5。本书后面提到的与此同名的书均是指的同一本。

图1-3-15　边界条件及分析设置定义

计算完成后，分别插入Temperature 2、Temperature 3、Temperature 4和Temperature 5（Scoping Method设置为Surface，Surface分别选择之前定义的截面Surface、Surface 2、Surface 3、Surface 4，Display Time设置为100 s）四项，后处理的计算结果如图1-3-16所示。

图1-3-16　瞬态后处理的结果

由图可知，在100 s时刻，对应解析解的4个节点温度分别为80.962℃、63.101℃、47.122℃、32.996℃，其解析解温度为80.328℃、66.638℃、50.527℃、32.397℃，其结果相差无几，而中间两点误差较大的原因是解析解在90～100 s步长设置过大。

4个节点的瞬态温度分析结果如图1-3-17所示。对比图1-3-11，其趋势均表现为：先逐渐升温，然后保持稳定，而图1-3-11所示的解析计算过程中出现的低于设定温度和计算曲线有拐点的现象在软件计算结果中均未出现。说明软件计算具有良好的鲁棒性。

图1-3-17　棒材中间4个节点瞬态温度分析结果

3．Steady-State Thermal模块计算

对该模型采用稳态分析。采用与瞬态分析一致的前处理和边界条件设置，因为稳态分析不考虑时间效应，所以仅将Analysis Settings中的Step End Time修改为1 s。具体原理参见《ANSYS Workbench有限元分析实例详解（静力学）》和《ANSYS Workbench有限元分析实例详解（动力学）》。

计算完成后，依然查看基于截面Surface、Surface 2、Surface 3、Surface 4的Temperature，计算结果如图1-3-18所示。

图1-3-18　稳态后处理的结果

由图可知，在稳态分析中对应 4 个节点的温度分别为84.001℃、68.001℃、52.001℃、36.001℃，与瞬态分析计算温度80.962℃、63.101℃、47.122℃、32.996℃相比，稳态分析计算结果为温度（一定过程时间后）不随时间变化的稳定状态，整体温度表现未必均匀一致；而瞬态分析计算结果为系统稳定过程中的整体表现（实时状态）。