Coupled Navier–Stokes/molecular dynamics simulations in nonperiodic domains based on particle forcing

This paper focuses on coupling methods for hybrid Navier–Stokes/molecular dynamics (MD) simulations. The computational domain is split in a continuum flow region, where a finite‐volume discretisation of the Navier–Stokes equations is used, and one or more particle domains, where molecular level modelling of the flow is employed. The domains are defined with a partial overlap, in which the flow states are coupled through an exchange of the velocity components. For the steady flows considered, an under‐relaxed Newton iteration method is used to drive the coupled system to convergence. The main focus of the present work is on methods to impose nonperiodic boundary conditions on the particle domain(s). A particle forcing is applied in the direction normal to the particle domain boundary to impose the boundary normal velocity component. A novel aspect of the present work is the extension of this method to more general nonplanar particle domain boundaries. The main contribution of the paper is the development of a particle forcing method in the direction tangential to the domain boundary, which is based on the equivalent continuum‐flow boundary shear stresses along with an iterative forcing strength adjustment based on the extrapolated particle boundary velocity. Furthermore, an adaptation scheme is presented, which uses the finite‐volume flux residuals of the particle bin averaged velocity field as a truncation criterion for the iterative force‐update scheme. It is demonstrated that by comparing the residual reduction for the momentum equation in the nonhomogeneous directions during the molecular dynamics simulations with that for a homogeneous direction, the forcing iteration at which the statistical noise in the velocity field dominates the uncertainty in the forcing strength can be determined. At this point the iteration can be truncated. It is shown that with adaptive schemes of this type, the total number of MD evaluations required in a coupled Navier–Stokes/MD simulation can be reduced relative to a hybrid scheme with a fixed number of forcing‐strength updates. Copyright © 2011 John Wiley & Sons, Ltd.     

An atomistic–continuum hybrid scheme for numerical simulation of droplet spreading on a solid surface

We present an atomistic–continuum hybrid method to investigate spreading dynamics of drops on solid surfaces. The Navier–Stokes equations are solved by the finite-volume method in a continuum domain comprised of the main body of the drop, and atomistic molecular dynamics simulations are used in a particle domain in the vicinity of the contact line. The spatial coupling between the continuum and particle domains is achieved through constrained dynamics of flux continuities in an overlap domain.     

Granular vortices: Identification,characterization and conditions for the localization of deformation

We relate the micromechanics of vortex evolution to that of force chain buckling and, on this basis, formulate the conditions for strain localization in a continuum model of dense granular media. Using the traditional bifurcation analysis of shear bands, we show that kinematic vortex fields are in fact solutions to the boundary value problem satisfying null boundary conditions. To establish an empirical basis for our study, we first develop a method to identify the location of the core and boundary of each vortex from a given displacement field in two dimensions. We then employ this method to characterize the residual deformation field (i.e., the deviation of particle motions from the continuum deformation) in a physical experiment and a discrete element simulation of dense granular samples submitted to biaxial compression. Vortices in the failure regime are essentially confined to the shear band. Primary vortices, the clear majority, rotate in the same direction as the shear band; secondary vortices, the so-called wakes, rotate in the opposite direction. Primary vortices align in spatial succession along the central axis of the band; wakes form next to the band boundaries, in between and beside two adjacent primary vortices. Force chain buckling, the governing mechanism for shear bands, is responsible for vortex formation in the failure regime. Vortex dynamics are consistent with stick-slip dynamics. From quiescent conditions of jamming or stick, vortical motions arise from force chain buckling and associated relative particle rotations and sliding; these in turn precipitate intermittent periods of unjamming or slip, evident in the attendant drops in stress ratio and bursts in both kinetic energy and local nonaffine deformation. A kinematic vortex field inside shear bands is proposed that is consistent with the equations of continuum mechanics and the underlying instability of force chain buckling: such a field is periodic with a repeating unit cell comprising a primary vortex at the center of the band, with two trailing wakes close next to the band boundaries.     

A Lagrangian-based continuum homogenization approach applicable to molecular dynamics simulations

The continuum notions of effective mechanical quantities as well as the conditions that give meaningful deformation processes for homogenization problems with large deformations are reviewed. A continuum homogenization model is presented and recast as a Lagrangian-based approach for heterogeneous media that allows for an extension to discrete systems simulated via molecular dynamics (MD). A novel constitutive relation for the effective stress is derived so that the proposed Lagrangian-based approach can be used for the determination of the “stress–deformation” behavior of particle systems. The paper is concluded with a careful comparison between the proposed method and the Parrinello–Rahman approach to the determination of the “stress–deformation” behavior for MD systems. When compared with the Parrinello–Rahman method, the proposed approach clearly delineates under what conditions the Parrinello–Rahman scheme is valid.     

An atomistically validated continuum model for strain relaxation and misfit dislocation formation

In this paper, molecular dynamics (MD) calculations have been used to examine the physics behind continuum models of misfit dislocation formation and to assess the limitations and consequences of approximations made within these models. Without compromising the physics of misfit dislocations below a surface, our MD calculations consider arrays of dislocation dipoles constituting a mirror imaged “surface”. This allows use of periodic boundary conditions to create a direct correspondence between atomistic and continuum representations of dislocations, which would be difficult to achieve with free surfaces. Additionally, by using long-time averages of system properties, we have essentially reduced the errors of atomistic simulations of large systems to “zero”. This enables us to deterministically compare atomistic and continuum calculations. Our work results in a robust approach that uses atomistic simulation to accurately calculate dislocation core radius and energy without the continuum boundary conditions typically assumed in the past, and the novel insight that continuum misfit dislocation models can be inaccurate when incorrect definitions of dislocation spacing and Burgers vector in lattice-mismatched systems are used. We show that when these insights are properly incorporated into the continuum model, the resulting energy density expression of the lattice-mismatched systems is essentially indistinguishable from the MD results.     

Stability and vibration of a helical rod constrained by a cylinder

The stability and vibration of an elastic rod with a circular cross section under the constraint of a cylinder is discussed. The differential equations of dynamics of the constrained rod are established with Euler’s angles as variables describing the attitude of the cross section. The existence conditions of helical equilibrium under constraint are discussed as a special configuration of the rod. The stability of the helical equilibrium is discussed in the realms of statics and dynamics, respectively. Necessary conditions for the stability of helical rod are derived in space domain and time domain, and the difference and relationship between Lyapunov’s and Euler’s stability concepts are discussed. The free frequency of flexural vibration of the helical rod with cylinder constraint is obtained in analytical form. The project supported by the National Natural Science Foundation of China (10472067). The English text was polished by Yunming Chen.     

Navier–Stokes Equations with Shear Thinning Viscosity. Regularity up to the Boundary

In this paper we consider a class of stationary Navier–Stokes equations with shear dependent viscosity, in the shear thinning case p < 2, under a non-slip boundary condition. We are interested in global (i.e., up to the boundary) regularity results, in dimension n = 3, for the second order derivatives of the velocity and the first order derivatives of the pressure. As far as we know, there are no previous global regularity results for the second order derivatives of the solution to the above boundary value problem. We consider a cubic domain and impose the non-slip boundary condition only on two opposite faces. On the other faces we assume periodicity, as a device to avoid effective boundary conditions. This choice is made so that we work in a bounded domain Ω and simultaneously with a flat boundary. The extension to non-flat boundaries is done in the forthcoming paper [7], by following ideas introduced by the author, for the case p > 2, in reference [5]. The results also hold in the presence of the classical convective term, provided that p is sufficiently close to the value 2.      

Molecular-dynamics simulation-based cohesive zone representation of intergranular fracture processes in aluminum

A traction-displacement relationship that may be embedded into a cohesive zone model for microscale problems of intergranular fracture is extracted from atomistic molecular-dynamics (MD) simulations. An MD model for crack propagation under steady-state conditions is developed to analyze intergranular fracture along a flat Σ99 [1 1 0] symmetric tilt grain boundary in aluminum. Under hydrostatic tensile load, the simulation reveals asymmetric crack propagation in the two opposite directions along the grain boundary. In one direction, the crack propagates in a brittle manner by cleavage with very little or no dislocation emission, and in the other direction, the propagation is ductile through the mechanism of deformation twinning. This behavior is consistent with the Rice criterion for cleavage vs. dislocation blunting transition at the crack tip. The preference for twinning to dislocation slip is in agreement with the predictions of the Tadmor and Hai criterion. A comparison with finite element calculations shows that while the stress field around the brittle crack tip follows the expected elastic solution for the given boundary conditions of the model, the stress field around the twinning crack tip has a strong plastic contribution. Through the definition of a Cohesive-Zone-Volume-Element—an atomistic analog to a continuum cohesive zone model element—the results from the MD simulation are recast to obtain an average continuum traction-displacement relationship to represent cohesive zone interaction along a characteristic length of the grain boundary interface for the cases of ductile and brittle decohesion.     

A gradient theory of single-crystal viscoplasticity that accounts for geometrically necessary dislocations

This study develops a gradient theory of single-crystal plasticity that accounts for geometrically necessary dislocations. The theory is based on classical crystalline kinematics; classical macroforces; microforces for each slip system consistent with a microforce balance; a mechanical version of the second law that includes, via the microforces, work performed during slip; a rate-independent constitutive theory that includes dependences on a tensorial measure of geometrically necessary dislocations. The microforce balances are shown to be equivalent to nonlocal yield conditions for the individual slip systems. The field equations consist of the yield conditions coupled to the standard macroscopic force balance; these are supplemented by classical macroscopic boundary conditions in conjunction with nonstandard boundary conditions associated with slip. As an aid to solution, a weak (virtual power) formulation of the nonlocal yield conditions is derived. To make contact with classical dislocation theory, the microstresses are shown to represent counterparts of the Peach-Koehler force on a single dislocation.     

Feasibility study of symmetric solution of molecular argon flow inside microscale nozzles

The computational cost of numerical methods in microscopic-scales such as molecular dynamics (MD) is a deterrent factor that limits simulations with a large number of particles. Hence, it is desirable to decrease the computational cost and run time of simulations, especially for problems with a symmetrical domain. However, in microscopic-scales, implementation of symmetric boundary conditions is not straight-forward. Previously, the present authors have successfully used a symmetry boundary condition to solve molecular flows in constant-area channels. The results obtained with this approach agree well with the benchmark cases. Therefore, it has provided us with a sound ground to further explore feasibility of applying symmetric solutions of micro-fluid flows in other geometries such as variable-area ducts. Molecular flows are solved for the whole domain with and without the symmetric boundary condition. Good agreement has been reached between the results of the symmetric solution and the whole domain solution. To investigate robustness of the proposed method, simulations are conducted for different values of affecting parameters including an external force, a flow density, and a domain length. The results indicate that the symmetric solution is also applicable to variable-area ducts such as micro-nozzles.     

Study on micronozzle flow and propulsion performance using DSMC and continuum methods

In this paper, both DSMC and Navier–Stokes computational approaches were applied to study micronozzle flow. The effects of inlet condition, wall boundary condition, Reynolds number, micronozzle geometry and Knudsen number on the micronozzle flow field and propulsion performance were studied in detail. It is found that within the Knudsen number range under consideration, both the methods work to predict flow characteristics inside micronozzles. The continuum method with slip boundary conditions has shown good performance in simulating the formation of a boundary layer inside the nozzle. However, in the nozzle exit lip region, the DSMC method is better due to gas rapid expansion. It is found that with decreasing the inlet pressure, the difference between the continuum model and DSMC results increases due to the enhanced rarefaction effect. The coefficient of discharge and the thrust efficiency increase with increasing the Reynolds number. Thrust is almost proportional to the nozzle width. With dimension enlarged, the nozzle performance becomes better while the rarefaction effects would be somewhat weakened.The project supported by the National Natural Science Foundation of China (10372099). The English text was polished by Boyi Wang     

离散颗粒集合体-Cosserat连续体模型的连接尺度方法

基于针对分子动力学-Cauchy连续体模型提出的连接尺度方法(BSM)[1,2],发展了耦合细尺度上基于离散颗粒集合体模型的离散单元法(DEM)和粗尺度上基于Cosserat连续体模型的有限元法(FEM)的BSM。仅在有限局部区域内采用DEM以从细观层次模拟非连续破坏现象,而在全域则采用花费计算时间和存储空间较少的FEM。通过连接尺度位移(包括平移和转动)分解,和基于作用于Cosserat连续体有限元节点和颗粒集合体颗粒形心的离散系统虚功原理,得到了具有解耦特征的粗细尺度耦合系统运动方程。讨论和提出了在准静态载荷条件下粗细尺度域的界面条件,以及动态载荷条件下可以有效消除粗细尺度域界面上虚假反射波的非反射界面条件(NRBC)。本文二维数值算例结果说明了所提出的颗粒材料BSM的可应用性和优越性,及所实施界面条件对模拟颗粒材料动力学响应的有效性。     

A new model for dissipative particle dynamics boundary condition of walls with different wettabilities

The implementation of solid-fluid boundary condition has been a major challenge for dissipative particle dynamics(DPD)method.Current implementations of boundary conditions usually try to approach a uniform density distribution and a velocity profile close to analytical solution.The density oscillations and slip velocity are intentionally eliminated,and different wall properties disappear in the same analytical solution.This paper develops a new wall model that combines image and frozen particles and a new strategy to emphasize different wall properties especially wettabilities.The strategy first studies the realistic wall-fluid system by molecular dynamics(MD)simulation depending on physical parameters.Then,a DPD simulation is used to match the density and velocity profiles with the new wall model.The obtained DPD parameters can simulate the systems with the same wall and fluid materials.With this method,a simulation of the Poiseuille flow of liquid argon with copper walls is presented.Other walls with super-hydrophilic,hydrophilic,and hydrophobic wettabilities are also simulated.The limitations of the analytical solution and the effect of the wall-fluid interaction are discussed.The results show that the method suggested in this paper can simulate the mesoscale behavior of the microchannel flow related to realistic systems.     

Nanoscale Poiseuille flow and effects of modified Lennard–Jones potential function

Numerical simulation of Poiseuille flow of liquid Argon in a nanochannel using the non-equilibrium molecular dynamics simulation (NEMD) is performed. The nanochannel is a three-dimensional rectangular prism geometry where the concerned numbers of Argon atoms are 2,700, 2,550 and 2,400 at 102, 108 and 120 K. Poiseuille flow is simulated by embedding the fluid particles in a uniform force field. An external driving force, ranging from 1 to 11 PN (Pico Newton), is applied along the flow direction to inlet fluid particles during the simulation. To obtain a more uniform temperature distribution across the channel, local thermostating near the wall are used. Also, the effect of other mixing rules (Lorenthz–Berthelot and Waldman–Kugler rules) on the interface structure are examined by comparing the density profiles near the liquid/solid interfaces for wall temperatures 108 and 133 K for an external force of 7 PN. Using Kong and Waldman–Kugler rules, the molecules near the solid walls were more randomly distributed compared to Lorenthz–Berthelot rule. These mean that the attraction between solid–fluid atoms was weakened by using Kong rule and Waldman–Kugler rule rather than the Lorenthz–Berthelot rule. Also, results show that the mean axial velocity has symmetrical distribution near the channel centerline and an increase in external driving force can increase maximum and average velocity values of fluid. Furthermore, the slip length and slip velocity are functions of the driving forces and they show an arising trend with an increase in inlet driving force and no slip boundary condition is satisfied at very low external force (<1 PN).     

Molecular dynamics simulation of liquid argon flow at platinum surfaces

The micro Poiseuille flow for liquid argon flowing in a nanoscale channel formed by two solid walls was studied in the present paper. The solid wall material was selected as platinum, which has well established interaction potential. We consider the intermolecular force not only among the liquid argon molecules, but also between the liquid argon atoms and the solid wall particles, therefore three regions, i.e. the liquid argon computation domain, the top and bottom solid wall regions are included for the force interaction. The present MD (Molecular Dynamics) simulation was performed without any assumptions at the wall surface. The objective of the study is to find how the flow and the slip boundaries at the wall surface are affected by the applied gravity force, or the shear rate. The MD simulations are performed in a nondimensional unit system, with the periodic boundary conditions applied except in the channel height direction. Once the steady state is reached, the macroscopic parameters are evaluated using the statistical mechanics approach. For all the cases tested numerically in the present paper, slip boundaries occur, and such slip velocity at the stationary wall surface increases with increasing the applied gravity force, or the shear rate. The slip length, which is defined as the distance that the liquid particles shall travel beyond the wall surfaces to reach the same velocity as the wall surface, sharply decreases at small shear rate, then slightly decreases with increasing the applied shear rate. We observe that the liquid viscosity remains nearly constant at small shear rates, and the Newtonian flow occurs. However, with increasing the shear rate, the viscosity increases and the non-Newtonian flow appears.     

Vorticity and Regularity for Viscous Incompressible Flows under the Dirichlet Boundary Condition. Results and Related Open Problems

In reference [7] it is proved that the solution of the evolution Navier–Stokes equations in the whole of R ³ must be smooth if the direction of the vorticity is Lipschitz continuous with respect to the space variables. In reference [5] the authors improve the above result by showing that Lipschitz continuity may be replaced by 1/2-H?lder continuity. A central point in the proofs is to estimate the integral of the term (ω · ∇)u · ω, where u is the velocity and ω = ∇ × u is the vorticity. In reference [4] we extend the main estimates on the above integral term to solutions under the slip boundary condition in the half-space R ₊³. This allows an immediate extension to this problem of the 1/2-H?lder sufficient condition. The aim of these notes is to show that under the non-slip boundary condition the above integral term may be estimated as well in a similar, even simpler, way. Nevertheless, without further hypotheses, we are not able now to extend to the non slip (or adherence) boundary condition the 1/2-H?lder sufficient condition. This is not due to the “nonlinear" term (ω · ∇)u · ω but to a boundary integral which is due to the combination of viscosity and adherence to the boundary. On the other hand, by appealing to the properties of Green functions, we are able to consider here a regular, arbitrary open set Ω.      

Molecular dynamics simulation of peeling a DNA molecule on substrate

Molecular dynamics (MD) simulations are performed to study adhesion and peeling of a short fragment of single strand DNA (ssDNA) molecule from a graphite surface. The critical peel-off force is found to depend on both the peeling angle and the elasticity of ssDNA. For the short ssDNA strand under investigation, we show that the simulation results can be explained by a continuum model of an adhesive elastic band on substrate. The analysis suggests that it is often the peak value, rather than the mean value, of adhesion energy which determines the peeling of a nanoscale material.The project supported by the Distinguished Young Scholar Fund of NSFC (10225209) and key project from the Chinese Academy of Sciences (KJCX2-SW-L2)