传热流体数值计算.doc

传热与流体流动的数值计算 传热与流体流动的数值计算 1 傅立叶定律傅立叶定律是导热理论的基础。其向量表达式为: (2-1)式中:热流密度，是向量，；温度梯度，是向量，/m ；导热系数，又称热导率，； 式中的负号表示的方向始终与相反。2 导热系数（thermal conductivity）及其影响因素PzxyesxnbtwWENSTB导热系数（）是一个比例常数，在数值上等于每小时每平方米面积上，当物体内温度梯度为1/m时的导热量。 导热系数是指在稳定传热条件下，1m厚的材料，两侧表面的温差为1度（K，C）,在1秒内，通过1平方米面积传递的热量，用表示，单位为瓦/米度，w/mk（W/mK,此处的K可用代替）。导热系数为温度梯度1/m，单位时间通过每平方米等温面的热传导热流量。单位是：W/(mK)。 3热传导微分方程推导 在t时刻w界面的温度梯度为在t时刻e界面的温度梯度为单位时间内六面体在方向流入的热流量为：；单位时间内六面体在方向流出的热流量为：；单位时间内六面体在方向流入的净热量为：图3-1 微分单元体各面上进出流量示意图 同理，单位时间内六面体在y方向流入的净热量为：； 单位时间内六面体在y方向流入的净热量为：； 单位时间内流入六面体的总热量为： (3-1)六面体内介质的质量为:。单位时间六面体内热量的变化量（增加）为:根据热量守恒定律:， ， ， ， 称为热扩散率或热扩散系数（thermal diffusivity），单位m2/s. :导热系数，单位W/(mK)； ：密度，单位kg/m3 c：热容，单位J/(kgK).思考：如果单元体内有热源：单位体积单位时间的散热量是q 方程怎么变？4岩石的热扩散率（导温系数) thermal diffusion coefficient ;thermal diffusivity; thermal degradation岩石的热扩散率也叫或热扩散系数，表示岩石在加热或冷却时各部分温度趋于一致的能力。它反映岩石的热惯性特征，是一个综合性参数。热扩散率越大的岩石，热能传播温度趋于一致的速度越大，透入的深度也越大。热扩散系数一般是根据岩石的导热系数(ramuda)、和密度(rou)的测量数据计算得到的。如果单元体内有热源：单位体积单位时间的散热量是q 则在 左边增加热源散发的热量qdxdydz，方程变为：如果是流体除了热传导外，还有流体流动带入和带出微元体的热量在t时刻w界面流体速度为U，流体温度为T单位时间流入微元体的流体质量为： 带入微元体的热量为：e界面流体速度为,流体温度为 单位时间流出微元体的流体质量为：带出微元体的热量为: 如果不考虑x方向速度变化，略去高阶微量，则e界面带出微元体的热量为：单位时间内在方向流入六面体的净热流量为：；同理， y方向： z方向：将上述三项加入到热传导方程（左边，进入微元体的热量）： （能量方程）2.2巷壁与风流间的对流换热运动着的流体与所接触的固体壁面之间的热量传递过程称为对流换热，它是流体(液体或气体)由于宏观相对运动，从某一区域迁移到温度不同的另一区域时引起热量传递的现象。固体壁面与流体之间存在温度差将产生对流换热，由于实际流体的粘性和壁面摩擦的共同影响，近壁流体分层流动，尤其与壁面直接接触的几何面上，总有一层很薄的流体粘附于表面，该层流体处于静止状态，所以热流通过表面层的传递只能依靠导热。显然，在流体发生热对流的同时，由于流体中温度分布的不均匀，也将伴随产生导热现象。因此，对流换热过程实际上是热对流和热传导的综合作用过程。牛顿冷却公式对流换热过程是一个受很多因素影响的复杂过程，如流体的流动状况、流体的物理性质、壁的形状和大小、表面粗糙度等。一般情况下对流换热的计算可采用牛顿冷却公式。根据对流换热定律，可以计算出从壁面某处进入通风风流的显热热流密度： (3) 式中：Tw= 巷道壁面的温度； T= 巷道内风流的平均温度；= 巷道壁面的换热系数。在围岩与风流的热交换过程中，多半是井巷低温风流流经高温岩壁，井巷壁面向风流放热，所以矿内常把上式中的对流换热系数（）称为巷壁与风流的换热系数，简称为放热系数。TTw圆形巷道（柱体）围岩与风流换热控制方程地热通过围岩向风流的传热现象与围岩本身的热传导、巷道壁面向风流的对流换热以及壁面上的水分蒸发等因素有关。由于实际情况下围岩的散热是一个很复杂的过程，为了方便本论文的研究，对要研究的物理模型做了简化和假设：1) 巷道为圆形、无限扩展，围岩岩石均质、各向同性；2) 不考虑围岩壁面的热辐射作用。根据上述假设，可得到描述考虑壁面水分蒸发时围岩与风流热质传递的数学方程，如式（3-1）： (3-1)式中：R调热圈半径，m；其他符号的意义同前章所述。图3-1 巷道围岩内节点划分0RNRP1PNRPI+1PIPI-1RI+1RIRI-1R1R0根据简化的数学模型，可将巷道围岩划分为一系列等间距 ()的同心圆，取垂直于长轴的巷道断面角度为，如图3-1所示。第二章 物理现象的数学描述控制微分方程：把控制传热、流体流动等有关过程的规律表达成数学形式；详细、完整的推导，阅读标准的教科书；数值解法方程的形式和意义：我们这里所提到的所有方程都具有一个共同的形式；形式上的一致是构成一个通用解法的基础。21 控制微分方程211 微分方程的意义守恒原理：各个微分方程各自代表着一定的守恒原理每一个方程以一定的物理量作为它的因变量，方程本身则代表着那些影响该因变量的各个因素之间必定存在着的某种平衡这些微分方程的因变量通常具有“比”的性质，即以单位质量为基础来表示各因变量。这种因变量的例子有：质量分量、速度(即单位质量的动量)以及比焓 这一类微分方程的各项代表着以单位容积为基础的效果例： 设想J表示一个典型因变量的流量密度让我们考虑如图21所示的尺寸为dx、dy及dz的控制容积(J在x方向的分量)：进入面积为dydz的流量密度： 离开与这个面相对曲面上的流量密度 通过该面的整个面积上流出的净流量是： dxdydz：所讨论区域的容积用同样的方法考虑y与z方向的贡献，有： （2.1）由于我们的数值方法是通过对一个控制容积进行平衡构成的(就象我们将要在后面看到的那样)，上述divJ的表达方式对我们来说是特别有用的 以单位容积为基础来表达一项的另一个例子是变化速率项 如果是某个“比”性质（单位质量）而是密度，那么就代表在单位容积内所包含的相应广延性质的大小于是是单位容积内有关性质的变化率 一个微分方程是这样一些项的组合；其中每一项代表一个以单位容积为基础的效应；而所有的项合在一起则反映着某种平衡或是守恒。我们现在以几个标准的微分方程为例，并由此找到一个通用的形式212 化学组分的守恒令ml代表一种化学组分l的质量分量当存在有速度场u时，可把ml的守恒表示为： （2.2）: 单位容积内化学组分l的质量变化率； uml:对流流量密度(即一般化流场u所携带的流量密度) Jl:扩散流量密度，它通常是由ml的梯度引起的 两部分流量密度(对流与扩散)的散度构成微分方程的第二项 Rl:单位容积的化学组分l的生成率.化学组分的生成系由化学反应所致。当然，依反应实际上是产生还是消毁组分l而定，Rl可能是正的，也可能是负的对于不参与化学反应的组分，Rl为0如果用菲克(Fick)扩散定律来表示扩散流量密度Jl，我们可以写出： (2.3) 其中，是扩散系数将方程(2.3)代入方程(2.2)，可以导出： (2.4)213 能量方程最通用形式表示的能量方程式含有相当数量的各种不同影响因素因为我们主要关心的是方程的形式而不是其细节，所以考虑某些限定的情况就足够了。对于可忽略粘性耗散作用的稳态的低速流，能量方程可写成 (2.5)其中h是比焓 k是导热系数 T是温度 是容积发热率根据傅里叶(Fourier)热传导定律，项代表流体中的热传导作用对理想气体以及固体与液体，我们可以应用下列关系： (2.6) 其中c是定压比热将该式代入能量方程，就得到： （2.7）如果c是常数，则h T关系可以简化为： （2.8）结果导出： （2.9）在这种情况下，不管是焓还是温度都可以选作为因变量。取方程中的速度u为0，就得到了稳态的势传导方程： （2.10）214 动量方程对于牛顿流体，控制给定方向上的动量守恒的微分方程可以写成类似线性的形式但是由于必须同时考虑切应力和正应力，加之与流体流动有关的斯托克斯粘性定律要比菲克(Fick)定律或傅里叶定律复杂，因此动量方程就要复杂得多用u表示x方向的速度，我们可以把相应的动量方程写成下列的形式： （2.11） 其中是粘度，p是压力，Bx是沿着x方向的单位容积内的体积力，Vx代表除去以所表示的粘性力项之外的其它所有粘性力项。215 紊流的时间平均方程在实际应用中，常常会遇到紊流的问题就工程中所遇到的实际问题，人们通常关心的是这种流动状态的时间平均特性因此人们通过平均运算的方法把不稳态的层流方程转化为对紊流流动的时间平均方程在进行这种平均运算的时候，人们假设：紊流中存在有相对平均值的快速而随机的脉动由平均运算所产生的附加项是所谓的雷诺（Reynolds）应力，紊流热流密度，紊流扩散流量密度等等紊流模型的任务就是确定用流动的平均性质来表示这些附加量的方法许多紊流模型采用紊流粘度或紊流扩散系数的概念来表示紊流应力及流量密度结果，紊流的时间平均方程就具有了与层流流动方程完全相同的形式。但是，诸如粘度、扩散系数以及导热系数这样一些层流交换系数则需要用相应的有效(即层流加紊流)交换系数取代。从计算机计算的观点看，按照这种格式描述的紊流即相当于具有一个相当复杂的粘度表达式的层流。 2.16紊流动能方程现今普遍流行的紊流“双方程模型”把紊流脉动动能k的方程做为其中的方程之一，该方程可以写作： （2.12）其中是k的扩散系数，G是紊流能量的生成率，以及是动能的耗散率量是方程中的净源项对变量的微分方程与此相类似2.17通用微分方程 通过上面简单地列举的一些有关微分方程，可以看出：在这里，我们所感兴超的所有的因变量似乎都服从一个通用的守恒原理如果用表示固变量，通用的微分方程就是： （2.13） 其中是扩散系数，S是源项对于特定意义的，具有特定的量和S（实际上，我们本应采用符号的，然而这样做会导致在随后的文章中应用过多的下标）上述通用微分方程中的四项分别是不稳态项、对流项、扩散项以及源项因变量可以代表各种不同的物理量，如：化学组分的质量分量、焓或温度、速度分量、紊流动能或紊流的长度尺度。与此相应，对于这些变量中的每一个都必须相对应的扩散系数以及源项S赋以适当的意义把作为扩散项：并不是所有的扩散流量密度都是受制于有关变量的梯度的。但是把作为扩散项的做法并没有把通用的方程只限于那些以梯度驱动的扩散过程凡是不能归入名义的扩散项的因子或项总是可以表示成源项的一部分。如果需要的话，我们甚至于可以把扩散系数取为0由于大多数的因变量确实需要有一个突出梯度驱动扩散这一特性的扩散项，所以我们把具有这一特性的项以显式的形式写在通用的方程中了出现在方程(213)中的密度可以通过状态方程与温度和质量分量这样一些变量相关联这些变量与速度分布都服从通用微分方程此外，流场还应当满足一个附加的约束条件，即质量守恒或连续性方程这个方程是： （2.14）上面，我们已经用向量的形式写出了有关的通用微分方程(213)和连续性方程(214)这些方程的另一种有用的表达方式是直角坐标的张量形式： （2.15） （2.16）这里下标j可以取信1、2、3分别代表三个空间坐标如果在一项内下标j重复出现，就意味着要取三项之和，例如： （2.17） （2.18）用直角坐标张量形式表达的一个直接好处是，只要简单地把下标j抹掉，就可以由这种形式得到该方程的一维形式把任何特定的微分方程改写成通用形式(213)的过程，就是把有关因变量的不稳态项、对流项及扩散项转换成共同的标准形式于是，把扩散项内梯度的系数取为对的表达式；而把方程右端的其余各项之和定义为源项S 尽管到现在为止我们一直都是把所有的变量当成有因次的量来考虑，但是，以无因次的变量来进行研究则往往比较方便。可以认为任何一个用无因次变量表示的微分方程都具有通用的形式213， (2.13) (重点)这时代表无因次的因变量，而和S分别代表扩散系数和源项的无因次形式在许多情况下，的无因次值可以简化为1，而S可以取值0或1热、质传递，流体流动，紊流以及有关的一些现象的所有各有关微分方程都可以看成是通用的方程的一个特殊情况我们需要关心的仅仅是方程(213)的数值解在编制计算机程序时，只需要写出一个求解方程(2.13)的通用程序就足够了们可以对不同意义的，重复使用这个程序；当然需要对相应的和S分别赋以各自合适的表达式，同时也需要给出相应合适的初始条件与边界条件通用的方程的概念使我们能够列出一个通用数值方法的公式，并编制通用的计算机程序22坐标的性质到现在为止，我们都在注意因变量现在我们将回过来考虑自变量（x,y,z,t），并从计算的观点讨论一下它们的性质221 自变量一般说来，因变量是三个空间坐标与时间的函数于是： (2.19)其中x，y，z以及t都是自变量并不是所有的问题都需要考虑所有四个自变量所涉及的自变量数愈少，需要计算值的位置或网格结点)也愈少一维问题：当有关的物理量只与一个空间坐标有关； 二维问题：所考虑的问题与二个空间坐标有关； 三维问题：所考虑的问题与三个空间坐标有关； 稳态问题：当所考虑的问题与时间无关时； 不稳态问题：当所考虑的问题与时间有关时。 把空间与时间两个因素放在一起考虑，我们将把一种状态说成是不稳态的一维问题，稳态的三维流动等等等温面-等温迁移法象方程(219)所表示的那种独立坐标的选择方式并不是唯一的形式代替把稳态温度分布写成T(x，y，z)，我们可以用另外一种写法： (2.20)这里z变为因变量，它代表在位置(x，y)相应于温度T的等温面高度一种采用该表达形式的方法已经由迪克斯(Dix)和西切克(Cizek)(1970)以及克兰克(Crank)及其合作者克兰克和费尔(Phahle)(1973)；克兰克和格普塔(Gupta)(1975)；克兰克和克劳利（Crowley)（1978）研究出来了这种方法叫做等温迁移法但是这种方法只限于对坐标是单调函数的温度场；对于更为一般的温度场，一定的T、x和y可能对应有几个高度z；这样，从计算的目的看，z就不适合作为因变量了2.22 坐标的合适选择恰当明智地选择坐标系统有时可以减少所得要的自变量数目由于网格结点的数目一般都与自变量的数目有关，因而以较少的自变量进行研究可以大大节省计算时间 我们用几个特殊的例子来说明坐标的选择究竟是怎样影响自变量的数目的1飞机周围的流体流动：在一个静止的坐标系上看以恒定速度飞行的飞机周围的流体流动是不稳态的；但是相对于固定在飞机上的移动坐标系而言，流动则是稳态的2在一圆管内的轴对称流动：直角坐标系内是三维的，但是在r、z的圆柱极坐标系内则是二维的，这是因为： (2.21) 与无关。3坐标变换可能用来进一步减少自变量的数量，例如：a在一块平板上的二维层流边界层给出了速度仅与有关的相似性状态，其中： (2.22)式中c是一个有因次的常数结果二维的问题便简化为一维的问题了b在半无限固体内的不稳态导热问题具有x和t这样两个自变量但是，对于某些简单的边界条件，可以把温度描写成只与有关，其中 (2.23) 式中C代表一个适当的有因次常数4改变因变量可能导致自变量数目的减少例如：a在一充分发展的通道流中，温度T与流动方向的坐标x以及横向坐标y有关但是在具有均匀壁温的热发展区，我们有： (2.24) 式中： Tb是整体温度，它随x而变化。b平面自由射流是一种二维流，但是我们可以写成： (2.25) 式中：, (2.26)这里uc代表中心线上的流速，y是横向坐标，是射流的特征宽度uc和这两个变量都随流动方向的坐标x而变化为了提高计算效率，应当选择合适的坐标系统进行数值计算2.23 单向与双向坐标 (重点)【定义】双向的坐标：如果在一个坐标上的一个给定位置处的条件，只受该位置两侧条件变化的影响的话，那么这个坐标就是一个双向的坐标 单向的坐标：如果在一个坐标上的一个给定位置处的条件，只受该位置一侧条件变化的影响，这样的坐标就是一个单向的坐标【例】通常，空间坐标是双向坐标。在一根园棒内的一维稳态导热是个双向坐标的例子在棒内任意给定点的温度会受到两端温度变化的影响时间坐标则总是单向坐标。在一块固体的不稳态冷却期间，其某一给定瞬时的温度只能受在该瞬时以前所发生的那些条件变化的影响这是一个常识问题：昨天的事情影响着今天发生的事情，但是明天的条件对于今天会发生什么事情则没有影响空间也可以作为单向坐标：在流体流动的作用下，甚至某种空间坐标也能变成非常接近于单向的坐标如果在一个坐标方向上有很强的单向流动，那么各种重要的影响只能是从上游传播到下游于是在某一点上的状态主要受其上游条件的影响，而受其下游条件的影响很小空间坐标的单向特性是一种近似确实，对流是一种单向的过程，而扩散(这个过程总是存在的)则具有双向的效应但是，当流量很大时，对流作用远远超过扩散作用，因而空间坐标就近乎是单向的了。抛物型、椭圆型与双曲型单向与双向坐标 看起来，通常用来进行微分方程分类的抛物型与椭圆型这两个数学术语相当于我们这里所用的计算方法上的概念单向与双向坐标抛物型这一术语表示一种单向的状态，而椭圆型则表示双向的概念。那么前面提到过的第三种类型，即双曲型，又是什么意思呢？碰巧双曲状态不能准确归入计算方法上的分类双曲型问题具有一种单向的特性，但是这种单向的特性并不是沿着坐标的方向，而是沿着一些称之为特征线的特殊线的方向。有一些采用特征线的数值方法，但是这些方法只限于双曲问题。另一方面，本书所提出的数值方法并没有利用双曲问题所具有的这一特殊性质我们将把双曲问题当成是椭圆问题这一总类中的一个组成部分(即全部是双向坐标)计算方法上的含义 上述关于单向和双向坐标的讨论，其动机在于：如果可以用一个单向的坐标来规定一个给定的状态，那么就有可能大大节省计算机的存储量与时间现在让我们来讨论一个不稳态的二维热传导问题我们将在计算域内构成一个二维网格结点的阵列。在任一瞬间，将存在一个相应的温度场在计算机内将对时间上顺序推移的每一瞬时计算这样一个温度场，但是，由于时间是一个单向坐标，因而某一时刻的温度场不受未来时刻温度场的影响实际上，整个不稳态问题可以简化为按要求重复进行的一个基本的步骤，即给定t时刻的温度场，求得tt时刻的温度场于是计算机内存只要供这两个温度场使用即可；对于所有的各个不同时间步可以一遍又一遍地使用同样的储存空间在这样的情况下，从一给定的初始温度场开始，我们就能够“前进”到在时间上顺序推移的各个瞬时在任何时间步，需要同时处理的未知量只有一个二维的温度数组这个温度数组与所有未来的温度值无关，而影响这一温度数组的以前瞬时的值则是已知的于是我们需要求解的就是一个简单得多的方程组，从而可以节省计算机的时间离散化方法主要任务：我们所感兴趣的现象是受一些微分方程支配的。我们已经把这些微分方程用一个变量为的通用方程全部归纳了现在我们的主要任务就是推导求解这个通用方程的方法 (2.13) 假设：为了便于理解，在本章中我们将假设变量仅仅是一个自变量x的函数3.1数值方法的本质3.11 任务数值解系与实验数据：一个微分方程的数值解系由一组构成因变量的分布的数所组成。在这个意义上讲，数值方法有点类似于在实验室中进行的实验。做实验时，仪器的一组读数为我们构成在所研究的域内被测量的物理量的分布，数值分析与实验室实验两者所获得的数据都只能是一些有限数量的数值。未知量：各个不同位置上的值让我们假设用一个关于x的多项式来代表的变化 （3.1）并采用数值方法来求得有限数量的系数把x的值以及各个a的值代入方程(31)，我们就可以计算出在任何位置x的值但是如果我们最终的兴趣是得到在各个不同位置上的值，那么这种做法就有点不太方便了。各个a值本身是没有什么特别的意义的要知道所需要的值，还必须进行前面所述的代入计算过程想法：建立一个把一系列给定点上的值作为原始的未知量的方法.实际上求解微分方程的大多数方法都属于这一类，因此我们将把我们的注意力放在这样一类方法上因而数值方法就是把计算域内有限数量位置(叫做网格结点)上的因变量值当作为基本的未知量来处理该方法的任务是提供一组关于这些未知量的代数方程并规定求解这组方程的算法3.12 离散化的概念离散化方法：把注意力集中在网格结点处的值，用离散的值取代包含在微分方程精确解中的连续信息于是我们已经离散了的分布，把这一类数值方法叫做离散化方法 图3-2 一维问题的网格节点群x 离散化方程：关于所选取的网格结点上未知值的代数方程(我们现在就把它们叫做离散化方程)系由支配的微分方程推导而得。网格结点之间的分布：推导过程中，我们必须对网格结点之间如何变化作某种假设尽管可以选择在整个计算域内满足一个简单表达式的关系作为的这种“分布”。采用分段分布是更为实际的方法，即：一定的段仅仅用一个小区域的内部及其边界上的网格结点上的值来描述该区间内的变化；于是一般将计算域分成一定数量的子域或单元每一个子域可以有一个独立的分布假设离散化的概念：这样，连续的计算域已经被离散开了这种对空间和因变量所作的系统的离散化使得我们有可能用比较容易求解的简单的代数方程取代前面提到过的控制微分方程。 一个离散化方程是连接一组网格结点处值的代数关系式;这样的一个方程系由支配的微分方程推导而得;表示与该微分方程相同的物理信息.3.13 离散化方程的结构 的分布 一定的离散化方程只与少数的几个网格结点有关，这种情况是我们选取分段分布这一特性的结果因此，在一个网格结点处的值只影响与其紧相邻的一些点上的分布可以预料到，当网格结点的数目变得很大时，离散化方程的解将趋近于相应微分方程的精确解。这是出自于这样的考虑：当网格结点紧挨在一起时，在相邻点之间的变化就变得很小，因此有关分布假设的实际细节就不那么重要了 离散化方程不是唯一的：对于一个已知的微分方程，可能的离散化方程决不是唯一的；尽管在网格结点数非常大的极限条件下，预计所有这些可能类型的离散化方程将会给出相同的解离散化方程的不同形式起因于分布假设以及推导方法的不同 有限差分和有限元法：可以认为这两种方法是离散化方法的两种可供选择的型式有限差分法与有限元法之间的区别来自选择分布和推导离散化方程的方法不同在本书内集中注意的方法具有有限差分法的外形，但是它采用了许多属于典型的有限元方法论所具有的思想3.2推导离散化方程的方法对于一个已知的微分方程，可以用许多方法推导出所要求的离散化方法。这里我们将简单地介绍几种普通的方法，然后再指出我们所喜爱的一个方法3.21泰勒级数公式通过截断泰勒级数来近似表示微分方程的导数构成。对于图31中位于结点1与结点3之间中点的结点2 ()，在2周围展开的泰勒级数给出：在第三项之后截断级数，将二个方程相加及相减： 把这二个表达式代入微分方程就推得有限差分方程这个方法含有这样的假设：的变化多少有点象是x的一个多项式，从而高阶导数是不那么重要的。但是当存在(比方说)指数形式的变化时，这种假设就可能导致人们不希望要的那种公式泰勒级数公式的推导是比较直截了当的，但是缺乏弹性并且其中各项的物理意义难以理解。3.22变分公式 3.23加权余数法 3.24控制容积公式 (略 见课本P31-34)3.3 一个说明性的例子让我们来讨论受下列方程控制的一维稳态热传导问题： （3.10） 图3-2 一维问题的网格节点群x其中k是导热系数，T是温度，以及S是单位容积的发热率准备 为了推导离散化方程，使用图3.2中所示的网格结点群l 集中注意网格结点Pl 以网格结点E及W作为它的两个邻点（E表示东侧，即正的x方向，而W表示西侧，或是负的x方向）l 虚线表示控制容积面；就目前来说，它们的准确位置是无关紧要的字母e与w代表这些面对于所考虑的一维问题，我们将假设在y与z方向为单位厚度于是，图中所示的控制容积是x11。如果我们在整个控制容积内积分方程，我们就得到： (3.10) (3.11)分布曲线的假设 需要一个分布假设或是一个内插分式图3.3中表示了两种简单的分布假设阶梯式分布：假设在一个网格结点处的T值代表它周围整个控制容积内的T值（阶梯式分布），对于这样的分布，斜率在控制容积面(即在w或e)上是不确定的分段线性的分布(图3.3b)这时，在网格结点之间采用线性的内插函数离散化方程 如果我们用分段线性分布来计算方程(3.11)中的， (3.11)所得的方程将为： (3.12)其中为S在整个控制容积内的平均值离散化方程可缩写成下列形式： (3.13)其中： 说明1. 离散化方程的标准形式: 方程(313)方程的左边: 在中心网格结点上的温度Tp；右侧:由相邻结点上的温度和常数b构成，对二维与三维的情况相邻结点的数目增加 一般说来，比较方便的是把方程(3.13)看成具有如下的形式： （3.15）式中下标nb表示一个相邻结点，表示对所有的相邻结点求和2. 分布假设: 在我们推导方程(313)时，我们已经采用了使我们能够估计值的最简单的分布假设，当然选用许多其它形式的内插函数本来也是可以的3. 此外，我们没有必要对所有的量都采用同样的分布函数，明白这一点是重要的例如，既没有必要用网格结点之间线性变化的S来计算，也没有必要由kp和ke之间线性变化的k来计算ke。4. 即便对于一个确定的变量，也没有必要对方程中所有各项都采用同样的分布函数假设例如，倘若在方程(310)中含有一个单独包含T的附加项，或许对该项使用阶梯分布函数也是可以允许的，而不必坚持如在计算时所采用的分段线性分布指导原则 上述有关选择分布函数的自由度，最终导致不同变型的离散化方程形式事实是，当网格结点的数目增加时，可以预料所有这些不同形式的方程都会给出相同的解。要求：即便是采用很粗的网格，解也总应该满足：（1）物理上的真实性；（2）总的平衡。但是，我们将在此提出一个附加的要求，这个要求的加入将使我们能够大大减少我们可以接受的公式(方程)的数目物理上的真实性是容易理解的图34所示的变化说明了这一概念一个真实的变化应当具有与准确变化相同的定性倾向在无内热源的热传导问题中，内部没有一处的温度可以超出边界温度所确定的温度范围之外当一块热的固体为绕流的流体所冷却时；固体的温度不可能降低到比该流体的温度还低我们将总是用这样的试验来检验我们的离散化方程总平衡的要求意味着对整个计算域应当满足积分守恒我们将坚持要求热流密度、质量流量以及动量通量必须准确地同相应的源和汇建立平衡，这种平衡不应当只是限于网格结点的数目变得很大时的情况，而是对于任何数目的网格结点都应该得到满足。（1）控制容积公式会使这一总的平衡成为可能，但是如同我们很快就会知道的那样，（2）在计算控制容积界面的热流密度、质量流量以及动量通量时还需要小心处置 物理上的真实性以及总的平衡这两个约束条件将用来指导我们选择分布假设以及所采用的有关措施在这些约束条件的基础上，我们将建立起几条基本的法则，应用这些法则我们就可以对现有的一些公式进行鉴别，并进而发展出一些新的公式来通常需要由数学上的考虑才能作出的一些决断，现在可以直接由物理学上的原理来进行指导了。源项的处理 一般来说，源项是因变量T本身的函数，因而在构成离散化方程的过程中需要知道这种函数关系。由于离散化方程需要用线性代数的技术来求解，因而我们在形式上只能考虑一种线性的函数关系。在下一章中将讨论对一个已知ST关系进行“线性化”处理的方法。这里把平均值表示成下列形式是足够的。 (3.16)式中，Sc代表的常数部分，而Sp是Tp的系数（显然Sp不代表在结点P所计算的S）在方程(316)中Tp的出现表明，在表示平均值时，我们已经假设：Tp值代表整个控制容积内的值，换句话说，已经采用了图3.3中中所示的阶梯式分布(应当注意：在对项采用分段线性分布的同时，我们可以自由地对源项采用阶梯式分布)。应用线性化的源项表达式，离散化方程的样子看起来仍然象方程(3.13)，但是系数的定义方程(3.24)要有所改变新的方程组是： (3.17)其中： 上述导言性的讨论为我们提供了足够的预备知识，这样，我们就可以来推导我们的离散化方程所应当服从的一些基本法则公式，以确保所得到的解满足物理上真实以及总的平衡这两个要求这些看来简单的法则具有深远的含义，它们将指导我们推演贯穿本书的那些方法3.4 四项基本法则法则1：在控制容积面上的连续性 当一个面作为两个相邻控制容积的公共面时，在这两个控制容积的离散化方程内必须用相同的表达式来表示通过该面的热流密度、质量流量以及动量通量讨论：显然，通过一个特定的面离开一个控制容积的热流通过同一个面进入第二个控制容积的热流密度相同，不然的话，总的平衡就不会满足尽管这一要求易于理解，但是稍不小心，就可能在一些微妙的地方违反。对于图3.2所示的控制容积，我们或许会用通过TW、TP及TE的二次曲线来计算界面上的热流密度。对下一个控制容积采用同一类的公式意味着：公共界面上的梯度系由不同的、与正要考察的控制容积有关的分布曲线算得所得到的 (因此热流密度)的不连续性绘于图35中另一种做法也可能导致热流密度的不连续性这时假定在给定的控制容积的各个表面上，热流密度完全为控制容积中心结点的导热系数kp所控制这样在考虑P点周围的控制容积时，在界面e处的热流密度(见图32)将表示成。而在把E作为控制容积的中心结点时，热流密度应为为了避免出现这种不连续性，记住这一点是有好处的，即：必须把通过界面上的热流密度看成是属于界面本身，而不是属于一定的控制容积的法则2：正系数 大多数的实际问题应该是这样的，即某个网格结点处的因变量值只是通过对流以及扩散的过程才受到相邻网格结点上的值的影响在其它条件不变的情况下，在一个网格结点处该因变量值的增加应当导致相邻网格结点上该值的增加（而不是减少）在方程（3.13)中，如果TE的增加必然导致TP增加的话，这就必然要求系数aE与aP具有相同的符号。挽句话说，对于通用性方程(315)，中心结点的系数aP与各相邻结点的系数anb全都必须具有相同的符号当然，让我们可以把它们全部取为正值或者全部取为负值。让我们决定这样来写我们的离散化方程，使所有的系数均为正值，这样，法则2可以表述如下：所有的系数(aP以及各相邻给点系数anb)必须总是正的 (3.13) （3.15）说明 由方程(314)所给出的系数的定义表明，我们所用的离散化方程方程(313)遵守正系数法则但是，正如我们将在后面看到的那样，有许多公式经常违反这一法则结果往往是得到一个物理上不真实的解存在负相邻结点系数就会导致这样的情况，在这种情况下，一个边界温度的增加会引起相邻网格结点上的温度的降低我们将只接受那些确保在所有情况下系数均为正的公式 (3.14)法则3：源项的负斜率线性化 要是我们研究一下方程(3.18)中系数的定义，就可以发现，即便所有相邻结点的系数均为正，由于SP项的关系，中心结点的系数aP仍可能变为负值当然，只要SP不为正值，这一危险就完全可以避免 (3.18)于是我们把法则3写成：当源项线性化为，系数SP必须总是小于或是等于0注意 这一法则看起来好像有点任意，其实并不然大多数物理过程确实在源项与因变量之间具有负的斜率关系，实际上，要是SP为正的话，物理状态就可能会变得不稳定了一个正的SP意味着，当TP增加时，源项也随着增加；如果这时没有有效的散热机构；这可能会反过来导致TP的增加，如此反复进行下去，造成温度飞升的不稳定现象从计算方法上讲，保持负的SP，使之不致产生不稳定性以及物理上的不真实解，这是至关重要的在下一章中，我们还要进一步讨论源项线性化的问题这里要充分注意到：为了使计算成功，负SP的原则是必不可少的法则4：相邻结点系数之和 控制微分方程往往只包含有因变量的导数项 （3.10）于是，如果T代表因变量，则函数T与T+c(其中c是一个任意常数)两者均满足微分方程微分方程所具有的这一特性也必定要反映在与之相对应的离散化方程中因此，当TP以及所有的Tnb都增加同一常数值时，方程(315)应当仍然适合 （3.15） 由这个要求可以得出结论：aP必须等于所有相邻结点的系数之和因此，我们可以把法则4表述如下：为了使微分方程在因变量增加一个常数之后也仍然能得到满足，我们要求： (3.19)讨论 很容易理解，方程(3.13，当源项与T无关时)确实满足这一法则 (3.13)本法则意味着：中心结点值TP是各相邻结点值Tnb的一个加权平均值与方程(3.13)不同，方程(3.17，当源项与T有关时)中的系数不服从这一法则 (3.17)当源项与T有关时:不能认为这种情况是对法则4的违背，而应该是本法则对这种情况不适用当源项与T有关时，T与Tc两者不能同时满足微分方程就是在这样的情况下，也不应当把这个法则忘记掉，而应该通过设想方程(3.17)的一个特殊情况来应用这个法则例如，如果我们取方程(3.17)中的SP为0，那么本法则就变得可以应用了，并且确实是被遵守的。温度场的确定性：当T以及T+c两者均满足微分方程时，所要求的温度场并不会变成是多值或不确定的T的值可以由适当的边界条件唯一地确定遵照法则4就可以确保：例如，倘若边界温度增加一个常数值，那么所有的温度就会准确地增加同一常数值理解法则4作用的另一方法是：若是源项不存在，而且所有的相邻结点的温度Tnb都相等，则中心结点的温度TP必定等于它们在这样的条件下，只有一个蹩脚的离散化方程才会算出TPTnb的结果来第四章 热传导41本章的对象 (略 见课本P46)42 一维稳态热传导421 基本方程作为解释四项基本法则的工具，在33节中介绍说明性例子的过程中，我们已经推导了稳态一维导热问题的离散化方程我们来复习一下那一部分的主要内容对这个问题的控制微分方程是：由这个方程可以推导出离散化方程 (4.2)其中： (4.3)网格点P、E及W示于图32中，图中还标明了各个不同的距离控制容积面位于网格点P及其相应的相邻点之间可以认为这些面的准确位置是任意的它们的位置可以有许多不同的安排方案。其中的一些将在461节中讨论就目前来说，我们将简单地认为e和w相对于网格点P、E和W的位置是已知的量Sc及SP由源项的线性化规定，其形式是： (4.4) 至于分布假设，梯度已经由T对x的分段线性的变化算得，而对于线性化的源项，假设TP值代表整个控制容积内的值当然，应当记住：只要不违背四项基本法则，选择其它形式的分布曲线是可能的，也是允许的这里，选择分布曲线的基本原则是：在四项基本法则的约束范围内宁可采纳简单一些的分布曲线，而只是在那些必要的场合，才采用高级而复杂形式的分布曲线一维热传导问题的许多重要方面仍然有待于讨论，我们该转到这些论题上来了422 网格间距非均匀网格对于在图32中所示的网格点，距离没有必要相等事实上，人们常常希望采用不均匀的网格间距这样做可以使我们能够有效地扩大我们计算的功能一般说来，只有在网格相当细的时候，我们才可能得到精确的解但是在因变量随x变化相当慢的区域内没有必要采用细的网格相反，在Tx变化比较陡的区域内则需要细的网格不均匀网格的准确度要比均匀网格的差？现在在一些人们中间似乎流行着一种错误的概念，这些人以为不均匀网格的准确度要比均匀网格的差些这样的论断是没有合理依据的网格间距应当直接与因变量在计算域内的变化方法联系起来此外，还没有一个通用的法则可以说明相邻网格间距之间最大(或最小)比值究竟应该保持多大如何设计一个合适的非均匀网格？因为在问题解决之前，Tx的分布是不知道的，我们能够怎样来设计一个合适的非均匀网格呢？（1）首先，人们通常对于所要得到的解有着某些定性的预计，他们可以从这些预计得到某些指导（2）其次，可以用初步的粗网格的解来求得Tx 变化的形式，而后可以构成一个合适的非均匀网格这就是为什么我们坚持要求我们的方法，即便是在粗网格的条件下也应该给出具有物理意义的解的原因之一。如果这种方法只是对足够细的网格才能给出合理的解，那么一个试探性的由粗网格求得的解将是无用的达到给定精确度所需要的网格点数及分布？达到给定精确度所需要的网格点数以及这些网格点在计算域内应取的分布方式与所要求解的问题的特性有关采用仅仅几个网格点进行试探性计算，为弄清有关解的情况提供了一个方便的途径最后，这恰恰是人们在实验室中相当普遍地采用的一种做法在那里，人们先进行预备性的实验或是探索性的试验，然后再应用由这些实验(或试验)所得到的数据资料来确定在最终的实验中所应安装的探头的位置与数目。423 界面导热系数界面导热系数：在方程(43)中，已经把导热系数ke用来代表通过控制容积面e的k值，类似地，kw代表界面w的k值当导热系数k是x的函数时，我们往往只知道在网格点W、P、E上的k值于是我们需要有一个如何由这些网格点值来计算界面导热系数(如ke)的规定。当然，下列讨论并不是针对均匀导热系数的情况而言的不均匀导热系数：不均匀导热系数的状况可以由材料的不均匀性(如组合的材料扳)引起即便对于均匀的材料，由于温度分布的不均匀，加上导热系数随温度而变化的函数关系，也会导致导热系数发生变化在处理对的通用微分方程时，扩散系数将起着与导热系数相同的作用在实际问题中，经常会遇到发生明显变化的情况例如，在紊流的情况下，可以代表紊流粘度，也可以表示紊流导热系数因此人们迫切希望能有个适合于不均匀的k或的公式。求得界面导热系数ke的最直截了当的方法是假设k在P点和E点之间呈线性变化。于是， (4.5)其中插入因子，是用图4.l所示的距离所定义的一个比值： (4.6)导热系数的突然变化的情况：我们马上就要说明：在某些情况下这种简单化的做法会导致相当不准确的结果；此外，这种做法不可能精确地处理在组合材料中可能遇到的导热系数的突然变化在阐明这种替代的办法时，我们认为我们主要关心的不是导热系数在界面e上的局部值我们的主要目标是要得到一个通过下式描述的界面热流密度qe的良好的表达式： (4.7)事实上，在推导离散化方程(4-2)时已经采用了这一表达式，我们想要采用的ke的表达式应当是一个可以由它表达“正确的”qe的式子让我们来讨论这样一种情况，围绕着网格点P的控制容积由具有均匀导热系数kp的材料所填满，而围绕着E点的控制容积则由导热系数为kE的