Finite Difference Framework Research Articles

Local Discontinuous Galerkin Methods for One-Dimensional Second Order Fully Nonlinear Elliptic and Parabolic Equations

This paper is concerned with developing accurate and efficient nonstandard discontinuous Galerkin methods for fully nonlinear second order elliptic and parabolic partial differential equations (PDEs) in the case of one spatial dimension. The primary goal of the paper to develop a general framework for constructing high order local discontinuous Galerkin (LDG) methods for approximating viscosity solutions of these fully nonlinear PDEs which are merely continuous functions by definition. In order to capture discontinuities of the first order derivative $$u_x$$ u x of the solution $$u$$ u , two independent functions $$q^-$$ q - and $$q^+$$ q + are introduced to approximate one-sided derivatives of $$u$$ u . Similarly, to capture the discontinuities of the second order derivative $$u_{xx}$$ u x x , four independent functions $$p^{- -}, p^{- +}, p^{+ -}$$ p - - , p - + , p + - , and $$p^{+ +}$$ p + + are used to approximate one-sided derivatives of $$q^-$$ q - and $$q^+$$ q + . The proposed LDG framework, which is based on a nonstandard mixed formulation of the underlying PDE, embeds a given fully nonlinear problem into a mostly linear system of equations where the given nonlinear differential operator must be replaced by a numerical operator which allows multiple value inputs of the first and second order derivatives $$u_x$$ u x and $$u_{xx}$$ u x x . An easy to verify set of criteria for constructing good numerical operators is also proposed. It consists of consistency and generalized monotonicity. To ensure such a generalized monotonicity property, the crux of the construction is to introduce the numerical moment in the numerical operator, which plays a critical role in the proposed LDG framework. The generalized monotonicity gives the LDG methods the ability to select the viscosity solution among all possible solutions. The proposed framework extends a companion finite difference framework developed by Feng and Lewis (J Comp Appl Math 254:81---98, 2013) and allows for the approximation of fully nonlinear PDEs using high order polynomials and non-uniform meshes. Numerical experiments are also presented to demonstrate the accuracy, efficiency and utility of the proposed LDG methods.

Mimetic discretization of two-dimensional magnetic diffusion equations

In case of non-constant resistivity, cylindrical coordinates, and highly distorted polygonal meshes, a consistent discretization of the magnetic diffusion equations requires new discretization tools based on a discrete vector and tensor calculus. We developed a new discretization method using the mimetic finite difference framework. It is second-order accurate on arbitrary polygonal meshes and a consistent calculation of the Joule heating is intrinsic within it. The second-order convergence rates in L2 and L1 norms were verified with numerical experiments.

Convergent finite difference methods for one-dimensional fully nonlinear second order partial differential equations

This paper develops a new framework for designing and analyzing convergent finite difference methods for approximating both classical and viscosity solutions of second order fully nonlinear partial differential equations (PDEs) in 1-D. The goal of the paper is to extend the successful framework of monotone, consistent, and stable finite difference methods for first order fully nonlinear Hamilton–Jacobi equations to second order fully nonlinear PDEs such as Monge–Ampère and Bellman type equations. New concepts of consistency, generalized monotonicity, and stability are introduced; among them, the generalized monotonicity and consistency, which are easier to verify in practice, are natural extensions of the corresponding notions of finite difference methods for first order fully nonlinear Hamilton–Jacobi equations. The main component of the proposed framework is the concept of a “numerical operator”, and the main idea used to design consistent, generalized monotone and stable finite difference methods is the concept of a “numerical moment”. These two new concepts play the same roles the “numerical Hamiltonian” and the “numerical viscosity” play in the finite difference framework for first order fully nonlinear Hamilton–Jacobi equations. In the paper, two classes of consistent and monotone finite difference methods are proposed for second order fully nonlinear PDEs. The first class contains Lax–Friedrichs-like methods which also are proved to be stable, and the second class contains Godunov-like methods. Numerical results are also presented to gauge the performance of the proposed finite difference methods and to validate the theoretical results of the paper.

Estimation of operating blast furnace reactor invisible interior surface using Differential Evolution

A blast furnace reactor that transforms processed minerals into liquid iron is the most important link in a steel making process chain. Estimating the condition of refractory shield in a running furnace is important for production campaign planning and hazardous leakage prevention. A method coupling Differential Evolution with heat transfer simulations is developed for estimating refractory inner surface within process real times. The fitness value that drives this evolutionary algorithm is calculated from the differences between temperatures measured from thermocouples embedded in the refractory, and those obtained from numerical solution of the heat transfer equations. A population of different refractory inner surfaces constitutes the candidate solution set. These evolving surfaces tend to become craggy which impact variably on accuracy of numerical solution with significant consequence on the final synthesized solution. These inter-candidate accuracy variations are neutralized by transforming each physical refractory grid configuration onto identical rectangular computational domains using a set of continuity-ensuring elliptic partial differential equations. The heat transfer equation is correspondingly transformed into the computational space and solved within a finite difference framework. A set of practical constraints are imposed on the Differential Evolution algorithm to guide evolution toward the optimal surface that best represents the measured temperatures. Rigorous verification of the various stages of development facilitating industrial acceptance of the method are presented. This simulation-embedded evolutionary inverse-design approach can be used generically for online estimation of unobservable critical parameters of infrastructure, equipment and processes in diverse domains.

A dispersively accurate compact finite difference method for the Degasperis–Procesi equation

In this paper we are aimed to solve the non-dissipative Degasperis–Procesi equation based on the u-P formulation. To resolve the computational difficulty at the wave crest where the first-order derivative may diverge and the shockpeakon solution may form, the first-order spatial derivative term in the two-step equations will be approximated in a conservative form. The resulting equations will be approximated by the symplecticity-preserving time-stepping scheme and the spatial discretization scheme that can optimize the numerical wavenumber for the first-order spatial derivative term. This scheme will be developed in a three-point grid stencil with the accuracy order of seventh within the combined compact finite difference framework. Besides the validation of numerical accuracy, we will in particular address the discrete conservation of Hamiltonians even when peakon collides with antipeakon and generates, as a result, a shockpeakon. We will also demonstrate the capability of applying the proposed numerical method to sharply resolve some important features of the third-order dispersive DP equation.

Sequential real rainbow options

We develop two models to value European sequential rainbow options. The first model is a sequential option on the better of two stochastic assets, where these assets follow correlated geometric Brownian motion processes. The second model is a sequential option on the mean-reverting spread between two assets, which is applicable if the assets are co-integrated. We provide numerical solutions in the form of finite difference frameworks and compare these with Monte Carlo simulations. For the sequential option on a mean-reverting spread, we also provide a closed-form solution. Sensitivity analysis provides the interesting results that in particular circumstances, the sequential rainbow option value is negatively correlated with the volatility of one of the two assets, and that the sequential option on the spread does not necessarily increase in value with a longer time to maturity. With given maturity dates, it is preferable to have less time until expiry of the sequential option if the current spread level is way above the long-run mean.

An Efficient 3D Model of Heterogeneous Materials for Elastic Contact Applications Using Multigrid Methods

A 3D graded coating/substrate model based on multigrid techniques within a finite difference frame work is presented. Localized refinement is implemented to optimize memory requirement and computing time. Validation of the solver is performed through a comparison with analytical results for (i) a homogeneous material and (ii) a graded material. The algorithm performance is analyzed through a parametric study describing the influence of layer thickness (0.01 < t/a < 10) and mechanical properties (0.005 < Ec/Es < 10) of the coating on the contact parameters (Ph, a). Three-dimensional examples are then presented to illustrate the efficiency and the large range of possibilities of the model. The influence of different gradations of Young’s modulus, constant, linear and sinusoidal, through the coating thickness on the maximum tensile stress is analyzed, showing that the sinusoidal gradation best accommodates the property mismatch of two successive layers. A final case is designed to show that full 3D spatial property variations can be accounted for. Two spherical inclusions of different size made from elastic solids with Young’s modulus and Poisson’s ratio are embedded within an elastically mismatched finite domain and the stress field is computed.

Detailed modelling of a moving heat source using multigrid methods

The actual contact between solid surfaces is generally rough and time dependent. This roughness can only be correctly described using a 3D model. Moreover, severe thermo-mechanical loading induces high stress and temperature gradients within a thin subsurface layer. Even higher gradients are obtained in the case of surface coatings. Consequently, the problem is strongly multiscale: from contact to coating to roughness, the characteristic dimensions range from millimeter to nanometer. A straightforward discretization of this multiscale problem would exceed the memory and CPU capabilities of current computers. The aim of this work is to propose a more efficient numerical model able to deal with this multiscale problem: using 10 9 points and 10 3 timesteps. This paper traces the history of the numerical solutions of the heat equation from the pioneering work of Carslaw and Jaeger, to the current era. The proposed model is based on multigrid techniques within a finite difference frame work. Localized refinement is implemented to optimize memory and computing time costs. The numerical performance of the solver is presented through a comparison with analytical results using different types of boundary conditions. A multisource contact is studied as a first approximation to real asperity interaction.

Simulating Three-Dimensional Free Surface Viscoelastic Flows using Moving Finite Difference Schemes

An efficient finite difference framework based on moving meshes methods is developed for the three-dimensional free surface viscoelastic flows. The basic model equations are based on the incompressible Navier-Stokes equations and the Oldroyd-B constitutive model for viscoelastic flows is adopted. A logical domain semi-Lagrangian scheme is designed for moving-mesh solution interpolation and convection. Numerical results show that harmonic map based moving mesh methods can achieve better ac- curacy for viscoelastic flows with free boundaries while using much less memory and computational time compared to the uniform mesh simulations. AMS subject classifications: 65M20, 65N22

Hybrid finite-difference integral equation solver for 3D frequency domain anisotropic electromagnetic problems

The modeling of the controlled-source electromagnetic (CSEM) and single-well and crosswell electromagnetic (EM) configurations requires fine gridding to take into account the 3D nature of the geometries encountered in these applications that include geological structures with complicated shapes and exhibiting large variations in conductivities such as the seafloor bathymetry, the land topography, and targets with complex geometries and large contrasts in conductivities. Such problems significantly increase the computational cost of the conventional finite-difference (FD) approaches mainly due to the large condition numbers of the corresponding linear systems. To handle these problems, we employ a volume integral equation (IE) approach to arrive at an effective preconditioning operator for our FD solver. We refer to this new hybrid algorithm as the finite-difference integral equation method (FDIE). This FDIE preconditioning operator is divergence free and is based on a magnetic field formulation. Similar to the Lippman-Schwinger IE method, this scheme allows us to use a background elimination approach to reduce the computational domain, resulting in a smaller size stiffness matrix. Furthermore, it yields a linear system whose condition number is close to that of the conventional Lippman-Schwinger IE approach, significantly reducing the condition number of the stiffness matrix of the FD solver. Moreover, the FD framework allows us to substitute convolution operations by the inversion of banded matrices, which significantly reduces the computational cost per iteration of the hybrid method compared to the standard IE approaches. Also, well-established FD homogenization and optimal gridding algorithms make the FDIE more appropriate for the discretization of strongly inhomogeneous media. Some numerical studies are presented to illustrate the accuracy and effectiveness of the presented solver for CSEM, single-well, and crosswell EM applications.

A third-order finite-volume residual-based scheme for the 2D Euler equations on unstructured grids

A residual-based (RB) scheme relies on the vanishing of residual at the steady-state to design a transient first-order dissipation, which becomes high-order at steady-state. Initially designed within a finite-difference framework for computations of compressible flows on structured grids, the RB schemes displayed good convergence, accuracy and shock-capturing properties which motivated their extension to unstructured grids using a finite volume (FV) method. A second-order formulation of the FV–RB scheme for compressible flows on general unstructured grids was presented in a previous paper. The present paper describes the derivation of a third-order FV–RB scheme and its application to hyperbolic model problems as well as subsonic, transonic and supersonic internal and external inviscid flows.

A numerical study of saline contamination of a coastal aquifer

A three-dimensional mathematical model of variable density groundwater flow and miscible solute transport is developed to investigate the saltwater contamination of a coastal aquifer in west Japan. The selected area is an agricultural region located on a coastal aquifer where the groundwater is being consumed continuously as the main water source for drinking, greenhouses and industries. The salinity of the pumped water has been analysed periodically and it has been found that the threat of saltwater intrusion is emerging due to the continuous exploitation of groundwater. Owing to the distributed nature of pumping-well locations and non-uniform geological conditions, a three-dimensional model is needed to understand the saline contamination of the coastal aquifer. A numerical model was built by using the finite difference framework, taking into account the development of a transition zone and the variation of fluid density within it. The model simulates the saturated and unsaturated zone groundwater flow under the influence of a variable density situation in the coastal aquifer. Moreover, it can be used as a numerical tool to simulate the salinity contamination of the freshwater in this coastal aquifer.

Conservative high order semi-Lagrangian finite difference WENO methods for advection in incompressible flow

In this paper, we propose a semi-Lagrangian finite difference formulation for approximating conservative form of advection equations with general variable coefficients. Compared with the traditional semi-Lagrangian finite difference schemes [5,25], which approximate the advective form of the equation via direct characteristics tracing, the scheme proposed in this paper approximates the conservative form of the equation. This essential difference makes the proposed scheme naturally conservative for equations with general variable coefficients. The proposed conservative semi-Lagrangian finite difference framework is coupled with high order essentially non-oscillatory (ENO) or weighted ENO (WENO) reconstructions to achieve high order accuracy in smooth parts of the solution and to capture sharp interfaces without introducing spurious oscillations. The scheme is extended to high dimensional problems by Strang splitting. The performance of the proposed schemes is demonstrated by linear advection, rigid body rotation, swirling deformation, and two dimensional incompressible flow simulation in the vorticity stream-function formulation. As the information is propagating along characteristics, the proposed scheme does not have the CFL time step restriction of the Eulerian method, allowing for a more efficient numerical realization for many application problems.

Proof of concept of regional scale hydrologic simulations at hydrologic resolution utilizing massively parallel computer resources

We present the results of a unique, parallel scaling study using a 3‐D variably saturated flow problem including land surface processes that ranges from a single processor to a maximum number of 16,384 processors. In the applied finite difference framework and for a fixed problem size per processor, this results in a maximum number of approximately 8 × 109 grid cells (unknowns). Detailed timing information shows that the applied simulation platform ParFlow exhibits excellent parallel efficiency. This study demonstrates that regional scale hydrologic simulations on the order of 103 km2 are feasible at hydrologic resolution (∼100–101 m laterally, 10−2–10−1 m vertically) with reasonable computation times, which has been previously assumed to be an intractable computational problem.

Multiscale modeling of the influence of Fe content in a Al–Si–Cu alloy on the size distribution of intermetallic phases and micropores

A multiscale model was developed to simulate the formation of Fe-rich intermetallics and pores in quaternary Al–Si–Cu–Fe alloys. At the microscale, the multicomponent diffusion equations were solved for multiphase (liquid-solid-gas) materials via a finite difference framework to predict microstructure formation. A fast and robust decentered plate algorithm was developed to simulate the strong anisotropy of the solid/liquid interfacial energy for the Fe-rich intermetallic phase. The growth of porosity was controlled by local pressure drop due to solidification and interactions with surrounding solid phases, in addition to hydrogen diffusion. The microscale model was implemented as a subroutine in a commercial finite element package, producing a coupled multiscale model. This allows the influence of varying casting conditions on the Fe-rich intermetallics, the pores, and their interactions to be predicted. Synchrotron x-ray tomography experiments were performed to validate the model by comparing the three-dimensional morphology and size distribution of Fe-rich intermetallics as a function of Fe content. Large platelike Fe-rich β intermetallics were successfully simulated by the multiscale model and their influence on pore size distribution in shape castings was predicted as a function of casting conditions.

A tensor artificial viscosity using a finite element approach

We derive a tensor artificial viscosity suitable for use in a 2D or 3D unstructured arbitrary Lagrangian–Eulerian (ALE) hydrodynamics code. This work is similar in nature to that of Campbell and Shashkov [1]; however, our approach is based on a finite element discretization that is fundamentally different from the mimetic finite difference framework. The finite element point of view leads to novel insights as well as improved numerical results. We begin with a generalized tensor version of the Von Neumann–Richtmyer artificial viscosity, then convert it to a variational formulation and apply a Galerkin discretization process using high order Gaussian quadrature to obtain a generalized nodal force term and corresponding zonal heating (or shock entropy) term. This technique is modular and is therefore suitable for coupling to a traditional staggered grid discretization of the momentum and energy conservation laws; however, we motivate the use of such finite element approaches for discretizing each term in the Euler equations. We review the key properties that any artificial viscosity must possess and use these to formulate specific constraints on the total artificial viscosity force term as well as the artificial viscosity coefficient. We also show, that under certain simplifying assumptions, the two-dimensional scheme from [1] can be viewed as an under-integrated version of our finite element method. This equivalence holds on general distorted quadrilateral grids. Finally, we present computational results on some standard shock hydro test problems, as well as some more challenging problems, indicating the advantages of the new approach with respect to symmetry preservation for shock wave propagation over general grids.

Numerical modelling of an atmospheric pressure dielectric barrier discharge in nitrogen: electrical and kinetic description

An atmospheric pressure molecular nitrogen dielectric barrier discharge, under relatively high frequency (130 kHz) sinusoidal excitation, is investigated by means of a coupled electro-dynamic and kinetic model. The spatio temporal evolution of charged particles and electric field along the discharge axis is analyzed in a 1.5D finite difference framework which consists of the solution of 1D continuity and 2D Poisson equations, respectively. A 0D kinetic model accounts for the evolution of N2 neutral particles considering the principal triplet electronic states , , N2 (B 3 Πg) and N2 (C 3 Πu) as well as the singlet metastable state . It is shown that the evolution of the discharge takes place in three distinct phases, each one characterized by different current waveforms. The discharge's efficiency on active species production is examined, revealing the dominant character of metastable species and . Comparison of the numerical model's results with experimental ones is performed in the framework of electrical and optical diagnostics.

A new element-by-element method for trajectory calculations with tetrahedral finite element meshes

A new method to track massless particles in a three-dimensional flow field is presented. The method is based on an element-by-element approach coupled with a predictor–corrector shooting scheme and does not use any time step. By analogy with time-dependent schemes, the number of shootings is related to an equivalent number of time steps. The method has been implemented in a finite element framework using unstructured tetrahedral finite element meshes. However, it is general enough so that it can be implemented in finite difference and finite volume frameworks as well. It has been tested on a variety of flow systems namely: the Poiseuille flow in an empty circular pipe, the rotating flow in a stirred tank, the shear flow in a square tank and the flow through a static mixer. Accuracy has been found to depend on the accuracy of the velocity computation, the number of points per element and the level of mesh refinement. Copyright © 2006 John Wiley & Sons, Ltd.

Interfacial operator approach to computing modes of surface plasmon polaritons for periodic structures

It has been difficult to compute the band structures for photonic crystals with metallic components included in the periodic units. The existence of modes of surface plasmon polariton presents the major difficulty not only because of the localized nature of the modes but also of the apparent necessity of handling a nonlinear eigenvalue problem. Here we show that by introducing an interfacial operator within the finite-difference framework, we are able to formulate the problem for computing modes of surface plasmon polariton in the format of standard eigenvalue problems. Results are uncovered by applying the method to periodic structures with corrugated interfaces between metals and dielectric materials, as well as other classes of interfaces.

Mass-conservative finite volume methods on 2-D unstructured grids for the Richards’ equation

The solution to the 2-D time-dependent unsaturated flow equation is numerically approximated by a second-order accurate cell-centered finite-volume discretization on unstructured grids. The approximation method is based on a vertex-centered Least Squares linear reconstruction of the solution gradients at mesh edges. A Taylor series development in time of the water content dependent variable in a finite-difference framework guarantees that the proposed finite volume method is mass conservative. A Picard iterative scheme solves at each time step the resulting non-linear algebraic problem. The performance of the method is assessed on five different test cases and implementing four distinct soil constitutive relationships. The first test case deals with a column infiltration problem. It shows the capability of providing a mass-conservative behavior. The second test case verifies the numerical approximation by comparison with an analytical mixed saturated–unsaturated solution. In this case, the water drains from a fully saturated portion of a 1-D column. The third and fourth test cases illustrate the performance of the approximation scheme on sharp soil heterogeneities on 1-D and 2-D multi-layered infiltration problems. The 2-D case shows the passage of an abrupt infiltration front across a curved interface between two layers. Finally, the fifth test case compares the numerical results with an analytical solution that is developed for a 2-D heterogeneous soil with a source term representing plant roots. This last test case illustrates the formal second-order accuracy of the method in the numerical approximation of the pressure head.