第三章-传热学数值计算方法.ppt

本章任务 前两章已阐述了获得物理现象的理论预测有一些明显的好处 且一些有趣的现象是受一些微分方程支配的 这些方程已用通用方程归纳了 接下来的任务就是推导求解这个通用方程的方法 即所谓的离散化方法 第三章离散化方法 推导前的假设 设通用变量 仅仅是一个变量的函数 一维问题 3 1数值方法的本质 3 1 1任务1 什么是微分方程的数值解 它是由一组可以构成因变量 分布的数组成的集合 即用一组数字表示待定变量在定义域内的分布 类似于在实验室中进行的实验 仪器的读数构成了所研究区域内被测物理量的分布 有限个离散点的值的集合 2 实验数据的处理 拟合关系式通常采用的方法是利用实验结果 拟合出遵循的方程 然后利用方程即可计算出任意位置点上因变量的值 eg 设 随x的变化遵循高次多项式的形式 利用实验数据 并采用数值方法求得有限数量的多项式中的各项系数值 进而得到相应的 与x遵循的方程式 由此即可求得任意点处的 值了 若最终的兴趣是得到不同位置上的 值 则此方法有些不便 因为各个系数a本身没有什么特别的意义 但要求得 还必须进行代入过程 3 数值方法及任务 建立一个把一系列给定点上的 值作为原始未知量的方程 实际上求解微分方程的多数方法均属此类 该方法的任务是提供一组关于这些未知量的代数方程并规定求解这组方程的算法 3 1 2离散化的概念 1 离散化方法及基本思想 把注意力集中在网格节点处的值 用离散的值取代包含在微分方程精确解中的连续信息 这样就离散了 的分布 这类数值方法叫离散化方法 根据实际研究对象 把定义域分为若干个有限的区域 在定义域内连续变化的待求变量场 由有限区域上的若干个点的待求变量值来表示 这就是离散化的基本思想 2 离散化方程 所取网格节点上未知因变量 值的代数方程 此方程由支配 的微分方程推导而得 在推导过程中 需对网格节点之间 如何变化作某种假设 变量 在节点间的分布形式不同 推导离散化方程的方法也就不同 另外可以选择在整个计算域内满足一个简单表达式的分布 更为实际的方法还是采用分段分布 即将计算区域分布一定数量的子域或单元 每个子域可以有一个独立的分布假设 3 1 3离散化方程的结构 1 离散化方程的结构 一个离散化方程是连接一组网格节点处 值的代数关系式 由支配 的微分方程推导而得 并表示与该微分方程相同的物理信息 当节点数很多时 离散方程的解接近于相应微分方程的精确解 相邻点之间 变化很小 有关 分段分布的细节就不那么重要了 相应于一个已知的微分方程 离散化方程的形式决不是唯一的 这起因于分布假设以及推导方法的不同 网格节点数非常多的极限条件下 所有可能类型的离散化方程将会给出相同的解 2 离散化方法 常见的方法主要有 有限差分法和有限元法 两种方法的区别来自于选择分布和推导离散化方程的方法不同 本书主要关注的方法具有有限差分的外形 但它采用了典型的有限元方法所具有的思想 把此方法叫有限差分法可能在于它坚持遵守习惯的有限差分法做法 3 2推导离散化方程的方法 对于一个已知的微分方程 可以用许多方法推导出所要求的离散化方程 3 2 1泰勒级数公式 1 定义 在有限差分法中 通过把控制方程中的各阶导数用相应的差分表达式来代替而形成离散方程 各阶导数的差分表达式可由泰勒级数展开而得 把这种建立离散方程的方法称为泰勒级数展开法 2 差分方程式的建立 节点i两侧分别有i 2 i 1 i 1 i 2 各节点间距都为h 用泰勒级数展开有 取左端及右端的前三项 并进行相加或相减 便可得中心差分的近似式 剩余项的最低阶导数前系数的次数 用同样的方法可以得到略去截断误差O h 的差分计算式 为了提高精度 可以得到截断误差更高阶的差分表达式 3 几点说明 差分表达式分子项系数的代数和为零 各阶导数差分表达式的量纲必须与导数的量纲一致 因而 一阶导数各个差分表达式的分母为 x 二阶为 x 2 给出一个差分表达式时 必须指明是对哪个点建立的 同样的节点数 不同的建格式的点 导致不同的截断误差 如 对i点只有一阶截差 但对i 1点则是二阶导数具有二阶截差的表达式 4 优缺点 推导比较直截了当 但其中各项的物理意义难以理解 3 2 2用多项式拟合法建立导数的差分表达式 导数的差分表达式也可以通过多项式的拟合来获得 相当于对未知函数的局部变化型线采用多项式来逼近 1 线性拟合 假设函数 x 在节点 i n 附近对x的变化关系近似为线性 则有 于是有 2 采用二次曲线拟合 可得到具有二阶精度的空间导数 于是有 主要用来处理对流项的高阶格式及边界条件 eg 如图所示 已知区域内部与边界节点的温度 物体的导热系数 const 试用多项式拟合法确定穿过壁面的热流密度 解 设壁面附近温度T按线性关系变化 则 如果取温度分布为二次曲线 则有 由上面三式可解得 3 2 3控制容积公式 1 控制容积法的基本思想 把计算区域分成许多互不重叠的控制容积 并使每一个网格节点都由一个控制容积所包围 对每一个控制容积积分微分方程 应用表示网格节点之间 变化的分段分布关系来计算所要求的积分 这样就得到了包含一组网格节点处 值的离散方程 2 控制容积法诱人的特征 所得结果将意味着任何一组的控制容积内 也是整个容积 诸如质量 动量以及能量的积分守恒都可以精确地得到满足 也就是说 不论网格划分的疏密情况如何 它的解都能满足控制容积的积分平衡 这个特点提供了在不失去物理上真实性的条件下 选择控制容积尺寸有更大的自由度 有限差分法 仅当网格极其细密时 离散方程才满足积分守恒 而有限体积法 即使在粗网格情况下 也显示出准确的积分守恒 有限差分法 只考虑网格点上的数值而不考虑值在网格点之间如何变化 有限体积法在寻求控制体积的积分时 必须假定值在网格点之间的分布 插值函数只用于计算控制体积的积分 得出离散方程之后 便可忘掉插值函数 上述观点 可使得对微分方程中不同的项采取不同的插值函数有完全的自由 3 积分法实施的步骤 将守恒型的控制方程在任一控制容积及时间间隔内对空间及时间积分 需要积分 须先设定待求变量在区域内的变化规律 即假定变量的分布函数 将其分布代入控制方程 在控制容积上积分 并整理成关于节点上未知值的代数方程 在控制容积积分前 须设定变量的分布规律 但得到离散方程后 节点间变量的分布规律就不再有什么意义了 因此对于不同变量可以采用不同的分布 对微分方程中的不同项可以采用不同的分布假设进行积分 3 3一个说明性的例子 一 方程的离散化 以一维稳态导热为例 一维稳态导热问题的控制方程 一维问题的网格节点群 给出网格节点群 划定控制容积 y z方向为单位长度 控制容积体积为 x 在整个控制容积内积分方程 选定未知函数及其导数的局部分布函数 分布函数 通常有阶梯式分布和分段线性分布两种 阶梯式分布 阶梯式分布 一个节点处的 值代表它周围整个控制容积的 值 它虽然简单 但不能用来计算变量在控制容积界面处的梯度值 故一般只用于源项 物性参数和变量在时域上的分布 P E W x 分段线性分布 分段线性分布 变量在网格节点间呈线性分布 可以用来计算变量的梯度 有时也用于计算变量在时域上的分布 离散化方程采用分段线性分布来计算积分 在整个控制容积内的积分平均值 上式可整理成如下形式 式中 说明 方程的推广 二维 三维的情况均适用 在推导公式时 采用了能够估算导数的最简单的分布假设 当然选用其它形式的内插函数也可以 没有必要对所有的量都采用同样的分布函数 即没有必要用网格节点之间线性变化的S来计算 也没有必要由kP和kE之间线性变化的k计算ke 对于一个确定的变量 没有必要对方程中所有各项都采用同样的分布函数假设 二 指导原则 分布函数的自由性将会导致不同变型的离散方程的形式 事实上 网格节点数的增加 所有这些不同形式的方程都会给出相同的解 附加要求 即使是采用很粗的网格 解也应该满足 物理上真实的性状和总的平衡 物理上的真实性 一个真实的变化应当具有与准确变化相同的定向性倾向 如 无内热源的热传导问题 热固体被绕流流体冷却 可以用此真实性来检验离散化方程的准确性 物理上真实与不真实的性状 总平衡的要求 对整个计算域应该满足积分守恒 要求q qm及动量通量 必须准确地同相应的源和汇建立平衡 这种平衡对任何数目的网格节点都应当得到满足 三 源项的处理 通常来讲 源项是因变量本身的函数 构成离散方程的过程中 需要知道这种函数关系 由于离散化方程需要用线性代数的技术来求解 所以 形式上只能考虑一种线性的函数关系 即 SP为TP的系数 不代表在节点P计算式的结果 TP代表整个控制容积的值 采用了阶梯式分布 应用线性化的源项表达式 离散化方程的形式一样 但系数有所改变 四 不同离散方法的比较 Taylor展开法与多项式拟合法偏重于从数学角度进行推导 把方程中的各阶导数用相应的差分式来表示 而控制容积法和平衡法则侧重于从物理观点来分析 每个离散方程都是有限大小容积上某种物理量守恒的表示式 Taylor展开法与多项式拟合法优点 易于对离散方程进行其数学特性的分析 缺点 变步长网格的离散方程形式比较复杂 导出过程的物理概念也不清晰 且不能保证所得方程具有守恒特性 控制容积法 平衡法 优点 导出过程物理概念清晰 离散系数具有一定的物理意义 并可以保证所得方程具有守恒特性 缺点 不便于对方程进行数学特性的分析 这两种方法分别展示了有限差分法与有限元法这两种数值解法的基本特点 有限容积法更具有吸引力 Taylor展开法 控制容积法 3 4四项基本法则 离散化方程应当服从的这些法则 可以确保所得的解满足物理上的真实性以及总的平衡这两个要求 法则1 在控制容积面上的连续性 在同一个界面上各物理量 及有关物性 及其一阶导数是连续的 所谓连续是指从界面两侧的两个控制容积写出的该界面上的值是相等的 即 e 如图所示 在P E两个控制容积的公共界面e上 离开P控制容积穿过e界面的q qm及动量通量 应各自等于穿过e界面进入E控制容积的相应的量 若公共界面e上的型线选择不妥 可能导致界面上连续性受到破坏如上图所示 界面上采用了二次曲线 由于从P控制容积及E控制容积来确定的二次曲线拟合点不完全相同 在e界面上 无论是 值还是其导数从P E两侧控制容积确定的值均不相等 使得格式失去守恒性 说明界面连续性的示意 若公共界面e上的导热系数选择不妥 也可能使界面上连续性被破坏 如 在给定的控制容积的各个表面上 热流密度完全为控制容积中心节点的导数系数kP或kE所控制 这样在考虑P点周围的控制容积时 在界面e处的热流密度将表示成 而在把E作为控制容积的中心节点时 界面e处的热流密度将表示成 从物理意义上看 若从界面的两侧计算所得通量不能相互抵消 相当于在界面上存在一个由计算而造成的源或汇 使总体计算误差增加 为避免出现这种不连续性 注意 必须把界面上的热流看成是属于界面本身 而不是属于一定的控制容积的 即物性参数用界面上的值 法则2 正系数 所有的系数 ap以及各相邻节点系数anb 必须总是正的 在一个网格节点处因变量值的增加 应当导致相邻网格节点上该值的增加 而不是减少 即方程 要求aE aW均与aP有相同的正负号 皆为正或负 规定取正 若相邻节点系数有正 有负 则往往不能确保得到物理上真实的解 以后的讨论只接受那些确保在所有情况下系数均为正的公式 法则3 源项的负斜率线性化 当源项线性化为时 系数SP必须满足 因为线性化后 系数aE aW没有变化 而ap变化为 若Sp 0 则中心节点系数ap有可能变为负值 若Sp 0 则可以保证中心节点系数ap 0 实际问题也表明 Sp不能为正 若如果这时没有有效的散热机构 可能会反过来使如此反复会造成温度飞升的不稳定现象