To use all functions of this page, please activate cookies in your browser.
With an accout for my.chemeurope.com you can always see everything at a glance – and you can configure your own website and individual newsletter.
- My watch list
- My saved searches
- My saved topics
- My newsletter
MUSCL stands for Monotone Upstream-centered Schemes for Conservation Laws, and the term was introduced in a seminal paper by Bram van Leer (van Leer, 1979). In this paper he constructed the first high-order, total variation diminishing (TVD) scheme where he obtained second order spatial accuracy. It is a finite volume method that provides high accuracy numerical solutions to partial differential equations which can involve solutions that exhibit shocks, discontinuities or large gradients.
The idea is to replace the piecewise constant approximation of Godunov's scheme by reconstructed states, derived from cell-averaged states obtained from the previous time-step. For each cell, slope limited, reconstructed left and right states are obtained and used to calculate fluxes at the cell boundaries (edges). These fluxes can, in turn, be used as input to a Riemann solver, following which the solutions are averaged and used to advance the solution in time. Alternatively, the fluxes can be used in Riemann-free-solver schemes, such as the Kurganov and Tadmor scheme outlined below.
Additional recommended knowledge
We will consider the fundamentals of the MUSCL scheme by considering the following simple first-order, scalar, 1D system, which is assumed to have a wave propagating in the positive direction,
Where represents a state variable and represents a flux variable.
The basic scheme of Godunov uses piecewise constant approximations for each cell, and results in a first-order upwind discretisation of the above problem with cell centres indexed as i. A semi-discrete scheme can be defined as follows,
This basic scheme is not able to handle shocks or sharp discontinuities as they tend to become smeared. An example of this effect is shown in the diagram opposite, which illustrates a 1D advective equation with a step wave propagating to the right. The simulation was carried out with a mesh of 200 cells and used a 4th order Runge-Kutta time integrator (RK-4).
To provide higher resolution of discontinuities, Godunov's scheme can be extended to use piecewise linear approximations of each cell, which results in a central difference scheme that is second-order accurate in space. The piecewise linear approximations are obtained from
Thus, evaluating fluxes at the cell edges we get the following semi-discrete scheme
where and are the piecewise approximate values of cell edge variables, i.e.
Whilst the above second-order scheme provides greater accuracy for smooth solutions, it is not a total variation diminishing (TVD) scheme and introduces spurious oscillations into the solution where discontinuities or shocks are present. An example of this effect is shown in the diagram opposite, which illustrates a 1D advective equation ut + ux = 0, with a step wave propagating to the right. This loss of accuracy is to be expected due to Godunov's theorem. The simulation was carried out with a mesh of 200 cells and used RK-4 for time integration.
MUSCL based numerical schemes extend the idea of using a linear piecewise approximation to each cell by using slope limited left and right extrapolated states. This results in the following high resolution, TVD discretisation scheme,
Which, alternatively, can be written in the more succinct form,
The numerical fluxes correspond to a nonlinear combination of first and second-order approximations to the continuous flux function.
The symbols and represent scheme dependent functions (of the limited extrapolated cell edge variables), i.e.
The function is a limiter function that limits the slope of the piecewise approximations to ensure the solution is TVD, thereby avoiding the spurious oscillations that would otherwise occur around discontinuities or shocks - see Flux limiter section. The limiter is equal to zero when and is equal to unity when . Thus, the accuracy of a TVD discretization degrades to first order at local extrema, but tends to second order over smooth parts of the domain.
The algorithm is straight forward to implement. Once a suitable scheme for has been chosen, such as the Kurganov and Tadmor scheme (see below), the solution can proceed using standard integration techniques such as Method of lines (MOL).
Kurganov and Tadmor central scheme
The Nessyahu and Tadmor central scheme (Nessyahu and Tadmor, 1990) is a Riemann-solver-free, second-order, high-resolution scheme that uses MUSCL reconstruction. It is straight forward to implement and can be used on scalar and vector problems. The algorithm is based upon central differences with comparable performance to Riemann type solvers when used to obtain solutions for PDE’s describing systems that exhibit high-gradient phenomena. The Kurganov and Tadmor scheme can be implemented as a fully-discrete or semi-discrete scheme. Here we consider the semi-discrete scheme later proposed by Kurganov and Tadmor (2000)
The calculation is shown below:
Where the local propagation speed, , is the maximum absolute value of the eigenvalue of the Jacobian of over cells given by
where represents the spectral radius of .
Beyond these CFL related speeds, no characteristic information is required.
The above flux calculation is sometimes referred to as local Lax-Friedrichs flux or Rusanov flux (Lax, 1954; Rusanov, 1961; Toro, 1999; Kurganov and Tadmor, 2000; Leveque, 2002).
An example of the effectiveness of using a high resolution scheme is shown in the diagram opposite, which illustrates the 1D advective equation , with a step wave propagating to the right. The simulation was carried out on a mesh of 200 cells, using the Kurganov and Tadmor central scheme with Superbee limiter and used RK-4 for time integration. This simulation result contrasts extremely well against the above first-order upwind and second-order central difference results shown above. This scheme also provides good results when applied to sets of equations - see results below for this scheme applied to the Euler equations. However, care has to be taken in choosing an appropriate limiter because, for example, the Superbee limiter can cause unrealistic sharpening for some smooth waves.
The scheme can readily include diffusion terms, if they are present. For example, if the above 1D scalar problem is extended to include a diffusion term, we get
for which Kurganov and Tadmor propose the following central difference approximation,
Full details of the algorithm (full and semi-discrete versions) and its derivation can be found in the original paper (Kurganov and Tadmor, 1999), along with a number of 1D and 2D examples. Additional information is also available in an earlier related paper (Nessyahu and Tadmor, 1990).
Note: Whilst this scheme was originally presented by Kurganov and Tadmor as a 2nd order scheme based upon linear extrapolation, a later paper (Kurganov and Levy, 2000) demonstrates that it can also form the basis of a third order scheme. A 1D advective example and a Euler equation example of their scheme, using parabolic reconstruction (3rd order), are shown in the parabolic reconstruction and Euler equation sections below.
Piecewise parabolic reconstruction
It is possible to extend the idea of linear-extrapolation to higher order reconstruction, and an example is shown in the diagram opposite. However, for this case the left and right states are estimated by interpolation of a second-order, upwind biased, difference equation. This results in a parabolic reconstruction scheme that is third-order accurate in space.
We follow the approach of Kermani (Kermani, et al, 2003), and present a third-order upwind biased scheme, where the symbols and again represent scheme dependent functions (of the limited reconstructed cell edge variables). But for this case they are based upon parabolically-reconstructed states, i.e.
Where = 1/3 and,
and the limiter function , is the same as above.
Parabolic reconstruction is straight forward to implement and can be used with the Kurganov and Tadmor scheme in lieu of the linear extrapolation shown above. This has the effect of raising the spatial solution of the KT scheme to 3rd order. It performs well when solving the Euler equations, see below. However, whilst this increase in spatial order has certain advantages over 2nd order schemes for smooth solutions, for shocks it is more disipative - compare diagram opposite with above solution obtained using the KT algorithm with linear extrapolation and Superbee limiter. This simulation was carried out on a mesh of 200 cells using the same KT algorithm but with parabolic reconstruction. Time integratration was by RK-4, and the alternative form of van Albada limiter, , was used to avoid spurious oscillations.
Euler equations example
For simplicity we consider the 1D case without heat transfer and without body force. Therefore, in conservation vector form, the general Euler equations reduce to
and where U is a vector of states and F is a vector of fluxes.
The equations above represent conservation of mass, momentum, and energy. There are thus three equations and four unknowns, ρ (density) u (fluid velocity), p (pressure) and E (total energy). The total energy is given by,
where represents specific internal energy.
In order to close the system an equation of state is required. One that suits our purpose is
where is equal to the ratio of specific heats for the fluid.
We can now proceed, as shown above in the simple 1D example, by obtaining the left and right extrapolated states for each state variable. Thus, for density we obtain
Having obtained the limited extrapolated states, we then proceed to construct the edge fluxes using these values. With the edge fluxes known, we can now construct the semi-discrete scheme, i.e.
The solution can now proceed by integration using standard numerical techniques.
The above illustrates the basic idea of the MUSCL scheme. However, for a practical solution to the Euler equations, a suitable scheme has to be chosen which will also define the function .
The diagram opposite shows a 2nd order solution to G A Sod's shock tube problem (Sod, 1978) using the above high resolution Kurganov and Tadmor Central Scheme (KT) with Linear Extrapolation and Ospre limiter. This illustrates clearly the effectiveness of the MUSCL approach to solving the Euler equations. The simulation was carried out on a mesh of 200 cells using Matlab code (Wesseling, 2001), adapted to use the KT algorithm and Ospre limiter. Time integration was performed by a 4th order SHK (equivalent performance to RK-4) integrator. The following initial conditions (SI units) were used:
The diagram opposite shows a 3rd order solution to G A Sod's shock tube problem (Sod, 1978) using the above high resolution Kurganov and Tadmor Central Scheme (KT) but with parabolic reconstruction and van Albada limiter. This again illustrates the effectiveness of the MUSCL approach to solving the Euler equations. The simulation was carried out on a mesh of 200 cells using Matlab code (Wesseling, 2001), adapted to use the KT algorithm with Parabolic Extrapolation and van Albada limiter. The alternative form of van Albada limiter, , was used to avoid spurious oscillations. Time integration was performed by a 4th order SHK integrator. The same initial conditions were used.
More information on these and other methods can be found in the references below.
|This article is licensed under the GNU Free Documentation License. It uses material from the Wikipedia article "MUSCL_scheme". A list of authors is available in Wikipedia.|