weno中文版

1迎风型WENO格式[1]JiangGuan-ShanandShuChi-Wang（舒其望）.EfficientImplementationofWeightedENOSchemes.JournalofComputationalPhysics,126,pp202-228,1996[2]LiuXu-Dong,S.OsherandT.Chan.WeightedEssentiallyNon-OscillatorySchemes.JournalofComputationalPhysics,115,pp200-212,1994[3]ShuChi-WangandS.Osher.EfficientImplementationofEssentiallyNon-OscillatoryShockCapturingSchemes.JournalofComputationalPhysics,77,pp439-471,1988[4]ShuChi-WangandS.Osher.EfficientImplementationofEssentiallyNon-OscillatoryShockCapturingSchemesII.JournalofComputationalPhysics,83,pp32-78,1989以上4篇文章是这篇关于WENO短文的参考文献，在阅读本文之前请先阅读上述参考文献，尤其是第一片。本文的目的在于提供一个关于WENO格式的中文参考，尤其是给出了当特征值为负数时数值通量表达式的形式，这点在参考文献[1]中并未明确给出。考虑一维标量双曲守恒型方程()0fuutx∂∂+=∂∂。在此，我们将只考虑它的空间离散，对于时间项稍后讨论。它的半离散守恒型格式如下()()()11221ˆˆjjjjxxfuuLufftxx+−=∂∂⎛⎞==−≈−−⎜⎟∂∂∆⎝⎠(1)其中L是空间离散算子。下面我们将通量()fu进行矢通量分裂为：()()()fufufu+−=+(2)在此的分裂规则是：()()0,0dfudfududu+−≥，可以采用多种分裂形式，例如采用Steger-Warming的分裂形式：()22(),2fuuλλελλ±∂±+==∂，其中ε是一个小量，使得f在声速点和驻点连续可微。这种分裂的精度高，并且f±满足多阶（大于2WENO格式的精度）连续可微。111111222222ˆˆˆˆˆˆ,jjjjjjffffff+−+−+++−−−=+=+(3)对于时间项，我们可以使用Runge-Kutta法来提高精度。三阶TVD型Runge-Kutta法（简写为RK-3）为：()()()()()()()()()1211221311444122333nnnnnuutLuuuutLuuuutLu+=+∆=++∆=++∆四阶非TVD型Runge-Kutta法（简写为RK-4）为：()()()()()()()()()()()()()()1213212331121211236nnnnnnuutLuuutLuuutLuuuuuutLu+=+∆=+∆=+∆=−++++∆下面我们假定()0fu′≥，讨论12ˆjf++的计算，并且略去上标“+”；当()0fu′时，12ˆjf−+与12ˆjf++关于12jx+对称。简要介绍一下迎风型WENO格式的重构，详细过程请参考Jiang和Shu的论文。其基本思想通过线性组合低阶通量得到高阶通量，这些系数叫做线性权重。WENO格式是基于ENO格式构造出来的，有关ENO格式的构造这里就不再赘述，可以参考Shu和Osher的文章。根据r阶精度的ENO格式的选取模板（即插值区域）思想，考虑迎风效应，选取初始模板为{}10jSx=，经过r-1次扩充，可能的模板有r个：{}()12,,,,0,1,1rkjkrjkrjkSxxxkr+−++−++==−LLr阶ENO格式是从上述模板中选取一个最光滑的模板rkS，在其上构造插值多项式()rkqx来逼近f在12jxx+=附近的数值通量，其数学表达式为：()()()()211122212221221ˆ,,rlrkrrrlkjkrjkkkljjjlljPxfqffqxPxaxx⎡⎤⎣⎦+−+++++=+⎛⎞∂===+∆⎜⎟∂⎝⎠∑L(4)3其中()rkPx是一个插值多项式，它满足下面的关系式：()()(){}121,,rrkjPxfxOxxxxxx++=+∆∀∈−∆(5)a2l为已知常数，可以从Shu和Osher的文章查到，例如：2411,245760aa=−=。现以三阶的ENO格式为例说明之。此时被选模板有三个：{}{}{}333021111212,,,,,,jjjjjjjjjSxxxSxxxSxxx−−−+++===它们对应的数值通量分别为：()()()333021111212,,,,,,,,jjjjjjjjjqfffqfffqfff−−−+++。首先我们来确定()()()12330210,,jjjjqfffqx−−+=，根据式(4)可以得到()()()2303320022PxqxPxaxx⎛⎞∂=+∆⎜⎟∂⎝⎠其中()30Px是在30S上构造的二次插值多项式，用Newton插值方式表示如下：()()()2212301112222jjjjjjjjfffffPxfxxxxxx−−−−−−−−+=+−+−∆∆于是我们可以得到()()()()2212301111222122224jjjjjjjjjjjfffffqxfxxxxfffxx−−−−−−−−−−+=+−+−−−+∆∆所以()()12330210211711,,366jjjjjjjqfffqxfff−−−−+==−+(6)同理可得()31qx和()32qx()()()()21111311122122224jjjjjjjjjjjfffffqxfxxxxfffxx+−+−+−−−+=+−+−−−+∆∆所以()()1233111111151,,663jjjjjjjqfffqxfff−+−++==−++(7)()()()()2221321112122122224jjjjjjjjjjjfffffqxfxxxxfffxx++++++++−−+=+−+−−−+∆∆所以()()1233212212151,,366jjjjjjjqfffqxfff+++++==+−(8)ENO的思想就是在若干个模板里选出一个最光滑的，并由它求出数值通量。在包含间断的模板里，我们非常希望采用这种方法，因为它能排除包含间断模板中的不是很精确甚至完全不精确的信息。然而在光滑的模板里面，这种做法却不4是我们所希望的，因为此时所有可能的插值模板都提供着同样精确的信息，从而浪费了其他r－1个模板。事实上对于r阶ENO格式而言，有r个被选模板，如果将它们合并在一起将形成一个包含2r－1个格点的插值区域，按照式子(4)可以给出网格界面处2r－1阶的近似：()()1122212111ˆ,,rrjrjrjjfqffqx−−−++−++==L(9)下面，我们依旧以三阶ENO为例，它的三个模板一共包含五个点，所以可以构造一个四次插值多项式：()()()()()21122112522342112211234881630161224224641224jjjjjjjjjjjjjjjjjjjjjjjfffffffffPxfxxxxfffffffffxxxxxxxx++−−++−−++−−++−−−+−−+−+=−−−×∆∆−+−−+−+−+−+−∆∆则数值通量函数可以构造为：()()()()254555242424PxPxqxPxaxaxxx⎛⎞⎛⎞∂∂=+∆+∆⎜⎟⎜⎟∂∂⎝⎠⎝⎠求得界面处的数值通量：()()11225521121134791ˆ3060602020jjjjjjjfqfqxfffff−−++++===−++−(10)它正好是五阶上游中心格式的数值通量，但是它不能排除包含间断的模板所带来的不精确信息，用它所计算的通量在间断附近会产生非物理的振荡。为此，构造WENO格式以解决此问题。WENO格式的数值通量表达式为：()12110ˆ,,rrkkjkrjkjkfqffω−+−+++==∑L(11)其中()1,,rkjkrjkqff+−++L由(4)确定，kω是赋予相应模板的权重，问题的关键是如何选择kω。为了获得基本无振荡的性质，应该使kω与相应模板的光滑性联系起来，使得包含间断的模板被赋予几乎为零的权重。同时要求在光滑的模板中权重的分布将使得式(11)与式(9)的通量相接近。其实可以组合r－1阶数值通量来求得2r－1阶数值通量：()()1211110,,,,rrrrkkjkrjkjrjrkCqffqff−−+−++−++−==∑LL(12)5通过简单的代数运算可以确定这些系数，并且它们满足101,0rrrkkkCC−==∑。对于三阶ENO格式而言，把式(6)(7)(8)和(10)带入(12)就可以确定这些系数，它们的值为：333012133,,10510CCC===(13)由式(11)和(12)等号两边相加，可得()()()121211110ˆ,,,,rrrrjrjrkkkjkrjkjkfqffCqffω−−−++−+−+++==+−∑LL因为101rrkkC−==∑，如果要求101rkkω−==∑，那么上式就可以改写为：()()()()11221211110ˆ,,,,rrrrjrjrkkkjkrjkjjkfqffCqfffuω−−−++−+−++++=⎡⎤=+−−⎢⎥⎣⎦∑LL(14)上式右端第一项是()12jfu+的2r-1阶近似，其精度为()21rOx−∆。因为式(14)第二项的求和符号内的两项的乘积的第二项的精度是()rOx∆，所以如果再要求权重系数满足()()1,0,1,1rrkkCOxkrω−=+∆=−L，那么式(14)右端第二项的精度在光滑区域将达到()21rOx−∆。也就是说如果kω的选择符合下述三个条件，那么式(11)的精度在光滑区域将有()21rOx−∆；同时此数值通量将有ENO基本无振荡的性质，此时kω称为最佳权重。①组合系数符合101,0rkkkωω−==∑；②在光滑区域内组合系数满足()1rrkkCOxω−=+∆；③任何包含间断的模板都被赋予几乎为零的权重，即0kω≈。Jiang和Shu选取如下符合上面三个条件的系数：()()011,,0,1,1rkkkkprkCkrISαωααααε−===−++++LL(15)其中，ε是一个小量，一般取10-6，用来避免分母为0；p是大于等于2的正整6数，它的作用是放大模板的不光滑性，使得包含间断的模板的权重kω更小；ISk是第k个模板的光滑性度量。Jiang和Shu采用如下方法度量模板的光滑性：()121221211jjlrrxklklxldqxISxdxdx+−−−=⎛⎞=∆⎜⎟⎝⎠∑∫(16)对于三阶ENO而言，将前面的()()3,0,1,2kqxk=带入上式，有()()()()()()2202121221111122212121312431241312124131234124jjjjjjjjjjjjjjjjjISffffffISfffffISffffff−−−−−+−+++++=−++−+=−++−=−++−+(17)当r=3时，我们验证一下kω是否符合最佳权重条件。首先，由上述各式构造的kω符合201,0kkkωω==∑。其次，分析kω的精度，将式(17)进行Taylor展开，我们可以得到：()()()()()()222350222361222352131221243131121243131221243jjjjjjjjjISxfxfxfOxISxfxfxfOxISxfxfxfOx⎛⎞′′′′′′=∆+∆−∆+∆⎜⎟⎝⎠⎛⎞′′′′′′=∆+∆+∆+∆⎜⎟⎝⎠⎛⎞′′′′′′=∆+∆−∆+∆⎜⎟⎝⎠如果0f′≠，有()()()()2210,1,2kISxfOxk′=∆+∆=；将kα进行Taylor展开，有()()()2321pkjkxfCOxα−′=∆+∆，把它带入kω的表达式，并利用Taylor展开和2301k