数值传热学第四章-数值计算-课件

*太太 原原 理理 工工 大大 学学1/72 ThermalThermal第四章第四章 热传导热传导主要包含以下内容：4.1本章的对象4.2 一维稳态热传导 4.3 不稳态一维热传导 4.4 二维与三维问题4.5超松弛与欠松弛 4.6某些几何上的考虑*太太 原原 理理 工工 大大 学学2/72 ThermalThermal4.1 本章的本章的对对象象一、本章研究一、本章研究对对象象 本本章章以以导导热热问问题题为为代代表表，介介绍绍扩扩散散方方程程的的数数值值求求解解法法。将将通通用用微微分分方方程程中中的的对对流流项项略略去去，整整个个方方法的介法的介绍绍将在第五章完成。将在第五章完成。*太太 原原 理理 工工 大大 学学3/72 ThermalThermal二、以二、以导热问题导热问题的数的数值值解作解作为为学学习习起点的原因起点的原因热传导作为物理过程易于理解，而且在数学上的复杂性最小，计算方法也比较成熟；工程流动与换热过程中的不少现象，其控制方程类似于热传导方程。如二维位势流动；常物性流体在直管内的充分发展对流换热；质扩散过程；轴承的润滑流动；某些通过多孔介质的流动。导热问题数值解过程中所采用的一些方法与技巧对于对流问题的数值解也适用。如边界条件的处理、源项的线性化及代数方程组的求解方法等。*太太 原原 理理 工工 大大 学学4/72 ThermalThermal.把热传导用作流体流动计算方案的基本组成部分 的做法有助于理解动量传递与热量传递之间的类 似性（用某种方法把速度与温度相比拟）。本章内容将是流本章内容将是流动动与与换热换热数数值值解的基解的基础础*太太 原原 理理 工工 大大 学学5/72 ThermalThermal4.2 一一维稳态热传导维稳态热传导4.2-1 基本方程基本方程一一维稳态导热问题的控制方程：的控制方程：式中：式中：相相应应的离散化方程：的离散化方程：*太太 原原 理理 工工 大大 学学6/72 ThermalThermal分布假分布假设设：由由T 对对 x 的的分段分段线线性的性的变变化化算得；算得；源源项项的的线线性化性化TP代表整个控制容代表整个控制容积积内的内的值值，即，即采用采用阶阶梯梯性分布性分布进进行行计计算的。算的。当然，不当然，不违违背四背四项项基本法基本法则则，选择选择其它形式的分布曲其它形式的分布曲线线也也是可以的，但尽可能采用是可以的，但尽可能采用简单简单一些的分布曲一些的分布曲线线。以下各以下各节节将将对对离散方程中的各离散方程中的各项给项给予予说说明明*太太 原原 理理 工工 大大 学学7/72 ThermalThermal4.2-2 网格网格间间距距1.1.采用不均匀的网格采用不均匀的网格间间距距WwPeEx x(x)w(x)e可以有效地可以有效地扩扩大大计计算功能。算功能。在温度在温度T 随随x 变变化化剧剧烈的区域上采用烈的区域上采用细细网格，而在网格，而在变变化化缓缓慢的区慢的区域采用域采用较较疏（粗）的网格。疏（粗）的网格。2.2.怎怎样设计样设计一个合适的非均匀网格一个合适的非均匀网格因因为为在在问题问题求解之前，求解之前，T x 的分布是不知道的，那么如何的分布是不知道的，那么如何设设计计网格呢网格呢？*太太 原原 理理 工工 大大 学学8/72 ThermalThermal.对对所要得到的解所要得到的解进进行某些定性的行某些定性的预计预计，使，使设计设计得到某些指得到某些指导导；.采用粗网格采用粗网格进进行行试试算算,求得求得Tx的的变变化形式，再化形式，再对对温度温度变变化化急急剧剧的区域加密，最后构成一个合适的非均匀网格。的区域加密，最后构成一个合适的非均匀网格。3.3.先疏后密的网格划分是有前提的先疏后密的网格划分是有前提的采用粗网格得到的数采用粗网格得到的数值值结结果必果必须须符合物理上的真符合物理上的真实实性性，要做到，要做到这这一点，就一点，就应该应该确保离散方程同确保离散方程同时满时满足四个基本法足四个基本法则则。达到达到给给定精度所需要的网格点数定精度所需要的网格点数,以及以及这这些网格点在些网格点在计计算域内算域内应应采取的分布方式与所求采取的分布方式与所求问题问题的特性有关。的特性有关。4.4.采用采用仅仅几个网格点几个网格点进进行行试试探性探性计计算，算，为为弄清有关解弄清有关解 的情况提供了一个方便的途径。的情况提供了一个方便的途径。也可来指也可来指导实验导实验。*太太 原原 理理 工工 大大 学学9/72 ThermalThermal4.2-3 界面界面导热导热系数系数1.1.问题问题的提出的提出 通用离散方程式通用离散方程式aE、aW分分别别是是节节点点E 与与 P 和和节节点点W 与与P 间间的的热导热导，热导热导的大小的大小反映了周反映了周围节围节点点对节对节点点P的影响程度的影响程度。系数。系数aE、aW中分中分别别含有交含有交界面界面导热导热系数系数ke与与kw。当。当k 是是x 的函数的函数时时，只知道，只知道kP、kE、kW，无法知道无法知道ke与与kw的的值值，而，而ke 与与 kw是决定交界面是决定交界面热热流量的关流量的关键键量。量。因此，因此，计计算算 ke 与与 kw的方法是否合理就的方法是否合理就显显得非常重要了。得非常重要了。式中：式中：*太太 原原 理理 工工 大大 学学10/72 ThermalThermal k 值值不均匀性不均匀性产产生的原因生的原因由材料的不均匀性引起（如由材料的不均匀性引起（如组组合材料板）；合材料板）；材料均匀，材料均匀，T分布的不均匀性也会分布的不均匀性也会导导致致k 的不均匀。的不均匀。2.2.求解方法求解方法.算算术术平均法平均法xPeE(x)e(x)e+(x)e-如如图图所示，所示，P、E之之间间，k与与x 呈呈线线性性关系，关系，则则由由P、E两点上的两点上的kP、kE 确确定定ke 的关系式的关系式为为：*太太 原原 理理 工工 大大 学学11/72 ThermalThermal显显然，然，这这相当于相当于线线性插性插值值。当界面。当界面e e位于两个位于两个节节点之点之间间的中点的中点时时，fe=0.5,此此时时*太太 原原 理理 工工 大大 学学12/72 ThermalThermal.调调和平均法和平均法利用利用传热传热学基本公式可以学基本公式可以导导出界面上出界面上当量当量导热导热系数的系数的调调和平和平均均公式。公式。据界面上据界面上热热流密度流密度连续连续的原的原则则，写出下式：写出下式：另一方面，按界面上当量另一方面，按界面上当量导热导热系数的含系数的含义义，应应有：有：比比较较两式可得：两式可得：可看成是串可看成是串联过联过程程热热阻叠加原阻叠加原则则的反映的反映(x)ePeE(x)e+(x)e-x*太太 原原 理理 工工 大大 学学13/72 ThermalThermal.两种方法的比两种方法的比较较算算术术平均法平均法简单简单方便，但在方便，但在处处理理导热导热性能相差很大的性能相差很大的组组合材料合材料导热时导热时存在明存在明显显缺陷。下面缺陷。下面讨论讨论两种极限情况：两种极限情况：kE0,即即设设想交界面想交界面e 是是k 相差很大的两种材料的分界面，相差很大的两种材料的分界面，节节点点E的控制容的控制容积积是是绝热绝热材料，材料，这时节这时节点点E、P之之间间的的导热导热量量应应该该小到接近于零，即两点小到接近于零，即两点间间的的热热阻阻应应接近于接近于。但用算。但用算术术平均平均法法计计算算,这时这时ke 与与kE无关，无关，仅仅与与kP 有关，有关，不符合物理不符合物理规规律。律。(x)ePeE(x)e+(x)e-x*太太 原原 理理 工工 大大 学学14/72 ThermalThermal P P控制容控制容积积是良是良导热导热体体 ,按算按算术术平均法平均法计计算算,当网格均分当网格均分时时,即即P、E两点两点间间的的热热阻阻为为 ，表，表明此明此时时P、E间间的的热热阻主要由阻主要由k 大的物体所决定，大的物体所决定，显显然不符合然不符合传热传热原理。原理。实际实际上，此上，此时时控制体控制体E 构成了构成了热热阻的主要部分。阻的主要部分。(x)ePeE(x)e+(x)e-x*太太 原原 理理 工工 大大 学学15/72 ThermalThermal调调和平均法和平均法可以合理求解上述两种极限情况可以合理求解上述两种极限情况 当当 时时表明表明ke 完全与完全与kP无关，无关，这这个个结结果是可以果是可以预预料的。因料的。因为围绕为围绕 P 点的材料点的材料k 值值高，其高，其热热阻与阻与围绕围绕E 点的材料相比可以忽略。点的材料相比可以忽略。当当 kE0 时时，由由计计算公式得算公式得 ke0，qe0。意味着在一个。意味着在一个绝绝 热层热层的表面上的表面上热热流密度流密度为为零，与零，与实际实际情况相符。情况相符。*太太 原原 理理 工工 大大 学学16/72 ThermalThermal调调和平均公式和平均公式的推的推导导是是对对于于稳态稳态、无内、无内热热源、源、k在相在相邻邻的两个控制容的两个控制容积积之之间发间发生生阶阶梯式梯式变变化化导导出的。从出的。从定性上，串定性上，串联热联热阻叠加的适阻叠加的适应应性不受上述条件的限性不受上述条件的限制。制。采用两种方法采用两种方法进进行数行数值计值计算的算的结结果表明，果表明，即使即使对对有内有内热热源或源或k 呈呈连续变连续变化的化的场场合，合，调调和平均也比算和平均也比算术术平均平均更好一些。更好一些。*太太 原原 理理 工工 大大 学学17/72 ThermalThermal4.2-4 非非 线线 性性若离散化方程是一个若离散化方程是一个线线性的代数方程，式中的各性的代数方程，式中的各项项系数均系数均为为已已知数，知数，联联立求解代数方程立求解代数方程组组可得到温度可得到温度场场。但但实际问题实际问题中，中，Kp、KE和和Ke或或线线性化源性化源项项的系数的系数SC、SP是温度是温度T的函数，的函数，这样这样离散化方程的系数离散化方程的系数aE、aW、aP本身也成本身也成为为温度的温度的函数，方程函数，方程非非线线性性采用采用拟线拟线性化的方法性化的方法求解（求解（迭代迭代）。）。具体求解步具体求解步骤骤：.先先设设定域内全部定域内全部节节点的温度点的温度值值（给给定一个初定一个初场场）；）；.由由这这些些设设定的定的T T值计值计算出离散化方程中系数的算出离散化方程中系数的试试探探值值；.解名解名义义上的代数方程上的代数方程组组，得到一，得到一组组新的新的T T值值；.以以这这些新的些新的T T值值作作为较为较好的估好的估计值计值，返回到，返回到，并重复此，并重复此过过程，直到重复程，直到重复计计算不再引起算不再引起T T值值任何意任何意义义上的上的变变化化为为止。止。收收敛敛*太太 原原 理理 工工 大大 学学18/72 ThermalThermal收收敛敛：计计算达到算达到最最终终不不变变的状的状态态称之称之为为迭代的收迭代的收敛敛。发发散：散：一次次的迭代永一次次的迭代永远远不会收不会收敛敛到一个解。到一个解。T T的的值值可能可能稳稳定地定地飘飘移或是以一个不断增大的振幅振移或是以一个不断增大的振幅振荡荡,这这种与收种与收敛敛相相对对的的过过程称之程称之为发为发散。散。一种好的数一种好的数值值方法方法应应当使当使发发散的可能性减散的可能性减为为最小。最小。*太太 原原 理理 工工 大大 学学19/72 ThermalThermal4.2-5 源源项项的的 线线 性化性化如果源如果源项项是常数，是常数，则则在离散方程的建立在离散方程的建立过过程中不会程中不会带带来任何困来任何困难难；当源；当源项项是所求是所求变变量的函数量的函数时时，源，源项项的数的数值处值处理十分重要，有理十分重要，有时时甚至是数甚至是数值值求解的关求解的关键键所在。所在。应应用用较为较为广泛的一种广泛的一种处处理方法是把源理方法是把源项项局部局部线线性化性化 SC常数，常数，SP 是是S 随随T 而而变变化的曲化的曲线线在在P点的斜率。点的斜率。表示在表示在TP的附近以直的附近以直线线代替曲代替曲线线。*太太 原原 理理 工工 大大 学学20/72 ThermalThermal几点几点说说明：明：.当源当源项为项为未知量的函数未知量的函数时时，线线性化的性化的处处理比假定源理比假定源项为项为常数更常数更为为合理；合理；因因为为S=f(T)，把，把S 作作为为常数常数处处理就是以上次迭代理就是以上次迭代计计算所得之算所得之T*来来计计算算 S，这样这样源源项项相相对对于于T永永远远有一有一个滞后。而个滞后。而线线性化性化处处理后，理后，TP 是迭代的当前是迭代的当前值值，这样这样使使 S能更快地跟上能更快地跟上TP 的的变变化。化。.线线性化性化处处理又是建立理又是建立线线性代数方程性代数方程所必需的；所必需的；.为为了保了保证证代数方程迭代求解的收代数方程迭代求解的收敛敛，要求，要求 ；*太太 原原 理理 工工 大大 学学21/72 ThermalThermal离散化方程的一般式：离散化方程的一般式：线线性代数方程迭代求解收性代数方程迭代求解收敛敛的一个充分条件是的一个充分条件是主主对对角占角占优优，即，即这这就要求就要求:.由代数方程迭代求解公式：由代数方程迭代求解公式：可可见见，SP绝对值绝对值的大小影响迭代的大小影响迭代过过程中温度的程中温度的变变化速化速度，度，SP绝对值绝对值越大，系越大，系统统的的惯惯性越大，相性越大，相邻邻两次迭代两次迭代之之间间TP 的的变变化越小，收化越小，收敛敛速度下降；但有利于克服迭速度下降；但有利于克服迭代代过过程的程的发发散。散。*太太 原原 理理 工工 大大 学学22/72 ThermalThermal源源项线项线性化方法性化方法举举例例：eg1.已知：已知：S=5-4T,可能的可能的线线性化形式有：性化形式有：相当于相当于设设S 为为常数，当常数，当S 的表达式很复的表达式很复杂时杂时，这样这样做或做或许许是唯一的一种是唯一的一种选择选择。这给这给出了比出了比实际实际S-T关系关系更陡的曲更陡的曲线线，其，其结结果使迭代果使迭代收收敛敛的速度减慢了。的速度减慢了。若在所研究的若在所研究的问题问题中中还还存在着其存在着其它的非它的非线线性性项时项时，这这种减慢的做法是受种减慢的做法是受欢欢迎的。迎的。*太太 原原 理理 工工 大大 学学23/72 ThermalThermaleg2.已知：已知：S=3+7 T,可能的可能的线线性化形式有：性化形式有：eg3.已知：已知：S=4-5 T3,可能的可能的线线性化形式有：性化形式有：已知的已知的S-T曲曲线线要比要比这这一关系所反映的曲一关系所反映的曲线线陡。陡。*太太 原原 理理 工工 大大 学学24/72 ThermalThermal(3)推荐的方法推荐的方法于是：于是：这这一一线线性化表示，在性化表示，在 点所点所选择选择的直的直线线与与S-T曲曲线线相切。相切。这这一一线线性化比已知的性化比已知的S-T曲曲线为线为陡，它使收陡，它使收敛敛速度速度减慢。减慢。*太太 原原 理理 工工 大大 学学25/72 ThermalThermal结论结论：在所有在所有负负斜率的直斜率的直线线中，与已知曲中，与已知曲线线相切的相切的直直线线通常通常为为最佳最佳(图图中中2线线)；较较陡的直陡的直线线是可以接受的是可以接受的(图图中中3线线)；欠陡的直；欠陡的直线线是不希望采是不希望采纳纳的的(图图中中1线线)，它不能体它不能体现现已知的已知的S 随随T 的下降速度。的下降速度。S ST T已知曲已知曲线线1 12 23 3*太太 原原 理理 工工 大大 学学26/72 ThermalThermal4.2-6 边边 界界 条条 件件1.1.问题问题的提出的提出前面所推得的离散化方程适用于前面所推得的离散化方程适用于稳态导热问题稳态导热问题的任何内部的任何内部节节点，点，为计为计算一个具体算一个具体问题问题，应应把把边边界条件也用离散方程表示。界条件也用离散方程表示。因因为为只有离散化方程的个数与待求只有离散化方程的个数与待求节节点点变变量的数目相等量的数目相等时时，代数方程代数方程组组才能封才能封闭闭。2.2.网格划分的两种方法网格划分的两种方法先定先定节节点位置，后定控制容点位置，后定控制容积积 A A类类网格网格 BiIeE x这这种网格划分方法，在种网格划分方法，在边边界上将出界上将出现现半控制容半控制容积积。当网格划分不均匀当网格划分不均匀时时，节节点位置并不落在控制容点位置并不落在控制容积积的的几何中心位置。几何中心位置。*太太 原原 理理 工工 大大 学学27/72 ThermalThermal先定控制容先定控制容积积，后定，后定节节点位置点位置 B B类类网网格格*这这种方法是在种方法是在边边界上附上一界上附上一层层厚度厚度为为零的控制容零的控制容积积，代表代表这这个控制容个控制容积积的的边边界界节节点恰好落在点恰好落在边边界上。无界上。无论论网格如何划分都不会出网格如何划分都不会出现现半控制容半控制容积积，并且所有，并且所有节节点点都位于控制容都位于控制容积积的几何中心。的几何中心。iIjJB由于在由于在边边界上将出界上将出现现不同的不同的控制容控制容积积，所以将根据不同，所以将根据不同的方法来离散的方法来离散边边界条件。界条件。三三类类不同的不同的边边界条件界条件*太太 原原 理理 工工 大大 学学28/72 ThermalThermal3.3.构建构建边边界界节节点的附加方程点的附加方程 A A类类网格网格 BiIeE x(x)i.第一第一类边类边界条件界条件 TB已知已知不必不必额额外增加外增加边边界界节节点方程，把点方程，把 TB代入代入邻邻近近节节点的代数方程即可。点的代数方程即可。.第二第二类边类边界条件界条件 qB已知已知可在可在边边界半控制容界半控制容积积内内对对微分方程微分方程积积分建立附加方程，分建立附加方程，热热平平衡式衡式为为：qB注：注：规规定以定以进进入入计计算区域的算区域的热热量量为为正。正。交界面交界面i 处处的的热热流流为为：*太太 原原 理理 工工 大大 学学29/72 ThermalThermal将将qi 代入上式得：代入上式得：qB已知已知时时，方程可整理成：，方程可整理成：式中：式中：.第三第三类边界条件界条件 已知已知 h、tf将将边边界界热热流流 代入上式的代入上式的B点方程，点方程，形式同上：形式同上：式中：式中：*太太 原原 理理 工工 大大 学学30/72 ThermalThermal3.3.构建构建边边界界节节点的附加方程点的附加方程 B B类类网格网格*.第一第一类边类边界条件界条件 TB已知已知不需要附加方程不需要附加方程.第二第二类边类边界条件界条件 qB已知已知iIjJB(x)bB点方程点方程为为：式中：式中：也可以看作是第一种情况，在离散区域也可以看作是第一种情况，在离散区域 x0时时的极限。的极限。BiIeE x(x)i*太太 原原 理理 工工 大大 学学31/72 ThermalThermal.第三第三类边界条件界条件 已知已知 h、tfB点方程形式同上点方程形式同上式中：式中：至此，已构成至此，已构成对对所有未知温度的足所有未知温度的足够够数量的方程数量的方程*太太 原原 理理 工工 大大 学学32/72 ThermalThermal4.2-7 线线性代数方程的求解性代数方程的求解1.1.求解方法求解方法线线性代数方程性代数方程组组的解法通常有的解法通常有迭代法和直接解法迭代法和直接解法两种。两种。由于一由于一维导热维导热数数值值求解的离散方程中，待求温度求解的离散方程中，待求温度仅仅与左与左右两个右两个节节点的温度有关，点的温度有关，这样这样形成的代数方程形成的代数方程组组的系数的系数矩矩阵阵将是将是三三对对角矩角矩阵阵，采用追赶法（，采用追赶法（TDMA）.2.2.TDMA算法要点算法要点 TDMA是一种是一种简单、方便、高效率的、方便、高效率的计算方法，求解算方法，求解过程包括程包括消去消去变变量量求系数的求系数的正正过过程和回代程和回代求温度求温度场场的的逆逆过过程程。*太太 原原 理理 工工 大大 学学33/72 ThermalThermal正正过过程：程：消元的目的是将每个包含消元的目的是将每个包含三个待求三个待求变变量量的的方程，消去一个方程，消去一个变变量，使之量，使之成成为为两个待求两个待求变变量量的方的方程，程，继继而把所有而把所有这这些两个些两个变变量的方程新系数量的方程新系数值计值计算算出来，直到出来，直到边边界上界上原来只有两个待求原来只有两个待求变变量的方程，量的方程，经经消元后消元后变变成了成了单变单变量方程量方程，并能直接求出，并能直接求出边边界界变变量量值为值为止。止。回代回代过过程：程：由已知的由已知的边边界界变变量量值值逐一按次序逐一按次序回代到回代到已已求出系数的求出系数的二二变变量方程中量方程中便可求出全部待求便可求出全部待求变变量。量。经过经过消元与回代便完成消元与回代便完成TDMA的全的全过过程。程。*太太 原原 理理 工工 大大 学学34/72 ThermalThermalTDMA算法的数学表述算法的数学表述21Nii+1i-13i节节点的离散方程可写点的离散方程可写为为：考考虑虑到到边边界界节节点点的特殊形式：的特殊形式：如果如果T1已知，已知，则则有：有：这这些条件意味着些条件意味着T1 是是T2的函数，而的函数，而T2与与T1、T3，T2 T3。这这一代入一代入过过程可一直程可一直继续继续到到TN TN+1的关的关系式。由于系式。由于TN+1并不存在，所以到此并不存在，所以到此实际实际上也就得到上也就得到了了TN的的值值。由此可以回代了，由。由此可以回代了，由TN TN-1 TN-2 T2 T1。*太太 原原 理理 工工 大大 学学35/72 ThermalThermal3.3.TDMA的数学推演的数学推演 假假设设在代入在代入过过程得到：程得到：试图试图找到一个关系：找到一个关系：将将 代入离散方程的一般表达式中，将代入离散方程的一般表达式中，将Ti-1用用Ti代即可代即可*太太 原原 理理 工工 大大 学学36/72 ThermalThermal 对对比比与与即可得到即可得到Pi、Qi的的递递推公式：推公式：特例：特例：求得求得TN后，可以回代了！后，可以回代了！*太太 原原 理理 工工 大大 学学37/72 ThermalThermal4.4.TDMA的的计计算机程序算机程序说说明明.由由P1、Q1的表达式的表达式计计算出算出P1、Q1的的值值；.由由P、Q的的递递推式推式计计算得出算得出P2、Q2；Pi、Qi.QN.令令 TN=QN.对对N-1,N-2,3，2，1应应用回代式用回代式 可求得可求得TN-1、TN-2 T3、T2、T1 。只要代数方程可以表达成只要代数方程可以表达成 的形式，的形式，三三对对角矩角矩阵阵的算法就是一个非常有效而又方便的求解的算法就是一个非常有效而又方便的求解方法。方法。TDMA法需要的法需要的计计算机算机贮贮存量及存量及计计算算时间时间与与N成正比，而非成正比，而非N2或或N3。它是一种直接解法。它是一种直接解法。*太太 原原 理理 工工 大大 学学38/72 ThermalThermal4.3 非非稳态稳态一一维热传导维热传导 通用控制方程中的通用控制方程中的扩扩散散项项和源和源项项的数的数值处值处理已理已经经解解决，下面决，下面讨论讨论非非稳态项稳态项的的处处理理(暂暂不考不考虑虑源源项项)。1.1.一一维维非非稳态稳态无源无源项项的的热传导热传导方程方程 任任务务：因因为时间为时间是一个是一个单单向坐向坐标标，即由一已知的初，即由一已知的初始温度开始，沿始温度开始，沿时间时间坐坐标标逐步向前求解。故逐步向前求解。故计计算任算任务务是已知是已知 时时刻刻T 在网格上的在网格上的值值，求得，求得 +时时刻刻T 在相在相应应网格点上的网格点上的值值。4.3-1 通用的离散化方程通用的离散化方程*太太 原原 理理 工工 大大 学学39/72 ThermalThermal2.2.对对整个控制容整个控制容积积积积分上述微分方程分上述微分方程WwPeEx x(x)w(x)e 需要需要给给出一个关于出一个关于TP、TE 和和TW 如何随如何随时间变时间变化的化的关系的假关系的假设设，可以采取的假，可以采取的假设设有好几种。有好几种。*太太 原原 理理 工工 大大 学学40/72 ThermalThermal 一些假一些假设设可可归纳为归纳为一般化的式子：一般化的式子：由此原由此原积积分整理分整理为为：*太太 原原 理理 工工 大大 学学41/72 ThermalThermal经经整理得：整理得：（将上（将上标标1去掉）去掉）式中：式中：*太太 原原 理理 工工 大大 学学42/72 ThermalThermal4.3-2 显显式、克式、克兰兰克克-尼科尼科尔尔森模式以及全森模式以及全隐隐式模式式模式上上节节离散方程中，离散方程中，对应对应特定的特定的f 值值，可得到不同的，可得到不同的格式。格式。f=0 显显式；式；f=0.5 克克兰兰克克-尼克尼克尔尔森模式森模式(C-N)；f=1 全全隐隐式。式。1.1.显显式格式式格式 f=0 方程方程变为变为：*太太 原原 理理 工工 大大 学学43/72 ThermalThermal +显显式格式式格式全全隐隐格式格式C-N格格式式T TP P TP 的的变变化曲化曲线见线见右右图图，实质实质上是假上是假设设上一上一时时刻刻 代表除了代表除了时时刻刻+之外的整之外的整个个时间间时间间隔上的隔上的值值。TP与其它未知量与其它未知量(TE、TW)无关，而直接可由上一无关，而直接可由上一时时刻刻已知温度求得已知温度求得“显显式式”*太太 原原 理理 工工 大大 学学44/72 ThermalThermal显显式格式：式格式：前的系数有可能前的系数有可能变为负值变为负值(可看可看成是成是 在在时间时间方向上的一个相方向上的一个相邻邻点点)，从而，从而违违背正背正系数法系数法则则。实际实际上，要使系数上，要使系数为为正，正，时间时间步步长长必必须须足足够够小，以使小，以使对对于于k 均匀分布，且均匀分布，且即即 的的选选取受上述条件的限制。取受上述条件的限制。*太太 原原 理理 工工 大大 学学45/72 ThermalThermal2.2.全全隐隐式格式式格式 f=1方程方程变为变为：+显显式格式式格式全全隐隐格式格式C-N格格式式T TP P方程中各方程中各项项系数恒大系数恒大于零，无条件于零，无条件稳稳定。定。TP 的的变变化曲化曲线见线见右右图图，在在时时刻刻，T TP P的的值值突然由突然由降到降到 ，而在随后的整个，而在随后的整个时时间间间间隔上保持隔上保持 。于是于是在整个在整个时时步内温度由步内温度由新新值值 所确定。所确定。与其它未与其它未知量知量(TE、TW)有关，必有关，必须须解解一一组联组联立方程立方程“隐隐式式”*太太 原原 理理 工工 大大 学学46/72 ThermalThermal3.3.克克兰兰克克-尼克尼克尔尔森模式森模式(C-N)f=0.5 这这种格式也称种格式也称半半隐隐格式格式,方程方程变为变为：表示求解表示求解TP可以用上一可以用上一时时刻刻 和当地和当地时时刻刻(+)的的值值，TP与与TE、TW有关，同有关，同样样需要需要联联立方程立方程组组求解。求解。但但 前的系数有可能出前的系数有可能出现负值现负值，所以也是，所以也是条件条件稳稳定定的，其条件是：的，其条件是：*太太 原原 理理 工工 大大 学学47/72 ThermalThermal两点两点说说明：明：.数学上的数学上的“稳定性定性”只能确保解的振只能确保解的振荡最最终将会消将会消失，它并不能确保一定会失，它并不能确保一定会给出物理上似乎合理的解。出物理上似乎合理的解。.对对于于小的小的时间时间步步，全，全隐隐格式不如格式不如C-N格式的解精格式的解精确。因确。因为对为对于小的于小的时间间时间间隔，温度隔，温度时间时间曲曲线线是接近是接近于于线线性的。性的。探索一种探索一种结结合两种模式合两种模式优优点而不分担它点而不分担它们们各自缺点的模式是各自缺点的模式是诱诱人的人的“指数格式指数格式”。.全全隐隐格式，表达形式格式，表达形式简单简单而且物理上又能达到而且物理上又能达到满满意要求，因此在后面的介意要求，因此在后面的介绍绍中将采用全中将采用全隐隐格式。格式。*太太 原原 理理 工工 大大 学学48/72 ThermalThermal4.3-3 全全隐隐式离散化方程式离散化方程1.1.全全隐隐格式的格式的优优点点 f=1 离散化方程无条件离散化方程无条件稳稳定。定。2.2.源源项项引入后的离散化方程形式引入后的离散化方程形式可可见见，当，当时时，此方程，此方程简简化化为稳态为稳态的离散方程。的离散方程。时间时间上的相上的相邻邻点点*太太 原原 理理 工工 大大 学学49/72 ThermalThermal3.全全隐隐格式的主要原格式的主要原则则:TP的新的新值值代表整个代表整个时间时间步上的步上的值值。因此如果因此如果导导热热系数系数kP是温度的函数，就是温度的函数，就应应当反复由当反复由TP的迭代算的迭代算得新的得新的kP。稳态稳态程序中的程序中的边边界条件、源界条件、源项项的的线线性化性化处处理以理以及及TDMA求解方程求解方程组组的算法全部适合非的算法全部适合非稳态问题稳态问题。*太太 原原 理理 工工 大大 学学50/72 ThermalThermal4.4 二二维维与三与三维问题维问题4.4-1 二二维问题维问题的离散化方程的离散化方程1.1.二二维维控制容控制容积积及及标标号号说说明明y y xWENSnwsex控制容控制容积积图图4-6 4-6 对对于二于二维维情况的控制容情况的控制容积积控制容控制容积积面相面相对对于网格点的于网格点的实际实际位置是自位置是自由的，但把控由的，但把控制容制容积积面布置面布置在两个网格点在两个网格点的中点是一种的中点是一种很自然的做法很自然的做法*太太 原原 理理 工工 大大 学学51/72 ThermalThermalP P、E E两点两点间间e e控制面上的控制面上的热热流密度流密度整个整个e e 控制面上的控制面上的热热流量流量2.2.二二维维非非稳态稳态、有内、有内热热源的控制方程源的控制方程其它其它 控制面上控制面上热热流的流的计计算与之算与之类类同。同。相相应应离散方程的一般形式：离散方程的一般形式：*太太 原原 理理 工工 大大 学学52/72 ThermalThermal式中各系数式中各系数为为：*太太 原原 理理 工工 大大 学学53/72 ThermalThermal4.4-2 三三维问题维问题的离散化方程的离散化方程与二与二维问题维问题的推的推导类导类似，再加一个似，再加一个顶顶部、底部方部、底部方向的相向的相邻邻点点T T点及点及B B点。点。系数的意系数的意义义：aE、aW.代表代表P点与相点与相邻邻点之点之间间的的热导热导表示表示 时时刻控制容刻控制容积积内部所包含的内能（除内部所包含的内能（除以以 ）表示表示内能内能项项与与SC所造成的在控制容所造成的在控制容积积内的内的发热发热率之和。率之和。表示表示相相邻邻空空间间与与时间时间系数之和，并包含一系数之和，并包含一项线项线性的源性的源项项所做的所做的贡贡献。献。*太太 原原 理理 工工 大大 学学54/72 ThermalThermal4.4-3 代数方程的求解代数方程的求解可以采用任何一种合适的求解方法，可以采用任何一种合适的求解方法，选择选择方法方法时时两者两者无相互影响；无相互影响；在在计计算机程序中，可以很方便地把两个步算机程序中，可以很方便地把两个步骤骤分成两个分成两个独立的独立的阶阶段，其中任何一段都可以段，其中任何一段都可以单单独独进进行修改。行修改。1.1.方程的离散化与求解方法分成两个步方程的离散化与求解方法分成两个步骤骤的好的好处处2.2.求解方法求解方法直接解法直接解法无需迭代的方法。但需要大量的无需迭代的方法。但需要大量的计计算机算机贮贮存量及存量及计计算算时间时间。对对于一个于一个只需一次求解的只需一次求解的线线性性问题问题，直接解法是可以接受的。直接解法是可以接受的。*太太 原原 理理 工工 大大 学学55/72 ThermalThermal迭代法：迭代法：对对于非于非线线性性问题问题，方程必方程必须须用最新的系数用最新的系数来重复求解，采用直接解法是很不来重复求解，采用直接解法是很不经济经济的。的。3.迭代解法的基本原理迭代解法的基本原理给给定定一个因一个因变变量量T T 的的初初场场，然后利用某种形式的代数，然后利用某种形式的代数方程方程求得一个改求得一个改进进的的场场，重复，重复进进行行这这一算法一算法过过程，直程，直到求得一个充分接近代数方程精确解的解到求得一个充分接近代数方程精确解的解为为止。止。在迭代在迭代计计算的算的过过程中，要程中，要注意注意计计算系数所付的代价与算系数所付的代价与求解方程所需花求解方程所需花费费的的时间时间的某种平衡。的某种平衡。在一在一组组固定的固定的系数系数值值下，把代数方程的求解下，把代数方程的求解过过程一直程一直进进行到最行到最终终收收敛敛是不必要、也是不明智的。在新的一是不必要、也是不明智的。在新的一组组系数系数值给值给定定之前之前应应用已用已经给经给定的系数定的系数值值迭代几次就已足迭代几次就已足够够了。了。*太太 原原 理理 工工 大大 学学56/72 ThermalThermal4.4.迭代求解的两种方法迭代求解的两种方法:高斯高斯-赛赛德德尔尔逐点逐点计计算法和逐行法算法和逐行法高斯高斯-赛赛德德尔尔逐点迭代法逐点迭代法（G-S）：迭代原理：迭代原理：按一定按一定顺顺序逐个序逐个访问访问每个网每个网格点，以格点，以计计算那里的算那里的变变量量值值。在。在计计算机内只需要算机内只需要储储存一存一组组T T值值，且存，且存贮贮器中相器中相应应的的T T值值交替改交替改变变如如下：下：离散化方程离散化方程为为：*太太 原原 理理 工工 大大 学学57/72 ThermalThermal计计算机存算机存贮贮器中所存在相器中所存在相邻邻点的温度点的温度值值。本次本次迭代迭代过过程中已程中已经经被被访问过访问过的相的相邻邻点，点，是新是新值值，而，而对对于那些尚待于那些尚待访问访问的相的相邻邻点点 是上次迭代得到的是上次迭代得到的值值。都是相都是相邻邻温度的最新温度的最新值值。当所有网格点都当所有网格点都访问访问过过 一次后，算是完成了一次高斯一次后，算是完成了一次高斯赛赛德德尔尔迭代。迭代。迭代迭代计计算收算收敛敛的准的准则则斯卡巴勒准斯卡巴勒准则则高斯高斯赛赛德德尔尔迭代收迭代收敛敛的的充分条件充分条件是：系数是：系数满满足足*太太 原原 理理 工工 大大 学学58/72 ThermalThermal这这个准个准则则是是充分条件充分条件，并非必要条件。即有，并非必要条件。即有时虽时虽然然违违背了背了这这个条件，但迭代仍可能收个条件，但迭代仍可能收敛敛。对对于于满满足正系数法足正系数法则则的方程，如非的方程，如非稳态稳态全全隐隐格式格式 ,由于各系数均大于零，且由于各系数均大于零，且几点几点说说明：明：对对于于稳态稳态有内有内热热源的源的问题问题，只要保，只要保证证 即可即可使得使得 成立（成立（当当SP=0时时，不等号由，不等号由边边界界节节点点获获得）。得）。*太太 原原 理理 工工 大大 学学59/72 ThermalThermal.对对于于稳态稳态无内无内热热源源的的问题问题，对对内内节节点只有点只有然而可以从然而可以从边边界界邻邻点点的离散方程中找出的离散方程中找出使得不等号成使得不等号成立立的条件。的条件。WEPSN如：假如：假设边设边界有一温度已知，如二界有一温度已知，如二维维网格中的网格中的W边边界，界，则对则对P点的方程点的方程为为：已知，已知，转转入到常数入到常数项项中，中，这时这时，计计算算 时时，只只对对那些未知那些未知的相的相邻邻系数求和，而系数求和，而计计算算aP时则时则是包含是包含边边界点系数在内的系界点系数在内的系数之和，数之和，所以所以 。对对第三第三类类BCs也成也成立立*太太 原原 理理 工工 大大 学学60/72 ThermalThermal高斯高斯赛赛德德尔尔点迭代点迭代计计算的缺点算的缺点收收敛敛速度慢，特速度慢，特别对别对于网格于网格节节点数多的情况尤点数多的情况尤为为明明显显。因因为为它是以一次迭代它是以一次迭代传递传递一个网格一个网格间间距的速率来距的速率来传递边传递边界界所所给给予的信息的予的信息的(仅仅将将扫扫描起始描起始边边的影响的影响传传到整个区域到整个区域)。高斯高斯-赛赛德德尔尔逐行迭代法逐行迭代法迭代原理：迭代原理：因点迭代不能很快地把因点迭代不能很快地把边边界信息界信息传传播到区播到区域内部，因此就域内部，因此就产产生了一种生了一种试图试图充分利用直接消元法和充分利用直接消元法和迭代法迭代法优优点的点的块块迭代法，即按照行或列划分成迭代法，即按照行或列划分成块块，在在块块内用内用TDMA，块块之之间间采用采用G-S法，法，这这就形成了一种高效就形成了一种高效的求解方法。的求解方法。因同一因同一块块内内各各节节点的点的值值是以是以隐隐含的方式相含的方式相互互联联系的，故又叫系的，故又叫隐隐式迭代法。式迭代法。*太太 原原 理理 工工 大大 学学61/72 ThermalThermal迭代方法迭代方法简简介：介：设设二二维导热维导热的离散化方程的离散化方程为为WPENSi 选线选线选块选块yx采用逐采用逐线线迭代迭代时时，将每，将每线线（行或列）（行或列）作作为为一一块块，线线内采用内采用TDMA求解，求解，线间线间采用采用G-S法求解。如法求解。如图图所示，所示，进进行逐行逐线线迭代迭代时时，i 列中列中P点的左右点的左右邻邻点点TW、TE都取已知的当前都取已知的当前值值，只有只有TP、TN、TS是未知是未知值值，这样这样同同一一维导热问题维导热问题，可以用，可以用TDMA法求法求解解线线内每个内每个节节点的温度，而点的温度，而线线两端两端的的边边界信息很快界信息很快传给线传给线内的所内的所有有节节点。点。x方向的迭代用方向的迭代用G-S法，从左向右法，从左向右扫扫描描时时，TW是是刚刚刚刚得到的当前得到的当前值值。*太太 原原 理理 工工 大大 学学62/72 ThermalThermal实际实际上，采用逐上，采用逐线线迭代迭代时时，二，二维导热维导热的离散化方程可的离散化方程可改写改写为为：逐逐线线迭代法能加速收迭代法能加速收敛敛由于同一由于同一线线上的上的所有所有节节点点都能在一次都能在一次TDMA求解中求解中获获得得边边界信息，而界信息，而扫扫描方向两端的描方向两端的边边界信息，界信息，虽虽然然还还要要依靠点迭代来依靠点迭代来传传播，但由于播，但由于线线数数毕毕竟要比竟要比节节点数少得点数少得多；再加上采用了逐点更新的多；再加上采用了逐点更新的G-S法，所以法，所以计计算速度算速度将会明将会明显显加快。加快。交替改交替改变应变应用逐行迭代的方向用逐行迭代的方向*太太 原原 理理 工工 大大 学学63/72 ThermalThermal逐逐线线迭代方向迭代方向选择选择的原的原则则（扫扫描方向）描方向）收收敛敛速度的快慢关速度的快慢关键键在于在于边边界信息的界信息的传传播速度。播速度。因此因此扫扫描的描的起始起始线线上上应该应该有有确定的确定的信息。信息。这样这样，绝热绝热或梯度或梯度为为零的零的边边界不界不应应作作为扫为扫描的起始描的起始边边。T1T2T3?当有当有对对流存在流存在时时，扫扫描方向尤描方向尤为为重重要。由上游向下游要。由上游向下游扫扫描的速度比相描的速度比相反方向反方向扫扫描的收描的收敛敛速度快得多。速度快得多。求解狭求解狭长长条区域的条区域的导热问题导热问题，应应以以长边长边的的扫扫描描为为主。（主。（热热阻大的方向）阻大的方向）当当y方向的系数方向的系数远远大于大于x 方向的系数方向的系数时时，对对 y 方向方向应应用用TDMA，沿沿x方向方向进进行行线线迭代收迭代收敛敛速度快。（速度快。（热热阻大的方向信息多）阻大的方向信息多）*太太 原原 理理 工工 大大 学学64/72 ThermalThermal4.5 超松弛与欠松弛超松弛与欠松弛1.1.定定义义：在代数方程迭代求解的在代数方程迭代求解的过过程中，往往希程中，往往希望望加快或是减慢前后两次迭代加快或是减慢前后两次迭代间间因因变变量的量的变变化。化。加加快快变变化化为为超松弛超松弛，减慢，减慢变变化化为为欠松弛。欠松弛。超松弛常用于和超松弛常用于和G-S法相法相结结合，合，这这种种组组合起来的模式合起来的模式是所是所谓谓的持的持续续超松弛。超松弛很少与逐行或逐列迭代超松弛。超松弛很少与逐行或逐列迭代法法结结合使用。合使用。在在强强烈的非烈的非线线性方程性方程组组的迭代求解的迭代求解过过程中，程中，采用欠松采用欠松弛避免弛避免发发散。散。*太太 原原 理理 工工 大大 学学65/72 ThermalThermal2.2.方程式的写法方程式的写法本次迭代所本次迭代所产产生的生的TP变变化化 为为松弛因子。松弛因子。1为为逐次超松弛逐次超松弛(SOR)，0欠松弛，欠松弛，i 1,超松弛。超松弛。左左侧侧系数比原来小了，右系数比原来小了，右侧侧多了一个多了一个负负值项值项；在；在惯惯性松弛中，性松弛中，i 0 的的结结果与上述相果与上述相对应对应，故，故也也为为超松弛。超松弛。0 的的结结果与上述相果与上述相对应对应，故也故也为亚为亚松弛。松弛。*太太 原原 理理 工工 大大 学学68/72 ThermalThermal5.i 值值的的选选取取不存在不存在选选取最佳取最佳i 值值的一般法的一般法则则，对对特定的特定的问题问题，可根据可根据实际经验实际经验确定。确定。i 值应值应当与当与aP值值不相上下，同不相上下，同时时i 的模的模值值越大，松弛越大，松弛作用也愈作用也愈剧剧烈。烈。6.6.稳稳定定问题问题与非与非稳稳定定问题问题求解方法上的求解方法上的类类似似 将将稳稳定定问题问题作作为为非非稳稳定定问题问题来求解是一种有效的方法。来求解是一种有效的方法。这时问题这时问题的的初始条件初始条件就是就是迭代迭代过过程的初始程的初始试试探探值值，每，每一一时时刻的解就相当于刻的解就相当于问题问题的一次迭代解，非的一次迭代解，非稳态问题稳态问题中的中的 相相应应表示上次迭代表示上次迭代值值 ，与与 i 起着相同的作起着相同的作用。因用。因 恒正，相当于一种特殊的欠松弛。恒正，相当于一种特殊的欠松弛。*太太 原原 理理 工工 大大 学学69/72 ThermalThermal 越小，越小，越大，欠松弛作用愈越大，欠松弛作用愈强强；当当 很大很大时时，很小，不再有欠松弛作用，而直很小，不再有欠松弛作用，而直接接转转化化为稳态问题为稳态问题的求解。的求解。7.7.加速加速边边界条件的影响界条件的影响传传入入计计算区域的方法算区域的方法Jacobi点迭代：完成一点迭代：完成一轮轮后，后，边边界的影响界的影响只能只能传传入到入到与与边边界相界相邻邻的一批的一批节节点上，即点上，即传传入一个网格入一个网格，且，且扫扫描的方向与收描的方向与收敛敛快慢无关。快慢无关。G-S点迭代：从左到右完成一点迭代：从左到右完成一轮扫轮扫描后，描后，左左边边界的影界的影响就响就传传到整个区域到整个区域。但右、上、下。但右、上、下边边界的影响也只界的影响也只能能传传入一个网格，迭代速度比入一个网格，迭代速度比Jacobi法快，且受迭代法快，且受迭代扫扫描方向的影响。描方向的影响。*太太 原原 理理 工工 大大 学学70/72 ThermalThermal线线迭代：自左向右迭代：自左向右扫扫描描(G-S法法)，不，不仅仅左左边边界的影界的影响逐步响逐步传传入，且在每列的直接求解中，上下端点入，且在每列的直接求解中，上下端点的影响全部的影响全部传传入到入到该该列的各个点上。因此，每完列的各个点上。因此，每完成一成一轮轮迭代，迭代，左、上、下左、上、下边边界的影响全部界的影响全部传传入，入，但但右右边边界的影响界的影响则仅则仅前前进进一个网格。一个网格。交替方向交替方向隐隐式式扫扫描法（描法（ADI法）：每一法）：每一轮轮的迭代的迭代包括了逐行与逐列的包括了逐行与逐列的扫扫描，描，因而在每一因而在每一轮轮迭代中，迭代中，所有所有边边界的影响均已界的影响均已传传入区域内部入区域内部，从而加快了，从而加快了收收敛敛速度。速度。*太太 原原 理理 工工 大大 学学71/72 ThermalThermal8.8.边边界条件界条件对对区域内部的影响程度区域内部的影响程度第一第一类边类边界条件：具有最界条件：具有最强强烈的影响，最有利于迭代烈的影响，最有利于迭代计计算算过过程的收程的收敛敛；第三第三类边类边界条件：界条件：规规定了参考点温度以及与壁面上温定了参考点温度以及与壁面上温度梯度成正比的系数度梯度成正比的系数h，对边对边界温度界温度规规定的明确程度定的明确程度较较第一第一类边类边界条件弱，界条件弱，对计对计算区域内部的影响算区域内部的影响较较弱；弱；第二第二类边类边界条件：界条件：仅规仅规定了定了边边界上温度界上温度变变化的斜率。化的斜率。因因为为在相同的斜率下，壁面温度可以在任意范在相同的斜率下，壁面温度可以在任意范围围内内变变动动，对对确定内部确定内部节节点温度提供的信息最少，不利于迭点温度提供的信息最少，不利于迭代代过过程的收程的收敛敛。*太太 原原 理理 工工 大大 学学72/72 ThermalThermal4.6-2 其它坐其它坐标标系系以上所提到的方法不只限于直角坐以上所提到的方法不只限于直角坐标标系，它系，它还还可以可以用于任意一种正交坐用于任意一种正交坐标标系。系。由一新的坐由一新的坐标标系引入的系引入的补补充特性主要是几何上的特充特性主要是几何上的特性，性，只要所要求的只要所要求的长长度、面度、面积积以及体以及体积计积计算得合适，算得合适，那就不需要什么新的原那就不需要什么新的原则则。于是，在任何正交坐于是，在任何正交坐标标系中的离散化方程可以按照同系中的离散化方程可以按照同样样的方法得到。的方法得到。*太太 原原 理理 工工 大大 学学73/72 ThermalThermal在极坐在极坐标标系中的控制容系中的控制容积积二二维维极坐极坐标问题标问题 (r,）二二维维极坐极坐标标系下的通用方程系下的通用方程(r,）坐坐标标系中的控制容系中的控制容积见积见右右图图，设设在在z方向厚度方向厚度为为1。控制容。控制容积积的体的体积为积为0.5（rn+rs）r。将微分方程在控制容将微分方程在控制容积积内内积积分，分，得到离散化方程：得到离散化方程：其形式与直角坐其形式与直角坐标标系下的相同，只是系数的表达式不同系下的相同，只是系数的表达式不同*太太 原原 理理 工工 大大 学学74/72 ThermalThermal式中各式中各项项系数系数为为：要求正交坐要求正交坐标标系的必要性！系的必要性！热热流流计计算算对应对应面面积积的要求的要求*太太 原原 理理 工工 大大 学学75/72 ThermalThermal4.7 结结 语语 本章本章迈迈开了开了发发展通用微分方程数展通用微分方程数值值方法的第一步