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

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