Search

We formulate hydrodynamic equations and spectrally accurate numerical methodsfor investigating the role of geometry in flows within two-dimensional fluidinterfaces. To achieve numerical approximations having high precision and levelof symmetry for radial manifold shapes, we develop spectral Galerkin methodsbased on hyperinterpolation with Lebedev quadratures for $L^2$-projection tospherical harmonics. We demonstrate our methods by investigating hydrodynamicresponses as the surface geometry is varied. Relative to the case of a sphere,we find significant changes can occur in the observed hydrodynamic flowresponses as exhibited by quantitative and topological transitions in thestructure of the flow. We present numerical results based on theRayleigh-Dissipation principle to gain further insights into these flowresponses. We investigate the roles played by the geometry especiallyconcerning the positive and negative Gaussian curvature of the interface. Weprovide general approaches for taking geometric effects into account forinvestigations of hydrodynamic phenomena within curved fluid interfaces.

We present an exact and efficient algorithm for reaction-diffusion-nucleation processes on a 2D-lattice. The algorithm makes use of first passage time (FPT) to replace the computationally intensive simulation of diffusion hops in KMC by larger jumps when particles are far away from step-edges or other particles. Our approach computes exact probability distributions of jump times and target locations in a closed-form formula, based on the eigenvectors and eigenvalues of the corresponding 1D transition matrix, maintaining atomic-scale resolution of resulting shapes of deposit islands. We have applied our method to three different test cases of electrodeposition: pure diffusional aggregation for large ranges of diffusivity rates and for simulation domain sizes of up to 4096 × 4096 sites, the effect of diffusivity on island shapes and sizes in combination with a KMC edge diffusion, and the calculation of an exclusion zone in front of a step-edge, confirming statistical equivalence to standard KMC simulations. The algorithm achieves significant speedup compared to standard KMC for cases where particles diffuse over long distances before nucleating with other particles or being captured by larger islands. © 2013 Elsevier Inc.

In this paper we present two unconditionally energy stable finite difference schemes for the modified phase field crystal (MPFC) equation, a sixth-order nonlinear damped wave equation, of which the purely parabolic phase field crystal (PFC) model can be viewed as a special case. The first is a convex splitting scheme based on an appropriate decomposition of the discrete energy and is first order accurate in time and second order accurate in space. The second is a new, fully second-order scheme that also respects the convex splitting of the energy. Both schemes are nonlinear but may be formulated from the gradients of strictly convex, coercive functionals. Thus, both are uniquely solvable regardless of the time and space step sizes. The schemes are solved by efficient nonlinear multigrid methods. Numerical results are presented demonstrating the accuracy, energy stability, efficiency, and practical utility of the schemes. In particular, we show that our multigrid solvers enjoy optimal, or nearly optimal complexity in the solution of the nonlinear schemes. © 2013 Elsevier Inc.

We analyze a new mathematical and numerical framework, the " Voronoi Implicit Interface Method" (" VIIM" ), first introduced in Saye and Sethian (2011) [R.I. Saye, J.A. Sethian, The Voronoi Implicit Interface Method for computing multiphase physics, PNAS 108 (49) (2011) 19498-19503] for tracking multiple interacting and evolving regions (" phases" ) whose motion is determined by complex physics (fluids, mechanics, elasticity, etc.). From a mathematical point of view, the method provides a theoretical framework for moving interface problems that involve multiple junctions, defining the motion as the formal limit of a sequence of related problems. Discretizing this theoretical framework provides a numerical methodolology which automatically handles multiple junctions, triple points and quadruple points in two dimensions, as well as triple lines, etc. in higher dimensions. Topological changes in the system occur naturally, with no surgery required. In this paper, we present the method in detail, and demonstrate several new extensions of the method to different physical phenomena, including curvature flow with surface energy densities defined on a per-phase basis, as well as multiphase fluid flow in which density, viscosity and surface tension can be defined on a per-phase basis.We test this method in a variety of ways. We perform rigorous analysis and demonstrate convergence in both two and three dimensions for a variety of evolving interface problems, including verification of von Neumann-Mullins' law in two dimensions (and its analog in three dimensions), as well as normal driven flow and curvature flow with and without constraints, demonstrating topological change and the effects of different boundary conditions. We couple the method to a second order projection method solver for incompressible fluid flow, and study the effects of membrane permeability and impermeability, large shearing torsional forces, and the effects of varying density, viscosity and surface tension on a per-phase basis. Finally, we demonstrate convergence in both space and time of a topological change in a multiphase foam. © 2012 Elsevier Inc..

The Ewald method is applied to accelerate the evaluation of the Green's function (GF) of an infinite equispaced linear array of point sources with linear phasing. Only a few terms are needed to evaluate Ewald sums, which are cast in terms of error functions and exponential integrals, to high accuracy. It is shown analytically that the choice of the standard "optimal" Ewald splitting parameter E0 causes overflow errors at high frequencies (period large compared to the wavelength), and convergence rates are analyzed. A recipe for selecting the Ewald splitting parameter is provided. © 2006 Elsevier Inc. All rights reserved.

A Radial Basis Function Generated Finite-Differences (RBF-FD) inspired technique for evaluating definite integrals over bounded volumes that have smooth boundaries in three dimensions is described. A key aspect of this approach is that it allows the user to approximate the value of the integral without explicit knowledge of an expression for the boundary surface. Instead, a tesselation of the node set is utilized to inform the algorithm of the domain geometry. Further, the method applies to node sets featuring spatially varying density, facilitating its use in Applied Mathematics, Mathematical Physics and myriad other application areas, where the locations of the nodes might be fixed by experiment or previous simulation. By using the RBF-FD-like approach, the proposed algorithm computes quadrature weights for N arbitrarily scattered nodes in only O (N log ⁡ N) operations with tunable orders of accuracy. • A fast and accurate method for quadrature over 3D volumes with smooth boundaries. • High orders of accuracy are achieved on node sets featuring no particular uniformity. • Can be distributed to any number of available processors. • Knowledge of an expression for the boundary surface is unnecessary. [ABSTRACT FROM AUTHOR]

We present three new families of fast algorithms for classical potential theory, based on Ewald summation and fast transforms of Gaussians and Fourier series. Ewald summation separates the Green function for a cube into a high-frequency localized part and a rapidly-converging Fourier series. Each part can then be evaluated efficiently with appropriate fast transform algorithms. Our algorithms are naturally suited to the use of Green functions with boundary conditions imposed on the boundary of a cube, rather than free-space Green functions. Our first algorithm evaluates classical layer potentials on the boundary of a d-dimensional domain, with d equal to two or three. The quadrature error is O(hm) +ε, where h is the mesh size on the boundary and m is the order of quadrature used. The algorithm evaluates the discretized potential using N elements at O(N) points in O(N log N) arithmetic operations. The constant in O(N log N) depends logarithmically on the desired error tolerance. Our second scheme evaluates a layer potential on the domain itself, with the same accuracy. It produces Md values using N boundary elements in O((N + Md) log M) arithmetic operations. Our third method evaluates a discrete sum of values of the Green function, of the type which occur in particle methods. It attains error ε{lunate} at a cost O(Na log N), where a = 2 (1 + D d) and D is the Hausdorff dimension of the set where the sources concentrate in the limit N √ ∞. Thus it is O(N log N) when the sources do not cluster too much and close to O(N log N) in the important practical case when the points are uniformly distributed over a hypersurface. We also sketch an O(N log N) algorithm based on special functions. Two-dimensional numerical results are presented for all three algorithms. Layer potentials are evaluated to second-order accuracy, in times which exhibit considerable speedups even over a reasonably sophisticated direct calculation. Discrete sum calculations are speeded up astronomically; our algorithm reduces the CPU time required for a calculation with 40,000 points from six months to one hour. © 1992.

The capability to rapidly and robustly simulate the scattering of linear waves by periodic, multiply layered media in two and three dimensions is crucial in many engineering applications. In this regard, we present a High-Order Perturbation of Surfaces method for linear wave scattering in a multiply layered periodic medium to find an accurate numerical solution of the governing Helmholtz equations. For this we truncate the bi-infinite computational domain to a finite one with artificial boundaries, above and below the structure, and enforce transparent boundary conditions there via Dirichlet–Neumann Operators. This is followed by a Transformed Field Expansion resulting in a Fourier collocation, Legendre–Galerkin, Taylor series method for solving the problem in a transformed set of coordinates. Assorted numerical simulations display the spectral convergence of the proposed algorithm. [ABSTRACT FROM AUTHOR]

The direct differentiation and improved least squares shadowing methods are both developed for accurately and efficiently calculating the sensitivity coefficients of time averaged quantities for chaotic dynamical systems. The key idea is to recast the time averaged integration term in the form of differential equation before applying the sensitivity analysis method. An additional constraint-based equation which forms the augmented equations of motion is proposed to calculate the time averaged integration variable and the sensitivity coefficients are obtained as a result of solving the augmented differential equations. The application of the least squares shadowing formulation to the augmented equations results in an explicit expression for the sensitivity coefficient which is dependent on the final state of the Lagrange multipliers. The LU factorization technique to calculate the Lagrange multipliers leads to a better performance for the convergence problem and the computational expense. Numerical experiments on a set of problems selected from the literature are presented to illustrate the developed methods. The numerical results demonstrate the correctness and effectiveness of the present approaches and some short impulsive sensitivity coefficients are observed by using the direct differentiation sensitivity analysis method. [ABSTRACT FROM AUTHOR]

We present a practical algorithm based on symplectic splitting methods intended for the numerical integration in time of the Schrödinger equation when the Hamiltonian operator is either time-independent or changes slowly with time. In the later case, the evolution operator can be effectively approximated in a step-by-step manner: first divide the time integration interval in sufficiently short subintervals, and then successively solve a Schrödinger equation with a different time-independent Hamiltonian operator in each of these subintervals. When discretized in space, the Schrödinger equation with the time-independent Hamiltonian operator obtained for each time subinterval can be recast as a classical linear autonomous Hamiltonian system corresponding to a system of coupled harmonic oscillators. The particular structure of this linear system allows us to construct a set of highly efficient schemes optimized for different precision requirements and time intervals. Sharp local error bounds are obtained for the solution of the linear autonomous Hamiltonian system considered in each time subinterval. Our schemes can be considered, in this setting, as polynomial approximations to the matrix exponential in a similar way as methods based on Chebyshev and Taylor polynomials. The theoretical analysis, supported by numerical experiments performed for different time-independent Hamiltonians, indicates that the new methods are more efficient than schemes based on Chebyshev polynomials for all tolerances and time interval lengths. The algorithm we present automatically selects, for each time subinterval, the most efficient splitting scheme (among several new optimized splitting methods) for a prescribed error tolerance and given estimates of the upper and lower bounds of the eigenvalues of the discretized version of the Hamiltonian operator. [ABSTRACT FROM AUTHOR]

Abstract: We present a mathematical framework for constructing and analyzing parallel algorithms for lattice kinetic Monte Carlo (KMC) simulations. The resulting algorithms have the capacity to simulate a wide range of spatio-temporal scales in spatially distributed, non-equilibrium physiochemical processes with complex chemistry and transport micro-mechanisms. Rather than focusing on constructing exactly the stochastic trajectories, our approach relies on approximating the evolution of observables, such as density, coverage, correlations and so on. More specifically, we develop a spatial domain decomposition of the Markov operator (generator) that describes the evolution of all observables according to the kinetic Monte Carlo algorithm. This domain decomposition corresponds to a decomposition of the Markov generator into a hierarchy of operators and can be tailored to specific hierarchical parallel architectures such as multi-core processors or clusters of Graphical Processing Units (GPUs). Based on this operator decomposition, we formulate parallel Fractional step kinetic Monte Carlo algorithms by employing the Trotter Theorem and its randomized variants; these schemes, (a) are partially asynchronous on each fractional step time-window, and (b) are characterized by their communication schedule between processors. The proposed mathematical framework allows us to rigorously justify the numerical and statistical consistency of the proposed algorithms, showing the convergence of our approximating schemes to the original serial KMC. The approach also provides a systematic evaluation of different processor communicating schedules. We carry out a detailed benchmarking of the parallel KMC schemes using available exact solutions, for example, in Ising-type systems and we demonstrate the capabilities of the method to simulate complex spatially distributed reactions at very large scales on GPUs. Finally, we discuss work load balancing between processors and propose a re-balancing scheme based on probabilistic mass transport methods. [Copyright &y& Elsevier]

Abstract: This paper gives a first numerical comparison of the Discontinuous-Galerkin-(DG) and Spectral-Difference-Method (SD) under almost equivalent conditions for smooth test cases. [Copyright &y& Elsevier]