# Wave-propagation algorithms¶

## One space dimension¶

Clawpack can be used to solve a system of equations of the form

$\kappa(x)q_t+f(q)_x = \psi(q,x,t),$

where $$q=q(x,t)\in{\cal R}^m$$. The standard case of a homogeneous conservation law has $$\kappa\equiv 1$$ and $$\psi\equiv 0$$,

(1)$q_t+f(q)_x=0.$

The flux function $$f(q)$$ can also depend explicitly on $$x$$ and $$t$$ as well as on $$q$$. Hyperbolic systems that are not in conservation form, e.g.,

$q_t+A(x,t)q_x=0,$

can also be solved. See [LeVeque-FVMHP] for more details about the wave-propagation algorithms that are briefly described here.

The basic requirement on the homogeneous system is that it be hyperbolic in the sense that a Riemann solver can be specified that, for any two states $$q_{i-1}$$ and $$Q_i$$, returns a set of $$M_w$$ waves $${\cal W}^p_{i-1/2}$$ and speeds $$s^p_{i-1/2}$$ satisfying

$\sum_{p=1}^{M_w} {\cal W}^p_{i-1/2} = Q_i-Q_{i-1} \equiv \Delta Q_{i-1/2}.$

The Riemann solver must also return a left-going fluctuation $${\cal A}^-\Delta Q_{i-1/2}$$ and a right-going fluctuation $${\cal A}^+\Delta Q_{i-1/2}$$. In the standard conservative case (1) these should satisfy

(2)${\cal A}^-\Delta Q_{i-1/2}+{\cal A}^+\Delta Q_{i-1/2} = f(Q_i)-f(Q_{i-1})$

and the fluctuations then define a “flux-difference splitting”.

${\cal A}^-\Delta Q_{i-1/2} = \sum_p (s^p_{i-1/2})^- {\cal W}^p_{i-1/2}, \qquad {\cal A}^+\Delta Q_{i-1/2} = \sum_p (s^p_{i-1/2})^+ {\cal W}^p_{i-1/2},$

where $$s^- = \min(s,0)$$ and $$s^+ = \max(s,0)$$. In the nonconservative case eqn{claw_1dnoncon}, there is no “flux function” $$f(q)$$, and the constraint (2) need not be satisfied.

## Godunov’s method¶

Only the fluctuations are used for the first-order Godunov method, which is implemented in the form

$Q_i^{n+1} = Q_i^n - \frac{\Delta t}{\Delta x}\left[{\cal A}^+\Delta Q_{i-1/2} + {\cal A}^-\Delta Q_{i+1/2}\right],$

assuming $$\kappa \equiv 1$$.

The Riemann solver must be supplied by the user in the form of a subroutine rp1, as described in Specifying the Riemann solver.

Typically the Riemann solver first computes waves and speeds and then uses these to compute $${\cal A}^+Q_{i-1/2}$$ and $${\cal A}^-Q_{i-1/2}$$ internally in the Riemann solver.

## High-resolution methods¶

The waves and speeds must also be returned by the Riemann solver in order to use the high-resolution methods described in [LeVeque-FVMHP], which reduce to a second-order accurate Lax-Wendroff type method when limiters are not used. By introducing wave limiters, non-physical oscillations near discontinuities or steep gradients in the solution can be avoided. The limiters are based on the theory of TVD methods for nonlinear scalar equations and extended in a natural way to systems of equations.

These methods take the form

$Q_i^{n+1} = Q_i^n - \frac{\Delta t}{\Delta x}\left[{\cal A}^+Q_{i-1/2} + {\cal A}^-Q_{i+1/2}\right] - \frac{\Delta t}{\Delta x}(\tilde F_{i+1/2} - \tilde F_{i-1/2})$

where

$\tilde F_{i-1/2} = \frac 1 2 \sum_{p=1}^{M_w} |s^p_{i-1/2}| \left( 1-\frac{\Delta t}{\Delta x} |s^p_{i-1/2}|\right) \tilde{\cal W}_{i-1/2}^p.$

Here $$\tilde{\cal W}_{i-1/2}^p$$ represents a limited version of the wave $${\cal W}_{i-1/2}^p$$, obtained by comparing $${\cal W}_{i-1/2}^p$$ to $${\cal W}_{i-3/2}^p$$ if $$s^p>0$$ or to $${\cal W}_{i+1/2}^p$$ if $$s^p<0$$.

## f-wave formulation¶

For equations with spatially-varying flux functions, such as

(3)$q_t+f(q,x)_x=0.$

it is often convenient to use the f-wave formulation of the algorithm, as proposed in [BaleLevMitRoss02]. Rather than decomposing the jump $$Q_i-Q_{i-1}$$ into waves, the idea of the f-wave approach is to decompose the flux difference $$f(Q_i,x_i) - f(Q_{i-1},x_{i-1})$$ into f-waves using appropriate eigenvectors of the Jacobian matrices to either side of the interface. See f-wave Riemann solvers for more details.

## Capacity functions¶

When a capacity function $$\kappa(x)$$ is present, the Godunov method becomes

$Q_i^{n+1} = Q_i^n - \frac{\Delta t}{\kappa_i \Delta x} \left[{\cal A}^+Q_{i-1/2} + {\cal A}^-Q_{i+1/2}\right],$

See [LeVeque-FVMHP] for discussion of this algorithm and its extension to the high-resolution method. Capacity functions are useful in particular for solving equations on a mapped grid.

## Source terms¶

If the equation has a source term, a routine src1 must also be supplied that solves the source term equation $$q_t=\psi$$ over a time step. A fractional step method is used to couple this with the homogeneous solution, as described in Using src for source terms.

## Boundary conditions¶

Boundary conditions are imposed by setting values in ghost cells each time step, as described in Boundary conditions.

TODO: Continue description – 2d and 3d, transverse solvers.