WENO3¶

WENO3 Link¶

code¶

WENO (Weighted Essentially Non-Oscillatory) schemes.

PyWENO¶

PyWENO contains two modules to help authors generate WENO codes

PyWENO.

WENO（Weighted Essentially Non-Oscillatory）格式是一种高阶精度的数值格式，广泛应用于计算流体动力学（CFD）中，特别是在求解含有激波或高梯度区域的偏微分方程时。WENO 格式通过自适应地选择光滑的子模板来构造高阶精度的数值通量，从而避免非物理振荡。

WENO 格式简介 WENO 格式的核心思想是通过加权多个低阶模板的数值通量来构造高阶精度的通量。对于三阶 WENO 格式（WENO3），通常使用两个二阶模板来构造三阶精度的数值通量。

WENO3 公式假设我们求解一维守恒律方程：

\[ \cfrac{\partial u}{\partial t}+\cfrac{\partial f}{\partial x}=0 \]

其中$f(u)$是通量函数。

通量分裂¶

首先对通量进行分裂（如 Lax-Friedrichs 分裂）：

\[ f(u)=f^{+}(u)+f^{-}(u) \]

其中：

\[ f^{+}(u)=\cfrac{1}{2}(f(u)+\alpha{u}), \quad f^{-}(u)=\cfrac{1}{2}(f(u)-\alpha{u}) \]

$\alpha$是通量的最大特征值。

模板选择¶

对于 WENO3，使用两个二阶模板：

模板 1：$\{u_{i-1},u_{i}\}$
模板 2：$\{u_{i},u_{i+1}\}$

数值通量计算¶

对于正通量$f^{+}$，使用左侧模板；对于负通量$f^{-}$，使用右侧模板。

正通量$f^{+}$:¶

\[ f^{+}_{i+1/2}=\omega_{1}f^{1}_{i+1/2}+\omega_{2}f^{2}_{i+1/2} \]

其中：

$ f_{i+1/2}^{(1)} = -\frac{1}{2} f_{i-1} + \frac{3}{2} f_i $
$ f_{i+1/2}^{(2)} = \frac{1}{2} f_i + \frac{1}{2} f_{i+1} $
权重 $ w_1 $ 和 $ w_2 $ 通过光滑度指标计算。

负通量 $ f^- $：¶

\[ f_{i+1/2}^- = w_1 f_{i+1/2}^{(1)} + w_2 f_{i+1/2}^{(2)} \]

其中：

$ f_{i+1/2}^{(1)} = \frac{1}{2} f_i + \frac{1}{2} f_{i+1} $
$ f_{i+1/2}^{(2)} = \frac{3}{2} f_{i+1} - \frac{1}{2} f_{i+2} $

4. 光滑度指标¶

光滑度指标用于计算权重 $ w_1 $ 和 $ w_2 $：

\[ \beta_1 = (u_i - u_{i-1})^2, \quad \beta_2 = (u_{i+1} - u_i)^2 \]

权重计算：

\[ w_k = \frac{\alpha_k}{\alpha_1 + \alpha_2}, \quad \alpha_k = \frac{C_k}{(\beta_k + \epsilon)^2} \]

其中 $C_1=\frac{2}{3}$, $C_2=\frac{1}{3}$ 是模板的系数，$ \epsilon $ 是一个小常数（通常取 $ 10^{-6} $）。

one-dimensional scalar conservation laws¶

ESSENTIALLY NON-OSCILLATORY AND WEIGHTED ESSENTIALLY NON-OSCILLATORY SCHEMES FOR HYPERBOLIC CONSERVATION LAWS

CHI-WANG SHU

One Space Dimension.

Reconstruction and Approximation in 1D. In this section we concentrate on the problems of interpolation and approximation in one space dimension.

Given a grid

Text Only
  [   1  |  ...    N-1    |    N    ]
  1/2   3/2             N-1/2      N+1/2
  a                                  b

\[ a=x_{\frac{1}{2}}<x_{\frac{3}{2}}<\cdots<x_{N-\frac{1}{2}}<x_{N+\frac{1}{2}}=b \]

We define cells, cell centers, and cell sizes by

\[ \begin{array}{l} I_{i}\equiv [x_{i-\frac{1}{2}}, x_{i+\frac{1}{2}}],\quad x_{i}\equiv\cfrac{1}{2}(x_{i-\frac{1}{2}}+x_{i+\frac{1}{2}})\\ \Delta \equiv x_{i+\frac{1}{2}}-x_{i-\frac{1}{2}},i=1,2,\cdots,N\\ \end{array} \]

We denote the maximum cell size by $$ \Delta x \equiv \underset{1\le i\le N}{\text{max}}\Delta x_{i} $$

Reconstruction from cell averages¶

One dimensional reconstruction

Given the cell averages of a function v(x):

\[ \bar{v}_{i}\equiv \cfrac{1}{\Delta x_{i}}\int_{x_{i-\frac{1}{2}}}^{x_{i+\frac{1}{2}}}v(\xi)d\xi, i=1,2,...,N\\ \]

find a polynomial $pi(x)$, of degree at most k − 1, for each cell $I_{i}$, such that it is a k-th order accurate approximation to the function $v(x)$ inside $I_{i}$:

\[ p_{i}(x)=v(x)+O(\Delta x^{k}),\quad x \in I_{i},\quad i=1,\cdots,N. \]

In particular, this gives approximations to the function v(x) at the cell boundaries

\[ {v}^{-}_{i+\frac{1}{2}}=p_{i}(x_{i+\frac{1}{2}}),\quad {v}^{+}_{i-\frac{1}{2}}=p_{i}(x_{i-\frac{1}{2}}) ,\quad i=1,\cdots,N. \]

which are k-th order accurate:

\[ {v}^{-}_{i+\frac{1}{2}}=v(x_{i+\frac{1}{2}})+O(\Delta x^{k}),\quad {v}^{+}_{i-\frac{1}{2}}=v(x_{i-\frac{1}{2}})+O(\Delta x^{k}) ,\quad i=1,\cdots,N. \]

Given the location $I_{i}$ and the order of accuracy $k$, we first choose a “stencil”, based on $r$ cells to the left, $s$ cells to the right, and $I_{i}$ itself if $r, s ≥ 0$, with $r + s +1= k$:

\[ S(i)\equiv \{I_{i-r},\cdots,I_{i+s}\} \]

There is a unique polynomial of degree at most $k−1 = r+s$, denoted by $p(x)$ (we will drop the subscript i when it does not cause confusion), whose cell average in each of the cells in $S(i)$ agrees with that of $v(x)$:

\[ \cfrac{1}{\Delta x_{j}}\int_{x_{j-\frac{1}{2}}}^{x_{j+\frac{1}{2}}}p(\xi)d\xi=\bar{v}_{j},\quad j=i-r,\cdots,i+s \]

\[ {v}^{-}_{i+\frac{1}{2}}=\sum_{j=0}^{k-1}c_{rj}{v}^{-}_{i-r+j},\quad {v}^{+}_{i-\frac{1}{2}}=\sum_{j=0}^{k-1}\tilde c_{rj}{v}^{-}_{i-r+j} \]

\[ \tilde{c}_{rj}=c_{r-1,j} \]

We summarize this as follows: given the $k$ cell averages

\[ \bar{v}_{i-r},\cdots,\bar{v}_{i-r+k-1} \]

there are constants $c_{rj}$ such that the reconstructed value at the cell boundary ${x}_{i+\frac{1}{2}}$

\[ {v}_{i+\frac{1}{2}}=\sum_{j=0}^{k-1}c_{rj}\bar{v}_{i-r+j} \]

is $k-th$ order accurate:

\[ {v}_{i+\frac{1}{2}}=v(x_{i+\frac{1}{2}})+O(\Delta x^{k}) \]

use the Lagrange form of the interpolation polynomial:

\[ {P}(x)=\sum_{m=0}^{k}V(x_{i-r+m-\frac{1}{2}})\prod_{l=0,l\ne m}^{k} \cfrac{x-x_{i-r+l-\frac{1}{2}}}{x_{i-r+m-\frac{1}{2}}-x_{i-r+l-\frac{1}{2}}} \\ \]

\[ c_{rj}=\sum^{k}_{m=j+1}\cfrac{\sum^{k}_{l=0,l\ne m}\prod^{k}_{q=0,q\ne m,l}(r-q+1)}{\prod^{k}_{l=0,l\ne m}(m-l)} \]

From Table 2.1, we would know, for example, that

\[ v_{i+\frac{1}{2}}=-\cfrac{1}{6}\bar{v}_{i-1}+\cfrac{5}{6}\bar{v}_{i}+\cfrac{1}{3}\bar{v}_{i+1}+O(\Delta x^{3}).\\ \]

calc crj

2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54

Python

import math class="kn">from fractions import Fraction class="k">def calculate_crj(r, j, k): result = Fraction(0, 1) # 外层求和：m 从 j+1 到 k for m in range(j + 1, k + 1): numerator = 0 # 内层求和：l 从 0 到 k，且 l ≠ m for l in range(0, k + 1): if l == m: continue  # 跳过 l = m 的情况 product = 1 # 内层连乘：q 从 0 到 k，且 q ≠ m, l for q in range(0, k + 1): if q == m or q == l: continue  # 跳过 q = m 或 q = l 的情况 product *= (r - q + 1) numerator += product denominator = 1 # 分母连乘：l 从 0 到 k，且 l ≠ m for l in range(0, k + 1): if l == m: continue  # 跳过 l = m 的情况 denominator *= (m - l) result += Fraction(numerator, denominator) return result class="k">for k in range(1,8): print(f"=== k = {k} ===") # 计算矩阵并存储到列表中 mat = [] for r in range(k): row = [] for j in range(k): c_rj = calculate_crj(r, j, k) row.append(c_rj) mat.append(row) # 计算每个分数的最大宽度 max_width = max(len(str(item)) for row in mat for item in row) for row in mat: for item in row: print(f"{str(item):^{max_width}}", end=' ') print() print() >