Central-upwind Scheme Research Articles

Conservative Third-Order Central-Upwind Schemes for Option Pricing Problems

In this paper, we propose the application of third-order semi-discrete central-upwind conservative schemes to option pricing partial differential equations (PDEs). Our method is a high-order extension of the recent efficient second-order “Black-Box” schemes that successfully priced several option pricing problems. We consider the Kurganov–Levy scheme and its extensions, namely the Kurganov–Noelle–Petrova and the Kolb schemes. These “Black-Box” solvers ensure non-oscillatory property and achieve desired accuracy using a third-order central weighted essentially non-oscillatory (CWENO) reconstruction. We compare the schemes using a European test case and observe that the Kolb scheme performs better. We apply the Kolb scheme to one-dimensional butterfly, barrier, American and non-linear options under the Black–Scholes model. Further, we extend the Kurganov–Levy scheme to solve two-dimensional convection-dominated Asian PDE. We also price American options under the constant elasticity of variance (CEV) model, which treats volatility as a stochastic instead of a constant as in Black–Scholes model. Numerical experiments achieve third-order, non-oscillatory and high-resolution solutions.

A positivity-preserving central-upwind scheme for isentropic two-phase flows through deviated pipes

Directional drilling in oil and gas extraction can encounter difficulties such as accumulation of solids in deviated pipes. Motivated by such phenomenon, we consider a model for isentropic two-phase flows through deviated pipes. The system of partial differential equations is aimed at simulating the dynamics between a particle bed and a gas phase. The pipe can be either horizontal or vertically deviated where the effects of gravity are incorporated. Furthermore, the acceleration or deceleration due to friction between phases is investigated and spectral properties of the hyperbolic system of balance laws are described. The existence and characterization of steady states under appropriate conditions is analyzed. A new type of steady states arises when a balance between gas and solid phases results in a non-uniform solid particle bed and vanishing solid velocity. This state corresponds to an accumulation of sedimented solids. A central-upwind scheme that preserves the positivity of the gas and solid densities and volume fractions is presented. Including an application of the model to an analysis of accumulation of solids, a variety of numerical tests is presented to show the merits of the scheme.

One-Dimensional/Two-Dimensional Coupling Approach with Quadrilateral Confluence Region for Modeling River Systems

We study shallow water flows in river systems. An accurate description of such flows can be obtained using the two-dimensional (2-D) shallow water equations, which can be numerically solved by a shock-capturing finite-volume method. This approach can, however, be inefficient and computationally unaffordable when a large river system with many tributaries and complex geometry is to be modeled. A popular simplified approach is to model flow in each uninterrupted section of the river (called a reach) as one-dimensional (1-D) and connect the reaches at the river junctions. The flow in every reach can then be modeled using the 1-D shallow water equations, whose numerical solution is dramatically less computationally expensive compared with solving its 2-D counterpart. Even though several point-junction models are available, most of them prove to be sufficiently accurate only in the case of a smooth flow though the junction. We propose a new 1-D/2-D river junction model, in which each reach of the river is modeled by the 1-D shallow water equations, while the confluence region, where the mixing of flows from the different directions occurs, is modeled by the 2-D ones. We define the confluence region to be a trapezoid with parallel vertical sides. This allows us to take into account both the average width of each reach and the angle between the directions of flow of the tributary and the principal river at the confluence. We choose a trapezoidal confluence region as it is consistent with the 1-D model of the river. We implement well-balanced positivity preserving second-order semi-discrete central-upwind schemes developed in Kurganov and Petrova (Commun Math Sci 5:133–160, 2007) for the 1-D shallow water equations and in Shirkhani et al. (Comput Fluids 126: 25–40, 2016) for the 2-D shallow water equations using quadrilateral grids. For the 2-D junction simulations in the confluence region we choose a very coarse 2-D mesh as the goal of our model is not to resolve the fine details of complex 2-D vortices that form around the junction, but to efficiently compute average water depth and velocity in the connected 1-D reaches. A special ghost cell technique is developed for coupling the reaches to the confluence region, which is one of the most important parts of a good 1-D/2-D coupling method. The proposed approach leads to very significant computational savings compared to numerically solving the full 2-D problem. We perform several numerical experiments to demonstrate plausibility of the proposed 1-D/2-D coupling model.

Debris flows with pore pressure and intergranular friction on rugged topography

The dynamic behavior of debris flows features the interplay of a non-hydrostatic pore-fluid pressure with the non-linear deformational behavior of the granular skeleton and the internal contact stress between grains. This complex physical background is considered by amending the classical depth-integrated modeling for granular-fluid flows by two additional fields, an extra pore-fluid pressure and a hypoplastic intergranular stress. A scaled and depth-integrated model is developed and transferred into a system of terrain-following coordinates, enabling the application on rugged topography. With this model, numerical investigations are carried out, using a non-oscillatory, shock-capturing central-upwind scheme. Parameter studies show the general impact of the additional fields, completed by comparison to the experimental results of a dam break scenario. Furthermore, application to the landslide event at the village of Hsiaolin in Taiwan, 2009, show the capability of the model to cope with large scale scenarios. The results show that the model and its implementation provide insights in the flow dynamics and the possibility to application on complex topography, considering an enhanced approach to the physics of debris flows.

Path-conservative central-upwind schemes for nonconservative hyperbolic systems

We develop path-conservative central-upwind schemes for nonconservative one-dimensional hyperbolic systems of nonlinear partial differential equations. Such systems arise in a variety of applications and the most challenging part of their numerical discretization is a robust treatment of nonconservative product terms. Godunov-type central-upwind schemes were developed as an efficient, highly accurate and robust ``black-box’’ solver for hyperbolic systems of conservation and balance laws. They were successfully applied to a large number of hyperbolic systems including several nonconservative ones. To overcome the difficulties related to the presence of nonconservative product terms, several special techniques were proposed. However, none of these techniques was sufficiently robust and thus the applicability of the original central-upwind schemes was rather limited. In this paper, we rewrite the central-upwind schemes in the form of path-conservative schemes. This helps us (i) to show that the main drawback of the original central-upwind approach was the fact that the jump of the nonconservative product terms across cell interfaces has never been taken into account and (ii) to understand how the nonconservative products should be discretized so that their influence on the numerical solution is accurately taken into account. The resulting path-conservative central-upwind scheme is a new robust tool for both conservative and nonconservative hyperbolic systems. We apply the new scheme to the Saint-Venant system with discontinuous bottom topography and two-layer shallow water system. Our numerical results illustrate the good performance of the new path-conservative central-upwind scheme, its robustness and ability to achieve very high resolution.

An asymptotic preserving scheme for the two-dimensional shallow water equations with Coriolis forces

We consider the two-dimensional Saint-Venant system of shallow water equations with Coriolis forces. We focus on the case of a low Froude number, in which the system is stiff and conventional explicit numerical methods are extremely inefficient and often impractical. Our goal is to design an asymptotic preserving (AP) scheme, which is uniformly asymptotically consistent and stable for a broad range of (low) Froude numbers. The goal is achieved using the flux splitting proposed in [Haack et al., Commun. Comput. Phys., 12 (2012), pp. 955–980] in the context of isentropic Euler and Navier-Stokes equations. We split the flux into the stiff and nonstiff parts and then use an implicit-explicit approach: apply an explicit hyperbolic solver (we use the second-order central-upwind scheme) to the nonstiff part of the system while treating the stiff part of it implicitly. Moreover, the stiff part of the flux is linear and therefore we reduce the implicit stage of the proposed method to solving a Poisson-type elliptic equation, which is discretized using a standard second-order central difference scheme.We conduct a series of numerical experiments, which demonstrate that the developed AP scheme achieves the theoretical second-order rate of convergence and the time-step stability restriction is independent of the Froude number. This makes the proposed AP scheme an efficient and robust alternative to fully explicit numerical methods.

A New Approach for Designing Moving-Water Equilibria Preserving Schemes for the Shallow Water Equations

We construct a new second-order moving-water equilibria preserving central-upwind scheme for the one-dimensional Saint-Venant system of shallow water equations. The idea is based on a reformulation of the source terms as integral in the flux function. Reconstruction of the flux variable yields then a third order equation that can be solved exactly. This procedure does not require any further modification of existing schemes. Several numerical tests are performed to verify the ability of the proposed scheme to accurately capture small perturbations of steady states.

Comparison of Shallow Water Solvers: Applications for Dam-Break and Tsunami Cases with Reordering Strategy for Efficient Vectorization on Modern Hardware

We investigate in this paper the behaviors of the Riemann solvers (Roe and Harten-Lax-van Leer-Contact (HLLC) schemes) and the Riemann-solver-free method (central-upwind scheme) regarding their accuracy and efficiency for solving the 2D shallow water equations. Our model was devised to be spatially second-order accurate with the Monotonic Upwind Scheme for Conservation Laws (MUSCL) reconstruction for a cell-centered finite volume scheme—and be temporally fourth-order accurate using the Runge–Kutta fourth-order method. Four benchmark cases of dam-break and tsunami events dealing with highly-discontinuous flows and wet–dry problems were simulated. To this end, we applied a reordering strategy for the data structures in our code supporting efficient vectorization and memory access alignment for boosting the performance. Two main features are pointed out here. Firstly, the reordering strategy employed has enabled highly-efficient vectorization for the three solvers investigated on three modern hardware (AVX, AVX2, and AVX-512), where speed-ups of 4.5–6.5× were obtained on the AVX/AVX2 machines for eight data per vector while on the AVX-512 machine we achieved a speed-up of up to 16.7× for 16 data per vector, all with singe-core computation; with parallel simulations, speed-ups of up to 75.7–121.8× and 928.9× were obtained on AVX/AVX2 and AVX-512 machines, respectively. Secondly, we observed that the central-upwind scheme was able to outperform the HLLC and Roe schemes 1.4× and 1.25×, respectively, by exhibiting similar accuracies. This study would be useful for modelers who are interested in developing shallow water codes.

Multiresolution-based adaptive central high resolution schemes for modeling of nonlinear propagating fronts

A class of wavelet-based methods, which combine adaptive grids together with second-order central and central-upwind high-resolution schemes, is introduced for accurate solution of first order hyperbolic systems of conservation laws and related equations. It is shown that the proposed class recovers stability which otherwise may be lost when the underlying central schemes are utilized over non-uniform grids. For simulations on irregular grids, both full-discretized and semi-discretized forms are derived; in particular, the effect of using certain shifts of the center of the computational allow to use standard slope limiters on non-uniform cells. Thereafter, the questions of numerical entropy production, local truncation errors, and the Total Variation Diminishing (TVD) criterion for scalar equations, are investigated. Central schemes are sensitives for irregular grids. To cure this drawback in the present context of multiresolution, the adapted grid is locally modified around high-gradient zones by adding new neighboring points. For an adapted point of resolution scale j, the new points are locally added in both resolution j and all successive coarser resolutions. It is shown that this simple grid modification improves stability of numerical solutions; this, however, comes at the expense that cell-centered central and central-upwind high resolution schemes, may not satisfy the TVD property. Finally, we present numerical simulations for both scalar and systems of non-linear conservation laws which confirm the simplicity and efficiency of the proposed method. These simulations demonstrate the high accuracy, the entropy production of wavelet-based algorithms which can effectively detect high gradient zones of shock waves, rarefaction regions, and contact discontinuities.

A robust numerical model for shallow water governing solute transport with wet/dry interfaces

In this study, a robust numerical model is developed to guarantee that the C-property for both hydrodynamic and solute transport models and positivity of both water depth and solute concentration. Moreover, special attention is paid to satisfying such properties in the unstructured triangular cells near wet–dry interfaces. In the developed numerical model, the finite-volume central-upwind scheme is adopted to approximate the advective fluxes, and Green’s theorem is applied to estimate the diffusive fluxes. The C-property for hydrodynamic model near wet–dry area is achieved by introducing an auxiliary free surface variable used in the piecewise linear reconstruction and discretizing the bed slope source term in a special form. And the C-property for solute transport is guaranteed by reconstructing secondary variables instead of the primary variables within the frame of central-upwind method, and using the coupling approach for the governing system. The positivity of the solute concentration is satisfied using a proposed diluting time step restriction and the positivity of the water depth is guaranteed by a draining time technique. The desired properties of the proposed numerical model are verified by several numerical experiments.

Numerical method to simulate detonative combustion of hydrogen-air mixture in a containment

During an accident accompanied by core damage in a water-cooled nuclear reactor, a large amount of hydrogen is generated by fuel clad oxidation and released into the reactor containment. Hydrogen mitigation in the containment during a severe accident is important because the integrity of the structures important to safety could be threatened by an explosive combustion of the hydrogen. Even when equipped with a hydrogen mitigation system, if a highly concentrated hydrogen mixture cloud is developed in a compartment of the containment, the compartment integrity with a pressure load from a local detonation must be evaluated. In this study, a reliable numerical method to simulate the detonation of a hydrogen-air mixture was developed for the evaluation of a pressure load from a detonation. The numerical method has the capability to capture shocks robustly as well as resolve combustion phenomena stably. The Euler equations of mass, momentum, energy and gas species are employed as the governing equations of detonation. A central-upwind scheme was adopted as a shock capturing scheme and a reduced chemical kinetic model was applied for combustion. The numerical method was validated by solving experimental cases and applied for understanding detonation shock behaviour in a generic containment.

Efficient conservative second-order central-upwind schemes for option-pricing problems

The conservative Kurganov–Tadmor (KT) scheme has been successfully applied to option-pricing problems by Germán I. Ramírez-Espinoza and Matthias Ehrhardt. These included the valuation of European, Asian and nonlinear options as Black–Scholes partial differential equations, written in the conservative form, by simply updating fluxes in the black box approach. In this paper, we describe an improvement of this idea through a fully vectorized algorithm of nonoscillatory slope limiters and the efficient use of time solvers. We also propose the application of second-order extensions of KT to option-pricing problems. Our test problems solve one-dimensional benchmark and convection-dominated European options as well as digital and butterfly options. These demonstrate the robustness and flexibility of the pricing methods and set a basis for complex problems. Further, the computation of option Greeks ensures the reliability of these methods. Numerical experiments are performed on barrier options, early exercisable American options and two-dimensional fixed and floating strike Asian options. To the authors’ knowledge, this is the first time American options have been priced by applying the early exercise condition on the semi-discrete formulation of central-upwind schemes. Our results show second-order, nonoscillatory and high-resolution properties of the schemes as well as computational efficiency.

Central-upwind scheme for 2D turbulent shallow flows using high-resolution meshes with scalable wall functions

In this paper, the Reynolds-averaged Navier Stokes equations supplemented by the algebraic stress model are solved using the central-upwind scheme for simulating 2D turbulent shallow flows. The model is of spatially and temporally second-order accurate. To increase the accuracy, high-resolution meshes up to 3.4 million cells (6.8 million edges) are used. Consequently, a strategy combining the hydrostatic and topography reconstructions and the scalable wall functions – is proposed to accurately simulate wet–dry phenomena near the interfaces (moving boundary geometries) thus ensuring a proper calculation for the turbulence properties. This strategy has been proven to be accurate and to not deteriorate the results for such very fine meshes thus giving flexibility to users in generating meshes.

Analysis of nonequilibrium effects and flow instability in immiscible two-phase flow in porous media

Two-dimensional, high-resolution, numerical solutions for the classical formulation and two widely accepted nonequilibrium models of multiphase flow through porous media are generated and compared with experimental observations from literature. Flow equations for simultaneous flow of two immiscible phases through porous media are written in a vorticity stream-function form. In the resulting system of equations, the vorticity stream-function equation is solved using a spectral method and the transport equation is discretized in space using a central-upwind scheme. The system of equations is solved for a two-dimensional domain using a semi-implicit time-stepper. The solutions reveal behavior that is not apparent in one-dimensional solutions, namely that sharpening of the saturation front caused by the inclusion of a dynamic capillary pressure results in propagation of viscous fingers compared to the classical formulation. The inclusion of nonequilibrium effects in the constitutive relations, in the form of effective saturation, introduces dispersion and smears the otherwise highly resolved viscous fingers in the saturation front. Once developed, the length of the mixing zone in the numerical solutions remains constant with time regardless of the degree of instability. This is contrary to the evolution of the mixing zone observed in unstable flow experiments where, unlike the numerical solutions, the propagation speed of the leading edge of the front appears to increase with time.

On the role of pore-fluid pressure evolution and hypoplasticity in debris flows

Abstract For debris flows, two characteristics play a crucial role in their behaviors: (i) the non-linear deformational behavior that stems from the internal contact stress between grains, and (ii) the pore-pressure feedback, which reduces the intergranular friction and, therefore, enhances the mobility of the whole mixture. In a previous work, the constitutive equations of a general grain–fluid multiphase mixture are deduced by exploiting the entropy principle based on the formulation of Muller and Liu. In the work at hand, this thermodynamically consistent and general model is applied to a two-phase debris flow. For this, it is scaled, depth-integrated, compared with other known models and employed for numerical investigations. The extra pore-fluid pressure and the hypoplastic stress are described as internal variables that are, respectively, governed by a pressure diffusion equation and a transport equation related to the hypoplastic material. Numerical investigation of a grain–fluid material, sliding down a slope into a runout, is carried out using a non-oscillatory, shock-capturing central-upwind scheme with the total variation diminishing property. Comparison of the obtained results with those from classical debris flow models shows that the proposed thermodynamic model provides a phenomenological insight into the influence of the pore-pressure feedback and intergranular friction in the flow dynamics. The extra pore-fluid pressure enhances the flow dynamics, while the effect of the intergranular stress is more significant during the settling of the flow and affects the shape of the granular bulk.

Well-balanced positivity preserving central-upwind scheme with a novel wet/dry reconstruction on triangular grids for the Saint-Venant system

In this paper, we construct an improved well-balanced positivity preserving central-upwind scheme for the two-dimensional Saint-Venant system of shallow water equations. As in Bryson et al. (2011) [7], our scheme is based on a continuous piecewise linear discretization of the bottom topography over an unstructured triangular grid. The main new technique is a special reconstruction of the water surface in partially flooded cells. This reconstruction is an extension of the one-dimensional wet/dry reconstruction from Bollermann et al. (2013) [3]. The positivity of the computed water depth is enforced using the “draining” time-step technique introduced in Bollermann et al. (2011) [4]. The performance of the proposed central-upwind scheme is tested on a number of numerical experiments.

An Alternative Formulation of Discontinous Galerkin Schemes for Solving Hamilton–Jacobi Equations

The aim of this paper is to develop an alternative formulation of discontinuous Galerkin (DG) schemes for approximating the viscosity solutions to nonlinear Hamilton–Jacobi (HJ) equations. The main difficulty in designing DG schemes lies in the inherent non-divergence form of HJ equations. One effective approach is to explore the elegant relationship between HJ equations and hyperbolic conservation laws: the standard DG scheme is applied to solve a conservation law system satisfied by the derivatives of the solution of the HJ equation. In this paper, we consider an alternative approach to directly solving the HJ equations, motivated by a class of successful direct DG schemes by Cheng and Shu (J Comput Phys 223(1):398–415, 2007), Cheng and Wang (J Comput Phys 268:134–153, 2014). The proposed scheme is derived based on the idea from the central-upwind scheme by Kurganov et al. (SIAM J Sci Comput 23(3):707–740, 2001). In particular, we make use of precise information of the local speeds of propagation at discontinuous element interface with the goal of adding adequate numerical viscosity and hence naturally capturing the viscosity solutions. A collection of numerical experiments is presented to demonstrate the performance of the method for solving general HJ equations with linear, nonlinear, smooth, non-smooth, convex, or non-convex Hamiltonians.

A well-balanced positivity preserving two-dimensional shallow flow model with wetting and drying fronts over irregular topography

This paper presents an improved well-balanced Godunov-type 2-D finite volume model with structured grids to simulate shallow flows with wetting and drying fronts over an irregular topography. The intercell flux is computed using a central upwind scheme, which is a Riemann-problem-solver-free method for hyperbolic conservation laws. The nonnegative reconstruction method for the water depth is implemented to resolve the stationary or wet/dry fronts. The bed slope source term is discretized using a central difference method to capture the static flow state over the irregular topography. Second-order accuracy in space is achieved by using the slope limited linear reconstruction method. With the proposed method, the model can avoid the partially wetting/drying cell problem and maintain the mass conservation. The proposed model is tested and verified against three theoretical benchmark tests and two experimental dam break flows. Further, the model is applied to predict the maximum water level and the flood arrival time at different gauge points for the Malpasset dam break event. The predictions agree well with the numerical results and the measurement data published in literature, which demonstrates that with the present model, a well-balanced state can be achieved and the water depth can be nonnegative when the Courant number is kept less than 0.25.

Numerical simulation of Argon fuelled self-field magnetoplasmadynamic thrusters using the central-upwind scheme flux interpolations

In the present paper, we describe the implementation of a density-based scheme including a divergence cleaning method for the simulation of the magnetohydrodynamic (MHD) equations under a finite volume formulation. This new algorithm was developed for both ideal and resistive MHD equations and makes use of the central-upwind schemes of Kurganov and Tadmor for flux calculation. The developed density-based code is validated through comparisons of the obtained results with experimental data and simulations of previous works. The obtained results are in good agreement with the experimental data and previous simulation results.

A new sequential method for three-phase immiscible flow in poroelastic media

A new computational model is developed to solve three phase immiscible compressible flow in strongly heterogeneous poroelastic media. Within the context of the iterative coupled formulation based on the fixed stress split algorithm the governing equations are decomposed into three subsystems associated with the geomechanics, hydrodynamics and transport problems. The hydrodynamics involves an overall compressibility which plays essential role in the magnitude of additional source in the pressure equation associated the time derivative of the total mean stress. We show that the additional coefficients in the source can be decomposed into a rock/grain and fluid components involving the formation volume factors of each phase. The flow equations are written in a proper conservative form for application of mixed-hybrid formulation whereas the nonlinear coupling between flow and transport is handled by a proper sequential iterative scheme. A new post-processing approach is proposed to update the Lagrangian porosity based on the freezing of the total mean stress which reproduces the one-way formulation in a straightforward fashion under the stationarity of the total mean stress. The geomechanical subsystem is also discretized by a mixed formulation based on the decomposition of the effective stress into spherical and deviatoric components in order to handle the nearly-incompressible undrained limit. The system of conservation laws for the liquid saturations is rephrased in an alternative form for an appropriate operator splitting method, leading to the appearance of an extra source term arising from the transient behavior of the Lagrangian porosity induced by the geomechanical coupling. Within the framework of a predictor–corrector scheme, the predictor step is discretized by an extended version of a higher order central-upwind scheme whereas the corrector, which captures the influence of rock skeleton deformation upon transport, is captured by a Godunov splitting. The potential of the new computational model is illustrated in numerical simulations of water-flooding and of WAG (Water Alternating Gas) injection problems in poroelastic media. Numerical results illustrate precisely the role of the transient nature of the total stress upon subsidence, rock compaction, oil and gas production and finger grow in primary and secondary recovery regimes.