Accurate computation of Stokes flow driven by an open immersed interface

We present numerical methods for computing two-dimensional Stokes flow driven by forces singularly supported along an open, immersed interface. Two second-order accurate methods are developed: one for accurately evaluating boundary integral solutions at a point, and another for computing Stokes solution values on a rectangular mesh. We first describe a method for computing singular or nearly singular integrals, such as a double layer potential due to sources on a curve in the plane, evaluated at a point on or near the curve. To improve accuracy of the numerical quadrature, we add corrections for the errors arising from discretization, which are found by asymptotic analysis. When used to solve the Stokes equations with sources on an open, immersed interface, the method generates second-order approximations, for both the pressure and the velocity, and preserves the jumps in the solutions and their derivatives across the boundary. We then combine the method with a mesh-based solver to yield a hybrid method for computing Stokes solutions at N² grid points on a rectangular grid. Numerical results are presented which exhibit second-order accuracy. To demonstrate the applicability of the method, we use the method to simulate fluid dynamics induced by the beating motion of a cilium. The method preserves the sharp jumps in the Stokes solution and their derivatives across the immersed boundary. Model results illustrate the distinct hydrodynamic effects generated by the effective stroke and by the recovery stroke of the ciliary beat cycle.     

A New Method to Measure Small Amounts of Solute Atoms on Planar Defects and Application to Inversion Domain Boundaries in Doped Zinc Oxide

We demonstrate the application of a new method of analytical transmission electron microscopy for measuring very accurately small amounts of solute atoms within a well-defined planar defect such as a stacking fault, grain boundary or an interface. The method is based on acquiring several spectra with different electron beam diameters from the same position centred on the defect. It can be applied to energy-dispersive X-ray microanalysis (EDXS) or electron energy-loss spectroscopy (EELS) and does not necessitate a scanning unit. The accuracy has been tested numerically under different conditions using simulations for a specific geometry and has been found to be substantially better than that of any other current standard technique. Our calculations suggest an extremely high accuracy theoretically achievable in the determination of e.g. the Gibbsian solute excess or the doping level of a grain boundary down to about ±1% of an effective monolayer, i.e. ±0.1 atoms/nm² under typical experimental conditions. The method has been applied to zinc oxide, which forms inversion domain boundaries (IDBs) when doped with different transition metal oxides such as SnO₂ or Sb₂O₃. We obtained an experimental precision of ±0.4 atoms/nm², which has allowed us to solve the atomic structure of the IDBs.     

A second order virtual node method for elliptic problems with interfaces and irregular domains

We present a second order accurate, geometrically flexible and easy to implement method for solving the variable coefficient Poisson equation with interfacial discontinuities or on irregular domains, handling both cases with the same approach. We discretize the equations using an embedded approach on a uniform Cartesian grid employing virtual nodes at interfaces and boundaries. A variational method is used to define numerical stencils near these special virtual nodes and a Lagrange multiplier approach is used to enforce jump conditions and Dirichlet boundary conditions. Our combination of these two aspects yields a symmetric positive definite discretization. In the general case, we obtain the standard 5-point stencil away from the interface. For the specific case of interface problems with continuous coefficients, we present a discontinuity removal technique that admits use of the standard 5-point finite difference stencil everywhere in the domain. Numerical experiments indicate second order accuracy in L^∞.     

The Elastoplast Discontinuous Galerkin (EDG) method for the Navier–Stokes equations

The present work details the Elastoplast (this name is a translation from the French “sparadrap”, a concept first applied by Yves Morchoisne for Spectral methods [1]) Discontinuous Galerkin (EDG) method to solve the compressible Navier–Stokes equations. This method was first presented in 2009 at the ICOSAHOM congress with some Cartesian grid applications. We focus here on unstructured grid applications for which the EDG method seems very attractive. As in the Recovery method presented by van Leer and Nomura in 2005 for diffusion, jumps across element boundaries are locally eliminated by recovering the solution on an overlapping cell. In the case of Recovery, this cell is the union of the two neighboring cells and the Galerkin basis is twice as large as the basis used for one element. In our proposed method the solution is rebuilt through an L² projection of the discontinuous interface solution on a small rectangular overlapping interface element, named Elastoplast, with an orthogonal basis of the same order as the one in the neighboring cells. Comparisons on 1D and 2D scalar diffusion problems in terms of accuracy and stability with other viscous DG schemes are first given. Then, 2D results on acoustic problems, vortex problems and boundary layer problems both on Cartesian or unstructured triangular grids illustrate stability, precision and versatility of this method.     

Studies on interdiffusion in Pd/Mg/Si films: Towards improved cyclic stability in hydrogen storage

The paper reports the diffusion coefficients of grain boundary diffusion and grain boundary assisted lattice diffusion of Pd in Mg in Pd/Mg/Si system, a useful material for hydrogen storage, at 473 K in vacuum. The grain boundary diffusivity is measured by Whipple model and grain boundary assisted lattice diffusivity by plateau rise method using Pd depth profiles constructed by Rutherford backscattering spectrometry. It is established that grain boundary diffusivities are about six orders of magnitude faster than lattice diffusivities. Fine grained microstructure of Pd film, high abundance of defects in Mg film and higher stability associated with Pd-Mg intermetallics are responsible for the diffusion of Pd into grain boundaries and subsequently in the interiors of Mg. Besides the indiffusion of Pd, annealing also brings about an outdiffusion of Mg into Pd film. Examination by nuclear reaction analysis involving ²⁴Mg(p,p′γ)²⁴Mg resonance reaction shows the occurrence of Mg outdiffusion. Minimization of surface energy is presumably the driving force of the process. In addition to Pd/Mg interface, diffusion occurs across Mg/Si (substrate) interface as well on increasing the annealing temperature above 473 K. These studies show that dehydrogenation of films accomplished by vacuum annealing should be limited to temperatures less than 473 K to minimize the loss of surface Pd, the catalyst of the hydrogen absorption-desorption process and Mg, the hydrogen storing element, by way of interfacial reactions.     

Photoelastic investigation of interfacial fracture between orthotropic and isotropic materials

Static and dynamic fracture of interfaces between orthotropic and isotropic materials were studied using photoelasticity. In this study, a bi-material specimen made of PSM-1^® and Scotchply^® 1002, a unidirectional glass fiber reinforced epoxy composite, was used. Two fiber orientations, fibers parallel to the interface (α=0°) and fibers perpendicular to the interface (α=90°) were considered. Center crack bi-material specimens having different crack lengths were loaded quasi-statically and the full-field isochromatics were recorded using a digital camera. The complex stress intensity factor corresponding to each crack length was calculated from the isochromatics and the values were compared to that obtained from boundary collocation method. Dynamic interfacial fracture was studied with an edge crack bi-material geometry for the two different fiber orientations. The isochromatics around the propagating crack were recorded using a digital high-speed camera. The fracture parameters such as crack speed, complex stress intensity factor and energy release rate were extracted from the isochromatics using the asymptotic stress field equations. The complex stress intensity factor obtained from the static experiments was in close agreement with that calculated using the boundary collocation method. The results also indicated that the fiber orientation with respect to the interface influences the fracture parameters for stationary and propagating cracks.     

Boundary integral solutions of coupled Stokes and Darcy flows

An accurate computational method based on the boundary integral formulation is presented for solving boundary value problems for Stokes and Darcy flows. The method also applies to problems where the equations are coupled across an interface through appropriate boundary conditions. The adopted technique consists of first reformulating the singular integrals for the fluid quantities as single and double layer potentials. Then the layer potentials are regularized and discretized using standard quadratures. As a final step, the leading term in the regularization error is eliminated in order to gain one more order of accuracy. The numerical examples demonstrate the increase of the convergence rate from first to second order and show a decrease in magnitude of the error. The coupled problems require the computation of the gradient of the Stokes velocity at the common interface. This boundary condition is also written as a combination of single and double layer potentials so that the same approach can be used to compute it accurately. Extensive numerical examples show the increased accuracy gained by the correction terms.     

Accurate measurement of grain boundary misorientation by Large Angle Convergent Beam Electron Diffraction

A new method using Large Angle Convergent Beam Electron Diffraction (LACBED) patterns is proposed to measure accurately the grain boundary misorientation. The LACBED patterns which are obtained with a defocused convergent electron beam having a convergence semi-angle in the range 1 to 5^o contain very sharp deficiency lines. Due to the good quality of the LACBED patterns, these sharp deficiency lines can be used to measure with great accuracy the grain boundary misorientation. In addition, since the LACBED method is a defocus mode method, the patterns contain at the same time information on the reciprocal space (the deficiency lines typical of the crystal orientation of the two grains on each side of the grain boundary) and on the real space (the image of the grain boundary). We describe a method which allows the identification of the misorientation from these LACBED patterns. The main point to consider is the accuracy which is about 0.05^o. It is much better than the one obtained from other conventional methods used to measure this misorientation.     

Fundamental frequency of clamped plates with circularly periodic boundaries

A boundary perturbation method is developed to determine the fundamental frequency of vibrating plates. The method is then applied to wavy, star shape and polygonal plates with clamped boundary conditions. Approximate analytical solutions of the fundamental frequency are obtained with an accuracy of O(ε⁴), where ε is the deviation from the unit circle.     

A boundary condition capturing immersed interface method for 3D rigid objects in a flow

To simulate the flow around an object, we can replace the object with the fluid enclosed by a singular force. We can then simulate the flow on a fixed domain with a fluid–fluid interface supporting the singular force. In this paper, we present a boundary condition capturing approach to determine the singular force for a 3D rigid object. We apply a discontinuous body force to enforce the rigid motion of the fluid replacing the object and compute the singular force based on the kinematics of the object. Due to the singular force and the body force, the flow is not smooth across the interface. We solve the flow using the immersed interface method. Our boundary condition capturing immersed interface method is very efficient and stable, and its accuracy based on the infinity norm is near second order for the velocity and above first order for the pressure.     

A model for sonar interrogation of complex bottom and surface targets in shallow-water waveguides

Many problems of current interest in underwater acoustics involve low-frequency broadband sonar interrogation of objects near the sea surface or sea floor of a shallow-water environment. When the target is situated near the upper or lower boundary of the water column the acoustic interactions with the target objects are complicated by interactions with the nearby free surface or fluid-sediment interface, respectively. A practical numerical method to address such situations is presented. The model provides high levels of accuracy with the flexibility to handle complex, three-dimensional targets in range-independent environments. The model is demonstrated using several bottom target scenarios, with and without locally undulating seabeds. The impact of interface and boundary interactions is considered with an eye toward using the sonar return signal as the basis for acoustic imaging or spectral classification.     

A second order virtual node method for elliptic problems with interfaces and irregular domains in three dimensions

We present a numerical method for the variable coefficient Poisson equation in three-dimensional irregular domains and with interfacial discontinuities. The discretization embeds the domain and interface into a uniform Cartesian grid augmented with virtual degrees of freedom to provide accurate treatment of jump and boundary conditions. The matrix associated with the discretization is symmetric positive definite and equal to the standard 7-point finite difference Poisson stencil away from embedded interfaces and boundaries. Numerical evidence suggests second order accuracy in the L^∞-norm. Our approach improves the treatment of Dirichlet and jump constraints in the recent work of Bedrossian et al. [1] and introduces innovations necessary in three dimensions. Specifically, we construct new constraint-based Lagrange multiplier spaces that significantly improve the conditioning of the associated linear system of equations; we provide a method for cell-local polyhedral approximation to the zero isocontour surface of a level set needed for three-dimensional embedding; and we show that the new Lagrange multiplier spaces naturally lead to a class of easy-to-implement multigrid methods that achieve near optimal efficiency, as shown by numerical examples. For the specific case of a continuous Poisson coefficient in interface problems, we provide an expansive treatment of the construction of a particular solution that satisfies the value jump and flux jump constraints. As in [1], this is used in a discontinuity removal technique that yields the standard 7-point stencil across the interface and only requires a modification to the right-hand side of the linear system.     

A volume of fluid approach for crystal growth simulation

A new approach to simulating the dendritic growth of pure metals, based on a recent volume of fluid (VOF) method with PLIC (piecewise linear interface calculation) reconstruction of the interface, is presented. The energy equation is solved using a diffuse-interface method, which avoids the need to apply the thermal boundary conditions directly at the solid front. The thermal gradients at both sides of the interface, which are needed to obtain the front velocity, are calculated with the aid of a distance function to the reconstructed interface. The advection equation of a discretized solid fraction function is solved using the unsplit VOF advection method proposed by López et al. [J. Comput. Phys. 195 (2004) 718–742] (extended to three dimensions by Hernández et al. [Int. J. Numer. Methods Fluids 58 (2008) 897–921]), and the interface curvature is computed using an improved height function technique, which provides second-order accuracy. The proposed methodology is assessed by comparing the numerical results with analytical solutions and with results obtained by different authors for the formation of complex dendritic structures in two and three dimensions.     

Multiscale molecular dynamics using the matched interface and boundary method

The Poisson-Boltzmann (PB) equation is an established multiscale model for electrostatic analysis of biomolecules and other dielectric systems. PB based molecular dynamics (MD) approach has a potential to tackle large biological systems. Obstacles that hinder the current development of PB based MD methods are concerns in accuracy, stability, efficiency and reliability. The presence of complex solvent-solute interface, geometric singularities and charge singularities leads to challenges in the numerical solution of the PB equation and electrostatic force evaluation in PB based MD methods. Recently, the matched interface and boundary (MIB) method has been utilized to develop the first second order accurate PB solver that is numerically stable in dealing with discontinuous dielectric coefficients, complex geometric singularities and singular source charges. The present work develops the PB based MD approach using the MIB method. New formulation of electrostatic forces is derived to allow the use of sharp molecular surfaces. Accurate reaction field forces are obtained by directly differentiating the electrostatic potential. Dielectric boundary forces are evaluated at the solvent-solute interface using an accurate Cartesian-grid surface integration method. The electrostatic forces located at reentrant surfaces are appropriately assigned to related atoms. Extensive numerical tests are carried out to validate the accuracy and stability of the present electrostatic force calculation. The new PB based MD method is implemented in conjunction with the AMBER package. MIB based MD simulations of biomolecules are demonstrated via a few example systems.     

A “superlattice” model for a smooth GaAs/AlAs (001) heterointerface

The effect of a smooth interface potential on the electronic states in GaAs/AlAs (001) structures is investigated using the pseudopotential method. In this approach, the transition region between GaAs and AlAs is assumed to be a layer corresponding to a half-period of the (AlAs)₂(GaAs)₂ superlattice, with the potential of this layer being close to the real potential near the heterointerface. In this case, the intervalley mixing occurs at two boundaries and in the transition layer rather than at one boundary, as in the model with a sharply cut-off potential. It is shown that a smooth potential has an appreciable effect on electron tunneling in structures with thin layers. This effect is especially important in the case where short-wavelength X states are involved. For one GaAs/AlAs (001) boundary, the transition layer acts as a quantum well localizing the charge density of a mixed Γ-X state near the boundary. In structures with a layer thickness of less than 2 nm, the differences in the resonance energies obtained in the models with a smooth heterointerface and with a sharp heterointerface can be as high as ～0.1 eV. The envelopes of the wave functions associated with Γ ₁ ⁽¹⁾ , Γ ₁ ⁽²⁾ , and Γ ₃ ⁽¹⁾ superlattice valleys and with Γ₁, X ₁, and X ₃ valleys of GaAs and AlAs are analyzed. It is shown that the matching matrices for the envelope functions at the GaAs/(AlAs)₂(GaAs)₂ and (AlAs)₂(GaAs)₂/AlAs boundaries depend only weakly on the electron energy near the bottom of the conduction band and that the probability densities calculated using these functions agree with the results of many-band calculations. Therefore, these functions can be used to construct a model with a smooth interface potential in the framework of the effective-mass method.     

Fixed-grid method for the modelling of unsteady partial oxidation of a spherical coal particle

The main goal of this work is the development of a fixed-grid method to model unsteady partial oxidation of a solid with implicit tracking of the interface. As a first step diffusive oxidation of a spherical coal particle is considered. The energy and species conservation equations formulated in spherical coordinates are discretised using the finite-volume approach. The boundary conditions for the temperature and species mass fractions at the solid–gas interface are modelled via special source terms activated in the interface cells. The numerical model was validated against analytic one- and two-film models for coal combustion in a dry-air atmosphere. Very good agreement was obtained. Based on the model developed a numerical study was carried out on the influence of water vapour on the partial oxidation of a spherical coal particle. Numerous numerical simulations were performed for particle diameters in the range 200×10^?6 m to 2×10^?2 m. The ambient temperature was varied in the range between 700 and 3000 K. The analysis of results showed that the addition of H₂O has an influence on the solution convergence due to the catalytic effect of water in the coal monoxide oxidation reaction making the whole system stiffer. However, at the same time, it was found that if the ambient mass fraction of water vapour is below 1×10^?3, its influence on combustion rates is minimal. The results of numerical simulations obtained for higher H₂O concentration (>1×10^?3) are discussed.     

A ghost fluid,level set methodology for simulating multiphase electrohydrodynamic flows with application to liquid fuel injection

In this paper, we present the development of a sharp numerical scheme for multiphase electrohydrodynamic (EHD) flows for a high electric Reynolds number regime. The electric potential Poisson equation contains EHD interface boundary conditions, which are implemented using the ghost fluid method (GFM). The GFM is also used to solve the pressure Poisson equation. The methods detailed here are integrated with state-of-the-art interface transport techniques and coupled to a robust, high order fully conservative finite difference Navier–Stokes solver. Test cases with exact or approximate analytic solutions are used to assess the robustness and accuracy of the EHD numerical scheme. The method is then applied to simulate a charged liquid kerosene jet.     

A nonconforming combination for solving elliptic problems with interfaces

Nonconforming combinations are provided for solving interface problems of elliptic equations. In these approaches, the Ritz-Galerkin method with particular solutions is used for the part of a solution domain where there are interface singular points; and the conventional finite element method is used for the rest of the solution domain. In addition, admissible functions chosen are constrained to be continuous only at the element nodes on the common boundary of the subdomains. Error bounds are derived in the Sobolev norms, and numerical experiments are given for solving a model interface problem of the equation, −Δu + U = 0. Moreover, a significant coupling relation, L + 1 = O(|ln h|), is found for interface problems by using the nonconforming combinations, where (L + 1) is the total number of particular solutions used in the Ritz-Galerkin method, and h is the maximal boundary length of triangular elements in the finite element method.     

An accurate, stable and efficient domain-type meshless method for the solution of MHD flow problems

The aim of the present paper is the development of an efficient numerical algorithm for the solution of magnetohydrodynamics flow problems for regular and irregular geometries subject to Dirichlet, Neumann and Robin boundary conditions. Toward this, the meshless point collocation method (MPCM) is used for MHD flow problems in channels with fully insulating or partially insulating and partially conducting walls, having rectangular, circular, elliptical or even arbitrary cross sections. MPC is a truly meshless and computationally efficient method. The maximum principle for the discrete harmonic operator in the meshfree point collocation method has been proven very recently, and the convergence proof for the numerical solution of the Poisson problem with Dirichlet boundary conditions have been attained also. Additionally, in the present work convergence is attained for Neumann and Robin boundary conditions, accordingly. The shape functions are constructed using the Moving Least Squares (MLS) approximation. The refinement procedure with meshless methods is obtained with an easily handled and fully automated manner. We present results for Hartmann number up to 10⁵${10}^5$. The numerical evidences of the proposed meshless method demonstrate the accuracy of the solutions after comparing with the exact solution and the conventional FEM and BEM, for the Dirichlet, Neumann and Robin boundary conditions of interior problems with simple or complex boundaries.     

An MPI parallel level-set algorithm for propagating front curvature dependent detonation shock fronts in complex geometries

We present a parallel, two-dimensional, grid-based algorithm for solving a level-set function PDE that arises in Detonation Shock Dynamics (DSD). In the DSD limit, the detonation shock propagates at a speed that is a function of the curvature of the shock surface, subject to a set of boundary conditions applied along the boundaries of the detonating explosive. Our method solves for the full level-set function field, φ(x, y, t), that locates the detonation shock with a modified level-set function PDE that continuously renormalises the level-set function to a distance function based off of the locus of the shock surface, φ(x, y, t)=0. The boundary conditions are applied with ghost nodes that are sorted according to their connectivity to the interior explosive nodes. This allows the boundary conditions to be applied via a local, direct evaluation procedure. We give an extension of this boundary condition application method to three dimensions. Our parallel algorithm is based on a domain-decomposition model which uses the Message-Passing Interface (MPI) paradigm. The computational order of the full level-set algorithm, which is O(N ⁴), where N is the number of grid points along a coordinate line, makes an MPI-based algorithm an attractive alternative. This parallel model partitions the overall explosive domain into smaller sub-domains which in turn get mapped onto processors that are topologically arranged into a two-dimensional rectangular grid. A comparison of our numerical solution with an exact solution to the problem of a detonation rate stick shows that our numerical solution converges at better than first-order accuracy as measured by an L1-norm. This represents an improvement over the convergence properties of narrow-band level-set function solvers, whose convergence is limited to a floor set by the width of the narrow band. The efficiency of the narrow-band method is recovered by using our parallel model.