Symmetric Finite Difference Research Articles

A regularized representation of the fractional Laplacian in n dimensions and its relation to Weierstrass-Mandelbrot-type fractal functions

We demonstrate that the fractional Laplacian (FL) is the principal characteristic operator of harmonic systems with {\it self-similar} interparticle interactions. We show that the FL represents the "{\it fractional continuum limit}" of a discrete "self-similar Laplacian" which is obtained by Hamilton's variational principle from a discrete spring model. We deduce from generalized self-similar elastic potentials regular representations for the FL which involve convolutions of symmetric finite difference operators of even orders extending the standard representation of the FL. Further we deduce a regularized representation for the FL $-(-\Delta)^{\frac{\alpha}{2}}$ holding for $\alpha\in \R \geq 0$. We give an explicit proof that the regularized representation of the FL gives for integer powers $\frac{\alpha}{2} \in \N\_0$ a distributional representation of the standard Laplacian operator $\Delta$ including the trivial unity operator for $\alpha\rightarrow 0$. We demonstrate that self-similar {\it harmonic} systems are {\it all} governed in a distributional sense by this {\it regularized representation of the FL} which therefore can be conceived as characteristic footprint of self-similarity.

A New Hybrid Stochastic Approximation Algorithm

We introduce Secant-Tangents AveRaged (STAR) Stochastic Approximation (SA), a new SA algorithm that estimates the gradient using a hybrid estimator, which is a convex combination of a symmetric finite difference and an average of two direct gradient estimators. For the deterministic weight sequence that minimizes the variance of the STAR gradient, we prove that for quadratic functions, the mean squared error (MSE) of the STAR-SA algorithm using this weight sequence is strictly less than that of the classical SA methods of Robbins-Monro (RM) and Kiefer-Wolfowitz (KW). We also prove convergence of the STAR-SA algorithm for general concave functions. Furthermore, we illustrate its effectiveness through numerical experiments by comparing the MSE of the STAR-SA algorithm against RM and KW for simple quadratic functions with various steepness and noise levels.

Numerical Investigation of the Steady State of a Driven Thin Film Equation

A third-order ordinary differential equation with application in the flow of a thin liquid film is considered. The boundary conditions come from Tanner's problem for the surface tension driven flow of a thin film. Symmetric and nonsymmetric finite difference schemes are implemented in order to obtain steady state solutions. We show that a central difference approximation to the third derivative in the model equation produces a solution curve with oscillations. A difference scheme based on a combination of forward and backward differences produces a smooth accurate solution curve. The stability of these schemes is analysed through the use of a von Neumann stability analysis.

An interpolation matched interface and boundary method for elliptic interface problems

An interpolation matched interface and boundary (IMIB) method with second-order accuracy is developed for elliptic interface problems on Cartesian grids, based on original MIB method proposed by Zhou et al. [Y. Zhou, G. Wei, On the fictious-domain and interpolation formulations of the matched interface and boundary method, J. Comput. Phys. 219 (2006) 228–246]. Explicit and symmetric finite difference formulas at irregular grid points are derived by virtue of the level set function. The difference scheme using IMIB method is shown to satisfy the discrete maximum principle for a certain class of problems. Rigorous error analyses are given for the IMIB method applied to one-dimensional (1D) problems with piecewise constant coefficients and two-dimensional (2D) problems with singular sources. Comparison functions are constructed to obtain a sharp error bound for 1D approximate solutions. Furthermore, we compare the ghost fluid method (GFM), immersed interface method (IIM), MIB and IMIB methods for 1D problems. Finally, numerical examples are provided to show the efficiency and robustness of the proposed method.

OPTIMIZED ARTIFICIAL DISSIPATION TERMS FOR ENHANCED STABILITY LIMITS

Finite difference approximations for the convection equation are developed, which exhibit enhanced stability limits for explicit Runge–Kutta integration. Stability limits are increased by adding artificial dissipation terms, which are optimized to yield greatest stable time steps. For the artificial dissipation terms, symmetric finite difference approximations of even-order derivatives are used with differencing stencils equal to the convective stencils. The spatial discretization inclusive of the added dissipation term is shown to be consistent with a first derivative. The formal order of accuracy in space is decreased by one order, while the order of time integration is not affected. As a result, the time step limits of originally stable Runge–Kutta integration is increased, for some combinations of spatial discretization and time integration by a factor of two. Algorithms, which are unstable without damping are stabilized. The dispersion properties of the algorithms are not influenced by the proposed damping terms. Spectral analysis of the algorithms show very low dissipation error for dimensionless wave numbers k Δ x < 0.5. Stability conditions based on von Neumann stability analysis are given for the proposed schemes for explicit Runge–Kutta time integration of orders up to ten.

Stability of a symmetric finite-difference scheme with approximate transparent boundary conditions for the time-dependent Schrödinger equation

Stability of a symmetric finite-difference scheme with approximate transparent boundary conditions for the time-dependent Schrödinger equation

Numerical analysis of the ab initio computation of the effects of ionization on the nonlinear susceptibility coefficients of the hydrogen atom

This paper provides a numerical analysis of a procedure for determining the effects of ionization on the nonlinear susceptibility coefficients of the hydrogen atom. To solve the relevant system of Schrodinger-type equations, we have developed a multidomain pseudospectral code with high accuracy symmetric finite differences to update cell boundary points. Using a conservative time stepping, one is then able to resolve the oscillatory solutions to the underlying equations, compute the ionization probability, and accurately determine the polarization. To gain insight into the physical mechanisms involved, we have calculated the full susceptibility and one which depends solely on the electronic bound-bound transitions. Our analysis reveals that saturation of the susceptibility occurs without including bound-continuum transitions. We have also found that a linear extrapolation of the total susceptibility versus the ionization probability to obtain the instantaneous susceptibility is quantitatively unreliable.

Fourth Order Symmetric Finite Difference Schemes for the Acoustic Wave Equation

We present an explicit, symmetric finite difference scheme for the acoustic wave equation on a rectangle with Neumann and/or Dirichlet boundary conditions. The scheme is fourth order accurate both in time and space. It is obtained by mass lumping of a finite element scheme. The accuracy and the difference approximations at the boundary are analyzed in terms of local and global errors.

Evacuation and outgassing of vacuum glazing

The evacuation and outgassing of vacuum glazing is studied by modeling and experimental methods. A two-dimensional finite difference method of modeling the pumpdown process is validated by comparison with a one dimensional, cylindrically symmetric finite difference method, and with measurements on experimental samples. In the absence of outgassing, the pressure in the glazing decreases linearly with time in the viscous flow region; in the molecular flow region, the pressure decreases exponentially with time. The time taken to evacuate the sample from atmospheric pressure to the transition pressure is approximately equal to one time constant for the subsequent exponential pressure decrease. The time constant associated with the exponential pressure decrease in the molecular flow region can be written τ=V/C, where V is the volume of the sample, and C is the series combination of the conductance associated with the radially inwards gas flow towards the pumpout port, and the conductance of the pumpout tube itself. For very small gaps, τ varies inversely with gap; for large gaps, τ is proportional to gap. The time constant is least for gaps in the range 0.1–0.5 mm, depending on the dimensions of the sample and the pumpout tube. The outgassing rate at elevated temperatures from the internal surfaces of the glazing is shown to decrease slowly, relative to the time constant for pressure reduction in a typical sample. Over a wide pressure range, the outgassing rate is unaffected by the pressure in the sample. The time required to evacuate a sample of vacuum glazing is therefore determined mainly by outgassing, rather than by pumping considerations.

Direct Numerical Simulation of Turbulent Compression Ramp Flow

A numerical procedure for the direct numerical simulation of compressible turbulent flow and shock–turbulence interaction is detailed and analyzed. An upwind-biased finite-difference scheme with a compact centered stencil is used to discretize the convective part of the Navier–Stokes equations. The scheme has a uniformly high approximation order and allows for a spectral-like wave resolution while dissipating nonresolved wave numbers. When hybridized with an essentially nonoscillatory scheme near discontinuities, the scheme becomes shock–capturing and its resolution properties are preserved. Diffusive parts are discretized with symmetric compact finite differences and an explicit Runge–Kutta scheme is used for time-advancement. The peculiarities of efficient upwinding and coupling procedures are described and validation results are given. Using direct numerical simulation data, some aspects of turbulent supersonic compression ramp flow are studied to demonstrate the effectiveness of the simulation procedure.

Modelling one-dimensional unsteady flows in ducts: Symmetric finite difference schemes versus galerkin discontinuous finite element methods

The Euler equations for one-dimensional unsteady flows in ducts have been solved resorting to classical symmetric shock-capturing methods with second-order accuracy and to the recent discontinuous Galerkin finite-element method, with second- and third-order accuracy. In particular, the finite difference techniques adopted are the two-step Lax—Wendroff method and the MacCormack predictor—corrector method, with the addition of the flux corrected transport (FCT) or of the Davis nonupwind TVD scheme to suppress the spurious oscillations in the vicinity of discontinuous solutions. The finite-element method adopted is based on the weak formulation of the Euler equations, which are solved by introducing a discontinuous finite-element space discretization. A dissipative mechanism has been considered to supplement the FEM with a “discontinuity capturing” operator, adding a “viscous like” term to damp minor numerical overshoots arising in proximity of steep gradients of the solution. The numerical tests chosen to carry out a comparison between these schemes are the shock-tube problem and the shock—turbulence interaction problem. Both the test cases considered show the superiority of third-order FEM calculations, whereas the comparison between the computer run times points out the greater computational effort required.

Computations of nonlinear wave interaction in shock-wave focusing process using finite volume TVD schemes

An upwind and a symmetric finite difference TVD scheme have been extended to cell-centered finite volume methods. These newly extended schemes have been used to analyze a shock-wave focusing phenomenon that is dominated by complex nonlinear interaction of shock waves. Numerical solutions obtained from these two schemes successfully resolved all the waves evolving through the focusing process with reasonable accuracy and resolution. The phenomenon near the focal region was found to be strongly affected by the gas dynamic nonlinearities of shock waves. The computed solutions showed good agreement with experimental data and, thus, demonstrated the accuracy and reliability of the present numerical schemes as applied to the analysis of complex unsteady shock-shock interaction problems. The present numerical simulation has clearly revealed the mechanism of shock-shock/expansion wave interaction in a shock-wave focusing process, which determines the wave-front geometries at and beyond the focus and confines the pressure amplification at the focal region. In addition, some comparisons have been made for the performance of both schemes considering the effect of the limiter functions and those results are proposed.

Finite Difference Calculus Invariant Structure of a Class of Algorithms for the Nonlinear Klein–Gordon Equation

In a previous work, the authors have presented a formalism for deriving systematically invariant, symmetric finite difference algorithms for nonlinear evolution differential equations that admit conserved quantities. This formalism is herein cast in the context of exact finite difference calculus. The algorithms obtained from the proposed formalism are shown to derive exactly from discrete scalar potential functions using finite difference calculus, in the same sense as that of the corresponding differential equation being derivable from its associated energy function (a conserved quantity). A clear ramification of this result is that the derived algorithms preserve certain discrete invariant quantities, which are the consistent counterpart of the invariant quantities in the continuous case. Results on the nonlinear stability of a class of algorithms that are derived using the proposed formalism, and that preserve energy or linear momentum, are discussed in the context of finite difference calculus. Some numerical experiments are presented to illustrate the conservation property of the proposed algorithms.

Convergence Results for Harmonic Gradient Estimators

Sensitivity analysis of steady state simulation outputs typically involves estimating gradients. This paper presents convergence results for the harmonic gradient estimators. Sufficient conditions are formulated that validate the interchange of the derivative and the expectation operators for these estimators. The relationship between these estimators and finite differences gradient estimators is discussed. In particular, the harmonic estimators are shown to be variations of finite differences gradient estimators. Exploiting the orthogonal property of the trigonometric basis results in harmonic gradient estimation procedures requiring two simulation runs. Computational results with various queueing system simulation models are included to compare and illustrate the different estimators. These results suggest that the harmonic gradient estimation procedures requiring two simulation runs may be an alternative to forward (symmetric) finite differences gradient estimation procedures with common random number streams requiring p + 1 (2p) simulation runs. INFORMS Journal on Computing, ISSN 1091-9856, was published as ORSA Journal on Computing from 1989 to 1995 under ISSN 0899-1499.

A Tensor product generalized ADI method for the method of planes

AbstractWe consider solving separable, second order, linear elliptic prtial differential equations in three independent variables. If the partial differential opertor separates into two terms, one depending on x and y, and one depending on z, then we use the method of planes to obtain a discrete problem, which we write in tensor product from as We apply a new interative method, the tensor product generalized alternating direction implicit method, to solve the discrete problem. We study a specific implementation that uses Hermite bicubic collocation in the xy direction and symmetric finite differences in the z direction. We demostrate that this method is a fast and accurate way to solve the large linear systems arising from three‐dimensional elliptic problems.

A tensor product generalized ADI method for elliptic problems on cylindrical domains with holes

We consider solving second order linear elliptic partial differential equations together with Dirichlet boundary conditions in three dimensions on cylindrical domains (nonrectangular in x and y) with holes. We approximate the partial differential operators by standard partial difference operators. If the partial differential operator separates into two terms, one depending on x and y, and one depending on z, then the discrete elliptic problem may be written in tensor product form as ( T z ⊗ I + I⊗ A xy ) U= F. We consider a specific implementation which uses a Method of Planes approach with unequally spaced finite differences in the xy direction and symmetric finite difference in the z direction. We establish the convergence of the Tensor Product Generalized Alternating Direction Implicit iterative method applied to such discrete problems. We show that this method gives a fast and memory efficient scheme for solving a large class of elliptic problems.

New symmetric finite difference methods for computing eigenvalues of a boundary‐value problem

AbstractWe develop a new finite difference method of order two and its refinement for computing the eigenvalues of a two‐point boundary‐value problem involving a fourth‐order linear ordinary differential equation. Our methods leads to generalized symmetric five‐ and seven‐band matrix eigenvalue problems. To demonstrate the effectiveness of the proposed methods, the numerical evidence is also included.

Collocation for Singular Perturbation Problems III: Nonlinear Problems without Turning Points

A class of nonlinear singularly perturbed boundary value problems is considered, with restrictions which allow only well-posed problems with possible boundary layers, but no turning points. For the numerical solution of these problems, a close look is taken at a class of general purpose, symmetric finite difference schemes arising from collocation. It is shown that if locally refined meshes, whose practical construction is discussed, are employed, then high order uniform convergence of the numerical solution is obtained. Nontrivial examples are used to demonstrate that highly accurate solutions to problems with extremely thin boundary layers can be obtained in this way at a very reasonable cost.

A new symmetric five-diagonal finite difference method for computing eigenvalues of fourth-order two-point boundary value problems

We present a new symmetric five-diagonal finite difference method for computing eigenvalues of two-point boundary value problems involving a fourth-order differential equation.

Difference operator for variable stiffness plates

AbstractA symmetric central finite difference operator is derived for the small deflection analysis of transversely loaded, ‘thin’, elastic plates, whose stiffnesses vary as functions of the in‐plane co‐ordinates. The operator is written here in Cartesian co‐ordinates for a square mesh, but the same technique can also be used to construct an operator for a rectangular mesh.