热物理过程的数值模拟-计算传热学2

02 数值传热学的基本原理
控制方程
控制方程
数值传热学的核心是求解控制方 程，这些方程描述了热量传递过 程中的物理规律。
偏微分方程
控制方程通常以偏微分方程的形 式给出，包含了温度、时间、空 间等变量的变化关系。
初始条件和边界条
件
为了求解控制方程，需要给出初 始条件和边界条件，这些条件限 定了问题的解的范围。
详细描述
传热过程模拟是数值传热学的另一重要应用，通过建立传热过程的数学模型，可以模拟物体内部的温 度分布和热量传递过程。这对于能源、化工、电子等领域中的热工设备设计和优化具有重要意义。
04 数值传热学面临的挑战与 解决方案
计算精度与稳定性问题
总结词
计算精度和稳定性是数值传热学中的核心问题，直接关系到模拟结果的准确性和可靠性。
详细描述
多尺度问题要求数值方法能够捕捉到不同尺度的物理现象，并准确地将它们联系起来。 这需要发展具有多尺度分辨率的数值方法，如多重网格法、谱方法和自适应网格法等。
非线性问题
总结词
非线性问题在传热过程中广泛存在，如 流动、相变和化学反应等，给数值模拟 带来很大难度。
VS
详细描述
非线性问题需要数值方法能够处理高度非 线性的物理方程，并能够准确地捕捉到非 线性现象。这需要发展高效的数值算法， 如有限元法和有限体积法等，同时还需要 考虑非线性问题的特殊性质，如初始条件 和边界条件等。
02
它涉及传热学的基本原理、数学 建模、数值计算和计算机技术等 多个领域，是计算流体动力学和 计算传热学的重要组成部分。
数值传热学的重要性
随着科技的发展，传热问题在能源、 环境、航空航天、化工等领域越来越 突出，数值传热学的应用也越来越广 泛。

热传导的基本原理与计算方法热传导是指热量从高温区向低温区传递的过程。

它是热力学的一种基本现象，广泛应用于物理学、化学、材料科学等领域。

热传导研究的是物质中热量的传导机制、热传导的速率和规律以及如何控制和改变热传导过程。

一、热传导的基本原理在物理学中，热量的传导可以用热传导定律来描述，即热传导的速率与热差成正比，与导热系数和传热面积成反比。

物质温度较高的区域传递给相邻温度较低的区域，热量的传导是靠原子、分子、电子等的热运动完成的。

这些粒子在物质内做无规则的振动、流动，高温区的热粒子向低温区运动，直到它们的热平衡达到。

热传导的基本原理可以用一维热传导方程来描述：$$\frac{\partial T}{\partial t}=\alpha\frac{\partial^2 T}{\partialx^2}.$$其中，T代表温度，x代表长度，t代表时间，α代表物质的导热系数。

方程的右侧表示温度梯度，表示热量的传递速度。

二、计算热传导的基本方法由于热传导过程的复杂性，通过简单的数学方程来计算热传导的速率是不可能的。

因此，人们开发了许多传热学模型和计算方法。

这些方法主要可以分为两种：一种是基于传热学原理和模型计算的解析解，另一种是基于数值方法求解的计算机模拟。

1. 解析解法解析解法是指根据物理模型和数学方程分析热传导的过程，得到解析解的方法。

这种方法的优点是计算结果精确，适用于简单的热传导问题，如一维热传导、恒定温差热传导等。

解析解法的缺点是只能用于特定情况下的计算，不适用于复杂的三维热传导问题。

2. 数值模拟法数值模拟法是指利用数字计算机来模拟热传导过程，在计算机上求解热传导方程。

这种方法的优点是可以模拟任意形状复杂的热传导问题，适用范围广，计算结果较为准确。

数值模拟法的缺点是需要高性能计算机进行计算，耗费时间和资源较多。

三、热传导应用范围热传导的应用范围非常广泛，涉及物理、化学、材料等多个领域。

在工程领域，热传导的应用与产品的保温、散热、冷却、加热等相关。

传热学数值模拟实例教程王志军编著邓权威河南理工大学二〇〇九年十二月前言一、实验说明导热问题实际上就是对导热微分方程（能量方程）在规定的定解条件下进行求解，而对流问题除了对能量方程进行求解外，往往还需对质量守恒方程以及动量方程进行求解。

对于少数几何形状以及边界条件简单的问题能获得分析解，但对于大多数工程技术中遇到的许多几何形状或边界条件复杂的导热对流问题，数学上还无法得除其分析解。

另一方面，在近几十年中，随着计算机技术的迅速发展，数值模拟技术得到了飞速的发展，其中CFD （计算流体力学）能解决流体流动，传热传质等很多工程问题，因而发展非常快。

Fluent 作为目前国际上最流行的商用CFD软件之一，在美国和中国的市场占有率都超过60%。

只要涉及到流体、热传递以及化学方法等问题都可以用Fluent进行求解。

它具有丰富的物理模型、先进的数值方法以及强大的前后处理功能，在航空航天、汽车设计、石油天然气、消防火灾、环境分析等方面都有着广泛的应用。

本模拟实例库主要是运用成熟的Fluent软件对传热学的一些简单问题进行数值求解，主要包括一维稳态导热问题的求解，二维多热源的稳态导热问题，二维方腔内自然对流和混合对流，管内强制对流换热问题的数值模拟。

模拟实验的目的在于是为同学们提供一个形象直观而又生动的工具，为本科传热学的学习提供一个新的视角，使传热学的学习从抽象的理论中解放出来，变得直接而有主动，增强他们学习的兴趣与动力，从枯燥的灌输中解放出来。

另一方面数值模拟还能加深学生对基本概念、基本规律的理解。

杨世铭说：“传热学课程的教学应当从以往的单纯地为后续专业课服务而转变到着重培养学生的素质与能力方面来。

通过将CFD数值模拟方法渗透到传热学的本科实验中，为培养学生的素质与能力提供一个强有力的工具，最终促进学生创新能力和应用能力的全面提升。

二、Fluent软件简介Fluent软件是美国Fluent公司开发的通用CFD流场计算分析软件，囊括了Fluent Dynamic International、比利时Polyflow和Fluent Dynamic International（FDI）的全部技术力量（前者是公认的粘弹性和聚合物流动模拟方面占领先地位的公司，而后者是基于有限元方法CFD软件方面领先的公司）。

t f (x, y, z)
t f (x, y, z, )
2. 等温线，等温面
1) 定义：同一瞬间温度相等的各点连成的线或面称为 等温线(Isotherm)或等温面(Isothermal surface)。
5/41
2）特点：
传热学 Heat Transfer 第5版
(1)等温线(面)不能相交(同一点不可能有两个温度)；
(1768-1830)
9/41
传热学 Heat Transfer 第5版
1. 导热基本定律的文字表达
在导热现象中，单位时间内通过给定截面的热量， 正比于垂直于该截面方向上的温度变化率和截面面 积，方向与温度梯度相反。
2. 导热基本定律的数学表达
q gradt t n
A
Φ
c
a c
称为热扩散率(Thermal diffusivity)
或导温系数，单位：m2/s，是物性参数；
2.λ=constant 并且t x 2

2t y 2

2t z 2
)
a2t
Laplace算子
28/41
传热学 Heat Transfer 第5版
4/41
传热学 Heat Transfer 第5版
按温度场随空间与时间的变化特性，可以区分为：
稳态温度场 t f (x, y, z) 非稳态温度场
t f (x, y, z, )
一维温度场 二维温度场 三维温度场
t f (x)
t f (x, )
t f (x, y)
t f (x, y, )
传热学 Heat Transfer 第5版
代入能量平衡式， （1）＋（2）＝（3） 得导热微分方程的基本形式

热物理过程的数值模拟（计算传热学）教学大纲（热工类各专业及机械类动力机械专业研究生适用）（30学时）重庆大学热工教研室二零零零年七月热物理过程的数值模拟（计算传热学）教学大纲第一章绪论1、热物理过程的预测方法及特点：预测的实质。

理论分析、实验研究与数值计算。

预测方法的选择。

2、计算传热学的发展；国内外计算传热的有关活动、研究内容、当前的发展方向。

第二章热物理过程的数学描述1、控制微分方程：微分方程的意义。

连续性方程、化学组分议程、动量议程、能量议程、湍流流动的时间平均方程、湍流动能方程、通用微分方程。

2、边界条件和初始条件：第一、二、三、四类边界条件、初始条件、拟非稳态的概念。

3、控制方程的简化：恰当坐标系、自变量和因变量的变换、无量纲化。

第三章离散化方法1、离散化的概念：数值方法的基本思想、因变量的离散化、区域的离散化。

2、空间区域的离散化方法：空间区域离散人及几何要素、内节点法和外节点法。

3、推导离散化方程的方法：泰勒级数展开法、变分法、加权残值法、控制容积积分法、控制容积平衡法。

4、二个指导原则和四项基本法则：解的物理真实性和总量平衡、控制容积界面上的相容性、正系数法则、源项的负斜率线性化、邻点系数和法则。

第四章热传导问题的数值解灶1、一维稳态热传导：基本方程、网络间距、界面导热系数、源项的线性化、非线性的处理、边界条件的处理、线性代数方程组的求解（TDMA）。

2、一维非稳态热传导：离散化方程的一般形式、显格式、全隐格式、C-N格式。

3、离散方程的性质：相容性、收敛性、稳定性、代数方程组的求解、稳定性与解的物理真实性。

4、多维稳态热传导：离散化方程、三种坐标系中系数的通用表达式、边界条件的处理、代数方程组的求解方法（点迭代法、块迭代法）、迭代法的收敛性及其改善。

5、多维非稳态热传导：离散化方程的形式、稳定性、交替方向迭代（ADI）。

第五章对流—扩散方程的差分格式1、一维稳态对流—扩散：中心差分、上风差风、指数格式、混合格式、乘方格式、通用化格式。

Tu-B管
4.3 管内流动沸腾传热分析
1. 管内受迫对流沸腾（两相流）特征 流型演化 泡状流：汽泡小而分散，并逐渐增多
块状流（栓塞流、弹状流）：小汽泡合形成大 汽泡
环状流：大汽泡进一步合并，在管中形成汽芯， 把液体排挤到壁面上。热量主要以对流方式通 过液膜，汽化过程主要发生在液汽界面上 单相流： 液层全部蒸发，蒸汽单相流 传热机理
成核条件
① 具有一定的过热度(twts) ② 壁面凹腔处，且r>Rmin
③ (twts)越大，则Rmin越小
Rmin 2 Ts 2 Ts r v t r v (tw ts )
壁面会形成更多的汽化核心N(r)
4.1 沸 腾 传 热
4. 汽泡动力学简介

研究汽泡的产生、长大、脱离过程
仅讨论蒸汽的膜状凝结
4.4 凝 结 传 热
2. 膜状凝结传热计算
（1）液膜流动状态及其对传热的影响 层流 紊流 紊流核心层
层流底层
影响表面传热系数的主要因素 液膜流态判定准则数：Re数 临界雷诺数Rec≈ 1600(1800) 若竖壁液膜已转变为紊流，则
h hl xc x ht 1 c L L
4.2 大容器沸腾传热计算
沸腾传热计算特点

牛顿冷却公式仍然适用
q h tw ts ht

影响因素多
h f [t , g ( 1 v ), r, , c p , ,, Cwl ]
与沸腾液体及表面材料有关的系数

计算公式分歧较大
4.2 大容器沸腾传热计算
1. 大容器饱和核态沸腾 米海耶夫计算式
hr
或
h hc 0.75hr

然 而 利 用 商 用 软 件 简化 这 一过 程 的 分 析 ，以下 是 简化 后 的过 程 分析 步骤 ：１）定 义单 元类 型 ；２）定义 所 测 材 料 的性 能 参 数 ；３）创 建 几 何 模 型 ，划 分 网格 ；４） 施 加相 应 的条 件 ，求解 ；５）查 看求 解 结果 ；６）利用 数 据 对所 求 问题 进行 分 析 。
科学技术 应用
数值模 拟在传 热学 中的运 用

徐 剑波 贵 州 省节 能监 测 中心 ，贵 州贵 阳 ５５０００１
摘 要 在传 热学 中 ，边 界 效应 一 直是 讨论 的重 要 问题 ，由于边 界 效应 的 不 确定 性 ，其 对 计 算和 实验 精 度 都会 存 在一定的影响。对影响边界效应因素的讨论和分析可以为计算和实验结果提供参考。本文在常功率二维非稳 态传热 的条件下 ，根据二维常功率平面热源法测量材料导热系数的基本原理 ，建立了考虑边界效应 的二 维传热模型 ，并进 行 数值 求 解 。讨论 了对 流 换 热 系数 、模 型 尺 寸大 小 、无量 纲数 毕 渥数对 边 界效 应 的影 响 。边界 效应 随 着对对 流换 热 系数的增 大而增强 ，随着模型尺寸的增大而减小。边界效应随着无量纲数毕渥数的变化分为两种情况 ：１）当导热 系数 和长度 不变 时 ，随毕 渥数 的增 大 而增 大 。２）当导 热 系数 和 对 流换 热 系数 不变 时 ，随毕 渥数 的增 大 而减 小。 关键 词 数 值模 拟 ；边界 效应 ；非稳 态 中图 分类 号 ＴＱＯ 文 献标 识码 Ａ 文章 编号 ２０９５—６ ３６３（２０１６）０５一Ｏ１ ３８一Ｏ１
伴 随 计 算 机 技 术 的 进 步 ，我们 之 前很 多遗 留 需 要 解 决 的传 热 问题 可 以用 数值 求解 的方 式 进 行 模 拟 解 决 。 数值 模 拟解 决 传热 学 问题 的基 本 的 思路 可 以大 概 总结 为 如 下 ：把 空间 、时 间坐 标 系 中 的温度 场 用根 据情 况 设 定 数量 的离散 点 上 的值 的集 合 来 替代 ，采 用计 算机 模 拟 求 解按 ～ 定方 程 建立 起 来 的这些 值 的循 环代 替方程 ，来 获 得 离散 点 上我 们 需要 的温 度场 的值 。这 一基 本 求解 思 路 描述 成 以 下几 点 ：１）根 据 实 际情 况 建 立 控 制 方 程及 定 解条 件 ；２）确 定研 究本 体 的相 关 节 点 ；３）建立 节 点物 理量 的相关代数方程 ；４）通过计算机设立温度场 的符 合研究对象的迭代初值 ；５）求解相关的代数方程组 ；６） 根据 计 算 出 的解 ，分析 我 们研 究所 要 达 到 的 目的 。这 ６ 个步 骤就 是 导热 问题 数值 求解 的基本 步骤 。

❖ 假设各层之间接触良好，可以近似地认 为接合面上各处的温度相等
t2
t3 t4
❖ 边界条件： x 0
n
x i i1
t t1 t tn1
❖ 热阻：
r1
1 1
,
, rn
n n
t1
t2
t3
t4
三层平壁的稳态导热
由热阻分析法：
q
t1 tn1
n
ri
i 1
t1 tn1
n i i1 i
问：现在已经知道了q，如何计算其中第 i 层的右侧壁温？
对上述方程(a)积分两次:
第一次积分
第二次积分
r
dt dr
c1
t c1 ln r c2
tw1 c1 ln r1 c2 ; tw2 c1 ln r2 c2
应用边界条件 获得两个系数
c1
tw2 tw1 ln( r2 r1)
;
c2
tw1
(tw2
tw1)
ln r1 ln( r2 r1)
t
对于变截面肋片来讲，由于从导热微分方程求得的肋片 散热量计算公式相当复杂，因此，人们仿照等截面直肋。 利用肋片效率曲线来计算方便多了，书中图2－14和2－15 分别给出了三角形直肋和矩形剖面环肋的效率曲线。
图 2－14
图 2－15
4. 通过接触面的导热
实际固体表面不是理想平整的，所以两固体表面直接接触的界 面容易出现点接触，或者只是部分的而不是完全的和平整的面 接触 —— 给导热带来额外的热阻 —— 接触热阻 (Thermal contact resistance）
h2
ql
1
tf1 tf 2 1 ln r2
1

数值传热学数值传热又称计算传热，是传热学与数值方法相结合的一门交叉学科，它采用数值方法描述流动和传热问题的控制方程，并用计算机求解。

数值换热，其基本思想是将原始坐标在空间和时间上连续的物理量场（如速度场、温度场和浓度场等），用一系列有限个离散点上的数值来代替，通过一定的原理建立离散点变量值之间的关系代数方程（称为离散方程）。

通过求解所建立的代数方程组，得到求解变量的近似值。

1简介数值传热学（numerical heat transfer）数值传热学，又称计算传热学，是指对描写流动与传热问题的控制方程采用数值方法，通过计算机求解的一门传热学与数值方法相结合的交叉学科。

数值传热学的基本思想是把原来在空间与时间坐标中连续的物理量的场（如速度场，温度场，浓度场等），用一系列有限个离散点上的值的集合来代替，通过一定的原则建立起这些离散点变量值之间关系的代数方程（称为离散方程）。

求解所建立起来的代数方程已获得求解变量的近似值。

2发展简史数值传热学，主要由20世纪中叶，S.V. Patankar和D.B.Spalding 等人在总结前人的研究基础上所提出。

E.M.Sparrow对数值传热学的发展也起到了一定的促进作用。

国内比较知名的学者是陶文铨教授。

陶文铨3研究方法数值传热学常用的数值方法1.有限差分法历史上最早采用的数值方法，对简单几何形状中的流动与换热问题最容易实施的数值方法。

其基本点是：将求解区域中用于坐标轴平行的一系列网格的交点所组成的点的集合来代替，在每个节点上，将控制方程中每一个导数用相应的差分表达式来代替，从而在每个节点上，形成一个代数方程，每个方程中包括了本节点及其附近一些节点上的未知值，求解这些代数方程就获得了所需的数值解。

2.有限容积法将所计算的区域划分成一系列控制容积划分为一系列控制容积，每个控制容积都有一个节点做代表。

通过将守恒型的控制方程对控制容积坐积分导出离散方程。

在导出过程中，需要对界面上的被求函数本身及其一阶导数的构成做出假定，是目前流动与换热问题的数值计算中应用最广的一种方法。

钢： q tw 1 tw 2 3 .4 6 3 1 0 0 1 .4 0 1 6 5 W 0 m 2 0 .05
铬砖：
q tw 1 tw 2 2 .3 3 2 1 0 0 9 .2 0 1 8 3 W 0 m 2 0 .05
硅藻土砖： q t w 1 t w 2 0 .2 3 4 1 0 2 0 9 .6 0 1 8 2 W 0 m 2 0 .05
空调房间的隔热设计， 墙体内表面温度保持恒定， 外表面周期变化
干热性气候区，白天在太阳辐射作用下， 墙体外表面温度高于墙体内表面温度，到太 阳下山直至夜间，又低于内表面温度。
在非稳态导热中，由于温度不稳定，围护机构不断吸收或释放热 量，即材料在导热的同时还伴随着蓄热量的变化，这是非稳态导 热区别与稳态导热的重要特点。非稳态导热计算极其繁琐，一般 可采取简化模型进行计算。
力最强。同时，黑体也能发射一切波长的辐射。
黑体辐射的控制方程： 斯蒂芬波尔兹曼定律
Stefan-Boltzmann 定律 AT4 Eb T4
真实物体则为：b AT4 ET4
式中： 表示黑体的辐射系数：5.68*10-
8w/m2·k4； 表示黑度，又称发射率，是实际物
体的辐射系数与黑体的辐射系数之比。T表示物 体表面的绝对温度，K。
(1)材料为铜，λ =375 w/(mK );
(2)材料为钢， λ =36.4 w/(mK );
(3)材料为铬砖， λ =2.32 w/(mK )；
(4)材料为铬藻土砖， λ =0.242 w/(mK )。
解：一维稳态导热公式有：
铜：
q tw 1 tw 2 3 7 3 5 0 1 0 1 .5 1 0 6 W 0 m 2 0 .05

四、非线笥问题迭代式解法的收敛性每一层次上满足迭代法求解的收敛条件+相邻次间代数方程的系数变化不太大（亦即未知量的变化不太大J多数情形下非线性问题迭代式解法是可以收敛的）。

使相邻两层次间未知量变化不太大的措施：1欠松弛迭代常用逐次欠弛线迭法（SLUR）：一组临时系数下逐线迭代求解+对所得的解施以欠松弛，再用欠松弛后的解去计算新的系数，常数，以进入下一层次的迭代。

实施：常把欠松弛处理纳入迭代过程，而不是在一个层次迭代完成后再行欠松弛。

.（n 1）川）.'a n bt n bt p =t p （t p ）ap（先）t p n1） = 7an b t n b b （1一•）屯t p n）co oa'p t p n 9 、a n bt n b b'a'p -a^ ■, b' = b （^ ）（a p ）t p n），用交替方向线迭代法求解这一方程，就实现了SLUR的迭代求解。

为一般化起见，上式中t n b上没有标以迭代层次的符号（J, GS时不相同）。

2、采用拟非稳态法前面已指出，稳态问题的迭代解法与非稳态问题的步进法十分相似。

对于非线性稳态问题，从代数方程的一组临时系数进入到另一组临时系数亦好象非稳态问题前进了一个时间层，非稳态问题的物理特性：系数热惯性越大（a； = PM v/也I ），温度变化越慢，仿此，对稳态非线性问题，可在离散方程中加入拟非稳态项，以减小未知量托两个层次间的变化，即由（=a n b -S p：V）t p n。

= Ua n bt n b b=（3a n b - S p：V a；）t p n。

=二a n bt n b b a p tfZa n bt n b - b - a；t p n）（n 1）t p oEa n b -S p心V +a p一直进行到t p,t n b收敛，虚拟时间步的大小通过计算实践确定。

3、采用Jacobi点迭代法中止迭代的判据（该层次迭代）除前述变化率判据外，还可以规定迭代的轮数，例如规定进行4-6次ADI线迭代就结束该层次上的计算。

第二章热传导一、名词解释1．温度场：某一瞬间物体内各点温度分布的总称。

一般来说，它是空间坐标和时间坐标的函数。

2．等温面(线)：由物体内温度相同的点所连成的面（或线）。

3．温度梯度：在等温面法线方向上最大温度变化率。

4．热导率：物性参数，热流密度矢量与温度降度的比值，数值上等于1 K／m的温度梯度作用下产生的热流密度。

热导率是材料固有的热物理性质，表示物质导热能力的大小。

5．导温系数：材料传播温度变化能力大小的指标。

6．稳态导热：物体中各点温度不随时间而改变的导热过程。

7．非稳态导热：物体中各点温度随时间而改变的导热过程。

8．傅里叶定律：在各向同性均质的导热物体中，通过某导热面积的热流密度正比于该导热面法向温度变化率。

9．保温(隔热)材料：λ≤0.12 W/(m·K)(平均温度不高于350℃时)的材料。

10．肋效率：肋片实际散热量与肋片最大可能散热量之比。

11．接触热阻：材料表面由于存在一定的粗糙度使相接触的表面之间存在间隙，给导热过程带来额外热阻。

12．定解条件(单值性条件)：使微分方程获得适合某一特定问题解的附加条件，包括初始条件和边界条件。

二、填空题1．导热基本定律是_____定律，可表述为。

（傅立叶，）2．非稳态导热时，物体内的_____场和热流量随_____而变化。

（温度，时间）3．导温系数的表达式为_____，单位是_____，其物理意义为_____。

（a=λ/cρ,m2/s,材料传播温度变化能力的指标）4．肋效率的定义为_______。

（肋片实际散热量与肋片最大可能散热量之比。

）5．按照导热机理，水的气、液、固三种状态中_______态下的导热系数最小。

（气）6．一般，材料的导热系数与_____和_____有关。

（种类，温度）7．保温材料是指_____的材料.（λ≤0.12 W/(m·K)(平均温度不高于350℃时)）8．已知材料的导热系数与温度的关系为λ＝λ0(1+bt)，当材料两侧壁温分别为t1、t2时，其平均导热系数可取下的导热系数。

《计算传热学》学期作业多尺度模拟方法概述摘要：本文简单介绍多尺度模拟的思想，应用及存在的问题。

关键词：数值模拟；多尺度模拟世界的本质是多尺度的，在不同的尺度下物质表现出不同的特征。

如流体在分子尺度下表现为离散的不确定的粒子，而在宏观尺度下表现为连续的确定性的介质。

在不同的时间和空间尺度下由于其尺度特性的不同，往往所采用的方法也不同,如图1[1]所示。

图1 各种空间时间尺度下适用的模拟方法文献[2]利用Kn数来鉴定何种特征尺度下流体流动适合用何种方法。

Kn数的物理意义是分子平均自由程与特征长度的比值。

Kn<10-3，流动符合连续介质假设，可用N-S方程；10-3<Kn<10-1，边界是滑移边界，速度和温度有跳跃，控制方程为N-S方程；10-1<Kn<10，过渡流动，N-S方程不再适用，可用格子Boltzmann方法；Kn>10，分子流动，可用分子动力学模拟方法。

模拟方法大致可分为宏观方法，介观方法，微观方法。

宏观方法即流动符合连续介质假设，传热的空间尺度和时间尺度符合傅立叶导热定律；微观方法是从分子运动碰撞理论来建立方程；介观方法是介于微观方法和宏观方法之间。

这三种方法各有优缺点。

宏观方法不能揭示微观的物理现象，但是方法成熟，应用方便。

微观或介观方法更适合描述极端尺度的物理现象，但是计算量巨大，方法不成熟，工程应用极少。

如果在采用宏观方法的过程中，可将微观尺度的信息带入，建立一种微观——宏观耦合的多尺度模拟方法可以结合两者的优点，又可以削弱两者的缺点。

多尺度问题表现[3]为: 已知一个模型的宏观描述, 但这种宏观描述在某些局部区域失效, 必须要用低尺度微观非线性描述代替。

模型的微观特性既受制于宏观上的作用因素, 又可能显著影响宏观性能。

但微观结构, 性能与状态何时、以怎样的途径去影响宏观性能并不清楚。

假定一个给定系统的微观行为可以使用微观模型变量u表示, 系统的宏观行为用宏观模型变量U表示, 那么宏观模型变量U与微观模型变量u可以通过压缩乘子Q或者重构算子R联系起来：U=Qu RU=u多尺度模拟的难度在于两种尺度的耦合，即如何建模。

图一
二、 计算原理
本次上机模拟实验选等温边界条件。墙角是中心对称的，所以取其 1/4 研究， 方便计算机计算。上机模拟选取网格划分方法同实际实验，可根据热平衡法列 出节点方程，各方向导入单元体的热量之和为零。该边界条件下共有四类节点，
Hale Waihona Puke 内节点、内边界点、外边界点和绝热边界点。
图二
四种节点的节点方程简化如下：
eps=1; temp=A[i][j]; A[i][j]=(A[i-1][j]+A[i+1][j]+A[i][j-1]+A[i][j+1])/4; eps=A[i][j]-temp;
} eps=1; temp=A[0][j]; A[0][j]=(A[0][j-1]+A[0][j+1]+2*A[1][j])/4; eps=A[i][15]-temp; } //计算墙体外表面导热量 q_out=0; for(i=1;i<11;i++) q_out=q_out+A[i][0]-A[i][1]; for(j=1;j<15;j++) q_out=q_out+A[11][j]-A[10][j]; q_out=q_out+(A[0][0]-A[0][1]+A[11][15]-A[10][15])/2; q_out=q_out*0.53; //计算墙体内表面导热量 q_in=0; for(i=1;i<7;i++) q_in=q_in+A[i][4]-A[i][5]; for(j=5;j<15;j++) q_in=q_in+A[7][j]-A[6][j]; q_in=q_in+(A[0][4]-A[0][5]+A[7][15]-A[6][15])/2; q_in=q_in*0.53; //计算平均导热量和相对误差 q=(q_in+q_out)/2; eps=abs(q_in-q_out); } //输出结果 for(i=11;i>5;i--) { for(j=0;j<16;j++) out<<setw(8)<<setprecision(2)<<A[i][j]<<" "; out<<endl; } for(i=5;i>=0;i--) { for(j=0;j<6;j++) out<<setw(8)<<setprecision(2)<<A[i][j]<<" "; out<<endl; } out<<"墙体内表面导热量q_in="<<q_in<<"\n"; out<<"墙体外表面导热量q_out="<<q_out<<"\n"; out<<"墙体平均导热量q="<<q<<"\n"; return 0; }

假设单管长度为l，圆筒壁的外半 径小于长度的1/10。 一维、稳态、无内热源、常物性：
d dt (r ) 0 dr dr
(a)
r r1时 t t w1 第一类边界条件： r r2 时 t t w 2
对上述方程(a)积分两次:
第一次积分
第二次积分 应用边界条件
dt r c1 t c1 ln r c2 dr
直接积分，得：
t t1
x
dt c1 t c1 x c2 dx
t2 t1 c 带入边界条件： 1 c2 t1
t2 o

t2 t1 t x t1 带入Fourier 定律 dt t2 t1 dx
t w1 c1 ln r1 c2 ; t w 2 c1 ln r2 c2
t w 2 t w1 ; c1 ln(r2 r1 )
获得两个系数
ln r1 c2 t w1 (t w 2 t w1 ) ln(r2 r1 )
将系数带入第二次积分结果
t 2 t1 t t1 ln(r r1 ) ln(r2 r1 )
a 几何条件：单层平板； b 物理条件：、c、 已知；无内热源 c 时间条件： 稳态导热 : t 0 d 边界条件：第一类
o

x
根据上面的条件可得：
t t c ( ) Φ x x
控制 方程
d 2t dx
2
0
边界 条件
x 0, t t w1 第一类边条： x , t t w2
通过球壳的导热自己推导
5 其它变面积或变导热系数问题 求解导热问题的主要途径分两步： (1) 求解导热微分方程，获得温度场； (2) 根据Fourier定律和已获得的温度场计算热流量； 对于稳态、无内热源、第一类边界条件下的一维导热 问题，可以不通过温度场而直接获得热流量。此时， 一维Fourier定律：

三三
三三三三三三三三三 三三
三三 三三
三三三三三三三
三三
三三三三三三三三
三三
三三三三三三三三三三三
18
第2章 稳态热传导
2.1 典型一维稳态导热问题的分析解
2.3.1 通过平壁的导热：
一维、稳态、常物性、无内热源情况，考察平板的导热情况。
c t


x
t x
t x
n
中，gradt表示空间某点的温度梯度，
n表示通过该点的等温线上的法向单位矢量，温度升高的方向。
利用等温线和热流线来定量且形 象地表述一个导热过程： 等温线表示热流梯度，而热流线 是与等温线处处垂直的一组曲线， 通过平面上任一点的热流线与该 点的热流密度相切。 相邻两条热流线之间所传递的热 流量处处相等，相当于构成了一 个热流通道。 该方法用于现代工程软件应用。
2.类似于非导电固体;(倾向于此观点)
2
第2章 稳态热传导
等温场(temperature field)：
温度场：物体中存在温度的场。 温度分布：各时刻物体中各点温度所组成的集合
分类：
稳态温度场：物体中各点温度不随时间而变。 t f x, y, z 瞬态温度场：物体中各点温度随时间变化。 t f x, y, z,
几何条件： 说明导热体的几何形状(平壁或圆筒壁)和大小(厚度、直径等)
物理条件:
说明导热体的物理特征如：物性参数λ、c 和 r 的数值，是否 随温度变化；有无内热源、大小和分布；是否各向同性 初始(时间)条件: 说明在时间上导热过程进行的特点 稳态导热过程不需要时间条件 — 与时间无关 对非稳态导热过程应给出过程开始时刻导热体内的温度分布
疏密可直观反映出不同区域温度热流密度的相对大小。

《传热学》上机实践大作业二维导热物体温度场的数值模拟 能动A02 赵凯 2010031134一、物理问题有一个用砖砌成的长方形截面的冷空气通道，其截面尺寸如下图所示，假设在垂直于纸面方向上冷空气及砖墙的温度变化很小，可以近似地予以忽略。

在下列两种情况下试计算：砖墙横截面上的温度分布；垂直于纸面方向的每米长度上通过砖墙的导热量。

第一种情况：内外壁分别均匀地维持在0C ︒及30C ︒； 第二种情况：内外表面均为第三类边界条件，且已知：Km K m W h C t Km W h C t •=•=︒=•=︒=∞∞/35.0/93.3,10/35.10,30222211λ砖墙导热系数二、数学描写1、控制方程该问题为无内热源的二维稳态导热问题，因此控制方程为导热微分方程：02222=∂∂+∂∂y t x t 2、边界条件该问题中，导热物体在x 方向上，y 方向上都是对称的，因此可以只取其中的四分之一部分作为研究对象，其他部分情况完全相同，如下图所示：对于上图所示各边界：边界1：由对称性可知：其为绝热边界，即0=w q 。

边界2：第一种情况：其为等温边界，满足第一类边界条件。

即： C t w ︒=0第二种情况：其为对流边界，满足第三类边界条件。

即:)()(2f w w w t t h ntq -=∂∂-=λ 边界3：第一种情况：其为等温边界，满足第一类边界条件。

即： C t w ︒=30 第二种情况：其为对流边界，满足第三类边界条件。

即：)()(1f w w w t t h ntq -=∂∂-=λ三、方程离散如下图所示，用一系列与坐标轴平行的间隔10cm 的网格线将求解区域划分成子区域。

可将上图所示各节点分成内节点与边界点两类。

分别利用热平衡法列各个节点的代数方程。

第一种情况（等温边界）： 边界点：边界1（绝热边界）：5~2),2(411,11,12,1,=++=+-m t t t t m m m m 11~8),2(411,161,16,15,16=++=+-n t t t t n n n n 边界2（内等温边界）： 7,16~7;7~1,6,0,=====n m n m t n m边界3（外等温边界）： 12,16~2;12~1,1,30,=====n m n m t n m内节点：11~8,15~6;11~2,5~2);(411,1,,1,1,====+++=-+-+n m n m t t t t t n m n m n m n m n m第二种情况（对流边界）： 边界点：边界1（绝热边界）：5~2),2(411,11,12,1,=++=+-m t t t t m m m m11~8),2(411,161,16,15,16=++=+-n t t t t n n n n边界2（内对流边界）：6~1,)2(222111,61,6,5,6=++++=∆∆-+n Bi t Bi t t t t n n n n16~7,)2(2221117,17,18,7,=++++=∆∆-+m Bi t Bi t t t t m m m m边界3（外对流边界）：11~1,)2(2222221,11,1,2,1=++++=∆∆-+n Bi t Bi t t t t n n n n16~2,)2(22222212,112,111,12,=++++=∆∆-+m Bi t Bi t t t t m m m m内角点： )3(22)(21116,67,78,67,57,6+++++=∆∆Bi t Bi t t t t t外角点： )1(222211,112,212,1+++=∆∆Bi t Bi t t t内节点：11~8,15~6;11~2,5~2);(411,1,,1,1,====+++=-+-+n m n m t t t t t n m n m n m n m n m（10,22121==∆=∞∆t t xh Bi λ；30,21212==∆=∞∆t t xh Bi λ）四、编程求解第一种情况（等温边界)：Fortran程序代码如下所示：Program denwengimplicit noneinteger::t1=0integer::t2=30integer m,nreal::t(16,12),ta(16,12),et(16,12)real::epslona=1realfainei,fainei1,fainei2,fainei3,fainei4,fainei5,fai nei6,fainei7realfaiwai,faiwai1,faiwai2,faiwai3,faiwai4,faiwai5 ,faiwai6,faiwai7real pianchado n=1,7t(6,n)=t1end dodo m=7,16t(m,7)=t1end dodo n=1,12t(1,n)=t2end dodo m=2,16t(m,12)=t2end dodo m=2,5do n=1,11t(m,n)=10end doend dodo m=6,16do n=8,11t(m,n)=10end doend doopen(01,file='dengwen.dat')do while(epslona>0.00000001)do m=2,5ta(m,1)=0.25*(2*t(m,2)+t(m-1,1)+t(m+1,1)) end dodo m=2,5do n=2,11ta(m,n)=0.25*(t(m+1,n)+t(m-1,n)+t(m,n+1)+t( m,n-1))end doend dodo m=6,15do n=8,11ta(m,n)=0.25*(t(m+1,n)+t(m-1,n)+t(m,n+1)+t( m,n-1))end doend dodo n=8,11ta(16,n)=0.25*(2*t(15,n)+t(16,n-1)+t(16,n+1)) end dodo n=1,7ta(6,n)=t1end dodo m=7,16ta(m,7)=t1end dodo n=1,12ta(1,n)=t2end dodo m=2,16ta(m,12)=t2end dodo m=1,16do n=1,12et(m,n)=abs(ta(m,n)-t(m,n))end doend doepslona=maxval(et(1:16,1:12))do m=1,16do n=1,12t(m,n)=ta(m,n)end doend doend dofainei1=0.5*lanbuda*t(5,1)fainei3=lanbuda*t(5,8)fainei5=0.5*lanbuda*t(16,8)fainei2=0do n=2,7fainei6=lanbuda*t(5,n)fainei2=fainei2+fainei6end dofainei4=0do m=6,15fainei7=lanbuda*t(m,8)fainei4=fainei4+fainei7end dofainei=4*(fainei1+fainei2+fainei3+fainei4+fai nei5)faiwai1=0.5*lanbuda*(30-t(2,1))faiwai3=lanbuda*(30-t(2,11))faiwai5=0.5*lanbuda*(30-t(16,11))faiwai2=0do n=2,10faiwai6=lanbuda*(30-t(2,n))faiwai2=faiwai2+faiwai6end dofaiwai4=0do m=3,15faiwai7=lanbuda*(30-t(m,11))faiwai4=faiwai4+faiwai7end dofaiwai=4*(faiwai1+faiwai2+faiwai3+faiwai4+ faiwai5)print*,' m n t 'do m=1,16do n=1,12print*, m,n,t(m,n)write(01,*) m,n, t(m,n)end doend dopiancha=abs(fainei-faiwai)/((fainei+faiwai)/2) print*,'内部热流量=',faineiprint*,'外部热流量=',faiwaiprint*,'热平衡偏差=',pianchaend program denweng运行结果如图所示：第二种情况(对流边界）: Fortran程序代码如下所示：program duiliuimplicit noneinteger::t1=10integer::t2=30integer m,nreal::t(16,12),ta(16,12),et(16,12)real::epslona=1real bi1,bi2realfainei,fainei1,fainei2,fainei3,fainei4,fainei5,fai nei6,fainei7realfaiwai,faiwai1,faiwai2,faiwai3,faiwai4,faiwai5 ,faiwai6,faiwai7real pianchabi1=h1*detax/lanbudabi2=h2*detax/lanbudado m=1,16do n=1,12t(m,n)=10end doend doopen(01,file='crs.dat')do while(epslona>0.000000001)do m=2,5ta(m,1)=0.25*(2*t(m,2)+t(m-1,1)+t(m+1,1)) end dodo n=8,11ta(16,n)=0.25*(2*t(15,n)+t(16,n-1)+t(16,n+1)) end dodo n=2,6 ta(6,n)=(2*t(5,n)+t(6,n+1)+t(6,n-1)+2*bi1*t1) /(2*bi1+4)end dodo m=7,15ta(m,7)=(2*t(m,8)+t(m+1,7)+t(m-1,7)+2*bi1* t1)/(2*bi1+4)end dodo n=2,11ta(1,n)=(2*t(2,n)+t(1,n+1)+t(1,n-1)+2*bi2*t2) /(2*bi2+4)end dodo m=2,15ta(m,12)=(2*t(m,11)+t(m+1,12)+t(m-1,12)+2 *bi2*t2)/(2*bi2+4)end dodo m=2,5do n=2,11ta(m,n)=0.25*(t(m+1,n)+t(m-1,n)+t(m,n+1)+t( m,n-1))end doend dodo m=6,15do n=8,11ta(m,n)=0.25*(t(m+1,n)+t(m-1,n)+t(m,n+1)+t( m,n-1))end doend dota(6,7)=(2*t(5,7)+2*t(6,8)+t(7,7)+t(6,6)+2*bi1*t1)/(2*bi1+6)ta(1,12)=(t(2,12)+t(1,11)+2*bi2*t2)/(2*bi2+2) ta(6,1)=(t(5,1)+t(6,2)+bi1*t1)/(bi1+2)ta(16,7)=(t(16,8)+t(15,7)+bi1*t1)/(bi1+2)ta(16,12)=(t(16,11)+t(15,12)+bi2*t2)/(bi2+2) ta(1,1)=( t(2,1)+t(1,2)+bi2*t2)/(bi2+2)do m=1,16do n=1,12et(m,n)=abs(ta(m,n)-t(m,n))end doend doepslona=maxval(et(1:16,1:12))do m=1,16do n=1,12t(m,n)=ta(m,n)end doend doend dofainei1=0.05*h1*(t(6,1)-10)fainei3=0.1*h1*(t(6,7)-10)fainei5=0.05*h1*(t(16,7)-10)fainei2=0do n=2,6fainei6=0.1*h1*(t(6,n)-10)fainei2=fainei2+fainei6end dofainei4=0do m=7,15fainei7=0.05*h1*(t(m,8)-10)fainei4=fainei4+fainei7end dofainei=4*(fainei1+fainei2+fainei3+fainei4+fai nei5)faiwai1=0.05*h2*(30-t(1,1))faiwai3=0.1*h2*(30-t(1,12))faiwai5=0.05*h2*(30-t(16,12))faiwai2=0do n=2,11 faiwai6=0.1*h2*(30-t(1,n))faiwai2=faiwai2+faiwai6end dofaiwai4=0do m=2,15faiwai7=0.1*h2*(30-t(m,12))faiwai4=faiwai4+faiwai7end dofaiwai=4*(faiwai1+faiwai2+faiwai3+faiwai4+ faiwai5)do n=1,12do m=1,16print*, m,n,t(m,n)write(01,*) m,n,t(m,n)end doend dopiancha=abs(fainei-faiwai)/((fainei+faiwai)/2) print*,'内部热流量=',faineiprint*,'外部热流量=',faiwaiprint*,'热平衡偏差=',pianchaclose(01)end program duiliuWORD完整版----可编辑----教育资料分享运行结果如图所示：----完整版学习资料分享----五、结果讨论1,、温度场分布图用以上数值模拟得到的各节点温度绘制温度场分布图。