Jan 18 2017

math.NA arXiv:1701.04798v1

This paper investigates oscillation-free stability conditions of numerical methods for linear parabolic partial differential equations with some example extrapolations to nonlinear equations. Not clearly understood, numerical oscillations can create infeasible results. Since oscillation-free behavior is not ensured by stability conditions, a more precise condition would be useful for accurate solutions. Using Von Neumann and spectral analyses, we find and explore oscillation-free conditions for several finite difference schemes. Further relationships between oscillatory behavior and eigenvalues is supported with numerical evidence and proof. Also, evidence suggests that the oscillation-free stability condition for a consistent linearization may be sufficient to provide oscillation-free stability of the nonlinear solution. These conditions are verified numerically for several example problems by visually comparing the analytical conditions to the behavior of the numerical solution for a wide range of mesh sizes.

Problems of stabilization of the unstable cycle of one-dimensional complex dynamical system are briefly discussed. These questions reduced to the problem of description of the ranges of polynomials $q(z) = q_1z + q_2z^2 +\dots + q_nz^n$ defined in the unit disk and normalized by the conditions $q(1) = 1 $ and this is the main subject of the present paper.

Jan 18 2017

math.NA arXiv:1701.04750v1

The main goal of this paper is to investigate the order reduction phenomenon that appears in the integral deferred correction (InDC) methods based on implicit-explicit (IMEX) Runge-Kutta (R-K) schemes when applied to a class of stiff problems characterized by a small positive parameter $\varepsilon$, called singular perturbation problems (SPPs). In particular, an error analysis is presented for these implicit-explicit InDC (InDC-IMEX) methods when applied to SPPs. In our error estimate, we expand the global error in powers of $\varepsilon$ and show that its coefficients are global errors of the corresponding method applied to a sequence of differential algebraic systems. A study of these errors in the expansion yields error bounds and it reveals the phenomenon of order reduction. In our analysis we assume uniform quadrature nodes excluding the left-most point in the InDC method and the globally stiffly accurate property for the IMEX R-K scheme. Numerical results for the Van der Pol equation and PDE applications are presented to illustrate our theoretical findings.

Jan 18 2017

math.NA arXiv:1701.04719v1

We develop a cut finite element method for the Darcy problem on surfaces. The cut finite element method is based on embedding the surface in a three dimensional finite element mesh and using finite element spaces defined on the three dimensional mesh as trial and test functions. Since we consider a partial differential equation on a surface, the resulting discrete weak problem might be severely ill conditioned. We propose a full gradient and a normal gradient based stabilization computed on the background mesh to render the proposed formulation stable and well conditioned irrespective of the surface positioning within the mesh. Our formulation extends and simplifies the Masud-Hughes stabilized primal mixed formulation of the Darcy surface problem proposed in [28] on fitted triangulated surfaces. The tangential condition on the velocity and the pressure gradient is enforced only weakly, avoiding the need for any tangential projection. The presented numerical analysis accounts for different polynomial orders for the velocity, pressure, and geometry approximation which are corroborated by numerical experiments. In particular, we demonstrate both theoretically and through numerical results that the normal gradient stabilized variant results in a high order scheme.

Jan 18 2017

math.NA arXiv:1701.04685v1

In this paper we derive a discretisation of the equation of quasi-static elasticity in homogenization in form of a variational formulation and the so-called Lippmann-Schwinger equation, in anisotropic spaces of translates of periodic functions. We unify and extend the truncated Fourier series approach, the constant finite element ansatz and the anisotropic lattice derivation. The resulting formulation of the Lippmann-Schwinger equation in anisotropic translation invariant spaces unifies and analyses for the first time both the Fourier methods and finite element approaches in a common mathematical framework. We further define and characterize the resulting periodised Green operator. This operator coincides in case of a Dirichlet kernel corresponding to a diagonal matrix with the operator derived for the Galerkin projection stemming from the truncated Fourier series approach and to the anisotropic lattice derivation for all other Dirichlet kernels. Additionally, we proof the boundedness of the periodised Green operator. The operator further constitutes a projection if and only if the space of translates is generated by a Dirichlet kernel. Numerical examples for both the de la Vallée Poussin means and Box splines illustrate the flexibility of this framework.

It is often claimed that error cancellation plays an essential role in quantum chemistry and first-principle simulation for condensed matter physics and materials science. Indeed, while the energy of a large, or even medium-size, molecular system cannot be estimated numerically within chemical accuracy (typically 1 kcal/mol or 1 mHa), it is considered that the energy difference between two configurations of the same system can be computed in practice within the desired accuracy. The purpose of this paper is to provide a quantitative study of discretization error cancellation. The latter is the error component due to the fact that the model used in the calculation (e.g. Kohn-Sham LDA) must be discretized in a finite basis set to be solved by a computer. We first report comprehensive numerical simulations performed with Abinit on two simple chemical systems, the hydrogen molecule on the one hand, and a system consisting of two oxygen atoms and four hydrogen atoms on the other hand. We observe that errors on energy differences are indeed significantly smaller than errors on energies, but that these two quantities asymptotically converge at the same rate when the energy cut-off goes to infinity. We then analyze a simple one-dimensional periodic Schrödinger equation with Dirac potentials, for which analytic solutions are available. This allows us to explain the discretization error cancellation phenomenon on this test case with quantitative mathematical arguments.

Saint-Venant equations can be generalized to account for a viscoelastic rheology in shallow flows. A Finite-Volume discretization for the 1D Saint-Venant system generalized to Upper-Convected Maxwell (UCM) fluids was proposed in [Bouchut \& Boyaval, 2013], which preserved a physically-natural stability property (i.e. free-energy dissipation) of the full system. It invoked a relaxation scheme of Suliciu type for the numerical computation of approximate solution to Riemann problems. Here, the approach is extended to the 1D Saint-Venant system generalized to the finitely-extensible nonlinear elastic fluids of Peterlin (FENE-P). We are currently not able to ensure all stability conditions a priori, but numerical simulations went smoothly in a practically useful range of parameters.

Jan 18 2017

math.NA arXiv:1701.04564v1

In the setting of continuum elasticity, phase transformations involving martensitic variants are modeled by a free energy density function that is non-convex in strain space. Here, we adopt an existing mathematical model in which we regularize the non-convex free energy density function by higher-order gradient terms at finite strain and derive boundary value problems via the standard variational argument applied to the corresponding total free energy, inspired by Toupin's theory of gradient elasticity. These gradient terms are to preclude existence of arbitrarily fine microstructures, while still allowing for existence of multiple solution branches corresponding to local minima of the total free energy; these are classified as metastable solution branches. The goal of this work is to solve the boundary value problem numerically in three dimensions, observe solution branches, and assess stability of each branch by numerically evaluating the second variation of the total free energy. We also study how these microstructures evolve as the length-scale parameter, the coefficient of the strain gradient terms in the free energy, approaches zero.

Jan 18 2017

math.NA arXiv:1701.04558v1

In this study, the numerical solutions of reaction-diffusion systems are investigated via the trigonometric quintic B-spline nite element collocation method. These equations appear in various disciplines in order to describe certain physical facts, such as pattern formation, autocatalytic chemical reactions and population dynamics. The Schnakenberg, Gray-Scott and Brusselator models are special cases of reaction-diffusion systems considered as numerical examples in this paper. For numerical purposes, Crank-Nicolson formulae are used for the time discretization and the resulting system is linearized by Taylor expansion. In the finite element method, a uniform partition of the solution domain is constructed for the space discretization. Over the mentioned mesh, dirac-delta function and trigonometric quintic B-spline functions are chosen as the weighted function and the bases functions, respectively. Thus, the reaction-diffusion system turns into an algebraic system which can be represented by a matrix equation so that the coeffcients are block matrices containing a certain number of non-zero elements in each row. The method is tested on different problems. To illustrate the accuracy, error norms are calculated in the linear problem whereas the relative error is given in other nonlinear problems. Subject to the character of the nonlinear problems, the occurring spatial patterns are formed by the trajectories of the dependent variables. The degree of the base polynomial allows the method to be used in high-order differential equation solutions. The algorithm produces accurate results even when the time increment is larger. Therefore, the proposed Trigonometric Quintic B-spline Collocation method is an effective method which produces acceptable results for the solutions of reaction-diffusion systems.

Jan 18 2017

math.NA arXiv:1701.04553v1

We prove in this paper the convergence of the Marker and cell (MAC) scheme for the dis-cretization of the steady-state and unsteady-state incompressible Navier-Stokes equations in primitive variables on non-uniform Cartesian grids, without any regularity assumption on the solution. A priori estimates on solutions to the scheme are proven ; they yield the existence of discrete solutions and the compactness of sequences of solutions obtained with family of meshes the space step of which tends to zero. We then establish that the limit is a weak solution to the continuous problem.

Jan 18 2017

math.NA arXiv:1701.04511v1

We consider the numerical approximations of a two-phase hydrodynamics coupled phase-field model that incorporates the variable densities, viscosities and moving contact line boundary conditions. The model is a nonlinear, coupled system that consists of incompressible Navier--Stokes equations with the generalized Navier boundary condition, and the Cahn--Hilliard equations with moving contact line boundary conditions. By some subtle explicit--implicit treatments to nonlinear terms, we develop two efficient, unconditionally energy stable numerical schemes, in particular, a linear decoupled energy stable scheme for the system with static contact line condition, and a nonlinear energy stable scheme for the system with dynamic contact line condition. An efficient spectral-Galerkin spatial discretization is implemented to verify the accuracy and efficiency of proposed schemes. Various numerical results show that the proposed schemes are efficient and accurate.

Jan 18 2017

math.NA arXiv:1701.04492v1

By viewing the nonuniform discrete Fourier transform (NUDFT) as a perturbed version of a uniform discrete Fourier transform, we propose a fast, stable, and simple algorithm for computing the NUDFT that costs $\mathcal{O}(N\log N\log(1/\epsilon)/\log\!\log(1/\epsilon))$ operations based on the fast Fourier transform, where $N$ is the size of the transform and $0<\epsilon <1$ is a working precision. Our key observation is that a NUDFT and DFT matrix divided entry-by-entry is often well-approximated by a low rank matrix, allowing us to express a NUDFT matrix as a sum of diagonally-scaled DFT matrices. Our algorithm is simple to implement, automatically adapts to any working precision, and is competitive with state-of-the-art algorithms. In the fully uniform case, our algorithm is essentially the FFT. We also describe quasi-optimal algorithms for the inverse NUDFT and two-dimensional NUDFTs.

Jan 18 2017

math.NA arXiv:1701.04476v1

In this paper, we propose a novel approach for coupling 2D/1D shallow water flow models. Efficiently coupling these models is vital for simulating the flow and flooding of open channels. Currently, existing methods couple the models either at the channel lateral boundaries (lateral methods) or at the location, along the channel flow direction, where the two sub-domains intersect (frontal methods). We classify these methods as horizontal methods since the coupling points are on the horizontal plane. The limitations of these methods include their inability to recover the 2D channel flow structure during flooding without losing efficiency, and their inability to switch between lateral and frontal methods, based on the problem. Here, we propose a new paradigm which is to think of the channel flow as a two layer flow. This leads to the vertical coupling method (VCM) which is able to recover 2D channel flow structure without losing efficiency, switch types depending on the flow and is a superset of some existing methods. The VCM is based on a user chosen elevation above which the user considers the channel to be full. From this elevation, all other flow quantities are defined including the two layers. The flows in the lower and upper layers are assumed to be 1D and 2D respectively, and the appropriate flow models with exchange terms are derived. The 2D shallow water flow model is retained for the floodplains. Finite volume methods (FVM) are formulated for the flood model; the FVM with the operator splitting approach are also formulated to solve and couple the two layer channel flow models. We show that the resulting method (i) is well-balanced (ii) preserves the "no-numerical" flooding property (iii) preserves conservation properties and (iv) adapts to the flow situation.