Non-periodic Boundary Conditions Research Articles

Toward the correction of effective electrostatic forces in explicit-solvent molecular dynamics simulations: restraints on solvent-generated electrostatic potential and solvent polarization.

Despite considerable advances in computing power, atomistic simulations under nonperiodic boundary conditions, with Coulombic electrostatic interactions and in systems large enough to reduce finite-size associated errors in thermodynamic quantities to within the thermal energy, are still not affordable. As a result, periodic boundary conditions, systems of microscopic size and effective electrostatic interaction functions are frequently resorted to. Ensuing artifacts in thermodynamic quantities are nowadays routinely corrected a posteriori, but the underlying configurational sampling still descends from spurious forces. The present study addresses this problem through the introduction of on-the-fly corrections to the physical forces during an atomistic molecular dynamics simulation. Two different approaches are suggested, where the force corrections are derived from special potential energy terms. In the first approach, the solvent-generated electrostatic potential sampled at a given atom site is restrained to a target value involving corrections for electrostatic artifacts. In the second approach, the long-range regime of the solvent polarization around a given atom site is restrained to the Born polarization, i.e., the solvent polarization corresponding to the ideal situation of a macroscopic system under nonperiodic boundary conditions and governed by Coulombic electrostatic interactions. The restraints are applied to the explicit-water simulation of a hydrated sodium ion, and the effect of the restraints on the structural and energetic properties of the solvent is illustrated. Furthermore, by means of the calculation of the charging free energy of a hydrated sodium ion, it is shown how the electrostatic potential restraint translates into the on-the-fly consideration of the corresponding free-energy correction terms. It is discussed how the restraints can be generalized to situations involving several solute particles. Although the present study considers a very simple system only, it is an important step toward the on-the-fly elimination of finite-size and approximate-electrostatic artifacts during atomistic molecular dynamics simulations.Electronic supplementary materialThe online version of this article (doi:10.1007/s00214-014-1600-8) contains supplementary material, which is available to authorized users.

Divergence-Free Wavelet Projection Method for Incompressible Viscous Flow on the Square

We present a wavelet numerical scheme for the discretization of two-dimensional Navier--Stokes equations with Dirichlet boundary condition on the square. This work is an extension to nonperiodic boundary conditions of the previous method of Deriaz and Perrier [E. Deriaz and V. Perrier, Multiscale Model. Simul., 7 (2008), pp. 1101--1129]. Here the temporal discretization is borrowed from the projection method. The projection operator is defined through a discrete Helmholtz--Hodge decomposition using divergence-free wavelet bases; this prevents the use of a Poisson solver as in usual methods, while improving the accuracy of the boundary condition. The stability and precision order of the new method are stated in the linear case of Stokes equations, confirmed by numerical experiments. Finally, the effectiveness, stability, and accuracy of the method are validated by simulations conducted on the benchmark problem of lid-driven cavity flow at Reynolds number $Re=1000$ and $Re=10000$.

A competitive lattice model Monte Carlo method for simulation competitive growth of different polytypes in close-packed crystals: 4H and 6H silicon carbide

A competitive lattice model Kinetic Monte Carlo method has been developed and implemented to simulate the SiC crystals growth by physical vapor transport method. The model is an advancement with respect to standard Monte Carlo algorithms, allowing to simulate polytypes competitive growth and defective evolution. The model can simulate the formation and evolution of threading screw dislocations in SiC as well. Kinetics processes with surface energy and bond energy, including deposition, evaporation and diffusion processes are considered in physical criterion in our model. Periodic boundary conditions and non-periodic boundary conditions are carried out. The competitive growth between 4H- and 6H-SiC polytypes suggests that retaining the step growth mode is an important factor to maintain a stable single 4H polytype during SiC growth. The threading screw dislocations simulation indicates that screw dislocations drive the growth of SiC. The evolution of surface morphologies demonstrates that fast surface diffusion rate and low growth rate may be one way to grow high quality, low-roughness 4H-SiC crystals.

Effect of Ferroelectric Polarization on Ionic Transport and Resistance Degradation in BaTiO3 by Phase‐Field Approach

We proposed a model to study the resistance degradation behavior of ferroelectric oxides in the presence of ferroelectric spontaneous polarization by combining the phase‐field model of ferroelectric domains and nonlinear diffusion equations for ionic/electronic transport. We took into account the nonperiodic boundary conditions for solving the electrochemical transport equations and Ginzburg–Landau equations using the Chebyshev collocation algorithm. We considered a single domain structure relative to a thin film BaTiO3 single crystal orientated to the normal of the electrode plates (Ni) in a single parallel plate capacitor configuration. The capacitor was subjected to a dc bias of 0.5 V either along the polarization direction or opposite to the polarization direction at 25°C. It is shown that the polarization bound charges at the metal/ferroelectric interface play an important role in charge carrier transport and leakage current evolution in BaTiO3 capacitor.

Hybrid simulation of whistler excitation by electron beams in two-dimensional non-periodic domains

We present a two-dimensional hybrid fluid-PIC scheme for the simulation of whistler wave excitation by relativistic electron beams. This scheme includes a number of features which are novel to simulations of this type, including non-periodic boundary conditions and fresh particle injection. Results from our model suggest that non-periodicity of the simulation domain results in the development of fundamentally different wave characteristics than are observed in periodic domains.

Determination of IVC breakpoint for josephson junction stack. Non-periodic boundary conditions with γ = 1

We prove that, in the case of non-periodic (with γ = 1) boundary conditions, the calculation of the current-voltage characteristic (IVC) for a stack of n intrinsic Josephson junctions reduces to solving a system of [(n + 1)/2] non-linear differential equations instead of the n original ones. The current voltage characteristic V(I) has the shape of a hysteresis loop. On the back branch of the loop V(I) decreases to zero rapidly near the breakpoint Ib. We succeeded to derive an algorithm determining the approximate breakpoint location and to improve simultaneously the mixed numerical-analytical algorithm of IVC calculation for a stack of Josephson junctions developed by us before. The efficiency of the improved algorithm is shown by the calculations of IVC for stacks consisting of various numbers of intrinsic Josephson junctions.

Exponential time differencing methods with Chebyshev collocation for polymers confined by interacting surfaces.

We present a fast and accurate numerical method for the self-consistent field theory calculations of confined polymer systems. It introduces an exponential time differencing method (ETDRK4) based on Chebyshev collocation, which exhibits fourth-order accuracy in temporal domain and spectral accuracy in spatial domain, to solve the modified diffusion equations. Similar to the approach proposed by Hur et al. [Macromolecules 45, 2905 (2012)], non-periodic boundary conditions are adopted to model the confining walls with or without preferential interactions with polymer species, avoiding the use of surface field terms and the mask technique in a conventional approach. The performance of ETDRK4 is examined in comparison with the operator splitting methods with either Fourier collocation or Chebyshev collocation. Numerical experiments show that our exponential time differencing method is more efficient than the operator splitting methods in high accuracy calculations. This method has been applied to diblock copolymers confined by two parallel flat surfaces.

General Model for an Entanglement-Enhanced Composed Quantum Game on a Two-Dimensional Lattice

We introduce a method of analyzing entanglement-enhanced quantum games on regular lattices of agents. Our method is valid for setups with periodic and non-periodic boundary conditions. To demonstrate our approach we study two different types games, namely the Prisoner's dilemma game and a cooperative Parrondo's game. In both cases we obtain results showing, that entanglement is a crucial resource necessary for the agents to achieve positive capital gain.

The Wave Equation in Non-classic Cases: Non-self Adjoint with Non-local and Non-periodic Boundary Conditions

In this paper the wave equation in some non-classic cases has been studied. In the rst case boundary conditions are non-local and non- periodic. At that case the associated spectral problem is a self-adjoint problem and consequently the eigenvalues are real. But in the second case the associated spectral problem is non-self-adjoint and consequently the eigenvalues are complex numbers, in which two cases, the solutions of the problem are constructed by the Fourier method. By compatibil- ity conditions and asymptotic expansions of the Fourier coecients, the convergence of series solutions are proved. Finally, series solutions are established and the uniqueness of the solution is proved by a special way which has not been used in classic texts.

Superconvergent error estimates for position-dependent smoothness-increasing accuracy-conserving (SIAC) post-processing of discontinuous Galerkin solutions

Superconvergence of discontinuous Galerkin methods is an area of increasing interest due to the ease with which higher order information can be extracted from the approximation. Cockburn, Luskin, Shu, and Suli showed that by applying a B-spline filter to the approximation at the final time, the order of accuracy can be improved from order k+1 to order 2k+1 in the L2-norm for linear hyperbolic equations with periodic boundary conditions (where k is the polynomial degree and h is the mesh element diameter) [Math. Comp. (2003)]. The applicability of this filter for linear hyperbolic problems with non-periodic boundary conditions was computationally extended and renamed a position-dependent smoothness-increasing accuracy-conserving (SIAC) filter by van Slingerland, Ryan, Vuik [SISC (2011)]. However, error estimates in the L2$-norm for this new position-dependent SIAC filter were never given. Furthermore, error estimates in the L-infinity-norm have not been established for the original kernel nor the position-dependent kernel. In this paper, for the first time we establish that it is possible to obtain order s, s=min{2k+1,2k + 2-\frac {d}{2}} accuracy in the L-infinity-norm for the position-dependent SIAC filter, where d is the dimension. Furthermore, we extend the error estimates given by Cockburn et al. so that they are applicable to the entire domain when implementing the position-dependent SIAC filter. We also computationally demonstrate the applicability of this filter for visualization of streamlines.

Combining the Fluctuating Charge Method, Non-Periodic Boundary Conditions and Meta-Dynamics: Aqua Ions as case studies.

We present the current status of development of our code for performing Molecular Dynamics (MD) simulations exploiting a polarizable force field based on the Fluctuating Charge (FQ) method and non-Periodic Boundary Conditions (NPBC). Continuing on the path set in a previous work, we increased the capabilities of the code by implementing a number of new features, including: a non-iterative algorithm for rigid trigonal molecule simulations; two additional temperature coupling schemes; a meta-dynamics based approach for effective free energy evaluations. Although these are well known algorithms, each present in one or more widely used MD packages, they have now been tested, for the first time, in the context of the FQ model coupled with NPBC. As case studies, we considered three aqueous ions of increasing charge, namely Na+, Ca2+ and La3+, at infinite dilution. In particular, by exploiting a computational approach recently proposed by our group and based on the metadynamics technique, we focused on the important role played by solvent polarization on ionic hydration structures, also investigating the free energy landscapes of ion coordination and the water exchange rates. Such an approach, previously tested with standard non-polarizable models, was applied here to evaluate the effects of explicit polarization on water exchange barriers between different solvent coordination structures. Moreover, we have analyzed and discussed in some detail non-linear electrostatic effects arising from solvent polarization while going from a mono- to a di- and trivalent ion.

Monomer-dimer problem on random planar honeycomb lattice

We consider the monomer-dimer (MD) problem on a random planar honeycomb lattice model, namely, the random multiple chain. This is a lattice system with non-periodic boundary condition, whose generating process is inspired by the growth of single walled zigzag carbon nanotubes. By applying algebraic and combinatorial techniques we establish a calculating expression of the MD partition function for bipartite graphs, which corresponds to the permanent of a matrix. Further, by using the transfer matrix argument we show that the computing problem of the permanent of high order matrix can be converted into some lower order matrices for this family of lattices, based on which we derive an explicit recurrence formula for evaluating the MD partition function of multiple chains and random multiple chains. Finally, we analyze the expectation of the number of monomer-dimer arrangements on a random multiple chain and the asymptotic behavior of the annealed MD entropy when the multiple chain becomes infinite in width and length, respectively.

A Local Pressure Boundary Condition Spectral Collocation Scheme for the Three-Dimensional Navier–Stokes Equations

A spectral collocation scheme for the three-dimensional incompressible $$({\varvec{u}},p)$$ formulation of the Navier–Stokes equations, in domains $$\varOmega $$ with a non-periodic boundary condition, is described. The key feature is the high order approximation, by means of a local Hermite interpolant, of a Neumann boundary condition for use in the numerical solution of the pressure Poisson system. The time updates of the velocity $${\varvec{u}}$$ and pressure $$p$$ are decoupled as a result of treating the pressure gradient in the momentum equation explicitly in time. The pressure update is computed from a pressure Poisson equation. Extension of the overall methodology to the Boussinesq system is also described. The uncoupling of the pressure and velocity time updates results in a highly efficient scheme that is simple to implement and well suited for simulating moderate to high Reynolds and Rayleigh number flows. Accuracy checks are presented, along with simulations of the lid-driven cavity flow and a differentially heated cavity flow, to demonstrate the scheme produces accurate three-dimensional results at a reasonable computational cost.

Numerical investigation on chaos in late boundary layer transition to turbulence

The mechanism of chaos in late boundary layer transition is a key issue of the laminar-turbulent transition process in a flat plate boundary layer. A careful study on the characteristic of chaos is carried out by high order direct numerical simulation (DNS). The process of flow chaos was originally considered as a result of large background noise and non-periodic spanwise boundary conditions. It was addressed that the large ring structures are affected by background noises first, and then the change of large ring structures affect the small scale vortices quickly, which directly lead to chaos and formation of turbulent flow. However, according to our DNS observation, the loss of symmetry starts from the middle level vortex rings while the top and bottom rings are still symmetric. The non-symmetric structure of second level vortex rings will influence the small scale vortices at the boundary layer bottom quickly. The loss of symmetry at the bottom of the boundary layer quickly spreads to upper level through ejections. This will lead to chaos of the whole flow field. Therefore, the internal instability of multiple level ring structures, especially the middle ring cycles, is the main reason for the process of flow chaos, but not the large background noise.

Dispersion Analysis of Complex Periodic Structures by Full-Wave Solution of Even-Odd-Mode Excitation Problems for Single Unit Cells

The Bloch modes of complex periodic structures are computed by superimposing results obtained from even-odd-mode full-wave driven simulations of individual unit cells. In order to emulate the periodic boundary conditions for the Bloch modes, one full-wave simulation is performed with magnetic and one with electric boundary conditions yielding the even- and odd-mode results, respectively. Therefore, the employed full-wave solver does not even need to support periodic boundary conditions. Since the non-periodic boundary conditions can be arbitrary, even open and radiating problems can be analyzed. The excitation of the structures is performed by discrete ports, which are located appropriately in order to excite the desired mode, typically the fundamental mode of the background structure. The found even- and odd-mode impedances deliver directly the complex propagation constant and the Bloch impedance of the periodic structure without any iterative search and under consideration of all electromagnetic interaction effects. This is in contrast to alternative solution methods, which require to solve computationally intensive eigenproblems or many excitation problems. Moreover, this approach allows to determine the field data corresponding to the Bloch mode, where arbitrary real and complex modes can be handled. Numerical results for different classes of problems including 1D and 2D microstrip structures as well as hollow waveguide based unit cells are presented.

Calculating the binding free energies of charged species based on explicit-solvent simulations employing lattice-sum methods: an accurate correction scheme for electrostatic finite-size effects.

The calculation of a protein-ligand binding free energy based on molecular dynamics (MD) simulations generally relies on a thermodynamic cycle in which the ligand is alchemically inserted into the system, both in the solvated protein and free in solution. The corresponding ligand-insertion free energies are typically calculated in nanoscale computational boxes simulated under periodic boundary conditions and considering electrostatic interactions defined by a periodic lattice-sum. This is distinct from the ideal bulk situation of a system of macroscopic size simulated under non-periodic boundary conditions with Coulombic electrostatic interactions. This discrepancy results in finite-size effects, which affect primarily the charging component of the insertion free energy, are dependent on the box size, and can be large when the ligand bears a net charge, especially if the protein is charged as well. This article investigates finite-size effects on calculated charging free energies using as a test case the binding of the ligand 2-amino-5-methylthiazole (net charge +1 e) to a mutant form of yeast cytochrome c peroxidase in water. Considering different charge isoforms of the protein (net charges -5, 0, +3, or +9 e), either in the absence or the presence of neutralizing counter-ions, and sizes of the cubic computational box (edges ranging from 7.42 to 11.02 nm), the potentially large magnitude of finite-size effects on the raw charging free energies (up to 17.1 kJ mol(-1)) is demonstrated. Two correction schemes are then proposed to eliminate these effects, a numerical and an analytical one. Both schemes are based on a continuum-electrostatics analysis and require performing Poisson-Boltzmann (PB) calculations on the protein-ligand system. While the numerical scheme requires PB calculations under both non-periodic and periodic boundary conditions, the latter at the box size considered in the MD simulations, the analytical scheme only requires three non-periodic PB calculations for a given system, its dependence on the box size being analytical. The latter scheme also provides insight into the physical origin of the finite-size effects. These two schemes also encompass a correction for discrete solvent effects that persists even in the limit of infinite box sizes. Application of either scheme essentially eliminates the size dependence of the corrected charging free energies (maximal deviation of 1.5 kJ mol(-1)). Because it is simple to apply, the analytical correction scheme offers a general solution to the problem of finite-size effects in free-energy calculations involving charged solutes, as encountered in calculations concerning, e.g., protein-ligand binding, biomolecular association, residue mutation, pKa and redox potential estimation, substrate transformation, solvation, and solvent-solvent partitioning.

Reverse Monte Carlo study of spherical sample under non-periodic boundary conditions: the structure of Ru nanoparticles based on x-ray diffraction data

A new method to fit experimental diffraction data with non-periodic structure models for spherical particles was implemented in the reverse Monte Carlo simulation code. The method was tested on x-ray diffraction data for ruthenium (Ru) nanoparticles approximately 5.6 nm in diameter. It was found that the atomic ordering in the ruthenium nanoparticles is quite distorted, barely resembling the hexagonal structure of bulk Ru. The average coordination number for the bulk decreased from 12 to 11.25. A similar lack of structural order has been observed with other nanoparticles (e.g. Petkov et al 2008 J. Phys. Chem. C 112 8907–11) indicating that atomic disorder is a widespread feature of nanoparticles less than 10 nm in diameter.

Effect of number of grains and boundary conditions on digital material representation deformation under plane strain

The main goal of this paper is to investigate the digital material representation of a single-phase polycrystalline unit cell. Particular attention is put on the amount of grains in the microstructure, which can be considered as the representative volume element of a sample subjected to plastic deformation conditions. Additionally, the influences of the periodic and non-periodic boundary conditions on deformation behavior of the unit cells are evaluated. Possibility of the application of periodic boundary conditions on the non-periodic unit cells using a buffer zone approach is also discussed. Obtained results of equivalent strains and reaction force on case studies are presented.

A study on boundary force model used in multiscale simulations with non-periodic boundary condition

This paper presents an investigation of the non-periodic boundary condition (NPBC) which is often used in multiscale atomistic–continuum simulations. The relationship between the boundary force exerted by the imaginary atoms outside the atomistic domain and the fluid state parameters including density and temperature at the boundary is studied. A fitting formula of the boundary force as a function of the fluid state has been proposed based on the relationship. The accuracy of the fitting formula is verified by the equilibrium molecular dynamics (MD) simulations. Poiseuille flow with viscous dissipation and unsteady heat transfer between two walls is then simulated using the proposed fitting formula. The elimination of density oscillation near the boundary of atomistic region and good agreement of velocity and temperature evolutions with time from pure MD and the multiscale simulations adopting NPBC further confirm the correctness of our fitting formula.

Graphics processing unit accelerated one‐dimensional blood flow computation in the human arterial tree

One-dimensional blood flow models have been used extensively for computing pressure and flow waveforms in the human arterial circulation. We propose an improved numerical implementation based on a graphics processing unit (GPU) for the acceleration of the execution time of one-dimensional model. A novel parallel hybrid CPU-GPU algorithm with compact copy operations (PHCGCC) and a parallel GPU only (PGO) algorithm are developed, which are compared against previously introduced PHCG versions, a single-threaded CPU only algorithm and a multi-threaded CPU only algorithm. Different second-order numerical schemes (Lax-Wendroff and Taylor series) are evaluated for the numerical solution of one-dimensional model, and the computational setups include physiologically motivated non-periodic (Windkessel) and periodic boundary conditions (BC) (structured tree) and elastic and viscoelastic wall laws. Both the PHCGCC and the PGO implementations improved the execution time significantly. The speed-up values over the single-threaded CPU only implementation range from 5.26 to 8.10 × , whereas the speed-up values over the multi-threaded CPU only implementation range from 1.84 to 4.02 × . The PHCGCC algorithm performs best for an elastic wall law with non-periodic BC and for viscoelastic wall laws, whereas the PGO algorithm performs best for an elastic wall law with periodic BC.