NASep 29, 2011
A Correction Function Method for Poisson Problems with Interface Jump ConditionsAlexandre Noll Marques, Jean-Christophe Nave, Rodolfo Ruben Rosales
In this paper we present a method to treat interface jump conditions for constant coefficients Poisson problems that allows the use of standard "black box" solvers, without compromising accuracy. The basic idea of the new approach is similar to the Ghost Fluid Method (GFM). The GFM relies on corrections applied on nodes located across the interface for discretization stencils that straddle the interface. If the corrections are solution-independent, they can be moved to the right-hand-side (RHS) of the equations, producing a problem with the same linear system as if there were no jumps, only with a different RHS. However, achieving high accuracy is very hard (if not impossible) with the "standard" approaches used to compute the GFM correction terms. In this paper we generalize the GFM correction terms to a correction function, defined on a band around the interface. This function is then shown to be characterized as the solution to a PDE, with appropriate boundary conditions. This PDE can, in principle, be solved to any desired order of accuracy. As an example, we apply this new method to devise a 4th order accurate scheme for the constant coefficients Poisson equation with discontinuities in 2D. This scheme is based on (i) the standard 9-point stencil discretization of the Poisson equation, (ii) a representation of the correction function in terms of bicubics, and (iii) a solution of the correction function PDE by a least squares minimization. Several applications of the method are presented to illustrate its robustness dealing with a variety of interface geometries, its capability to capture sharp discontinuities, and its high convergence rate.
COMP-PHDec 16, 2016
High order solution of Poisson problems with piecewise constant coefficients and interface jumpsAlexandre Noll Marques, Jean-Christophe Nave, Rodolfo Ruben Rosales
We present a fast and accurate algorithm to solve Poisson problems in complex geometries, using regular Cartesian grids. We consider a variety of configurations, including Poisson problems with interfaces across which the solution is discontinuous (of the type arising in multi-fluid flows). The algorithm is based on a combination of the Correction Function Method (CFM) and Boundary Integral Methods (BIM). Interface and boundary conditions can be treated in a fast and accurate manner using boundary integral equations, and the associated BIM. Unfortunately, BIM can be costly when the solution is needed everywhere in a grid, e.g. fluid flow problems. We use the CFM to circumvent this issue. The solution from the BIM is used to rewrite the problem as a series of Poisson problems in rectangular domains - which requires the BIM solution at interfaces/boundaries only. These Poisson problems involve discontinuities at interfaces, of the type that the CFM can handle. Hence we use the CFM to solve them (to high order of accuracy) with finite differences and a Fast Fourier Transform based fast Poisson solver. We present 2-D examples of the algorithm applied to Poisson problems involving complex geometries, including cases in which the solution is discontinuous. We show that the algorithm produces solutions that converge with either 3rd or 4th order of accuracy, depending on the type of boundary condition and solution discontinuity.
NANov 18, 2011
Jet schemes for advection problemsBenjamin Seibold, Jean-Christophe Nave, Rodolfo Ruben Rosales
We present a systematic methodology to develop high order accurate numerical approaches for linear advection problems. These methods are based on evolving parts of the jet of the solution in time, and are thus called jet schemes. Through the tracking of characteristics and the use of suitable Hermite interpolations, high order is achieved in an optimally local fashion, i.e. the update for the data at any grid point uses information from a single grid cell only. We show that jet schemes can be interpreted as advect-and-project processes in function spaces, where the projection step minimizes a stability functional. Furthermore, this function space framework makes it possible to systematically inherit update rules for the higher derivatives from the ODE solver for the characteristics. Jet schemes of orders up to five are applied in numerical benchmark tests, and systematically compared with classical WENO finite difference schemes. It is observed that jet schemes tend to possess a higher accuracy than WENO schemes of the same order.
NAApr 2, 2012
A comparative study of the efficiency of jet schemesPrince Chidyagwai, Jean-Christophe Nave, Rodolfo Ruben Rosales et al.
We present two versions of third order accurate jet schemes, which achieve high order accuracy by tracking derivative information of the solution along characteristic curves. For a benchmark linear advection problem, the efficiency of jet schemes is compared with WENO and Discontinuous Galerkin methods of the same order. Moreover, the performance of various schemes in tracking solution contours is investigated. It is demonstrated that jet schemes possess the simplicity and speed of WENO schemes, while showing several of the advantages as well as the accuracy of DG methods.
NAFeb 11, 2014
A Sharp-Interface Active Penalty Method for the Incompressible Navier-Stokes EquationsDavid Shirokoff, Jean-Christophe Nave
The volume penalty method provides a simple, efficient approach for solving the incompressible Navier-Stokes equations in domains with boundaries or in the presence of moving objects. Despite the simplicity, the method is typically limited to first order spatial accuracy. We demonstrate that one may achieve high order accuracy by introducing an active penalty term. One key difference from other works is that we use a sharp, unregularized mask function. We discuss how to construct the active penalty term, and provide numerical examples, in dimensions one and two. We demonstrate second and third order convergence for the heat equation, and second order convergence for the Navier-Stokes equations. In addition, we show that modifying the penalty term does not significantly alter the time step restriction from that of the conventional penalty method.
NADec 20, 2016
A Correction Function Method for the Wave Equation with Interface Jump ConditionsDavid S. Abraham, Alexandre Noll Marques, Jean-Christophe Nave
In this paper a novel method to solve the constant coefficient wave equation, subject to interface jump conditions, is presented. In general, such problems pose issues for standard finite difference solvers, as the inherent discontinuity in the solution results in erroneous derivative information wherever the stencils straddle the given interface. Here, however, the recently proposed Correction Function Method (CFM) is used, in which correction terms are computed from the interface conditions, and added to affected nodes to compensate for the discontinuity. In contrast to existing methods, these corrections are not simply defined at affected nodes, but rather generalized to a continuous function within a small region surrounding the interface. As a result, the correction function may be defined in terms of its own governing partial differential equation (PDE) which may be solved, in principle, to arbitrary order of accuracy. The resulting scheme is not only arbitrarily high order, but also robust, having already seen application to Poisson problems and the heat equation. By extending the CFM to this new class of PDEs, the treatment of wave interface discontinuities in homogeneous media becomes possible. This allows, for example, for the straightforward treatment of infinitesimal source terms and sharp boundaries, free of staircasing errors. Additionally, new modifications to the CFM are derived, allowing compatibility with explicit multi-step methods, such as Runge-Kutta (RK4), without a reduction in accuracy. These results are then verified through numerous numerical experiments in one and two spatial dimensions.
NADec 7, 2016
Conservative methods for dynamical systemsAndy T. S. Wan, Alexander Bihlo, Jean-Christophe Nave
We show a novel systematic way to construct conservative finite difference schemes for quasilinear first-order system of ordinary differential equations with conserved quantities. In particular, this includes both autonomous and non-autonomous dynamical systems with conserved quantities of arbitrary forms, such as time-dependent conserved quantities. Sufficient conditions to construct conservative schemes of arbitrary order are derived using the multiplier method. General formulas for first-order conservative schemes are constructed using divided difference calculus. New conservative schemes are found for various dynamical systems such as Euler's equation of rigid body rotation, Lotka-Volterra systems, the planar restricted three-body problem and the damped harmonic oscillator.
NAJan 25, 2013
Convecting reference frames and invariant numerical modelsAlexander Bihlo, Jean-Christophe Nave
In the recent paper by Bernardini et al. [J. Comput. Phys. 232 (2013), 1-6] the discrepancy in the performance of finite difference and spectral models for simulations of flows with a preferential direction of propagation was studied. In a simplified investigation carried out using the viscous Burgers equation the authors attributed the poorer numerical results of finite difference models to a violation of Galilean invariance in the discretization and propose to carry out the computations in a reference frame moving with the bulk velocity of the flow. Here we further discuss this problem and relate it to known results on invariant discretization schemes. Non-invariant and invariant finite difference discretizations of Burgers equation are proposed and compared with the discretization using the remedy proposed by Bernardini et al..
NANov 14, 2015
A Fourier penalty method for solving the time-dependent Maxwell's equations in domains with curved boundariesRyan Galagusz, David Shirokoff, Jean-Christophe Nave
We present a high order, Fourier penalty method for the Maxwell's equations in the vicinity of perfect electric conductor boundary conditions. The approach relies on extending the smooth non-periodic domain of the equations to a periodic domain by removing the exact boundary conditions and introducing an analytic forcing term in the extended domain. The forcing, or penalty term is chosen to systematically enforce the boundary conditions to high order in the penalty parameter, which then allows for higher order numerical methods. We present an efficient numerical method for constructing the penalty term, and discretize the resulting equations using a Fourier spectral method. We demonstrate convergence orders of up to 3.5 for the one-dimensional Maxwell's equations, and show that the numerical method does not suffer from dispersion (or pollution) errors. We also illustrate the approach in two dimensions and demonstrate convergence orders of 2.5 for transverse magnetic modes and 1.5 for the transverse electric modes. We conclude the paper with numerous test cases in dimensions two and three including waves traveling in an irregular waveguide, and scattering off of a windmill-like geometry.
NAMar 20, 2019
Area-Preserving Geometric Hermite InterpolationGeoffrey McGregor, Jean-Christophe Nave
In this paper we establish a framework for planar geometric interpolation with exact area preservation using cubic Bézier polynomials. We show there exists a family of such curves which are $5^{th}$ order accurate, one order higher than standard geometric cubic Hermite interpolation. We prove this result is valid when the curvature at the endpoints does not vanish, and in the case of vanishing curvature, the interpolation is $4^{th}$ order accurate. The method is computationally efficient and prescribes the parametrization speed at endpoints through an explicit formula based on the given data. Additional accuracy (i.e. same order but lower error constant) may be obtained through an iterative process to find optimal parametrization speeds which further reduces the error while still preserving the prescribed area exactly.
NAMay 25, 2016
A fast-marching algorithm for non-monotonically evolving frontsAlexandra Tcheng, Jean-Christophe Nave
The non-monotonic propagation of fronts is considered. When the speed function $F:\mathbb{R}^{n} \times [0,T]\rightarrow \mathbb{R}$ is prescribed, the non-linear advection equation $ϕ_{t}+F|\nabla ϕ|=0$ is a Hamilton-Jacobi equation known as the level-set equation. It is argued that a small enough neighbourhood of the zero-level-set $\mathcal{M}$ of the solution $ϕ: \mathbb{R}^{n} \times [0,T] \rightarrow \mathbb{R}$ is the graph of $ψ:\mathbb{R}^{n} \rightarrow \mathbb{R}$ where $ψ$ solves a Dirichlet problem of the form $H(\vec{u},ψ(\vec{u}),\nabla ψ(\vec{u}))=0$. A fast-marching algorithm is presented where each point is computed using a discretization of such a Dirichlet problem, with no restrictions on the sign of $F$. The output is a directed graph whose vertices evenly sample $\mathcal{M}$. The convergence, consistency and stability of the scheme are addressed. Bounds on the computational complexity are estimated, and experimentally shown to be on par with the Fast Marching Method. Examples are presented where the algorithm is shown to be globally first-order accurate. The complexities and accuracies observed are independent of the monotonicity of the evolution.
NAApr 3, 2017
A Parametric Interpolation Framework for 1D Scalar Conservation Laws using the Equal Area PrincipleGeoffrey McGregor, Jean-Christophe Nave
In this paper we develop a novel framework for numerically solving scalar conservation laws in one space dimension. Utilizing the method of characteristics in conjunction with the equal area principle we develop an approach where the weak solution is obtained purely as the solution of a parametric interpolation problem. As this framework hinges on the validity of the equal area principle, we provide a rigorous discussion of the equal area principle and show that, indeed, the equal area principle is equivalent to the Rankine-Hugoniot condition, within the specific context studied in this paper. Combining these results with properties of the characteristic equations yields the desired setting to define the equivalent parametric interpolation problem. We conclude by applying this framework to Burgers' equation and show how one obtains machine precision in the shock position when the initial condition can be represented exactly in the chosen space of parametric polynomials.
NAOct 11, 2018
On the arbitrarily long-term stability of conservative methodsAndy T. S. Wan, Jean-Christophe Nave
We show the arbitrarily long-term stability of conservative methods for autonomous ODEs. Given a system of autonomous ODEs with conserved quantities, if the preimage of the conserved quantities possesses a bounded locally nite neighborhood, then the global error of any conservative method with the uniformly bounded displacement property is bounded for all time, when the uniform time step is taken suciently small. On nite precision machines, the global error still remains bounded and independent of time until some arbitrarily large time determined by machine precision and tolerance. The main result is proved using elementary topological properties for discretized conserved quantities which are equicontinuous. In particular, long-term stability is also shown using an averaging identity when the discretized conserved quantities do not explicitly depend on time steps. Numerical results are presented to illustrate the long-term stability result.
NAJul 8, 2015
The multiplier method to construct conservative finite difference schemes for ordinary and partial differential equationsAndy T. S. Wan, Alexander Bihlo, Jean-Christophe Nave
We present the multiplier method of constructing conservative finite difference schemes for ordinary and partial differential equations. Given a system of differential equations possessing conservation laws, our approach is based on discretizing conservation law multipliers and their associated density and flux functions. We show that the proposed discretization is consistent for any order of accuracy when the discrete multiplier has a multiplicative inverse. Moreover, we show that by construction, discrete densities can be exactly conserved. In particular, the multiplier method does not require the system to possess a Hamiltonian or variational structure. Examples, including dissipative problems, are given to illustrate the method. In the case when the inverse of the discrete multiplier becomes singular, consistency of the method is also established for scalar ODEs provided the discrete multiplier and density are zero-compatible.
NAMay 27, 2015
A low complexity algorithm for non-monotonically evolving frontsAlexandra Tcheng, Jean-Christophe Nave
A new algorithm is proposed to describe the propagation of fronts advected in the normal direction with prescribed speed function F. The assumptions on F are that it does not depend on the front itself, but can depend on space and time. Moreover, it can vanish and change sign. To solve this problem the Level-Set Method [Osher, Sethian; 1988] is widely used, and the Generalized Fast Marching Method [Carlini et al.; 2008] has recently been introduced. The novelty of our method is that its overall computational complexity is predicted to be comparable to that of the Fast Marching Method [Sethian; 1996], [Vladimirsky; 2006] in most instances. This latter algorithm is O(N^n log N^n) if the computational domain comprises N^n points. Our strategy is to use it in regions where the speed is bounded away from zero -- and switch to a different formalism when F is approximately 0. To this end, a collection of so-called sideways partial differential equations is introduced. Their solutions locally describe the evolving front and depend on both space and time. The well-posedness of those equations, as well as their geometric properties are addressed. We then propose a convergent and stable discretization of those PDEs. Those alternative representations are used to augment the standard Fast Marching Method. The resulting algorithm is presented together with a thorough discussion of its features. The accuracy of the scheme is tested when F depends on both space and time. Each example yields an O(1/N) global truncation error. We conclude with a discussion of the advantages and limitations of our method.
MATH-PHAug 1, 2013
Invariant Discretization Schemes Using Evolution-Projection TechniquesAlexander Bihlo, Jean-Christophe Nave
Finite difference discretization schemes preserving a subgroup of the maximal Lie invariance group of the one-dimensional linear heat equation are determined. These invariant schemes are constructed using the invariantization procedure for non-invariant schemes of the heat equation in computational coordinates. We propose a new methodology for handling moving discretization grids which are generally indispensable for invariant numerical schemes. The idea is to use the invariant grid equation, which determines the locations of the grid point at the next time level only for a single integration step and then to project the obtained solution to the regular grid using invariant interpolation schemes. This guarantees that the scheme is invariant and allows one to work on the simpler stationary grids. The discretization errors of the invariant schemes are established and their convergence rates are estimated. Numerical tests are carried out to shed some light on the numerical properties of invariant discretization schemes using the proposed evolution-projection strategy.
NANov 26, 2009
A gradient-augmented level set method with an optimally local, coherent advection schemeJean-Christophe Nave, Rodolfo Ruben Rosales, Benjamin Seibold
The level set approach represents surfaces implicitly, and advects them by evolving a level set function, which is numerically defined on an Eulerian grid. Here we present an approach that augments the level set function values by gradient information, and evolves both quantities in a fully coupled fashion. This maintains the coherence between function values and derivatives, while exploiting the extra information carried by the derivatives. The method is of comparable quality to WENO schemes, but with optimally local stencils (performing updates in time by using information from only a single adjacent grid cell). In addition, structures smaller than the grid size can be located and tracked, and the extra derivative information can be employed to obtain simple and accurate approximations to the curvature. We analyze the accuracy and the stability of the new scheme, and perform benchmark tests.