[Paper Review] When WENO Smears Small Eddies — Fu (2019) Low-Dissipation TENO Finite-Volume Scheme

Under that same 10 percent shift, the normalized TENO weights stay near one quarter and the wide stencil stays on.

At Stanford, Lin Fu was running turbulence DNS and saw something off. Fifth-order WENO-JS caught shocks. It also blurred small eddies far from those shocks. The scheme was labeled fifth-order. The resolution in hand looked more like a third-order central difference.

Moving the $C_T$ slider from $10^{-5}$ to $10^{-3}$ or $10^{-7}$ tightens or loosens the cutoff, and the same bump is shaved differently in orange and in cyan.

N cells: 120

CFL: 0.40

TENO cutoff C_T: 1e-5

WENO-JSTENO-5L2 err — WENO: 0.0000 | TENO: 0.0000

After one period, check whether both square-pulse edges smear by only a cell or two. Also check whether a small ring appears beside the cyan pulse at $10^{-7}$ .

Paper Info and Background#

Author: Lin Fu (Stanford Turbulence Research Center)
Journal: Computer Physics Communications 235 (2019) 25–39
DOI: 10.1016/j.cpc.2018.10.009
Keywords: TENO, High-order schemes, Shock-capturing, Low-dissipation, Finite-volume method

The Hidden Cost of WENO-JS#

Classic WENO-JS (Jiang–Shu, 1996) measures a smoothness indicator $\beta_k$ on each small stencil $S_0, S_1, S_2$ and blends them with nonlinear weights.

w_k = \frac{\alpha_k}{\sum_j \alpha_j}, \qquad \alpha_k = \frac{d_k}{(\beta_k + \varepsilon)^2}

$d_k$ is the optimal linear weight that restores fifth order. $\varepsilon$ , near $10^{-6}$ , keeps the denominator off zero. Those weights are continuous. A tiny stencil wobble moves $w_k$ in a region that should have stayed at $w_k=d_k$ , and weight leaks onto a lower-order stencil.

A smooth wave loses amplitude and piles up phase error over tens of steps. In LES, where shocks and turbulent eddies share the domain, that dissipation is fatal.

Two Ideas Behind TENO: Strong Scale Separation and a Sharp Cutoff#

Fu changed two things.

The stencil set changed. Beside the three 3-point stencils $S_0, S_1, S_2$ there is one large 5-point stencil $S_3$ . If $S_3$ is clean, stop. The others are unused.

The weights were made binary. The smooth-or-not test is the cutoff $C_T$ .

\gamma_k = \left(1 + \frac{\tau_K}{\beta_k + \varepsilon}\right)^q, \qquad \chi_k = \frac{\gamma_k}{\sum_j \gamma_j}

\delta_k = \begin{cases} 0, & \chi_k < C_T \\ 1, & \text{otherwise} \end{cases}

$\tau_K=|\beta_0-\beta_2|$ is the global reference. $q=6$ stretches the ratio. $\delta_k$ is 0 or 1. A stencil is dropped completely or kept at its optimal weight. There is no continuous fade. Smooth regions keep the fifth-order linear scheme on $S_3$ . Only near a discontinuity does the scheme pick stencils the ENO way. Whether one cutoff can carry every wave is a separate question. The argument that each wave family needs its own reconstruction doubts the attempt.

From Continuous Weights to a 0/1 Verdict#

If $S_3$ is smooth, $\delta_3=1$ , the reconstruction is the fifth-order linear formula.

u_{i+1/2}^L = \frac{1}{60}\bigl(2\bar u_{i-2} - 13\bar u_{i-1} + 47\bar u_i + 27\bar u_{i+1} - 3\bar u_{i+2}\bigr)

$\bar u_j$ is the cell average. The formula is fully linear. When $C_T$ is near zero it is about as nondissipative as a central difference. If $\delta_3=0$ , only the surviving small stencils are blended with $d_k$ .

u_{i+1/2}^L = \sum_{k=0}^{2} w_k\, u_{k,i+1/2}^L, \qquad w_k = \frac{d_k\, \delta_k}{\sum_j d_j\, \delta_j}

That is the ENO property. A stencil that crosses a discontinuity has $\delta_k$ switched off, and the scheme moves to another stencil.

A 5th-Order Stencil Reconstruction in NumPy#

The code is short. From five points it builds three smoothness indicators and three candidate polynomials. TENO then cuts with $C_T$ . WENO-JS, the comparison, blends by inverse squares. This function itself is the TENO branch.

import numpy as np
 
def teno_beta_indicators(um2, um1, u0, up1, up2):
    # Jiang-Shu smoothness indicators for three 3-point stencils
    b0 = (13/12)*(um2 - 2*um1 + u0)**2 + 0.25*(um2 - 4*um1 + 3*u0)**2
    b1 = (13/12)*(um1 - 2*u0 + up1)**2 + 0.25*(um1 - up1)**2
    b2 = (13/12)*(u0 - 2*up1 + up2)**2 + 0.25*(3*u0 - 4*up1 + up2)**2
    return b0, b1, b2
 
def teno_weighted_reconstruction(u, i, CT=1e-5, q=6):
    N = len(u)
    um2, um1, u0 = u[(i-2) % N], u[(i-1) % N], u[i]
    up1, up2 = u[(i+1) % N], u[(i+2) % N]
    b0, b1, b2 = teno_beta_indicators(um2, um1, u0, up1, up2)
    b3 = max(b0, b1, b2)
    tau = abs(b0 - b2)
    eps = 1e-40
    g = np.array([(1 + tau/(b + eps))**q for b in (b0, b1, b2, b3)])
    chi = g / g.sum()
    delta = (chi >= CT).astype(float)
    if delta[3] == 1:  # S3 smooth -> use 5th-order linear
        return (2*um2 - 13*um1 + 47*u0 + 27*up1 - 3*up2) / 60
    d = np.array([1/10, 6/10, 3/10])
    alpha = d * delta[:3]
    if alpha.sum() < 1e-14:
        return u0  # last-resort donor cell
    w = alpha / alpha.sum()
    p0 = (2*um2 - 7*um1 + 11*u0) / 6
    p1 = (-um1 + 5*u0 + 2*up1) / 6
    p2 = (2*u0 + 5*up1 - up2) / 6
    return w[0]*p0 + w[1]*p1 + w[2]*p2

With an Euler step and periodic boundaries, a smooth Gaussian bulge comes back nearly intact under TENO and about 30 percent lower under WENO-JS. Same skeleton. Only the weights differ.

Substitute a small wobble into the weight formulas. Ideal linear weights are one tenth, six tenths, and three tenths. If all three β equal 1, WENO weights return there. Raise only the third by 10 percent, so β is 1, 1, 1.1. The third unnormalized weight is 0.3 over 1.21, about 0.248. Normalized weights are about 0.105, 0.633, and 0.262. The L1 gap from the ideal set is about 0.077.

The same numbers enter TENO. Tau is 0.1 and q is 6. The code sets the large-stencil indicator to the max of the three β values and normalizes four together. Chi sits near 0.256, 0.256, 0.244, 0.244. Against a cutoff of $10^{-5}$ that is four orders up. Every verdict stays 1, the wide stencil survives, and the fifth-order linear formula is used as written.

Epsilon in the displayed formula is near $10^{-6}$ . Inside the function it is $10^{-40}$ . Near β of 1 neither shows up. Only a stencil with β at 0 splits the denominator. The fifth-order coefficients 2, −13, 47, 27, and −3 sum to 60, so a constant field is unchanged while the wide stencil stays on.

What stays when C_T moves#

The advection above pushes one initial condition, a smooth sine-squared bump plus a square pulse. Orange is WENO-JS. Cyan is TENO-5. At $C_T=10^{-5}$ , after one period, both square-pulse edges smear by a cell or two. Only the orange sine bump drops in a way you can see. Raising $C_T$ to $10^{-3}$ makes the cut picky, small stencils win often near the pulse, and TENO dissipates more. Dropping it to $10^{-7}$ makes TENO insist on the fifth-order linear scheme, and a small ring shows beside the pulse. $C_T$ is the balance between stability and low dissipation.

Critical Notes — The Tuning Trap#

A few places the paper does not claim out loud.

$C_T$ sticks to the problem. Fu recommends $10^{-5}$ . With a strong shock and a very weak acoustic wave, that value lifts the acoustic wave a little, and $10^{-6}$ rings near the shock. Groups porting this into OpenFOAM or SU2 report retuning the value per case.

In finite volume, the smoothness $\beta_3$ of $S_3$ is expensive. Derivatives of the 5-point polynomial fill about six pages. A finite-difference form is much shorter. Keeping this paper's low-dissipation Riemann-flux switch wants the finite-volume layout.

The match with time integration is empty. SSP-RK3 is fine, buried under the spatial accuracy. Mixed with SSP-RK2 or an implicit method, the jump in $\delta_k$ is reported to bother Newton convergence.

In practice, run the grid study at one fixed $C_T$ . Check whether the spatial order actually comes back as 5. If it does not, $C_T$ is probably large enough that TENO stays in the nonlinear mode all the time.

The verification that breaks if C_T stays fixed#

Refine the mesh with $C_T$ frozen at $10^{-5}$ , and fifth order disappears the moment chi on a smooth cell falls under the cut. The 10 percent wobble above did not cross the cut. A shock sitting next to a very weak wave is a different story.

If an order check then comes back below 1 instead of 5, the cutoff ate the order, not the spatial operator. The verification where one code reported order 0.79 and order 2.00 at the same time asks what the case is actually measuring. Which of those two numbers does TENO with a frozen $C_T$ land on?