16901讲义笔记20023.13 主题:*有限体积格式 我们已经学习了偏微分方程组的有限差分格式。因此,我们今天将看一下有 限体积格式。首先让我们考虑1维传导方程: aTaT =0,v>0并且为常量 这个偏微分方程实际上是材料(即气体,固体等等)在某个体积上的守恒原 理的应用,然后再将这块体积“收缩”到无穷小。在多维情况,按照下面的方法 来做 固定体积V和表面 交=向外法线方向 假设我们使用的是能量守恒原理。也就是: 体积内能量随时间总的变化率=对体积做的功+对体积加的热量 假设我们处理的问题没有做功也没有加热,那么, 体积内能量随时间的总变化率=0 随时间的总变化率可以分为两部分: (1)能量局部的变化 (2)进入或离开体积的能量,即“流量” 如果设每单位体积的能量与温度成比例: 由假设→E≡pCT,其中p≡物质密度,c≡比热 d≡体积内随时间的变化率 乐E下戒△S=单位时间内从表面离开体积V的能量。 守恒式变为:
16.901 讲义笔记 2002.3.13 主题： * 有限体积格式 我们已经学习了偏微分方程组的有限差分格式。因此，我们今天将看一下有 限体积格式。首先让我们考虑 1 维传导方程： 0 T T v t x ∂ ∂ + = ∂ ∂ ， 并且为常量。 v > 0 这个偏微分方程实际上是材料（即气体，固体等等）在某个体积上的守恒原 理的应用，然后再将这块体积“收缩”到无穷小。在多维情况，按照下面的方法 来做。 假设我们使用的是能量守恒原理。也就是： 体积内能量随时间总的变化率=对体积做的功+对体积加的热量 假设我们处理的问题没有做功也没有加热，那么， 体积内能量随时间的总变化率=0 随时间的总变化率可以分为两部分： ⑴ 能量局部的变化 ⑵ 进入或离开体积的能量，即“流量” 如果设每单位体积的能量与温度成比例： 由假设→ E c ≡ ρ T ，其中 ρ ≡物质密度，c ≡ 比热。 V E dV t ∂ ≡ ∂ ∫∫∫ 体积内随时间的变化率 ， S ∫∫ Ev ndS ⋅ ≡ V K K w 单位时间内从表面离开体积 的能量。 守恒式变为：
2m+=0 我们可以用高斯定理(散度定理)作替换: Ev nidS V·(Ev)dl 其中(E)=2(uE)+(mE)+(oE),下=+v+ok ∫ (Ev) dv=0 VEF)=0因为是任 意的,故处处成立 将这些降为1维: (vE)=0 ot 由E≡pcT得: 最后,注意到质量为常数,可得: 因此,方程(1)变为 用T代替cT则得到我们的精确的模型方程: OT aT at ax 所以,让我们向前几步回到能量的积分方程 ∫。a+乐E戒S=0 1维的问题里,在从x到x1的“体积内”,我们有: ddds+中E元hhz=0
S 0 V E dV Ev ndS t ∂ + ⋅ = ∂ ∫∫∫ ∫∫ K K w 我们可以用高斯定理（散度定理）作替换： S ( ) V Ev ndS Ev dV ⋅ = ∇⋅ ∫∫ ∫∫∫ K K K w 其中 () () () ( Ev uE vE E xyz ω ) ∂∂∂ ∇ ⋅ ≡ + + ， ∂∂∂ K v ui vj k ≡ ++ω K K K K ( ) 0 [ ( )] V E Ev V t E Ev dV t ∂ +∇⋅ = ∂ 0 ∂ ⇒ +∇⋅ ∂ ∫∫∫ K = K 因为 是任 意的，故处处成立. 将这些降为 1 维： ( )0 E vE t x ∂ ∂ + = ∂ ∂ 由 E c ≡ ρ T 得： () ( ) cT v cT t x ρ ρ 0 ∂ ∂ + = ∂ ∂ （1） 最后，注意到质量为常数，可得： () ( ) 0 v t x ρ ρ ∂ ∂ + = ∂ ∂ 因此，方程（1）变为： ( ) ( )0 cT v cT t x ∂ ∂ + = ∂ ∂ 用 T 代替 cT 则得到我们的精确的模型方程： 0 T T v t x ∂ ∂ + = ∂ ∂ 所以，让我们向前几步回到能量的积分方程： S 0 V E dV Ev ndS t ∂ + ⋅ = ∂ ∫∫∫ ∫∫ K K w 1 维的问题里，在从 i x 到 i 1 x + 的“体积内”，我们有： Syz S 0 V E dxdydz Ev ndydz t ∂ + ⋅ = ∂ ∫∫∫ ∫∫ K K w
心下+观+E S是无关量,可以消去。 此外,=-7 由此得到 aE 假设我们给这一问题加上网格 并且假设未知数都是平均数 E edx d( eo dt dx dt 假设网格并不随时间变化 dE x5/=Ax1-,<到目前为止这是准确的 →Ax +(ED-(Ev)=0 dt 有两个问题使之变为离散方程: (1)如何随时间积分? (2)如何求(E)或(Ev)? 对于(1),我们可以用标准的时间积分(即常微分方程)方法。例如,向前欧拉 方法
i+1 1 i x x ( ) i i yz yz x x E dx S E v n E v n S t + ⎛ ⎞ ∂ + ⋅+⋅ = 0 ⎜ ⎟ ⎝ ⎠ ∂ ∫ K K KK Syz是无关量，可以消去。 此外， i x n i = − K K i 1 x n i + = K K 由此得到： i+1 i x 1 x () () 0 i i E dx Ev Ev t + ∂ + − = ∂ ∫ 假设我们给这一问题加上网格： 并且假设未知数都是平均数： i+1 i 1/2 x 1 x 2 1 i i E E x + + ≡ Δ ∫ dx i+1 i+1 i i x x 1 1 x x 2 2 ( ) i i Ed d dx Edx E x t dt dt + + ∂ ⇒ = =Δ ∂ ∫ ∫ 假设网格并不随时间变化: i+1 i 1 x 2 1 x 2 i i dE E dx x t dt + + ∂ = Δ ∂ ∫ ←到目前为止这是准确的 1 2 1 1 2 () () 0 i i i i dE x Ev Ev dt + + + ⇒Δ + − = ) Ev i+ 有两个问题使之变为离散方程： （1） 如何随时间积分？ （2） 如何求( 1 或 ？ ( ) Ev i 对于（1），我们可以用标准的时间积分（即常微分方程）方法。例如，向前欧拉 方法：
t+(Bn):1-(E)"=0 对于(2),问题有些微妙并且存在多个解(相同的wODE选项) 选项: (1)(文中没有给出。译者注) (2) 我们如何处理这个扩散问题呢? OT aT 即: ax ax 首先,扩散项是由于我们先前忽略了热量的增加才出现的: J。a+乐、印下S=热量增加 在表面热量才可以传递。典型的线性模型假设: ox(_jq=(qx+×)n kVT·n q →体积的热量增加=4AVT戒S aE d+元=事 kVTnids 注意:散度定理的应用导致了 dE +V(VE)=V·(kVT) 在1维中,由积分常数法则给出: x+(Ev)1-(Ev)=(k 现在,我们必须决定如何估计(k)和(k T T 个简单的方法就是用有限差分
1 1 2 2 1 1 1 2 () () 0 i i n n n n i i i E E x Ev E dt + + + + + − Δ +− v = 对于（2），问题有些微妙并且存在多个解（相同的 w|ODE 选项） 选项： （1） （文中没有给出。译者注） （2） 我们如何处理这个扩散问题呢？ 即： ( ) TT T v k t xx x ∂ ∂ ∂∂ + = ∂ ∂∂ ∂ 首先，扩散项是由于我们先前忽略了热量的增加才出现的： V S E dV Ev ndS t ∂ + ⋅= ∂ ∫∫∫ ∫∫ K K w 热量增加 在表面热量才可以传递。典型的线性模型假设： ( ) x y T q k x q qi q j n T q k Tn q k y ∂ ⎫ = ⎪ ∂ ⎪ ⎧ = ×+× ⋅ ⎬ ⎨ ⇒ ∂ ⎩ =∇⋅ = ⎪ ∂ ⎪⎭ K K K K ⇒体积的热量增加= S k T ndS ∇ ⋅ ∫∫ K w VSS E dV Ev ndS k T ndS t ∂ ⇒ + ⋅= ∇ ∂ ∫∫∫ ∫∫∫ ∫∫ ⋅ K K K x w 注意：散度定理的应用导致了： () ( ) E vE t ∂ + ∇⋅ =∇⋅ ∇ ∂ K k T 在 1 维中，由积分常数法则给出： i+1 i x 1 1 x () () ( ) ( ) ii i E T dx Ev Ev k k t x + + ∂ ∂ + −= − ∂ ∂ ∫ i T x ∂ ∂ 现在，我们必须决定如何估计( )i T k x ∂ ∂ 和 1 ( )i T k x + ∂ ∂ 一个简单的方法就是用有限差分：
aT (k)=k (△x,+Ax.) 若使用均匀的网格和常数k则得到: 且控制方程变为(假设k,Ax是常数): dx+(Ev4l-(Ev),=k △x dx+Ev 4-(Ev) Ax(3-27,+7 另一种办法是什么呢?
1 1 2 2 1 1 2 2 ( ) 1 ( ) 2 i i i i i i T T T k k x x x − + + − − ∂ = ∂ Δ +Δ 若使用均匀的网格和常数 k 则得到： 1 1 2 2 ( ) i i i T T T k k x x + − − ∂ = ∂ Δ 且控制方程变为（假设k ，Δx 是常数）： i+1 i 31 11 x 22 2 1 x () () ii ii i i TT TT E dx Ev Ev k k t x 2 x + + + + − − − ∂ + −= − ∂ Δ ∫ Δ i+1 i 2 2 2 x 1 3 1 x 2 2 2 () () ( 2 ) i i i ii T x x E k dx Ev Ev T T T t x + + +− ∂ Δ ∂ ∂ ⇒ + − = −+ ∂ Δ ∫ 1 2 阶中心差分 另一种办法是什么呢？