当前位置：文档之家 › 二维稳态热传导问题中的辛方法

二维稳态热传导问题中的辛方法

  • 格式：pdf
  • 大小：329.78 KB
  • 文档页数：7

第５卷第４期　２００６年ｉ２月　热　科　学　与　技　术　Ｊｏｕｒｎａｌ　ｏｆ　Ｔｈｅｒｍａｌ　Ｓｃｉｅｎｃｅ　ａｎｄ　Ｔｅｃｈｎｏｌｏｇｙ　ＶｏＩ．５　Ｎｏ．４　Ｄｅｃ．２００６　

文章编号：１６７１－８０９７（２００６）０４－０２８８－０７　

二维稳态热传导问题中的辛方法　

徐新生，　周震寰　

（大连理工大学工程力学系工业装备结构分析国家重点实验室，辽宁大连１１６０２４）　

摘要：通过　Ｊ入哈密顿体系，提出一种求解二维稳态热传导问题的辛方法。将热传导问题归结为哈密顿体系　

下的本征值和本征解问题。利用辛本征解空间的完备性，建立一套封闭的求解问题方法　这种辛方法可直接求　解各种边界条件问题，包括混合边界条件。研究结果表明，零本征值本征解描述了基本的均匀问题，而非零本　征值本征解则显示着端部效应影响特点；数值算例给出了辛本征值和本征解的一些规律和具体例子，这些数　

值例子说明了由于非均匀端部的温度和热流影响的衰减规律；这种方法也为研究其他问题提供了一条路径。　

关键词：温度场；热流密度；热传导；哈密顿体系；辛本征解；辛方法　

中图分类号：ＴＫ１２４　文献标识码：Ａ　

０　引　言　

热传导问题一直是热门问题，随着科学技术　

的发展更多的领域和方向与之紧密相关，此问题　

的研究也越来越受到关注。长期以来人们进行了　大量研究工作，提出了各种方法。ＺｕｂａｌｒＬ１　通过研　

究得到了半无限体热传导问题温度和热流的一个　

解析解。Ｏｃｈｉａｉ［２　采用了边界元方法讨论了含热　

源的稳态热传导问题，又用积分方程方法，即区域　

积分方法改进了边界元方法，并应用到了功能梯　度材料的热传导问题。段毅文［５　等在稳态热传导　

问题中给出了一种分形方法。段雅丽　］等针对二　

维热传导问题，提出了时间为三阶和空间为二阶　的无条件稳定的ＥＴＦ—ＦＤＳ—ＭＧ算法，分析了稳　

定性和收敛性。ＨｓｉｅｈＬ７　等用Ｆｏｕｒｉｅｒ变化研究了　

各向异性薄层介质的热传导问题，并将其等效地　归结为各向同性问题。Ｍａ［８　等将Ｈｓｉｅｈ的方法推　

广到了各向异性多层介质的热传导问题，提出了　

一个精确解。Ｍｏｎｔｅ［９。ｏ］在一维问题的启发下，提　

出一种问题的近似，将计算分段均匀体的热传导　

问题化为计算特征值问题，并将此方法推广到三　

维问题中。Ｙｅｎ　Ｌｌ妇等研究了瞬态热传导问题的初　

始边界，采用Ｇｒｅｅｎ函数和能量方法，分析了瞬态　

解的性质。Ｂｅｃｋ　Ｌｌ。　等用Ｌａｐｌａｃｅ变换和分离变量　方法给出两类Ｇｒｅｅｎ函数，讨论了瞬态热传导问　题的解析解。Ｍｉｌｏｓｅｒｉｃ　Ｌｌ副等通过分步求解偏微分　

方程，讨论了热脉冲引起热传导问题，并得到了一　个二维瞬态解。刘继军［１‘］等针对二维热传导方程　

问题构造了一种显式差分格式并进行了数值计　算。陈光南［１　５］等用变分原理在不规则结构网格上　

建立热流通量形式的差分格式，研究了二维热传　导方程的差分数值模拟。薛齐文Ｌｌ６］等采用一种时　

域精细算法来计算了二维双曲热传导问题。程新　

广Ｌ１　等在Ｂｉｏｔ［１８］变分原理的基础上将热传导问　

题与弹性力学比拟，提出一种变分方法。然而以上　

的方法是以温度为基本函数，属于一类变量的拉　

格朗日体系。对于热流等物理量很难直接给出，需　

通过微分或积分才能得到。本文将哈密顿体系引　

入到稳态热传导问题中。以温度为基本函数，求得　

对偶函数，该函数恰好就是问题的热流函数。从而　

建立了哈密顿正则方程组，将问题归结为辛本征　

值和本征解问题，这样形成一种辛方法。并给出了　

复杂边界条件的数值算例。　

１　稳态热传导问题的哈密顿体系　

考虑平面一长２ｚ宽２６的矩形区域。采用直角　

坐标（２，　）（见图１），温度Ｔ（ｚ，ｚ），由Ｆｏｕｒｉｅｒ定　

律，其热流密度为　

收稿日期　２００６－０７—１８｝修回日期；２００６．１０－１６．　基金项目；国家自然科学基金资助项目（１０２７２０２４）．　作者简介：徐新生（１９５７－），男，教授，博士，博士生导师，主要研究方向为计算固体力学．

　维普资讯 http://www.cqvip.com 第４期　徐新生等：二维稳态热传导问题中的辛方法　

图１　坐标系　

Ｆｉｇ．１　Ｃｏｏｒｄｉｎａｔｅ　ｓｙｓｔｅｍ　

（２，　）一一ｋ　Ｏ７＂三一矗　，　

，ｚ）一一ｋ　三一ｋＧＴ　

式中：ｋ为热传导系数；拉格朗日密度函数Ｌ（　）　

＝：＝一ｋＥ（Ｏ．Ｔ）　＋（　）　３／２＋ＴＳ，其中Ｓ为体热；　

哈密顿作用量可以写成Ｉ　Ｉ　Ｌ（Ｔ）ｄｘｄｚｔ＇ｔ　。由哈密　ｒ６　Ｊ—ｆＪ一６　顿原理得　

Ｌ（Ｔ）ｄ　ｄ　一０　（２）　

由变分方程（２）可得到通常情况下的稳态热传导　

方程，即ｋＶ２Ｔ＋Ｓ＝０。为导入哈密顿体系，以空　

间坐标ｚ为模拟时间，记（。）＝　（）。拉格朗日密度　

函数Ｌ（Ｔ，　）一一矗［（　丁）　＋Ｔ　３／２十　ｓ。温度　

函数的对偶变量为　

一　一一ｌ－ｅｊ１　（３）一＿＿－一一』　Ｌ　Ｊ　’　

可见　一　即为　向热流密度。哈密顿函数为　

Ｈ（Ｔ，　）一Ｔ　一Ｌ（Ｔ，Ｔ）一ｋ（　）。／２一　

／（２ｋ）一　。哈密顿正则方程为　

琶　一　或　ｌ　

一　：　＋ｓ　‘　

［　一　ｆ　１　，　一　

（４）　

可以看出，哈密顿正则方程（４）与稳态热传导方　

程是完全等价的。注意到在得到哈密顿正则方程　

的同时，边界条件也可以得到。不失一般性记为　

［户１Ｔ—Ｐ２　＋Ｐ２　７１］扪＝Ｐ３　（５）　

其中，ｌ一｛　，　）　为边界ａ　的外法向。当Ｐ　≠０　

和Ｐ　一０时为第一类热边条件；当Ｐ　＝０和Ｐｚ≠　

０时为第二类热边条件；当Ｐ　≠０和Ｐ　≠０时为　

第三类热边条件。考虑侧边条件（　＝士６）：［户　

十Ｐ２　］　一±６一　幻（　），其中　０和　＝１。　

若进行一个变换Ｔ＝一Ｔ＋［　∞　（６＋　）。＋　（６　

一　）　］／［４ｂ（ｐ　＋Ｐ　６）］和　＝９（Ｐ　≠一Ｐ。６），式　（４）和（５）可以归结为　

｛　｝一［　一　／矗］｛　｝＋｛　｝；ｃ户　＋　

Ｐ２ｏ．Ｔ３　：±６＝０　（６）　

式中：ｆ，＝一［　∞　（６＋　）。＋　卜∞（６一　

）。］／［－４ｂ（ｐ　＋Ｐ　ｂ）３，ｆ２一Ｓ＋ｋｅｇ‘”＋　

卜们］／Ｅ２ｂ（ｐ。＋Ｐ　６）］。记ｆ一｛ｆ　，ｆ。）　，这样问题　

转化为非齐次方程和齐次侧边条件问题。只需讨　

论齐次方程和齐次侧边条件的通解和求一个特解　

即可。　

２　辛共轭正交归一关系　

取全状态向量　一｛Ｔ，　，式（４）简写为　

ＨｈＵ。由广义分离变量法知其解为　—　』（　）　，　

此时特征方程Ｈ　一　Ｊ，其中　是本征值，　』　是本征向量，　为哈密顿算子矩阵。定义泛函　

ｒ　＜　，Ｊ，　＞：＝＝Ｉ　’．，　ｄｚ，其中　和　为两个　Ｊ一６　本征向量，．，为单位旋转矩阵（辛单位矩阵）［１　。可　

以证明Ｌｌ　：如果一是本征值，则一　也是本征　

值；本征值集合可以分成两类，即｛　：Ｒｅ（ｔ￣　）＜０　

或Ｒｅ（，ｕ　）＝０，Ｉｍ（，ｕ。）＜０）和｛　：　＝一　）；　存在辛共轭正交归一关系　

＜　：’，Ｊ，ｈｇｐ＞一　ｊ，＜　，，　。＞一一　ｆ，＜　：，，，　Ｊ　Ｊ　

＞一（　。，‘，，　＞＝０　（７）　

３　零本征值本征解　

以第二类边条件为例，即关于热流　的侧边　

绝热条件。只需讨论侧边绝热条件（以下均采用该　

条件，其他条件类似讨论）。讨论方程Ｈ　一０，可　

以得到基本的零本征值本征解　

∞一　∞＝＝＝｛１，０）　（８）　

这个本征解的物理意义是等温度场。考虑其约当　

型本征解，由Ｈ　５”一　５∞可得　５¨一｛０，一ｋ）　。　

约当型本征解应该写成　

”　５”＋ｚ　一｛ｚ，一ｋ）　（９）　

此本征解的物理意义是ｚ方向线性温度场或等热　

流场，可以证明不存在下一阶的约当型本征解。记　

＝ｊ７５ｏ），ｊ７５　＝一　”／（２ｂｋ），由式（７）知ｊ７　和　

ｊ７　满足辛共轭和归一关系。如果在平均意义下　

讨论问题，问题的解（零本征值本征解）可以写为　

一ｃ１ｊ７　＋Ｃ２ｊ７　（１Ｏ）

　维普资讯 http://www.cqvip.com 热　科　学　与　技　术　第５卷　

边界上　ｌ；。　；＝＝　士　，由辛共轭正交归一关系式　

（７）得Ｃ１　（Ｔｔ＋Ｔ　）／２三Ｔ　＋（　）伫三　

Ｔ一，一（１／ｋ）　，，Ｃ２＝２锄三２ｂ￣ｏ＿，，　三　，三　

（ｋ／２１）（Ｔｌ—Ｔ一　）。可见每个端部的温度和热流　

是不能同时给出的，应满足以上恒等式。可以说零　

本征值本征解可以描述等热流问题，或在非均匀　

问题中的在平均意义下问题的近似解（等效解）。　

因此零本征值本征解是不够的，应考虑非零本征　

值本征解。　

４　非零本征值本征解　

讨论非零本征值的特征方程（　～／ｄ）　

Ｏ（　≠Ｏ），其解为　

一｛一　Ｃｌｃｏｓ（，ｕｘ）－＋－Ｃ２ｓｉｎ（，ｕｘ）］　可以验证本征解式（１４）和（１５）满足辛共轭正交　

归一关系式（７），这样在辛几何空间问题的解总可　

以用零本征值本征解和非零本征值本征解展开，　

统一写成　

一∑４　（　）　＋５２ｂ　（　）　一　（１６）　月＝０　＝０　其中：　和　（　一０，１，２，…）为待定系数，并可由　

端部条件确定。由辛共轭正交归一关系式（７），待　

定系数可以表示为　

＝一＜　，－，　＞　，ａ　一＜　，Ｊ，　＞　一　

（　＝０，１，２，３，…）　（１７）　

一　５　非齐次方程的特解　

代入齐次侧边（绝热）条件式（５）知　

一Ｃ１ｓｉｎ（　）＋Ｃ２ＣＯＳ（　）一０；　

一Ｃ１ｓｉｎ（一ｚｂ）＋Ｃ２ＣＯＳ（一　）一０（１２）　

有非零解的条件是系数行列式为零　

ｆ－ｓｉｎ（￣ｈｏｓ，／Ａ　，ｆ—ｏ或ｓｉｎ（２　）　ｏ　ｌ　ｓｉｎ（ｚｂ）　ＣＯＳ（　）ｌ　一　’　

得到非零本征值　；　，ｃ／（２ｂ）（　＝＝＝０，士１，±２，　

士３，…）。口类本征解　

＝　｛＿　｝ｏｓ（　）＋　

可以看出本征解中由两部分组成；第一项为对称　

解（　为偶数）和第二项为反对称解（　为奇数）。卢　

类本征解也由对称解和反对称解组成：　

一　１　ｋｎｎ　

１　２ｂ　＋　由于辛本征解之间存在正交归一关系，讨论　

非齐次方程（４）时可将非齐次项按辛本征解展　

开：　ｆ＝∑Ｉ－ｄ　（　）　。　（　）＋　（ｚ）　（　）］　（１８）　

式中包括了零本征解，系数ｄ　（　）＝＜厂，－，，　

：　（ｚ）＞，ｇ　（　）　一＜／　，．，，　（　））。　

设非齐次方程特解为　

一∑［－Ｄ　（　）　（　）十Ｇ　（　）　（－ｚ）］（１９）　

可得　

Ｄ　（　）＝　Ｄ，。（　）＋ｄ　（　）；Ｇ　（　）一一　Ｇ　）＋　

ｇ　（　）　（２Ｏ）　

则　

Ｄ　（　）　一　Ｉ　ｄ　（　）　一。ｄ手；Ｇ　（　）　一　

Ｉ　ｇ　（　）Ｐ一　‘　ｄ手　（２１）　

有了特解，可得非齐次方程的解　＝　十　＝＝＝∑｛［Ｄ（　）＋ｉ２ｎｅ．ｕ．　］　。　（　）＋　

ＥＧ　（　）＋ｂ．ｅ一　］　（　）｝　（２２）　

６　端部条件的辛本征解表示　

端部（　：士　）可以简单给定温度条件或热　

流条件，也可以给出复杂的混合边界条件。这里先　

考虑混合边界条件，即每一个边部分给定温度条　

件，而另一部分给定热流条件。此时考虑一般性边　

界条件记为　

ｆ丁ｌ　一　＝Ｔ—　（ｚ）　（　≥），一　）　

【　ｌ　一一　一　（ｚ）　（　＜），一　）’　

ｆ丁ｌ　：　—Ｔｌ（ｚ）　（　≤　）　１　

：　一　（　）　（　＞　）　‘２３　¨　３　２　一　警　ｎ　，ｅ　兀一７一　一

　维普资讯 http://www.cqvip.com

合集下载

无单元 Galerkin 方法数值求解二维瞬态热传导问题

无单元 Galerkin 方法数值求解二维瞬态热传导问题

无单元 Galerkin 方法数值求解二维瞬态热传导问题

摘要：采用无单元 Galerkin 方法数值求解具有狄利克雷边界条件的二维瞬态热传导问题。首先离散该问题的时间变量，将该问题转化为与时间无关的边值问题；然后采用罚函数法处理 Dirichlet 边界条件，得到数值离散方程组，再利用 Matlab 软件求解给出的算例，结果表明该方法得到的数值结果与解析解吻合较好，该方法具有较高的计算精度和较好的收敛性。

关键词： 二维瞬态热传导问题； 无单元 Galerkin 方法； 罚函数法；

Matlab 软件中图分类号：G210.7 文献标识码：A 文章编号：20200168685

有限差分法（FDM）[1]、有限元法（FEM）[2]、边界元法（BEM） 得

[3]、无网格法[4]等数值方法是解决瞬态热传导问题的常用方 T k1

T k1

法。Burlayenko 等人[5]用梯度有限元法计算了梯度材料中的瞬T k1vd kvd k vd

y2 y

态温度。Sutradhar 和 Paulino[6]提出了一种简单的边界元法， 该法只考虑边界条件，适用于梯度材料中的三维非定常热传导问 T k1vd k1Q vd

k1vd

题。Sladek 等人[7]用无网格局部边界积分方程方法考虑梯度材料中的不稳定温度场。

无单元 Galerkin 方法是应用最为广泛的无网格方法之一。本文借鉴文献[8]应用无单元 Galerkin 方法数值求解具有狄利克雷边界条件的二维瞬态热传导问题。首先离散该问题的时间变

量，将该问题转化为与时间无关的边值问题；然后采用罚函数法

（7）

接下来，根据高斯公式化简（7）式，于是得与混合边值问

题（6）等价的变分问题：求Tk1 xH1 ，使得 ， 有

处理 Dirichlet 边界条件，得到数值离散方程组，再利用 Matlab 软件求解给出的算例，最后给出了一个数值算例来验证理论误差aTk1,v （8）Tk1vd

二维热传导方程有限元体积法

二维热传导方程有限元体积法

二维热传导方程有限元体积法

摘要：

一、二维热传导方程的概述

二、有限元体积法的基本原理

三、二维热传导方程的有限元体积法求解

四、结论与展望

正文：

一、二维热传导方程的概述

二维热传导方程是研究物体在二维空间中热传导现象的数学模型，它是基于物质的热流密度与温度梯度之间的关系建立的。在实际应用中，很多问题都涉及到二维热传导，例如电子器件的散热问题、建筑墙体的保温问题等。解决这类问题需要对二维热传导方程进行求解。

二、有限元体积法的基本原理

有限元体积法是一种求解偏微分方程的数值方法，它通过将求解域进行网格划分，将偏微分方程转化为有限元方程组，然后求解方程组得到数值解。这种方法可以有效地解决复杂区域的问题，具有较高的灵活性和适应性。

在有限元体积法中，需要构造适当的网格，并对每个网格节点分配一个体积单元。在每个体积单元内，偏微分方程被离散化为一个或多个代数方程。通过求解这些代数方程，可以得到每个体积单元内的解，进而得到整个求解域内的解。

三、二维热传导方程的有限元体积法求解 对于二维热传导方程，可以将其表示为以下形式：

T = -Q/k

其中，T 表示温度，Q 表示热流密度，k 表示热传导系数。

采用有限元体积法求解二维热传导方程，首先需要对求解域进行网格划分。在每个网格节点处，构造一个四边形单元，并在四边形内部划分多个三角形子单元。在每个三角形子单元内，根据热传导方程的离散形式，可以得到一个二维线性代数方程。

通过对所有子单元的代数方程进行求解，可以得到每个网格节点的温度值。接着，根据相邻节点的温度值，可以计算出每个网格节点的热流密度。最后，根据所有网格节点的热流密度，可以得到整个求解域内的热传导解。

四、结论与展望

有限元体积法是一种有效的求解二维热传导方程的数值方法，具有较高的灵活性和适应性。通过合理地选择网格划分和单元形状，可以提高求解精度和计算效率。在实际应用中，二维热传导方程广泛应用于各种工程领域，如电子器件散热、建筑墙体保温等。

一类二维稳态对流——扩散方程的有限差分法

一类二维稳态对流——扩散方程的有限差分法

一类二维稳态对流——扩散方程的有限差分法

一维稳态扩散方程描述了物质在一维空间中的扩散行为。然而，在某些情况下，我们需要研究物质在二维平面中的扩散行为，例如热传导、流体传输等。本文将介绍一类二维稳态对流-扩散方程的有限差分法。

二维稳态对流-扩散方程可以写作：

∇·(D∇u) + ∇·(cu) + fu = 0 —— (1)

其中，D是扩散系数，c是速度场，u是待求解的物理量，f是源项。在这个方程中，第一项表示物质的扩散项，第二项表示对流项，第三项表示源项。我们需要求解方程(1)，找到u的分布。

为了应用有限差分法来求解二维稳态对流-扩散方程，需要将二维空间离散化为一个网格。假设我们将x方向离散为Nx个等距的节点，y方向离散为Ny个等距的节点，那么我们可以得到一个(Nx+1)×(Ny+1)的网格。我们在网格节点上定义未知量u，然后将方程(1)对节点处的u进行离散化。

首先，我们对方程(1)的扩散项进行离散化。我们使用五点差分格式来近似二维Laplace算符∇·(D∇u)。对于网格节点(x,y)，我们可以得到以下差分格式：

(Dij(xi+1,yj)ui+1,j + Dij(xi-1,yj)ui-1,j +

Dij(xi,yj+1)ui,j+1 + Dij(xi,yj-1)ui,j-1 -4Dij(xi,yj)ui,j) /

∆x^2 + (Dij(xi,yj)ui,j) / ∆y^2

其中，∆x和∆y是网格步长，Dij是扩散系数。

接下来，我们对方程(1)的对流项进行离散化。我们使用中心差分格式来近似二维梯度算符∇·(cu)。对于网格节点(x,y)，我们可以得到以下差分格式：

(cxi+1/2,yj(ui+1,j - ui,j)) / ∆x + (cxi-1/2,yj(ui,j - ui-1,j)) / ∆x + (cyi,j+1/2(ui,j+1 - ui,j)) / ∆y + (cyi,j-1/2(ui,j - ui,j-1)) / ∆y 其中，cxi+1/2,yj、cxi-1/2,yj、cyi,j+1/2和cyi,j-1/2是速度场在节点(x,y)处的中心点处的x和y分量。

二维肋片稳态导热问题的数值计算

二维肋片稳态导热问题的数值计算

传热学：第四章 导热问题的数值求解

1 例题4-5 二维肋片稳态导热问题的数值计算

（1）自主编程，编程语言自定，最后提交源程序

（2）提交电子报告（word格式），包括：

（a）给出空间离散示意图（网格划分）

（b）节点离散方程

（c） 图示温度等值线（可以利用origin或matlab）

解：（a）空间离散示意图(由origin的graph作图)

（b）节点离散方程

由上图所示得各节点的节点离散方程

结点1：T(m,n)=0.25*(T(m+1,n)+T(m-1,n)+T(m,n+1)T(m,n-1))

结点2：T(M,n)= 1/(4+2*Bi)*(T(M,n-1)+T(M,n+1)+2*T(M-1,n))

结点3：T(m,N)=1/(4+2*Bi)*(T(m-1,N)+T(m+1,N)+2*T(m,N-1))

结点4：T(M,N)= 1/(2+2*Bi)*(T(M-1,1)+T(M,2))

结点5：T(m,1)=0.25*(T(m-1,1)+T(m+1,1)+2*T(m,2))

结点6：T(M,1)= 1/(2+2*Bi)*(T(M,N-1)+T(M,2))

（c） 温度等值线 传热学：第四章 导热问题的数值求解

2 等温线图，工况1（Bi=0.01）η= 0.9670

温度与y轴分布图（工况1） 传热学：第四章 导热问题的数值求解

3 工况2（Bi=1）η=0.1910

温度与y轴分布图（工况2） 传热学：第四章 导热问题的数值求解

4 图像分析：四幅图的显示来看，结果是可信的。要是网格划分过松，就会发现，在肋板顶端的绝热边界上温度的分布是有问题的，温度的最高值并不是在半肋板顶端边界n=1处，而是在n>1的不远处的离散点上，这是和我们的预期是相违背的，但是当网格划分到达一定的密度，就可以避免这个问题，虽然在图像上看不出来了，但是这个问题还是存在的，不过由于足够小的网格，是它可以忽略。