The thermal radiative transfer (TRT) equations form an integro-differential system that describes the propagation and collisional interactions of photons. Computing accurate and efficient numerical solutions TRT are challenging for several reasons, the first of which is that TRT is defined on a high-dimensional phase space that includes the independent variables of time, space, and velocity. In order to reduce the dimensionality of the phase space, classical approaches such as the P$_N$ (spherical harmonics) or the S$_N$ (discrete ordinates) ansatz are often used in the literature. In this work, we introduce a novel approach: the hybrid discrete (H$^T_N$) approximation to the radiative thermal transfer equations. This approach acquires desirable properties of both P$_N$ and S$_N$, and indeed reduces to each of these approximations in various limits: H$^1_N$ $\equiv$ P$_N$ and H$^T_0$ $\equiv$ S$_T$. We prove that H$^T_N$ results in a system of hyperbolic partial differential equations for all $T\ge 1$ and $N\ge 0$. Another challenge in solving the TRT system is the inherent stiffness due to the large timescale separation between propagation and collisions, especially in the diffusive (i.e., highly collisional) regime. This stiffness challenge can be partially overcome via implicit time integration, although fully implicit methods may become computationally expensivemore »
A self-adaptive theta scheme using discontinuity aware quadrature for solving conservation laws
Abstract We present a discontinuity aware quadrature (DAQ) rule and use it to develop implicit self-adaptive theta (SATh) schemes for the approximation of scalar hyperbolic conservation laws. Our SATh schemes require the solution of a system of two equations, one controlling the cell averages of the solution at the time levels, and the other controlling the space-time averages of the solution. These quantities are used within the DAQ rule to approximate the time integral of the hyperbolic flux function accurately, even when the solution may be discontinuous somewhere over the time interval. The result is a finite volume scheme using the theta time stepping method, with theta defined implicitly (or self-adaptively). Two schemes are developed, self-adaptive theta upstream weighted (SATh-up) for a monotone flux function using simple upstream stabilization, and self-adaptive theta Lax–Friedrichs (SATh-LF) using the Lax–Friedrichs numerical flux. We prove that DAQ is accurate to second order when there is a discontinuity in the solution and third order when it is smooth. We prove that SATh-up is unconditionally stable, provided that theta is set to be at least 1/2 (which means that SATh can be only first order accurate in general). We also prove that SATh-up satisfies the maximum more »
- Award ID(s):
- 1912735
- Publication Date:
- NSF-PAR ID:
- 10356572
- Journal Name:
- IMA Journal of Numerical Analysis
- ISSN:
- 0272-4979
- Sponsoring Org:
- National Science Foundation
More Like this
-
-
Simulation of flow and transport in petroleum reservoirs involves solving coupled systems of advection-diffusion-reaction equations with nonlinear flux functions, diffusion coefficients, and reactions/wells. It is important to develop numerical schemes that can approximate all three processes at once, and to high order, so that the physics can be well resolved. In this paper, we propose an approach based on high order, finite volume, implicit, Weighted Essentially NonOscillatory (iWENO) schemes. The resulting schemes are locally mass conservative and, being implicit, suited to systems of advection-diffusion-reaction equations. Moreover, our approach gives unconditionally L-stable schemes for smooth solutions to the linear advection-diffusion-reaction equation in the sense of a von Neumann stability analysis. To illustrate our approach, we develop a third order iWENO scheme for the saturation equation of two-phase flow in porous media in two space dimensions. The keys to high order accuracy are to use WENO reconstruction in space (which handles shocks and steep fronts) combined with a two-stage Radau-IIA Runge-Kutta time integrator. The saturation is approximated by its averages over the mesh elements at the current time level and at two future time levels; therefore, the scheme uses two unknowns per grid block per variable, independent of the spatial dimension. Thismore »
-
The Langevin Dynamics (LD) method (also known in the literature as Brownian Dynamics) is routinely used to simulate aerosol particle trajectories for transport rate constant calculations as well as to understand aerosol particle transport in internal and external fluid flows. This tutorial intends to explain the methodological details of setting up a LD simulation of a population of aerosol particles and to deduce rate constants from an ensemble of classical trajectories. We discuss the applicability and limitations of the translational Langevin equation to model the combined stochastic and deterministic motion of particles in fields of force or fluid flow. The drag force and stochastic “diffusion” force terms that appear in the Langevin equation are discussed elaborately, along with a summary of common forces relevant to aerosol systems (electrostatic, gravity, van der Waals, …); a commonly used first order and a fourth order Runge-Kutta time stepping schemes for linear stochastic ordinary differential equations are presented. A MATLAB® implementation of a LD code for simulating particle settling under gravity using the first order scheme is included for illustration. Scaling analysis of aerosol transport processes and the selection of timestep and domain size for trajectory simulations are demonstrated through two specific aerosol processes:more »
-
Seismic waves in earth media usually undergo attenuation, causing energy losses and phase distortions. In the regime of high-frequency asymptotics, a complex-valued eikonal is an essential ingredient for describing wave propagation in attenuating media, where the real and imaginary parts of the eikonal function capture dispersion effects and amplitude attenuation of seismic waves, respectively. Conventionally, such a complex-valued eikonal is mainly computed either by tracing rays exactly in complex space or by tracing rays approximately in real space so that the resulting eikonal is distributed irregularly in real space. However, seismic data processing methods, such as prestack depth migration and tomography, usually require uniformly distributed complex-valued eikonals. Therefore, we have developed a unified framework to Eulerianize several popular approximate real-space ray-tracing methods for complex-valued eikonals so that the real and imaginary parts of the eikonal function satisfy the classic real-space eikonal equation and a novel real-space advection equation, respectively, and we dub the resulting method the Eulerian partial-differential-equation method. We further develop highly efficient high-order methods to solve these two equations by using the factorization idea and the Lax-Friedrichs weighted essentially nonoscillatory schemes. Numerical examples demonstrate that our method yields highly accurate complex-valued eikonals, analogous to those from ray-tracing methods.more »
-
The YBJ equation (Young & Ben Jelloul, J. Marine Res. , vol. 55, 1997, pp. 735–766) provides a phase-averaged description of the propagation of near-inertial waves (NIWs) through a geostrophic flow. YBJ is obtained via an asymptotic expansion based on the limit $\mathit{Bu}\rightarrow 0$ , where $\mathit{Bu}$ is the Burger number of the NIWs. Here we develop an improved version, the YBJ + equation. In common with an earlier improvement proposed by Thomas, Smith & Bühler ( J. Fluid Mech. , vol. 817, 2017, pp. 406–438), YBJ + has a dispersion relation that is second-order accurate in $\mathit{Bu}$ . (YBJ is first-order accurate.) Thus both improvements have the same formal justification. But the dispersion relation of YBJ + is a Padé approximant to the exact dispersion relation and with $\mathit{Bu}$ of order unity this is significantly more accurate than the power-series approximation of Thomas et al. (2017). Moreover, in the limit of high horizontal wavenumber $k\rightarrow \infty$ , the wave frequency of YBJ + asymptotes to twice the inertial frequency $2f$ . This enables solution of YBJ + with explicit time-stepping schemes using a time step determined by stable integration of oscillations with frequency $2f$ . Other phase-averaged equations have dispersion relations with frequency increasing as $k^{2}$more »