向前差分格式求解二维热传导方程
二维热传导方程求解
二维热传导方程求解二维热传导方程是描述平面内物体温度分布随时间变化的数学模型，被广泛应用于工业制造、城市规划和环境模拟等领域。
本文将介绍二维热传导方程的求解方法及其应用。
一、二维热传导方程的基本形式二维热传导方程可以写成以下形式：∂u/∂t = α(∂^2u/∂x^2+∂^2u/∂y^2)其中，u表示温度分布，t表示时间，x和y分别表示平面内的水平和竖直坐标，α为热传导系数。
二、二维热传导方程的求解方法为了求解二维热传导方程，需要确定初始条件和边界条件。
初始条件指在t=0时刻温度分布的初始状态，边界条件指平面内边界的温度（或热流）分布。
常见的求解方法有有限差分法、有限元法和边界元法等。
这里以有限差分法为例。
有限差分法是将待解区域划分成一个个小网格，用数值方法近似代替微分方程，然后逐步迭代求解。
假设在(x_i,y_j,t_n)处的温度为u_(i,j,n)，则可以用以下式子近似代替热传导方程：u_(i,j,n+1) = u_(i,j,n) + αΔt/Δx^2(u_(i+1,j,n)+u_(i-1,j,n)-2u_(i,j,n))+ αΔt/Δy^2(u_(i,j+1,n)+u_(i,j-1,n)-2u_(i,j,n))其中，Δt为时间步长，Δx和Δy为空间步长。
通过迭代计算，即可得到平面内任意位置随时间的温度变化规律。
三、应用实例二维热传导方程的应用范围非常广泛。
在工业制造中，可以用来分析材料的热处理过程，优化生产工艺；在城市规划中，可以用来预测城市内部的热岛效应，为城市绿化提供科学依据；在环境模拟中，可以用来模拟地下水温度变化、河流水温变化等。
例如，在炼钢过程中，需要控制钢材的温度分布，以保证钢材的物理性能。
通过建立二维热传导方程模型，可以计算出钢材表面的温度分布，进而调整生产参数，达到最佳的钢材质量。
在城市规划中，针对不同的城市形态和环境条件，可以建立相应的二维热传导方程模型，预测城市内不同区域的温度分布情况，并提出合理建议。
6第六讲 典型模型方程-热传导方的差分格式程
ADE Methods
同一时间步，同时左右扫描
6、Hopscotch Method
Comments
2、Richardaon’s Method: CTCS
3、Simple Implicit(Laasonen) Method
算子表示：
其中：
未知，每一个时间步长需要求解三对角方程组
放大因子：
4、Crank-Nicolson Method: famous
1 N 1 1 uN uN j 1 2u j j 1
where
G 1
3-D: ADI Methods
假设
∆t左移，右端合并即为C-N格式
其中 rx x 2 , ry y 2
at
at
3-D: ADI Methods
4、Splitting or Fractional-step Methods
5、ADE Methods
1-D:
或
Simple Explicit Method
差分格式的放大因子与精确解的放大因子比较
精度高
Simple Explicit Method
Simple Explicit Method: Example
Simple Explicit Method: Example
没有相位误差，幅值误差：1.88%
2-D: Crank-Nicolson Scheme
通常采用迭代方 法，需要比三对 角更多的计算资 源
3、2-D: ADI Methods
time
1/ 2 n 1 / 2 1/ 2 uin , uin 1, j , ui , j 1, j
1 n 1 n 1 uin, j 1 , ui , j , ui , j 1
二维热传导方程的可视化计算
二维热传导方程的可视化计算二维热传导方程是描述二维物体热传导过程的数学模型。
在工程领域中，通过求解二维热传导方程，可以预测物体内部的温度分布，进而进行热设计和优化。
热传导是指物体内部由高温区向低温区传递热量的过程。
二维热传导方程是基于热传导定律和能量守恒定律建立的，它可以描述物体内部温度的时空变化。
二维热传导方程的一般形式如下：∂²T/∂x² + ∂²T/∂y² = α ∂T/∂t其中，T是温度，x和y是空间坐标，t是时间，α是热扩散系数。
为了求解二维热传导方程，需要给定边界条件和初始条件。
边界条件是指在物体表面的温度分布情况，而初始条件是指在初始时刻物体内部各点的温度分布。
通常情况下，我们采用数值方法来求解二维热传导方程，其中最常用的方法是有限差分法。
有限差分法将连续的空间和时间离散化，将二维热传导方程转化为一组离散的代数方程。
在计算机中，可以使用计算软件来实现二维热传导方程的可视化计算。
首先，需要将物体的几何形状离散化为一个个小区域，然后对每个小区域进行温度计算。
在计算过程中，可以使用迭代方法来逐步求解离散方程，直到达到收敛条件。
通过迭代计算，可以得到物体在不同时间点的温度分布情况。
在可视化计算中，可以将温度用不同的颜色表示，从而直观地显示物体内部的温度分布。
通过观察温度分布的变化，可以了解物体的热传导特性，并对其进行优化设计。
除了温度分布的可视化，还可以计算物体的热流量、热传导速率等热学参数。
这些参数对于热设计和工程优化非常重要，可以帮助工程师在设计过程中做出准确的决策。
二维热传导方程的可视化计算在工程领域中具有重要的应用价值。
通过求解二维热传导方程，可以预测物体内部的温度分布，为工程设计提供参考依据。
同时，可视化计算也为工程师提供了直观的数据展示方式，帮助他们更好地理解和分析热传导过程。
边界节点法计算二维瞬态热传导问题
Vol．35，No． 2，Feb．15，2014 2014 年2 月15 日出版文章编号: 1000-0887( 2014) 02-0111-10 应用数学和力学编委会，ISSN 1000-0887* 边界节点法计算二维瞬态热传导问题师晋红，傅卓佳，陈文( 1．河海大学力学与材料学院，南京210098;2．河海大学水文水资源与水利工程科学国家重点实验室，南京210098)( 我刊编委陈文来稿)摘要: 采用边界节点法( BKM) 结合双重互易法( DＲM) 求解二维瞬态热传导问题．采用差分格式处理时间变量，可将原瞬态热传导方程转化为一系列非齐次修正的Helmholtz 方程．随后，方程的解可分为特解和齐次解两部分计算，引入双重互易法在区域内部配点求解方程的特解，采用边界节点法仅需边界配点求解方程的齐次解．给出的数值算例显示该方法计算精度高，适用性好，具有很好的稳定性和收敛性，适合求解瞬态热传导问题．关键词: 瞬态热传导; 边界节点法; 双重互易法; 差分格式; 修正的Helmholtz 方程中图分类号: O241．82; O343．6文献标志码: Adoi: 10．3879 /j．issn． 1000-0887．2014．02． 001引言近年来，功能梯度材料、形状记忆合金等许多新兴材料被广泛用于高温涡轮发动机\［1 \］及伸缩翼系统\［2 \］的热障涂层等航天工程中．这些材料经常需要在高温环境中长期连续工作，所以很有必要了解这些材料及相关结构随时间演变的热传导性能，因此亟需发展计算瞬态热传导问题的数值算法，以便为航天工程材料的研发设计提供参考与指导．目前，很多数值算法\［3-9 \］已用于求解瞬态热传导问题，例如，有限差分法、有限元法、边界元方法等．有限差分法和有限元法在求解问题时要求在整个求解域上进行离散，在处理复杂几何区域等问题时存在计算量过大、网格生成困难等难题．边界元方法是一种边界型数值算法，将计算域降低了一维，然而边界元法选用含奇异性的基本解作为插值基函数，不可避免地需要处理费时费力的奇异积分计算问题．最近几十年发展起来的无网格方法\［10 \］，基于点的近似，克服网格依赖的缺陷，减少因网格畸变引起的计算困难，已广泛用于求解瞬态热传导问题\［5-7，11-12 \］．其中，基本解方法\［13 \］是* 收稿日期:基金项目:2013-10-07; 修订日期: 2013-11-19国家重点基础研究发展计划( 973 计划)( 2010CB832702 ) ; 国家杰出青年基金( 11125208) ; 国家自然科学基金( 11372097; 11302069 ) ; 高等学校学科创新引智计划( “111”计划) ( B12032)师晋红( 1989—) ，女，山西人，硕士生( E-mail: shijhhhu@ 163． com) ;作者简介:陈文( 1967—) ，男，江苏人，教授，博士，博士生导师( 通讯作者．edu． cn) ．111E-mail: chenwen@ hhu．112师 晋 红 傅 卓 佳 陈 文与边界元法相对应的一种无网格方法． 该方法无需积分，编程容易，不需要生成网格; 并只需在计算域边界上配点，且比边界元法的数值收敛速度快得多． 然而为了避免基本解的源点奇异 性，该法在物理边界外引入了虚假边界． 虚假边界的设置有相当大的随意性，对于复杂几何区 域和多连通区域，易造成计算不稳定．为了避免基本解方法布置虚假边界的问题，陈文等提出了边界节点法\［14-15 \］( boundaryknot method ，BKM) ． 该方法使用控制微分方程的非奇异径向基函数通解代替奇异基本解，在克 服基本解方法虚假边界问题的同时，保留了该法其他优点． 边界节点法计算二维和三维几何复 杂域各类物理问题的精度和稳定性都很高\［16 \］，有关论文总共已被他人引用 400 余次． 目前该方法已被广泛应用于声学和扩散问题\［17-18 \］．\［19 \］结合边界型数值算法被 在过去几十年中，双重互易法( dual reciprocity method ，D ＲM) 用来处理偏微分方程的非齐次项，并取得良好效果． 尤其是双重互易边界元法( dual reciprocity boundary element method ，D ＲBEM) 在边界元领域就深受大家欢迎，被广泛用于求解瞬态热传导问题\［19-21 \］． 双重互易法的基本想法是用一系列的近似特解来得到问题的一个特解，使用方便、高效、灵活，是求方程特解的一个有效方法．基于边界节点法和双重互易法近年来取得的众多成果，本文将采用边界节点法结合双重 互易法求解二维瞬态热传导问题． 文章首先采用差分格式处理热传导方程中的时间变量，将原 热传导方程转化为一系列非齐次修正的 Helmholtz 方程． 方程的解可分为特解和齐次解两部 分． 本文采用双重互易法求解方程的特解，采用边界节点法求解方程的齐次解． 文章第 1 节描 述热传导问题，并采用差分格式转化原热传导方程; 第 2 节介绍双重互易法和边界节点法; 第 3 节给出 3 个二维问题的数值算例并分析讨论; 第 4 节给出结论．1 问 题 描 述设热传导问题的考察区域为 Ω，Γ Ω 为其边界． 则问题的控制方程可表示如下:= ρc T ( x ，t ) ，k 2 T ( x ，t ) + Q ( x ，t ) = ( 1)x ∈ Ω．t 方程( 1) 受以下边界条件及初始条件的约束:－T ( x ，t ) q ( x ，t ) T ( x ，0)= T ( x ，t ) ， x ∈ ΓD ， ( 2) ( 3)( 4)为内热源．T ( x ，t ) = q －( x ，t ) ， ，= x ∈ Γ n ( x )N = T 0 ，x ∈ Ω， 其中 T ( x ，t ) 表示 t 时刻 x 点的温度，k 为导热系数，ρ 为密度，c 为比热容，Q ( x ，t ) －T ( x ，t ) ，q －( x ，t ) 分别表示边界上的温度和热流，T 表示初始温度，Γ 为 Dirichlet 边界，Γ 为0 D N Neumann 边界，Γ ΓD + ΓN ．= 采用差分格式处理原控制方程中的时间变量，对于时间间隔 ［t n ，t n +1］ ［0，t ］，物体中 任一点的温度 T ( x ，t ) 及其对时间的一阶导数和内热源可以近似为\［22 \］T ( x ，t ) = T n +1( x ) ，( 5) ( 6) ( 7)n +1 n T ( x ，t ) t T ( x ) － T ( x ) ， = τ Q ( x ，t ) = Qn +1 ( x ) ，边界节点法计算二维瞬态热传导问题 113其中，τ t n+1 － t n为迭代过程中的时间步长．= 将方程( 5) 和( 6) 代入到方程( 1) 中，整理后可得(2ρc)T n+1( x ) = － ρc T n( x ) － Q n+1 ( x ) ． 1 ( 8) － τk k τk方程( 8) 中的右端项包含前一个时间步 t = t n 的近似解 T n ( x ) ，本文选取已知的初始条件T ( x ， 0) = T 0 作为迭代的第一步．方程( 8) 可简化为( 2 － λ2 ) T ( x )= f ( x ) ，( 9)其中ρc ， ( 10)=λ 槡f ( x ) = －ρcT n( x ) 1 －Q n+1( x ) ．( 11)τk方 法 介 绍k2 方程( 9) 的解可写成如下特解与齐次解和的形式:T ( x ) = T p ( x ) + T h ( x ) ， 式中，T p ( x ) 和 T h ( x ) 分别表示方程的特解和齐次解．特解 T p ( x ) 满足( 12)( 2 － λ2 ) T ( x ) = f ( x ) ，( 13) x ∈ Ω．p 齐次解 T h ( x ) 满足 ( 2 － λ2 ) T ( x ) = 0，x ∈ Ω，x ∈ ΓD ， x ∈ ΓN ．( 14) ( 15) ( 16)h －T h ( x ) q h ( x ) = T ( x ) － T p ( x )= q ( x ) － q p ( x )= g ( x ) ， = h ( x ) ，－双重互易法采用双重互易法求解方程的特解 T p ( x ) ． 首先用如下线性组合来逼近右端项 f ( x )N If ( x ) =2． 1 \［23 \］:∑αi φi( x ) ， x ∈ Ω，( 17)i = 1 其中 αi 为待求系数，φi ( x ) = φ( ‖x － x i ‖) 确定 f ( x ) 后，特解 T p ( x ) 可写为N IT p ( x ) =其中 Ψi ( x ) 满足为径向基函数，x i 为区域 Ω 内N I 个不同的数据点． ∑αi Ψi ( x ) ，i = 1( 18)( 2－ λ2 ) Ψ ( x ) φ ( x ) ．( 19) = i i 如下\［24 \］:这里选取径向基函数 φ( x ) = r 2ln r ．( 20)φ 根据方程( 19) ，可得到 Ψ( x ) 的取值:2{4K0 ( λr ) = － 4 － 4 l n r r ln rλ2 r ≠ 0，Ψ －－ λ4 λ4 4 λ ( 21)4 + 4γ λ4 ( 2 )4 λ ln ， r = 0，= － +Ψ λ4 λ4114师 晋 红 傅 卓 佳 陈 文其中 γ ≈ 0． 577 215 664 901 532 8 为 Euler 常数，K 0 为零阶第二类修正的 Bessel 函数． 边界节点法采用边界节点法求解方程的齐次解 T h ( x ) ，T h ( x ) N ST h ( x ) =2． 2 \［14 \］ 由非奇异通解的线性组合来逼近 : ∑βi T i *( x ) ，x ∈ Ω，( 22)i = 1* * 2 其中 βi 为待求系数，T i( x ) ( ‖x － x i ‖) 为修正的 Helmholtz 算子 ( Δ － λ 的非奇异通) = T 解，x i 为源点集合，N S 为源点总数． 对二维问题，非奇异通解 T *( x ) 为 T* ( λr ) ， ( 23)= I 0 其中I 0 为零阶第一类修正的 Bessel 函数． 方程( 22) 中的 T h ( x ) 自动满足偏微分方程( 14) ，选取系数使得 T h ( x ) ( 15) 和( 16) ． 本文采用配置点法来达到该目的，可得到如下线性方程组:N S也满足边界条件 ∑βi T i *( x j ) = g ( x j ) ， = 1，2，…，n 1 ， ( 24) j i = 1N S T *( x j ) i ∑βi = h ( x j ) ， = n 1 + 1，n 1 + 2，…，n b ，( 25)j n i = 1其中，n 1 为 Dirichlet 边界上的配置点个数，n b 为边界上的配点个数． 写成矩阵形式则为 = b ， ( 26) TA β T其中 b = ( g ( x 1 ) ，g ( x 2 ) ，…，g ( x n ) ，h ( x n +1 ) ，h ( x n +2 ) ，…，h ( x n ) ) ，β ( β1 ，β2 ，…，βN ) ， = 1 1 1bSA = ( A ij ) 为插值矩阵:{T* ( x j ) ， = 1，2，…，N S ，j = 1，2，…，n 1 ，i i ( 27)A ij = T *( x j ) ， i = 1，2，…，N ，j = n + 1，n + 2，…，n i ．S 1 1 bn因此，方程( 9) 的解可表示为N IN ST ( x ) = T p ( x ) + T h ( x ) =∑αi Ψi ( x )i = 1*( x ) ， ( 28)+ ∑βi T i i = 1x ∈ Ω．数 值 算 例本节采用边界节点法求解 3 个数值算例，从而检验该方法在求解瞬态热传导问题中的适 用性及有效性． 前两个算例所考察问题不含内热源，第 3 个算例包含内热源．3 为了量化数值解的好坏，分别引入绝对误差 E err ( T ) 和平均相对误差Ｒerr ( T ) 如下: E err ( T ) ，( 29)= T i － T i N∑( T i － T i )i = 12Ｒerr ( T ) ， ( 30)=N∑T2槡 i i = 1其中T i 和 T i 分别表示 x i 点的精确解和数值解，N 为测试点总个数． 例 1 考虑一个热传导 问题的经典算 例，考察区域为正方形区域 Ω = { ( x 1 ，x 2 ) 0 ＜ x 1 ，x 2 ＜ 3 } ，各系数取值为 ρc = 1，k = 1． 25，Q = 0． 初始时刻整个区域内的温度均边界节点法计算二维瞬态热传导问题 115为 30 ℃ ，边界温度均为 0 ℃ ，则该问题的边界条件及初始条件可以表述为 T ( 0，y ，t ) = T ( x ，0，t ) = T ( 3，y ，t ) = T ( x ，3，t ) = 0，T 0 ( x ，y ) 该问题的解析解为该问题的几何区域及边界条件见图 1．= 30． ∞∞–k π2 ( i 2 + j 2) t T ( x ，y ，t ) = ∑∑A i sin sin exp [i πx j πy ]， 3 3 23 i = 1 j = 1其中［( － 1) i－ 1］［( － 1) j－ A i = 4 × 30 × ． 2ij π目前已有很多学者采用有限元法( FEM) \［3 \］、边界元法 ( BEM) \［8 \］、双重互易边界元法\［19 \］、Trefftz 有限元法( Trefftz FEM) \［7 \］等方法求解该算例． 这里将采用本文的方法 ( D ＲBEM) 求解该算例并与上述几种方法的数值结果进行比较．图1 例 1 的考察区域及边界条件 Fig ． 1 The domain and boundaryconditions in example 1图2 Fig ． 2 例 1 区域的边界及内部配点 The boundary and interpolationpoints in example 1在区域边界均匀配点 40 个，在区域内部均匀配点 81 个，见图 2． 在计算中，时间步长取 τ = 0． 05，求 t = 1． 2 时刻的数值解． 在区域中取 5 个测试点，采用边界节点法计算结果与其他方 法的数值结果进行比较，见表 1． 从表 1 中可以看出，本文方法所得数值解与解析解吻合得很 好，并且比其他几种方法的数值精度更高．表1 = 1．2 时刻各方法数值结果比较 t Table 1 C omparison of numerical results bewteen different methods at t = 1． 2point x y BKM FEM BEM D ＲBE M Trefftz FEM exact 1 2 3 452． 4 2． 4 1． 8 1． 8 1． 51． 5 2． 4 1． 5 1． 8 1． 51． 081 1 0． 631 2 1． 778 5 1． 690 8 1． 871 01． 139 0． 670 1． 843 1． 753 1． 9381． 114 0． 657 1． 798 1． 713 1． 8871． 099 0． 645 1． 784 1． 695 1． 8771． 103 0． 660 1． 797 1． 715 1． 8941． 065 0． 626 1． 723 1． 639 1． 812为了检验该方法的收敛性，在边界取不同个数的配点分别进行计算． 在整个区域均匀分布 441 个测试点，则 t 1． 2 时刻区域的平均相对误差见图 3． 从图 3 中可以看出，数值解的精度= 随边界配点个数的增加而提高，最后趋于平稳，可见，该方法对边界配点个数是收敛的．例 2 为了检验边界节点法在求解不光滑区域问题上的适用性，考虑 L 形区域见图 4，各116师 晋 红 傅 卓 佳 陈 文图3 t = 1．2 时不同边界配点个数的数值解精度 Fig ．3 Numerical accuracy variation with respect to the number of boundary points at t = 1． 2图4 例 2 中考察区域及 A ，B 测试点 Fig ． 4 The domain and 2 test pointsA ，B in example 2图5 Fig ． 5 例 2 区域的边界及内部配点 The boundary and interpolationpoints in example 2图6 例 2 中 A 点的数值解 Fig ． 6 Numerical solution at point A in example 2图 7 例2 中 B 点的数值解 Fig ． 7 Numerical solution at point B in example 2= 1，k = 1，Q = 0． 该问题的精确解为 T ( x ，y ，t ) = cos( 3x ) sinh( 2y ) e －5t． 系数取值为 ρc 该问题所考察区域的边界及内部配点见图 5，边界配点个数为 80，内部配点个数为 261． 在计算中，时间步长取 τ = 0． 02． 在区域中取 A ，B 两个测试点见图4，A ，B 两点随时间的温度变边界节点法计算二维瞬态热传导问题 117化分别见图 6、图 7． 从图中可以看出，计算所得的数值解与精确解吻合得非常好，精度很高． 可见，对于区域不光滑的问题，该方法同样适用．当取不同时间步长时，A ，B 两点随时间变化的数值解的绝对误差见图 8． 从图中可以看 出，时间步长取值越小，数值解的精度越高，同时所需的计算时间也越多．图8 例 2 中时间步长对数值解绝对误差的影响 Fig ． 8 Effect of time step length on the absolute numerical error in example 2例 3 考虑含内热源的瞬态热传导问题，考察区域见图 9． 各系数取值为 ρc = 1，k = 1， Q ( x ，y ，t ) = sin x sin 2y ( 5cos t － sin t ) ． 该问题的精确解为T ( x ，y ，= sin x sin 2y cos t ． 在区域边界均匀配点 108 个，在区域内部均匀配点 415 个，见图 10． 在计算中，时间步长取 τ = 0． 05． 在整个区域均匀分布1 643 个测试点，t =2 时刻整个区域的数值解及绝对误差分 布分别见图 11 和图 12． 其中，最大绝对误差为 8． 196 × 10 － 5，其相应的数值解约为 0． 19，相对误差为 4． 262 × 10 － 4． 因此，可以看出，该方法同样适用于求解含内热源的瞬态热传导问题，计算结果具有很高的精度．图9 例 3 的考察区域 Fig ． 9 The domain in example 3图10 例 3 区域的边界及内部配点 Fig ． 10 The boundary and interpolation points in example 3结 论4 本文将边界节点法结合双重互易技术应用于求解瞬态热传导问题． 其中，采用差分格式处 理热传导方程中的时间变量，将原热传导方程转化为一系列修正的 Helmholtz 方程，并采用双 重互易法求解方程的特解，采用边界节点法求解方程的齐次解． 文中给出 3 个数值算例，从不118师 晋 红 傅 卓 佳 陈 文图11 t = 2 时整个区域的数值解Fig ． 11 The numerical solution in the whole domain at t = 2图 12 t = 2 时整个区域数值解的绝对误差分布Fig ． 12 The absolute numerical error distribution in the whole domain at t = 2同边界形状、是否含内热源等方面检验了该方法的适用性及有效性，并与其他方法进行比较．数值结果表明，边界节点法无需积分，编程容易、不需要生成网格和划分虚假边界，计算精度 高、适用性好，且具有很好的稳定性和收敛性，适合求解瞬态热传导问题．此方法还有望高效稳定求解瞬态热传导的三维问题及非线性问题，这将是我们以后的研 究方向．参考文献( Ｒeferences) :［1］Burris K W ，Beardsley M B ，C huzho y L ． Component having a functionally graded material coating for improved performance \［P \］． U S Patent: 6087022，2000-7-11．Ｒ e nner E ． Thermal engine \［P \］． U S Patent: 3937019，1976-2-10．Bruch J C ，Z yvoloski G ． T ransient tw o -dimensional heat conduction problems solved by the finite ele-m ent m etho d \［J\］． I nt ernat ional Journal f or N umerical M et hods in Engineering ，1974，8( 3) : 481-494．B lobner J ，Bia ecki Ｒ A ，Kuhn G ． T ransient non-linear heat conduction-radiation pro blems —a bo und- ary elem ent fo rmulatio n \［J\］． I nt ernational Journal for N umeric al M et hods in Engineering ，1999，46( 11) : 1865-1882．［2］ ［3］［4］边界节点法计算二维瞬态热传导问题 119［5］Li Q H ，C h en S S ，Kou G X ． Tra ns ient hea t c onducti on a na l ys is us ing the M L PG method a nd modi-fied precise time step integration method \［J\］． Journal of Computational Physic s ，2011，230( 7) : 2736- 2750．Valtchev S S ，Ｒ oberty N C ． A time-marching M FS scheme for heat conduction problems\［J\］． Eng i- neering A nal y sis W i t h Boundary Elements ，2008，32( 6) : 480-493．Jirousek J ，Q in Q H ． Application of hybrid-T refftz element approach to transient heat conduction anal- ys is\［J\］． C omput ers ＆ St ruc t ures ，1996，58( 1) : 195-201．Brebbia C A ，T elles J C F ，W robel L C ． Boundary Element Techniques : T heory and A pplic ations in Eng i- neering \［M \］． Berlin: Springer-Verlag ，1984．欧阳华江． 广义热传导方程有限元算法的计算准则 \［J\］． 应用数学和力学，1992，13( 6) : 563- 571． ( O U YAN G H ua-jiang ． Criteria for finite element algorithm o f generalized heat conduction equa- tion \［J\］． A pplied M athematic s a nd M ec hanic s ，1992，13( 6) : 563-571． ( i n C hinese) )张雄，刘岩． 无网格法\［M \］． 北京: 清华大学出版社，2004． ( Z HAN G Xiong ，LIU Yan ． Mesh- less Methods \［M \］． Beijing: T singhua U n iversity Press ，2004． ( in Chinese) )Wang H ，Q in Q ，Kang Y ． A meshless model for transient heat conduction in functionally graded ma- teria ls\［J\］． C omputat ional M echanic s ，2006，38( 1) : 51-60．Fu Z J ，C hen W ，Q in Q H ． Three bounda ry meshless methods for heat conduction analysis in nonlin- ear FGM s w ith Kirchhoff and Laplace transformation \［J\］． A dvances in A pplied Mathematic s and Mechan- ic s ，2012，4( 5) : 519-542．Fairw eather G ，Karageorghis A ． The method of fundamental solutions fo r elliptic boundary value pro b- lems\［J\］． A dvances i n C omputational M at hematic s ，1998，9( 1 /2) : 69-95．Chen W ，T anaka M ． A meshless ，integratio n-free ，and boundary-only Ｒ BF technique \［J\］． Comput- er s ＆ M at hematic s W i t h A ppl i cations ，2002，43( 3) : 379-391．Fu Z J ，C hen W ，Q in Q H ． Boundary knot method for heat conduction in nonlinear functionally gra- ded ma teria l\［J\］． Engineering A nal y si s W i t h Boundary Elements ，2011，35( 5) : 729-734． Chen W ，H on Y C ． N umerical investigation on conve rgence of bounda ry knot method in the analysis of homogeneous Helmho ltz ，modified Helmholtz ，and convection-diffusion problems\［J \］． C omputer M et hods in A pplied M echanics and Engineering ，2003，192( 15) : 1859-1875．Hon Y C ，C hen W ． Boundary knot method fo r 2D and 3D Helmholtz and convection-diffusion prob-lems under complica ted geometry \［J\］． I nt ernat ional Journal for N umeric al M et hods in Engineering ，2003，56( 13) : 1931-1948． Jin B ，Z heng Y ． Boundary knot method for some inverse problems associated w ith the Helmholtz e-qua tion \［J\］． I nt ernat ional Journal f or N umerical M et hods in Engineering ，2005，62( 12) : 1636-1651．P artridge P W ，Brebbia C A ，Wro bel L C ． The Dual Ｒec iproc ity Boundary Element M ethod \［M \］．Southa mpton: Computa ti ona l M e c ha nics Publ ica ti ons ，1992．Bia ecki Ｒ A ，Jurga s＇ P ，Kuhn G ． D ual reciprocity BEM w ithout m atrix inversion fo r transient heat conduction \［J\］． Engineering A nal y si s W i t h Boundary Elements ，2002，26( 3) : 227-236．Bulgakov V ， arler B ，Kuhn G ． Iterative solution of systems of equations in the dual reciprocity bound- ary element method for the diffusion equation \［J \］． International J ournal for Numeric al Methods in Eng i- neering ，1998，43( 4) : 713-732． Chapko Ｒ ，Kress Ｒ ． Ｒ othe ＇s metho d for the heat equation and boundary integral equations\［J \］．Journal of I ntegral Equat ions a nd A ppli c ations ，1997，9( 1) : 47-69． Golberg M ． Ｒ e cent developments in the numerical evaluation of particular solutions in the boundary el-［6］［7］［8］［9］［10］［11］ ［12］［13］［14］ ［15］ ［16］ ［17］ ［18］ ［19］ ［20］［21］［22］［23］120师晋红傅卓佳陈文ement method \［J\］．A ppl i ed M at hematic s and C omputat ion，1996，75( 1) : 91-101．［24］ Chen C ，Ｒ ashed Y F． Evaluation o f thin plate spline based particular solutions fo r Helm holtz-type o p-era tors for the DＲ M \［J\］．M echanic s Ｒesearch C ommunicat ions，1998，25( 2) : 195-201．Boundary Knot Method for 2D Transient HeatConduction P roblem sSHI Jin-hong， FU Zhuo -jia， CHEN Wen( 1．C oll e ge of M echanic s and M at eri a l s，H o hai U ni v ersi t y，N anjing210098，P．Ｒ． C hina;2．Sta te Key Laboratory of H ydrolog y-W a ter Ｒesourc es and H ydraul i c Eng ineering，H ohaiU niversity，N a njing 210098，P．Ｒ． C hina)( Contributed by CHEN Wen，M ．AM M Editorial Boa rd)Abstract: The bounda ry knot method ( BKM ) in c onjunction w i th t he dual reciprocity method ( DＲM) w a s i ntroduce d to solve 2D tra nsi ent he at conduction proble ms．With the fini te diffe rence scheme appli ed to dea l w ith the time deri vative term ，the t ra nsie nt heat c onduction equation w as conver-ted to a se t of nonhomogeneous modifie d Helmholtz equati ons． Then the numerical sol ution to the non-homogeneous problems w a s divide d int o tw o pa rts: the partic ula r soluti on and the homoge neoussolution．The DＲM w i th few inner int erpol a tion nodes w a s empl oye d to get the pa rticul a r s ol uti on，a nd the BKM w i th bounda ry-onl y nodes us ed to obtain the homogeneous solution．N u m erical res ults show that the prese nt c ombi ned method has t he merits of high acc urac y，w ide a pplicabi lity，good sta bil it y and ra pi d conve rge nce，w h ic h w e re appe ali ng t o sol vi ng t ra nsie nt heat conduction proble ms．Key wor ds: transie nt heat c onduction; boundary knot method; dual reci procit y method; di ffere nce scheme; modi fied Helmholtz equati onFoundation item: The N ational Basic Ｒe sear ch Program of China ( 973 Program) ( 2010CB832702) ;The N ational Scie nce Fund for Distinguishe d Young Sc hola rs of China( 11125208) ;The N a tional N a tural Science Foundation of China( 11372097; 11302069)。
第九章_热传导方程的差分解法_郑大昉
类似地,其偏微分用差分近似为: 类似地 其偏微分用差分近似为 近似为
∂ui, j,k ui, j,k+1 − ui, j,k = ∂t τ 2 ∂ ui, j,k ui+1, j,k − 2ui, j,k + ui−1, j,k = 2 ∂x h2 ∂2ui, j,k ui, j+1,k − 2ui, j,k + ui, j−1,k = 2 ∂y h2
∂ui,k ∂x ∂ui,k − ∂x + h
(9-18)
二阶中心差商可近似为 二阶中心差商可近似为: 可近似为
∂2ui,k ∂x
即:
2
=
−
(9-19)
ui+1,k − 2ui,k + ui−1,k ∂2u = 2 2 ∂x i,k h
(9-20)
时间的一阶差商近似为 近似为: 另, 对时间的一阶差商近似为
(9-27)
u(x, y,0) = ϕ(x, y)
(9-28)
其边界条件留待后面给出 边界条件留待后面给出. 留待后面给出
差分方法 仍设空间步长 h 仍设空间步长: 空间步长 时间步长: 时间步长 空间为: 网格. 空间为 N× M 网格
τ
则:
Nh = l,
M =s h
t = kτ , k = 0,1 2,... , x = ih, i = 0,1 N ,..., y = jh, j = 0,1,..., M
∆t
∂u ∆Q = −K(x, y, z, t)∆t∆S ∂n
(9-1)
t1
t2 t1
和
Q =∫ 1
∂u dt ∫∫ K(x, y, z, t) dS ∂n (S)
二维热传导方程的差分格式
二维热传导方程的差分格式
二维热传导方程是描述热量在二维平面内传导的方程。
在数值计算中，我们通常采用差分法来求解二维热传导方程。
差分法的基本思想是将求解区域划分为若干个小区域，然后在每个小区域内近似求解热传导方程。
具体来说，我们可以采用有限差分法来离散化热传导方程，将偏导数转化为有限差分近似值，然后再用迭代方法求解离散化后的方程组。
常用的差分格式有显式差分格式、隐式差分格式、Crank-Nicolson差分格式等。
其中，显式差分格式计算简单，但是需要满足一定的时间步长条件才能保证稳定性；隐式差分格式较为稳定，但是计算量较大；Crank-Nicolson差分格式结合了显式和隐式的优点，是一种较为稳定且计算量较小的差分格式。
在实际应用中，我们可以根据具体问题的性质选择合适的差分格式来求解二维热传
导方程。
- 1 -。
【毕业设计(论文)】二维热传导方程有限差分法的MATLAB实现
第1章前言1.1问题背景在史策教授的《一维热传导方程有限差分法的MATLAB实现》和曹刚教授的《一维偏微分方程的基本解》中，对偏微分方程的解得MATLAB实现问题进行过研究，但只停留在一维中，而实际中二维和三维的应用更加广泛。
诸如粒子扩散或神经细胞的动作电位。
也可以作为某些金融现象的模型，诸如布莱克-斯科尔斯模型与Ornstein-uhlenbeck过程。
热方程及其非线性的推广形式也被应用与影响分析。
在科学和技术发展过程中，科学的理论和科学的实验一直是两种重要的科学方法和手段。
虽然这两种科学方法都有十分重要的作用，但是一些研究对象往往由于他们的特性（例如太大或太小，太快或太慢）不能精确的用理论描述或用实验手段来实现。
自从计算机出现和发展以来，模拟那些不容易观察到的现象，得到实际应用所需要的数值结果，解释各种现象的规律和基本性质。
科学计算在各门自然科学和技术科学与工程科学中其越来越大的作用，在很多重要领域中成为不可缺少的重要工具。
而科学与工程计算中最重要的内容就是求解科学研究和工程技术中出现的各种各样的偏微分方程或方程组。
解偏微分方程已经成为科学与工程计算的核心内容，包括一些大型的计算和很多已经成为常规的计算。
为什么它在当代能发挥这样大的作用呢？第一是计算机本身有了很大的发展；第二是数值求解方程的计算法有了很大的发展，这两者对人们计算能力的发展都是十分重要的。
1.2问题现状近三十年来，解偏微分方程的理论和方法有了很大的发展，而且在各个学科技术的领域中应用也愈来愈广泛，在我国，偏微分方程数值解法作为一门课程，不但在计算数学专业，而且也在其他理工科专业的研究生的大学生中开设。
同时，求解热传导方程的数值算法也取得巨大进展，特别是有限差分法方面，此算法的特点是在内边界处设计不同于整体的格式，将全局的隐式计算化为局部的分段隐式计算。
而且精度上更好。
目前，在欧美各国MATLAB的使用十分普及。
在大学的数学、工程和科学系科，MATLAB苏佳园：二维热传导方程有限差分法的MATLAB实现被用作许多课程的辅助教学手段，MATLAB也成为大学生们必不可少的计算工具，甚至是一项必须掌握的基本技能。
热传导方程的差分格式讲解
热传导方程的左分格式—上机卖习报告二零一gg年五月一维抛物方程的初边值问题分别用向前差分格式、向后差分格式、六点对称格式，求解下列问题:du d2u”(兀0) = sin兀X、0 <x <1w(0,O = z/(l，O = 0, r >0在f = 0.05,0.1和0.2时刻的数值解，并与解析解u^t) = e-7：l sm(^x)进行比较。
1差分格式形式设空间步长h = l/N，时间步长r>0, T=M T,网比r = r/h2.(1)向前差分格式向前差分格式，即Z = /C\) ‘“；=0 =心),必=吆=0其中，丿= 1,2,…,N —1,R = 1,2,…,M—l. ^r^at/h2表示网比。
(1)式可改写成如下:M*+1 = + (i-2r)Uj + rw*_! + tfj此格式为显格式。
其矩阵表达式如下：Q-2r r)r l-2r(j、r 1一2广rl吐7、厂1一2、用丿加(2)向后差分格式(1)向后差分格式，即=0=久形)上：=WN =a其中j = 12・・\N_l,k = H,M_L (2)式可改写成- rw ：［： + (l+2r )叶' -中；；=0 + 叭此种差分格式被称为隐格式。
其矩阵表达式如下：rl + 2r -r( j \ I”-r l + 2r-r l + 2r -rW.V-1-r 1 + 2广丿MJ< UN >(3) 六点对称格式六点差分格式：喟-0 _ a加:-2喟+唸;唏- 2”； +吃,—T2L戸 戸 J眄=0产久XJM=H ；=O.将(3)式改写成-g 唸;+ (1 + 时-1 昭=g 略 + (1 - 诃 * * 咯 + /其矩阵表达式如下：(1 + r -r/2<l-r r/2 ) ( j\ -r/2 l + rr/2 1-rui-r/2 l + r -r/2r/2 1-r r/2X-r 1+2匚M丿r/2 l-2r ；E >2利用MATLAB 求解问题的过程对每种差分格式依次取N = 40., r=l/1600, r=l/3200, el/6400,用 MATLAB 求解并图形比较数值解与精确解，用表格列出不同剖分时的Z?误差。