Conservative High Order Positivity-Preserving Discontinuous Galerkin Methods for Linear Hyperbolic and Radiative Transfer Equations | Journal of Scientific Computing Skip to main content
Log in

Conservative High Order Positivity-Preserving Discontinuous Galerkin Methods for Linear Hyperbolic and Radiative Transfer Equations

  • Published:
Journal of Scientific Computing Aims and scope Submit manuscript

Abstract

We further investigate the high order positivity-preserving discontinuous Galerkin (DG) methods for linear hyperbolic and radiative transfer equations developed in Yuan et al. (SIAM J Sci Comput 38:A2987–A3019, 2016). The DG methods in Yuan et al. (2016) can maintain positivity and high order accuracy, but they rely both on the scaling limiter in Zhang and Shu (J Comput Phys 229:8918–8934, 2010) and a rotational limiter, the latter may alter cell averages of the unmodulated DG scheme, thereby affecting conservation. Even though a Lax–Wendroff type theorem is proved in Yuan et al. (2016), guaranteeing convergence to weak solutions with correct shock speed when such rotational limiter is applied, it would still be desirable if a conservative DG method without changing the cell averages can be obtained which has both high order accuracy and positivity-preserving capability. In this paper, we develop and analyze such a DG method for both linear hyperbolic equations and radiative transfer equations. In the one-dimensional case, the method uses traditional DG space \(P^k\) of piecewise polynomials of degree at most k. A key result is proved that the unmodulated DG solver in this case can maintain positivity of the cell average if the inflow boundary value and the source term are both positive, therefore the positivity-preserving framework in Zhang and Shu (2010) can be used to obtain a high order conservative positivity-preserving DG scheme. Unfortunately, in two-dimensions this is no longer the case. We show that the unmodulated DG solver based either on \(P^k\) or \(Q^k\) spaces (piecewise kth degree polynomials or piecewise tensor-product kth degree polynomials) could generate negative cell averages. We augment the DG space with additional functions so that the positivity of cell averages from the unmodulated DG solver can be restored, thereby leading to high order conservative positivity-preserving DG scheme based on these augmented DG spaces following the framework in Zhang and Shu (2010). Computational results are provided to demonstrate the good performance of our DG schemes.

This is a preview of subscription content, log in via an institution to check access.

Access this article

Subscribe and save

Springer+ Basic
¥17,985 /Month
  • Get 10 units per month
  • Download Article/Chapter or eBook
  • 1 Unit = 1 Article or 1 Chapter
  • Cancel anytime
Subscribe now

Buy Now

Price includes VAT (Japan)

Instant access to the full article PDF.

Fig. 1
Fig. 2
Fig. 3
Fig. 4
Fig. 5
Fig. 6
Fig. 7
Fig. 8
Fig. 9
Fig. 10

Similar content being viewed by others

References

  1. Carlson, B.G., Lathrop, K.D.: Transport theory-the method of discrete ordinates. In: Greenspan, H., Kelber, C.N., Okrent, D. (eds.) Computing Methods in Reactor Physics, pp. 165–266. Gordon and Breach, New York (1968)

    Google Scholar 

  2. Cockburn, B., Hou, S., Shu, C.-W.: The Runge–Kutta local projection discontinuous Galerkin finite element method for conservation laws IV: the multidimensional case. Math. Comput. 54, 545–581 (1990)

    MathSciNet  MATH  Google Scholar 

  3. Cockburn, B., Lin, S.-Y., Shu, C.-W.: TVB Runge–Kutta local projection discontinuous Galerkin finite element method for conservation laws III: one dimensional systems. J. Comput. Phys. 84, 90–113 (1989)

    Article  MathSciNet  Google Scholar 

  4. Cockburn, B., Shu, C.-W.: TVB Runge-Kutta local projection discontinuous Galerkin finite element method for conservation laws II: General framework. Math. Comput. 52, 411–435 (1989)

    MathSciNet  MATH  Google Scholar 

  5. Cockburn, B., Shu, C.-W.: The Runge–Kutta discontinuous Galerkin method for conservation laws V: multidimensional systems. J. Comput. Phys. 141, 199–224 (1998)

    Article  MathSciNet  Google Scholar 

  6. Fiveland, W.A.: Discrete-ordinates solutions of the radiative transport equation for rectangular enclosures. J. Heat Transf. 106, 699–706 (1984)

    Article  Google Scholar 

  7. Lathrop, K.D., Carlson, B.G.: Discrete Ordinate Angular Quadrature of the Neutron Transport equation. Technical report LA-3186, Los Alamos Scientific Laboratory (1965)

  8. Lesaint, P., Raviart, P.A.: On a finite element method for solving the neutron transport equation. In: de Boor, C.A. (ed.) Mathematical Aspects of Finite Elements in Partial Differential Equations, pp. 89–145. Academic Press, New York (1974)

    Chapter  Google Scholar 

  9. Lewis, E.E., Miller Jr., W.F.: Computational Methods of Neutron Transport. Wiley, New York (1984)

    MATH  Google Scholar 

  10. Qin, T., Shu, C.-W.: Implicit positivity-preserving high order discontinuous Galerkin methods for conservation laws. SIAM J. Sci. Comput. 40, A81–A107 (2018)

    Article  MathSciNet  Google Scholar 

  11. Reed, W.H., Hill, T.R.: Triangular mesh methods for the neutron transport equation. Technical report LA-UR-73-479, Los Alamos Scientific Laboratory (1973)

  12. Shu, C.-W.: TVB uniformly high order schemes for conservation laws. Math. Comput. 49, 105–121 (1987)

    Article  MathSciNet  Google Scholar 

  13. Shu, C.-W., Osher, S.: Efficient implementation of essentially non-oscillatory shock-capturing schemes. J. Comput. Phys. 77, 439–471 (1988)

    Article  MathSciNet  Google Scholar 

  14. Yuan, D., Cheng, J., Shu, C.-W.: High order positivity-preserving discontinuous Galerkin methods for radiative transfer equations. SIAM J. Sci. Comput. 38, A2987–A3019 (2016)

    Article  MathSciNet  Google Scholar 

  15. Zhang, X., Shu, C.-W.: On positivity-preserving high order discontinuous Galerkin schemes for compressible Euler equations on rectangular meshes. J. Comput. Phys. 229, 8918–8934 (2010)

    Article  MathSciNet  Google Scholar 

Download references

Author information

Authors and Affiliations

Authors

Corresponding author

Correspondence to Juan Cheng.

Additional information

Juan Cheng: Research is supported in part by NSFC Grant 11471049 and U1630247 and Science Challenge Project No. TZ2016002. Chi-Wang Shu: Research is supported in part by ARO Grant W911NF-15-1-0226 and NSF Grant DMS-1719410.

Appendices

Appendix A: The Specific Test Function v for the Augmented DG Spaces

In this appendix we give explicit formulas for the specific test function v satisfying (4.1) for the augmented DG space \(R^1\) defined in (4.5). Since we just consider the case \(\alpha ,\beta >0\) in (2.7), we can simplify it as \(\alpha u_x+u_y+\gamma u=f\). Besides, for the DG scheme (2.8), we also take the size of the cell \(S_{i,j}\) as \(\Delta x_i=\Delta y_j\), which will only effect the ratio of \(\alpha \) and \(\beta \). The formula for v is then given as

$$\begin{aligned} \begin{aligned} v(\xi ,\eta )&=\frac{\Delta x_i}{\varLambda _{i,j}}\big (44{\widetilde{\gamma }}^3 +3{\widetilde{\gamma }}^2(\alpha +1)f_1(\xi ,\eta ) +{\widetilde{\gamma }}(f_2(\xi ,\eta )+2\alpha f_3(\xi ,\eta )+\alpha ^2f_4(\xi ,\eta ))\\&+\,\alpha (f_5(\xi ,\eta )+2\alpha f_6(\xi ,\eta )+786\alpha ^2(1-\xi ))+786(1-\eta ) \big ), \end{aligned} \end{aligned}$$
(A.1)

where \(\varLambda _{i,j}\) reads as

$$\begin{aligned} \begin{aligned} \varLambda _{i,j}&=4\big (393+393\alpha ^4+536 {\widetilde{\gamma }}+310\gamma ^2+88{\widetilde{\gamma }}^3 +11{\widetilde{\gamma }}^4\\&+\,2\alpha ^3(667+268{\widetilde{\gamma }})+\alpha ^2 (1882+1481{\widetilde{\gamma }}+310{\widetilde{\gamma }}^2)\\&+\,\alpha (1334+1481{\widetilde{\gamma }}+583{\widetilde{\gamma }}^2 +88{\widetilde{\gamma }}^3)\big ) \end{aligned} \end{aligned}$$

with \({\widetilde{\gamma }}=\gamma \Delta x_i\) and

$$\begin{aligned} \begin{aligned} f_1(\xi ,\eta )&=106-5\xi ^2-44\eta -5\eta ^2-20\xi \eta ,\\ f_2(\xi ,\eta )&=820-15\xi ^2-15\eta ^2-558\eta -60\xi -60\xi \eta ,\\ f_3(\xi ,\eta )&=848+30\xi ^2+30\eta ^2+120\xi \eta -441\xi -441\eta ,\\ f_4(\xi ,\eta )&=820-15\xi ^2-15\eta ^2-558\xi -60\eta -60\xi \eta ,\\ f_5(\xi ,\eta )&=2206+300\xi ^2+300\eta ^2+1200\xi \eta -1158\xi -1944\eta ,\\ f_6(\xi ,\eta )&=1103+150\xi ^2+150\eta ^2-972\xi +600\xi \eta -579\eta . \end{aligned} \end{aligned}$$
(A.2)

We can easily prove that the functions \(f_k(\xi ,\eta ),k=1,\dots ,6\) are all positive when \((\xi ,\eta )\in [-1,1]^2\). Hence \(v(\xi ,\eta )\) is positive on \(S_{i,j}\).

Appendix B

In this appendix we prove that, for the unmodulated DG scheme, all the upstream solution values at the interfaces \(u_{i+\frac{1}{2}}^-\) remain positive, when the physical inflow condition and the source term are both positive and when the mesh size h is small enough. Consider the case for \(\alpha >0\) in (2.1), then (2.1) can simplified as \(u_x+\gamma u=f\). By means of the similar technique utilized before, we would like to seek a special test function which satisfies

$$\begin{aligned} \int _{S_i}u\left( \gamma v-v'+v_{i+\frac{1}{2}}^{-}\delta _k(x)\right) =u_{i+\frac{1}{2}}^{-} \end{aligned}$$
(B.1)

with \(\delta _k(x)\) defined in (3.1). Thanks to Lemma 3.1, we obtain an ordinary differential equation for v as follows:

$$\begin{aligned} \gamma v-v'+v_{i+\frac{1}{2}}^{-}\delta _k(x)=\delta _k(x). \end{aligned}$$
(B.2)

Once again, the fact that v(x) is a polynomial with degree of k enables us to formulate v(x) explicitly by the form of

$$\begin{aligned} v(x)=\left\{ \begin{array}{ll} \frac{\sum \nolimits _{l=0}^{k}\left( \frac{1}{\gamma }\right) ^{l+1} \delta _k^{(l)}(x)}{1+\sum \nolimits _{l=0}^{k} \left( \frac{1}{\gamma }\right) ^{l+1}\delta _k^{(l)}\left( x_{i+\frac{1}{2}}\right) }, &{}\quad \text {if}~\gamma >0,\\ 1,&{}\quad \text {if}~\gamma =0. \end{array}\right. \end{aligned}$$
(B.3)

Obviously, v(x) is positive when \(\gamma =0\). Now let us analyze the case \(\gamma >0\). We can express \(\delta _k^{(l)}(x)\) more explicitly as follows

$$\begin{aligned} \delta _k^{(l)}(x)=\frac{1}{\Delta x_i}\sum \limits _{n=0}^k(2n+1)p_n^{(l)}(\xi (x)) =2^l\left( \frac{1}{\Delta x_i}\right) ^{l+1}\sum \limits _{n=l}^k(2n+1)p_n^{(l)}(\xi ). \end{aligned}$$
(B.4)

Meanwhile, we can reformulate v(x) as

$$\begin{aligned} \begin{aligned} v(x)&=\frac{\sum \nolimits _{l=0}^{k}\left( \frac{1}{\gamma }\right) ^{l+1} \delta _k^{(l)}(x)}{1+\sum \nolimits _{l=0}^{k} \left( \frac{1}{\gamma }\right) ^{l+1}\delta _k^{(l)}\left( x_{i+\frac{1}{2}}\right) } =\frac{\sum \nolimits _{l=0}^{k}2^l\left( \frac{1}{\gamma \Delta x_i}\right) ^{l+1}\sum \nolimits _{n=l}^k(2n+1)p_n^{(l)}(\xi )}{1+\sum \nolimits _{l=0}^{k}2^l\left( \frac{1}{\gamma \Delta x_i}\right) ^{l+1}\sum \nolimits _{n=l}^k(2n+1)p_n^{(l)}(1)}\\&=\frac{\sum \nolimits _{l=0}^{k}2^l(\gamma \Delta x_i)^{k-l} \sum \nolimits _{n=l}^k(2n+1)p_n^{(l)}(\xi )}{(\gamma \Delta x_i)^{k+1}+\sum \nolimits _{l=0}^{k}2^l (\gamma \Delta x_i)^{k-l}\sum \nolimits _{n=l}^k(2n+1)p_n^{(l)}(1)}\\&=\frac{2^k(2k+1)p_k^{(k)}(\xi )+\sum \nolimits _{l=0}^{k-1}2^l (\gamma \Delta x_i)^{k-l}\sum \nolimits _{n=l}^k(2n+1)p_n^{(l)} (\xi )}{2^k(2k+1)p_k^{(k)}(1)+(\gamma \Delta x_i)^{k+1} +\sum \nolimits _{l=0}^{k-1} 2^l(\gamma \Delta x_i)^{k-l}\sum \nolimits _{n=l}^k(2n+1)p_n^{(l)}(1)}. \end{aligned} \end{aligned}$$
(B.5)

Here we notice that \(p_k^{(k)}(\xi )=p_k^{(k)}(1)\) is a positive constant, which reads as

$$\begin{aligned} p_k^{(k)}(\xi )=p_k^{(k)}(1)=\frac{(2k)!}{2^k\cdot k!} \end{aligned}$$

since \(p_k(\xi )=\frac{1}{2^k\cdot k!}\frac{d^n}{d\xi ^n}(\xi ^2-1)^k\). Therefore, when the mesh size \(\Delta x_i\) is small enough, v(x) will be positive with a value close to 1.

For smaller k, we also list the range of \(\gamma \Delta x_i\) to ensure the positivity of v. Let us denote \({\widetilde{\gamma }}=\gamma \Delta x_i\) once again. We can write v out as a function of the variable \(\xi \) for the cases \(k=1,2,3,4\) as follows:

  1. (1)

    \(k=1\):

    $$\begin{aligned} v(\xi )=\frac{6+{\widetilde{\gamma }}+3{\widetilde{\gamma }}\xi }{6+4{\widetilde{\gamma }}+{\widetilde{\gamma }}^2}. \end{aligned}$$

    It turns out that \(v(\xi )\ge 0\) if \({\widetilde{\gamma }}\le 3\).

  2. (2)

    \(k=2\):

    $$\begin{aligned} v(\xi )=\frac{120+12{\widetilde{\gamma }}(1+5\xi ) +3{\widetilde{\gamma }}^2(-1+2\xi +5\xi ^2)}{120+72{\widetilde{\gamma }} +18{\widetilde{\gamma }}^2+2{\widetilde{\gamma }}^3}. \end{aligned}$$

    \(v(\xi )\) will be non-negative when \({\widetilde{\gamma }}\le 5\sqrt{\frac{2}{3}}\).

  3. (3)

    \(k=3\):

    $$\begin{aligned} v(\xi )=\frac{1680+120{\widetilde{\gamma }}(1+7\xi ) +30{\widetilde{\gamma }}^2(-1+2\xi +7\xi ^2) +{\widetilde{\gamma }}^3(-3-15\xi +15\xi ^2+35\xi ^3)}{1680+960{\widetilde{\gamma }}+240{\widetilde{\gamma }}^2 +32{\widetilde{\gamma }}^3+2{\widetilde{\gamma }}^4}. \end{aligned}$$

    When \({\widetilde{\gamma }}\le 5-\frac{5^{2/3}}{(1+\sqrt{6})^{1/3}}+\big (5(1+\sqrt{6})\big )^{1/3}\approx 5.6484\) we have \(v(\xi )\ge 0\).

  4. (4)

    \(k=4\):

    $$\begin{aligned}\begin{aligned} v(\xi )&=\frac{5}{\Lambda }\big (24192+1344{\widetilde{\gamma }} (1+9\xi )+336{\widetilde{\gamma }}^2(-1+2\xi +9\xi ^2)\\&+\,24{\widetilde{\gamma }}^3(-1-7\xi +7\xi ^2+21 \xi ^3) +{\widetilde{\gamma }}^4(3-12\xi -42\xi ^2+28\xi ^3+63\xi ^4)\big ) \end{aligned} \end{aligned}$$

    where

    $$\begin{aligned} \Lambda =8(15120+8400{\widetilde{\gamma }}+2100{\widetilde{\gamma }}^2 +300{\widetilde{\gamma }}^3+25{\widetilde{\gamma }}^4 +{\widetilde{\gamma }}^5). \end{aligned}$$

    It is easy to prove that \(v(\xi )\ge 0\) if

    $$\begin{aligned} {\widetilde{\gamma }}\le \frac{1}{2}\left( 7-\frac{11\cdot 7^{2/3}}{(3(81+4\sqrt{2157}))^{1/3}} +\frac{(7(81+4\sqrt{2157}))^{1/3}}{3^{2/3}}\right) \approx 4.2923. \end{aligned}$$

Rights and permissions

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Ling, D., Cheng, J. & Shu, CW. Conservative High Order Positivity-Preserving Discontinuous Galerkin Methods for Linear Hyperbolic and Radiative Transfer Equations. J Sci Comput 77, 1801–1831 (2018). https://doi.org/10.1007/s10915-018-0700-3

Download citation

  • Received:

  • Revised:

  • Accepted:

  • Published:

  • Issue Date:

  • DOI: https://doi.org/10.1007/s10915-018-0700-3

Keywords

Navigation