Numerical solution of a phase field model for incompressible two-phase flows based on artificial compressibility

We present an artificial compressibility based numerical method for a phase field model for simulating two-phase incompressible viscous flows. The phase model was proposed by Liu and Shen [Physica D. 179 (2003) 211–228], in which the interface between two fluids is represented by a thin transition region of fluid mixture that stores certain amount of mixing energy. The model consists of the Navier–Stokes equations coupled with the Allen–Cahn equation (phase field equation) through an extra stress term and a transport term. The extra stress in the momentum equations represents the phase-induced capillary effect for the mixture due to the surface tension. The coupled equations are cast into a conservative form suitable for implementation with the artificial compressibility method. The resulting hyperbolic system of equations are then discretized with weighted essentially non-oscillatory (WENO) finite difference scheme. The dual-time stepping technique is applied for obtaining time accuracy at each physical time step, and the approximate factorization algorithm is used to solve the discretized equations. The effectiveness of the numerical method is demonstrated in several two-phase flow problems with topological changes. Numerical results show the present method can be used to simulate incompressible two-phase flows with small interfacial width parameters and topological changes.     

Adaptive WENO Methods Based on Radial Basis Function Reconstruction

We explore the use of radial basis functions (RBF) in the weighted essentially non-oscillatory (WENO) reconstruction process used to solve hyperbolic conservation laws, resulting in a numerical method of arbitrarily high order to solve problems with discontinuous solutions. Thanks to the mesh-less property of the RBFs, the method is suitable for non-uniform grids and mesh adaptation. We focus on multiquadric radial basis functions and propose a simple strategy to choose the shape parameter to control the balance between achievable accuracy and the numerical stability. We also develop an original smoothness indicator which is independent of the RBF for the WENO reconstruction step. Moreover, we introduce type I and type II RBF-WENO methods by computing specific linear weights. The RBF-WENO method is used to solve linear and nonlinear problems for both scalar and systems of conservation laws, including Burgers equation, the Buckley–Leverett equation, and the Euler equations. Numerical results confirm the performance of the proposed method. We finally consider an effective conservative adaptive algorithm that captures moving shocks and rapidly varying solutions well. Numerical results on moving grids are presented for both Burgers equation and the more complex Euler equations.     

Adaptive finite differences and IMEX time-stepping to price options under Bates model

In this paper, we consider numerical pricing of European and American options under the Bates model, a model which gives rise to a partial-integro differential equation. This equation is discretized in space using adaptive finite differences while an IMEX scheme is employed in time. The sparse linear systems of equations in each time-step are solved using an LU-decomposition and an operator splitting technique is employed for the linear complementarity problems arising for American options. The integral part of the equation is treated explicitly in time which means that we have to perform matrix-vector multiplications each time-step with a matrix with dense blocks. These multiplications are accomplished through fast Fourier transforms. The great performance of the method is demonstrated through numerical experiments.     

Propagation of non-linear waves in a non-ideal relaxing gas

We have derived an evolution equation governing the far-field behaviour of small amplitude waves in a non-ideal relaxing gas for planar and converging flow. Asymptotic expansions of the flow variables for small amplitude waves have been used to derive the evolution equation. This equation turns out to be a generalized Burger's equation. The numerical solution of this equation is obtained by using the homotopy analysis method (HAM) proposed by Liao with two different initial conditions. Using the HAM, we have studied the effect of relaxation and nonlinearity. The convergence control parameter enables us to find a good approximate solution for such a complex flow problem. This method also confirms the capabilities and usefulness of convergence control parameter and HAM for complex and highly non-linear problems.     

A coarse-grid projection method for accelerating incompressible MHD flow simulations

Coarse grid projection (CGP) is a multiresolution technique for accelerating numerical calculations associated with a set of nonlinear evolutionary equations along with stiff Poisson’s equations. In this article, we use CGP for the first time to speed up incompressible magnetohydrodynamics (MHD) flow simulations. Accordingly, we solve the nonlinear advection–diffusion equation on a fine mesh, while we execute the electric potential Poisson equation on the corresponding coarsened mesh. Mapping operators connect two grids together. A pressure correction scheme is used to enforce the incompressibility constrain. The study of incompressible flow past a circular cylinder in the presence of Lorentz force is selected as a benchmark problem with a fixed Reynolds number but various Stuart numbers. We consider two different situations. First, we only apply CGP to the electric potential Poisson equation. Second, we apply CGP to the pressure Poisson equation as well. The maximum speed-up factors achieved here are approximately 3 and 23, respectively, for the first and second situations. For the both situations, we examine the accuracy of velocity and vorticity fields as well as the lift and drag coefficients. In general, the results obtained by CGP are in an excellent to reasonable range of accuracy. The CGP results are significantly more accurate compared to the numerical simulations of the advection–diffusion and electric potential Poisson equations on pure coarse scale grids.

     

Numerical simulation of contaminant biodegradation by higher order methods and adaptive time stepping

In this work we present and analyze a reliable and robust approximation scheme for biochemically reacting transport in the subsurface following Monod type kinetics. Water flow is modeled by the Richards equation. The proposed scheme is based on higher order finite element methods for the spatial discretization and the two step backward differentiation formula for the temporal one. The resulting nonlinear algebraic systems of equations are solved by a damped version of Newtons method. For the linear problems of the Newton iteration Krylov space methods are used. In computational experiments conducted for realistic subsurface (groundwater) contamination scenarios we show that the higher order approximation scheme significantly reduces the amount of inherent numerical diffusion compared to lower order ones. Thereby an artificial transverse mixing of the species leading to a strong overestimation of the biodegradation process is avoided. Finally, we present a robust adaptive time stepping technique for the coupled flow and transport problem which allows efficient long-term predictions of biodegradation processes.     

A level-set method for steady-state and transient natural convection problems

This paper introduces a topology optimization method for 2D and 3D, steady-state and transient heat transfer problems that are dominated by natural convection in the fluid phase and diffusion in the solid phase. The geometry of the fluid-solid interface is described by an explicit level set method which allows for both shape and topological changes in the optimization process. The heat transfer in the fluid is modeled by an advection-diffusion equation. The fluid velocity is described by the incompressible Navier-Stokes equations augmented by a Boussinesq approximation of the buoyancy forces. The temperature field in the solid is predicted by a linear diffusion model. The governing equations in both the fluid and solid phases are discretized in space by a generalized formulation of the extended finite element method which preserves the crisp geometry definition of the level set method. The interface conditions at the fluid-solid boundary are enforced by Nitsche’s method. The proposed method is studied for problems optimizing the geometry of cooling devices. The numerical results demonstrate the applicability of the proposed method for a wide spectrum of problems. As the flow may exhibit dynamic instabilities, transient phenomena need to be considered when optimizing the geometry. However, the computational burden increases significantly when the time evolution of the flow fields needs to be resolved.     

Numerical simulation of the transient flow behaviour in chemical reactors using a penalisation method

We present a numerical scheme to study the transient flow behaviour in complex geometries. The Navier–Stokes equations are discretized with a high-resolution Fourier pseudo-spectral discretization with adaptive time-stepping. Using a penalisation technique solid boundaries of arbitrary shape can be easily taken into account. As application we present different simulations of unsteady flows, typically encountered in chemical reactors. We study transitional flows in tube bundles (arrays of cylinders and squares) for different Reynolds numbers and angles of incidence and a channel flow with an obstacle.     

Fast and Stable Polynomial Equation Solving and Its Application to Computer Vision

This paper presents several new results on techniques for solving systems of polynomial equations in computer vision. Gröbner basis techniques for equation solving have been applied successfully to several geometric computer vision problems. However, in many cases these methods are plagued by numerical problems. In this paper we derive a generalization of the Gröbner basis method for polynomial equation solving, which improves overall numerical stability. We show how the action matrix can be computed in the general setting of an arbitrary linear basis for ?[x]/I. In particular, two improvements on the stability of the computations are made by studying how the linear basis for ?[x]/I should be selected. The first of these strategies utilizes QR factorization with column pivoting and the second is based on singular value decomposition (SVD). Moreover, it is shown how to improve stability further by an adaptive scheme for truncation of the Gröbner basis. These new techniques are studied on some of the latest reported uses of Gröbner basis methods in computer vision and we demonstrate dramatically improved numerical stability making it possible to solve a larger class of problems than previously possible.     

Finite element analysis of particle motion in steady inspiratory airflow

To accurately model the inhaled particle motion, equations governing particle trajectories in carrier flow are solved together with the Navier–Stokes equations. Under the relatively dilute particle condition in the mixture, equations for two phases are coupled through the interface drag shown in the solid-phase momentum equations. The present study investigates bifurcation flow in the human central airway using the finite element method. In the gas phase, we employ the biquadratic streamline upwind Petrov–Galerkin finite element model to simulate the incompressible air flow. To solve the equations of motion for the inhaled particles, we apply another biquadratic streamline upwind finite element model. A feature common to two models applied to each phase of equations is that both of them provide nodally exact solutions to the convection–diffusion and the convection–reaction equations, which are prototype equations for the gas-phase and the solid-phase equations, respectively. In two dimensions, both models have ability to introduce physically meaningful artificial damping terms solely in the streamline direction. With these terms added to the formulation, the discrete system is enhanced without compromising the numerical diffusion error. Tests on inspiratory problem were conducted, and the results are presented, with an emphasis on the discussion of particle motion.     

Mixed Volume Element-Characteristic Fractional Step Difference Method for Contamination from Nuclear Waste Disposal

A nonlinear system with boundary-initial value conditions of convection–diffusion partial differential equations is presented to describe incompressible nuclear waste disposal contamination in porous media. The flow pressure is determined by an elliptic equation, the concentrations of brine and radionuclide are formulated by convection–diffusion equations, and the transport of temperature is defined by a heat equation. The pressure appears in convection–diffusion equations and heat equation in the form of Darcy velocity and controls the physical processes. The fluid pressure and velocity are solved by the conservative mixed volume element and the computation accuracy of Darcy velocity is improved one order. A combination method of the mixed volume element and the approximation of characteristics is applied to solve the brine and heat, where the diffusion is discretized by a mixed volume element method and the convection is treated by the method of characteristics. The characteristics can confirm strong computation stability at sharp fronts and it can avoid numerical dispersion and nonphysical oscillation. Larger time-steps along the characteristics are shown to result in smaller time-truncation errors than those resulting from standard methods. The mixed volume element method has the property of conservation on each element and it can obtain numerical solutions of the brine and adjoint vectors. The radionuclide is solved by a coupled method of characteristics and fractional step difference. The computational work is reduced greatly by decomposing a three-dimensional problem into three successive one-dimensional problems and using the algorithm of speedup. Using numerical analysis of priori estimates of differential equations, we demonstrate an optimal second order estimate in \(l^2\) norm. Numerical data are appropriate with the scheme and it is shown that the method is a powerful tool to solve the well-known problems in porous media.     

A Spectral Stochastic Semi-Lagrangian Method for Convection-Diffusion Equations with Uncertainty

Real life convection-diffusion problems are characterized by their inherent or externally induced uncertainties in the design parameters. This paper presents a spectral stochastic finite element semi-Lagrangian method for numerical solution of convection-diffusion equations with uncertainty. Using the spectral decomposition, the stochastic variational problem is reformulated to a set of deterministic variational problems to be solved for each Wiener polynomial chaos. To obtain the chaos coefficients in the corresponding deterministic convection-diffusion equations, we implement a semi-Lagrangian method in the finite element framework. Once this representation is computed, statistics of the numerical solution can be easily evaluated. These numerical techniques associate the geometrical flexibility of the finite element method with the ability offered by the semi-Lagrangian method to solve convection-dominated problems using time steps larger than its Eulerian counterpart. Numerical results are shown for a convection-diffusion problem driven with stochastic velocity and for an incompressible viscous flow problem with a random force. In both examples, the proposed method demonstrates its ability to better maintain the shape of the solution in the presence of uncertainties and steep gradients.     

Solving inverse problems involving the Navier–Stokes equations discretized by a Lagrange–Galerkin method

In this article, we are investigating the numerical approximation of an inverse problem involving the evolution of a Newtonian viscous incompressible fluid described by the Navier–Stokes equations in 2D. This system is discretized using a low order finite element in space coupled with a Lagrange–Galerkin scheme for the nonlinear advection operator. We introduce a full discrete linearized scheme that is used to compute the gradient of a given cost function by ensuring its consistency. Using gradient based optimization algorithms, we are able to deal with two fluid flow inverse problems, the drag reduction around a moving cylinder and the identification of a far-field velocity using the knowledge of the fluid load on a rectangular bluff body, for both fixed and prescribed moving configurations.     

Exact Boundary Condition for Semi-discretized Schrödinger Equation and Heat Equation in a Rectangular Domain

An approach to solve finite time horizon suboptimal feedback control problems for partial differential equations is proposed by solving dynamic programming equations on adaptive sparse grids. A semi-discrete optimal control problem is introduced and the feedback control is derived from the corresponding value function. The value function can be characterized as the solution of an evolutionary Hamilton–Jacobi Bellman (HJB) equation which is defined over a state space whose dimension is equal to the dimension of the underlying semi-discrete system. Besides a low dimensional semi-discretization it is important to solve the HJB equation efficiently to address the curse of dimensionality. We propose to apply a semi-Lagrangian scheme using spatially adaptive sparse grids. Sparse grids allow the discretization of the value functions in (higher) space dimensions since the curse of dimensionality of full grid methods arises to a much smaller extent. For additional efficiency an adaptive grid refinement procedure is explored. The approach is illustrated for the wave equation and an extension to equations of Schrödinger type is indicated. We present several numerical examples studying the effect the parameters characterizing the sparse grid have on the accuracy of the value function and the optimal trajectory.     

Nichtlineare Stabilität und Phasenuntersuchung adaptiver Nyström-Runge-Kutta-Methoden

The stability of adaptive Nyström-Runge-Kutta procedures is studied for a wide class of nonlinear stiff systems of second order differential equations. We show that for a large class of semi-discrete hyperbolic and parabolic problems the restriction of the stepsize is not due to the stiffness of the differential equation. Furthermore we use the scalar test equation \(y'' = - \omega ^2 y + q \cdot e^{iv(t - t_0 )} \) to derive conditions which ensure that the numerical forced oscillation is in phase with the analytical forced oscillation. The order of adaptive Nyström-Runge-Kutta methods (with a stability-matrix based on a diagonal Padéapproximation) for which the forced oscillation is in phase with its analytical counterpart cannot be greater than two. This barrier of order is not true forr-stage implicit Nyström methods of orderp=2r.     

Refinement of flexible space–time finite element meshes and discontinuous Galerkin methods

In this paper we present an algorithm to refine space–time finite element meshes as needed for the numerical solution of parabolic initial boundary value problems. The approach is based on a decomposition of the space–time cylinder into finite elements, which also allows a rather general and flexible discretization in time. This also includes adaptive finite element meshes which move in time. For the handling of three-dimensional spatial domains, and therefore of a four-dimensional space–time cylinder, we describe a refinement strategy to decompose pentatopes into smaller ones. For the discretization of the initial boundary value problem we use an interior penalty Galerkin approach in space, and an upwind technique in time. A numerical example for the transient heat equation confirms the order of convergence as expected from the theory. First numerical results for the transient Navier–Stokes equations and for an adaptive mesh moving in time underline the applicability and flexibility of the presented approach.     

Two-Dimensional Force-Free Relaxation and Superdiffusion Governed by Fractional Kinetic Equations

Starting from the integral Chapman-Kolmogorov equation and the Langevin equations with a multivariate Levy noise, we consider possible generalization of a one-dimensional fractional kinetic equations. It is assumed that the Levy noise has independent components. As the result, we get fractional Fokker-Planck equation for the probability density function in a phase space and fractional Einstein-Smoluchowski equation for the probability density function in a real space. We consider force-free superdiffusion and relaxation described by fractional equations and carry out numerical simulation of these processes by solving the corresponding Langevin equations. Analytical and numerical results agree quantitatively.     

A Galerkin Method for the Simulation of the Transient 2-D/2-D and 3-D/3-D Linear Boltzmann Equation

Many production steps used in the manufacturing of integrated circuits involve the deposition of material from the gas phase onto wafers. Models for these processes should account for gaseous transport in a range of flow regimes, from continuum flow to free molecular or Knudsen flow, and for chemical reactions at the wafer surface. We develop a kinetic transport and reaction model whose mathematical representation is a system of transient linear Boltzmann equations. In addition to time, a deterministic numerical solution of this system of kinetic equations requires the discretization of both position and velocity spaces, each two-dimensional for 2-D/2-D or each three-dimensional for 3-D/3-D simulations. Discretizing the velocity space by a spectral Galerkin method approximates each Boltzmann equation by a system of transient linear hyperbolic conservation laws. The classical choice of basis functions based on Hermite polynomials leads to dense coefficient matrices in this system. We use a collocation basis instead that directly yields diagonal coefficient matrices, allowing for more convenient simulations in higher dimensions. The systems of conservation laws are solved using the discontinuous Galerkin finite element method. First, we simulate chemical vapor deposition in both two and three dimensions in typical micron scale features as application example. Second, stability and convergence of the numerical method are demonstrated numerically in two and three dimensions. Third, we present parallel performance results which indicate that the implementation of the method possesses very good scalability on a distributed-memory cluster with a high-performance Myrinet interconnect.     

An adaptive phase field method for the mixture of two incompressible fluids

This paper develops an adaptive moving mesh method to solve a phase field model for the mixture of two incompressible fluids. The projection method is implemented on a half-staggered, moving quadrilateral mesh to keep the velocity field divergence-free, and the conjugate gradient or multigrid method is employed to solve the discrete Poisson equations. The current algorithm is composed by two independent parts: evolution of the governing equations and mesh-redistribution. In the first part, the incompressible Navier-Stokes equations are solved on a fixed half-staggered mesh by the rotational incremental pressure-correction scheme, and the Allen-Cahn type of phase equation is approximated by a conservative, second-order accurate central difference scheme, where the Lagrangian multiplier is used to preserve the mass-conservation of the phase field. The second part is an iteration procedure. During the mesh redistribution, the phase field is remapped onto the newly resulting meshes by the high-resolution conservative interpolation, while the non-conservative interpolation algorithm is applied to the velocity field. The projection technique is used to obtain a divergence-free velocity field at the end of this part. The resultant numerical scheme is stable, mass conservative, highly efficient and fast, and capable of handling variable density and viscosity. Several numerical experiments are presented to demonstrate the efficiency and robustness of the proposed algorithm.