Elsevier

Journal of Computational Physics

Volume 234, 1 February 2013, Pages 452-471
Journal of Computational Physics

Lax–Friedrichs fast sweeping methods for steady state problems for hyperbolic conservation laws

https://doi.org/10.1016/j.jcp.2012.10.008Get rights and content

Abstract

Fast sweeping methods are efficient iterative numerical schemes originally designed for solving stationary Hamilton–Jacobi equations. Their efficiency relies on Gauss–Seidel type nonlinear iterations, and a finite number of sweeping directions. In this paper, we generalize the fast sweeping methods to hyperbolic conservation laws with source terms. The algorithm is obtained through finite difference discretization, with the numerical fluxes evaluated in WENO (Weighted Essentially Non-oscillatory) fashion, coupled with Gauss–Seidel iterations. In particular, we consider mainly the Lax-Friedrichs numerical fluxes. Extensive numerical examples in both scalar and system test problems in one and two dimensions demonstrate the efficiency, high order accuracy and the capability of resolving shocks of the proposed methods.

Introduction

Hyperbolic conservation laws and Hamilton–Jacobi equations are first order nonlinear partial differential equations which arise in many applications such as gas dynamics [19], [40], shallow water waves [41], geometrical optics [42], crystal growth [23], [25], etching, photolithography [33], computer vision [32], [24] and seismology [7]. The solutions of these equations can develop singularity, such as discontinuity in the solutions or their derivatives. Under these circumstances, the solutions do not satisfy the equations in the classical sense and one must resort to weak solutions. For hyperbolic conservation equations, “vanishing viscosity solution” and “entropy solution” are introduced to define the weak solution uniquely (see [20], [18] and references therein). In early 1980s, Crandall and Lions [10], [11] introduced “viscosity solution”, among weak solutions, to study of the existence, uniqueness, and stability properties of Hamilton–Jacobi equations. Since then, many numerical methods have been proposed to approximate the viscosity solutions [34], [18], [9], [38]. The challenges of the numerical schemes generally are to achieve high order accuracy and be able to resolve the shock or singularity well, while maintaining conservation. One key to the scheme design of hyperbolic conservation laws and Hamilton–Jacobi equations, in order to correctly capture the viscosity solution, is the use of consistent and conservative numerical fluxes and Hamiltonians.

For Hamilton–Jacobi equations,ϕt(x,t)+H(x,t,ϕ,ϕ(x,t))=0after Osher proved in [22] the link between static and time dependent Hamilton–Jacobi equations, in which the zero level set of the viscosity solution of the time dependent problem at a later time t is the t-level set of the static problem, fast sweeping methods for static Hamilton–Jacobi equations became popular due to their high efficiency. A fast sweeping method mainly consists of the following three essential ingredients: (1) an efficient local solver on a given Cartesian mesh [39], [46], [15], [17], [45], [26] or triangulation [26], [16] based on monotone numerical Hamiltonians, (2) systematic orderings of solution nodes according to some pre-determined information-flowing directions, and (3) Gauss–Seidel type iterations based on a given order of solution nodes. Among fast sweeping schemes, the methods based on upwind Hamiltonians are most efficient for convex Hamiltonians [39], [46], while the methods based on Lax–Friedrichs fluxes [15], [45] are most flexible to deal with general non-convex Hamiltonians.

Although for hyperbolic conservation laws with source terms,ut+·f(u)=s(u),in which the Jacobian matrix f(u) is diagonalizable with all the eigenvalues being real for any u, there is no such link between time dependent problems and static problems as in Hamilton–Jacobi equations, it is still favorable to have efficient numerical methods for steady state hyperbolic problems. A class of schemes, called “residual distribution schemes” [1], [2], [3], [12], [29], [37], were proposed for solving steady state problems with pseudo time stepping. The spirit of these methods is to distribute the residuals, defined through integrating the flux and source terms on triangular or quadrilinear cells in a conservative fashion, and march in pseudo time. Residual distribution schemes were later generalized to high order schemes by Abgrall and Roe [4] through a finite element based approach. Based on the same distribution principles, Chou and Shu [8] developed a finite difference based residual distribution scheme which works on curvilinear meshes, and their scheme achieves high order accuracy and low computational cost as in finite difference methods, but without the constraint of uniform meshes.

These residual distribution schemes, though more efficient by using only pseudo time stepping, are still greatly constrained by the CFL condition for stability, which can substantially limit the speed of the schemes. In this paper, we develop a Gauss–Seidel type iterative method to accelerate the speed to compute the steady state solutions of hyperbolic equations. Inspired by fast sweeping methods for time independent Hamilton–Jacobi equations, we propose methods which discretize the steady state hyperbolic conservation laws directly, by approximating the spatial derivatives with consistent and conservative numerical fluxes, and iterate with Gauss–Seidel type nonlinear method with a finite number of alternating sweeping directions. In particular, we use the Lax–Friedrichs fluxes evaluated in WENO (Weighted Essentially Non-oscillatory) fashion [35], [36], [21], [14], [34], to achieve high order accuracy as well as high resolution of shocks.

It is worth pointing out here that, while most steady state hyperbolic problems have unique steady states, with initial conditions reasonably perturbed from the solutions, there are some problems whose steady states are totally dependent on the initial conditions through mass conservation [30]. In those cases, Gauss–Seidel type sweeping may not conserve the mass through the iterations, and therefore an additional constraint needs to be imposed in order to select the particular steady state.

This paper is organized as follows. In Section 2, we present the high order Lax–Friedrichs sweeping methods for one-dimensional scalar and system problems. In Sections 3, the method is extended to two-dimensional equations. Section 4 describes an efficient accuracy-preserving stopping criterion for the fast sweeping iterative scheme. Section 5 contains extensive numerical simulations for one and two-dimensional scalar and system steady state problems to demonstrate high order accuracy, efficiency and robustness of our scheme. Conclusions are given in Section 6.

Section snippets

One-dimensional scalar problems

In this section, we consider the one-dimensional scalar steady state problemf(u)x=s(u,x),x[a,b]subject to an initial guess and appropriate boundary conditions.

To obtain a numerical scheme, the interval is first discretized uniformly into N cells, and the grid points are denoted by {xj}j=0N, where xj=a+jΔx,j=0,,N and Δx=(b-a)/N. The midpoint of a cell is defined as xj+12=(xj+xj+1)/2, j=0,,N-1. The numerical approximation of u on the grid points xj are denoted by uj,j=0,,N. A conservative

Two-dimensional problems

The sweeping method we described in the previous sections can be easily extended to the two-dimensional steady state problemf(u)x+g(u)y=s(u,x,y),(x,y)[a,b]×[c,d],Let {(xi,yj)},i=0,,Nx,j=0,,Ny denote the grid points of a uniform discretization of the computational domain, with Δx=(b-a)/Nx and Δy=(d-c)/Ny as the mesh sizes for x and y direction, respectively. We use ui,j to represent the numerical solution of u at grid point (xi,yj). A conservative finite difference discretization of (14) can

Stopping criteria for convergence

Iterative schemes always require stopping criteria by which the convergence of the numerical scheme is determined. Traditional Gauss–Seidel iterations adopt the stopping criteria in which the algorithm stops when a particular norm of the difference of successive iterations, called residue, is smaller than a number. This stopping criteria is effective when the residue of the scheme is monotonically decreasing, but not robust enough for schemes with oscillatory convergence history, such as high

Numerical results

Here we show the numerical results of the proposed Lax–Friedrichs WENO sweeping method for hyperbolic conservation laws with source terms in both scalar and system test problems in one and two dimensions. The efficiency and high order accuracy of the proposed scheme will be demonstrated.

All the spatial discretizations in our numerical results are uniform. The values of σ in Lax–Friedrichs fluxes and α in the flux splitting are updated at the beginning of each directional sweeping and relaxed by

Conclusion

In this paper, we proposed the first fast sweeping method for solving steady state hyperbolic equations with source terms. The method is based on the Lax–Friedrichs fluxes with WENO (Weighted Essentially Non-oscillatory) reconstruction to achieve high order accuracy. The alternating sweeping directions with Gauss–Seidel iterative updates are used to accelerate the speed of convergence. Furthermore, the modified stopping criteria is proposed to stop the algorithm when the residue of the scheme

Acknowledgement

We thank reviewers for their constructive comments and Professor Chi-Wang Shu for the valuable discussion.

References (46)

  • C.-W. Shu et al.

    Efficient implementation of essentially non-oscillatory shock-capturing schemes

    J. Comput. Phys.

    (1988)
  • C.-W. Shu et al.

    Efficient implementation of essentially non-oscillatory shock-capturing schemes, II

    J. Comput. Phys.

    (1989)
  • Nail K. Yamaleev et al.

    Third-order energy stable WENO scheme

    J. Comput. Phys.

    (2009)
  • R. Abgrall et al.

    Residual distribution schemes for conservation laws via adaptive quadrature

    SIAM J. Scient. Comput.

    (2002)
  • R. Abgrall et al.

    High order fluctuation scheme on triangular meshes

    J. Scient. Comput.

    (2003)
  • W. Cai et al.

    Essentially nonoscillatory spectral Fourier methods for shock wave calculations

    Math. Comput.

    (1989)
  • C. Chapman

    Fundamentals of Seismic Wave Propagation

    (2010)
  • B. Cockburn et al.

    Discontinuous Galerkin Methods: Theory, Computation and Applications, Lecture Notes in Computational Science and Engineering

    (2000)
  • M.G. Crandall et al.

    Viscosity solutions of Hamilton–Jacobi equations

    Trans. Amer. Math. Soc.

    (1983)
  • M.G. Crandall et al.

    Two approximations of solutions of Hamilton–Jacobi equations

    Math. Comput.

    (1984)
  • H. Deconinck et al.

    P

    (1993)
  • P. Embid et al.

    Multiple steady states for 1-d transonic flow

    SIAM J. Scient. Stat. Comput.

    (1984)
  • C.-Y. Kao et al.

    Fast sweeping methods for static Hamilton–Jacobi equations

    SIAM J. Numer. Anal.

    (2005)
  • Cited by (29)

    • Spatial convergence study of Lax-Friedrichs WENO fast sweeping method on the S<inf>N</inf> transport equation with nonsmoothness

      2022, Annals of Nuclear Energy
      Citation Excerpt :

      Weighted essentially non-oscillatory (WENO) is a classical scheme with wide stencil for solving hyperbolic equations, which achieves high order accuracy in smooth regions and resolves discontinuities in an essentially non-oscillatory fashion (Shu, 2020). Inspired by fast sweeping methods for time independent Hamilton-Jacobi equations, the Lax-Friedrichs fast sweeping method based on the classical WENO scheme of Jiang and Shu (WENO-JS) (Jiang and Shu, 1996) was applied for solving steady-state hyperbolic conservation laws (Chen et al., 2013). Then the Lax-Friedrichs WENO fast sweeping method was applied for the SN transport equation (Wang and Byambaakhuu, 2019).

    • Implicit fast sweeping method for hyperbolic systems of conservation laws

      2021, Journal of Computational Physics
      Citation Excerpt :

      Their paraxial formulation used shock matching to produce sharp shock surfaces, but their approach required a selection principle to determine the correct solution branch. A somewhat different FSM for steady state conservation laws was proposed by Chen et al. [18] and later extended in [19,20]. Their method relied on a direct discretization of the static form of Eq. (1) using Lax-Friedrichs fluxes that eliminated the need for inversion.

    • Fixed-point fast sweeping weighted essentially non-oscillatory method for multi-commodity continuum traffic equilibrium assignment problem

      2018, Applied Mathematical Modelling
      Citation Excerpt :

      However, the computational efficiency of this approach is restricted by a Courant–Friedrichs–Lewy (CFL) condition for the discrete time step size and it takes significant time to evolve the solution of the pseudo-time problem to steady state [23]. Recently, a fast sweeping method, which was originally designed for solving static Hamilton–Jacobi equations [24–27], has been applied to solve hyperbolic conservation laws with source terms [23,28]. This method is an efficient Gauss–Seidel iterative numerical scheme to the boundary value problems with solution information propagating along characteristics starting from the boundary.

    View all citing articles on Scopus
    1

    This author is supported by NSF DMS1020625.

    2

    This author is partially supported by NSF DMS1216742.

    View full text