Dear all,
I recently read a paper (Sivak A. David et al., The journal of physical chemistry, 2014) in which it is pointed out the importance of using integrator splitting procedures when dealing with Langevin dynamics in a MD simulation.
I am actually using LAMMPS to do exactly this (system of colloidal particles in a Langevin dynamics context) and i followed the iFLD (inertia-Fast Lubrication Dynamics) scheme proposed here (Bolintineanu S. Dan et al., Comp. Part. Mech., 2014) to build up my simulation setup.
My question is: in light of what is underlined in the paper of Sivak et al., is it the iFLD scheme consistent in reproducing the Langevin dynamics ? Is there a fix in lammps which implements some integrator splitting procedures as the ones describes in Sivak et al. ?
Cheers,
Gianluca
Hi @Gianluca,
The langevin implementation allows for the use of the Grønbech-Jensen-Farago (gjf) integrator. Please see the manual page as well as the reference article for more details.
As a side note, a recent systematic review from Grønbech-Jensen states that the Sivak timestep splitting procedure is a special case for the gjf procedure (see what is stated about reference 36). I can´t tell how this relates to the iFLD procedure though.
I hope this helps answering your questions.
Hi Germain,
thanks for the references, it helps !
By contrast, I am a bit puzzled on how to couple the langevin fix with the lubricate pair style (that I need to get the hydrodynamic mobility tensor with the lubrication contribution).
Since the langevin fix also adds the viscous damp (that is already taken into account in the lubricate pair_style) it kind of seems to me that a straightforward coupling of the langevin fix and the lubricate pair_style “double” the viscous damping contribution.
Maybe, I simply modify the lagevin fix suppressing the viscous damp contribution already in the lubricate pair_style ? Would it to do the job ?
thanks again,
Gianluca
I don’t think you need to use fix langevin at all. Mind you, it is only a thermostat, it does not implement the basic time integration. That would be done by fix nve or equivalent. Please see this note in the pair style lubricate docs:
When using the FLD terms, these pair styles are designed to be used with explicit time integration and a correspondingly small timestep. Thus either fix nve/sphere or fix nve/asphere should be used for time integration. To perform implicit FLD, see the pair_style lubricateU command.
and the text following it.
Hi Alex,
thanks for the note.
Yes, I am currently using the classical fix nve/sphere integrator but I see a bizarre behavior for a simple setup. Basically, I have a simple pair of particles and I let the system evolve until the pair of particles forms a stable contact (I made some modifications in order to allow for particles’ contact and then use the granular pair_style). I keep track of the time evolution of the interparticle distance (h vs t) and, starting from a relatively large timestep, I reduce it. The point is that I don’t see any convergence of the h vs t curves reducing the timestep…
I was, then, wondering if the velocity-verlet algorithm used in the fix nve/sphere is a proper choice for the time integration of the Langevin dynamics (and that’s how I fall on the Sivak et al. paper).
Gianluca
I am sorry, but this is far outside of my area of expertise and beyond that not really a LAMMPS issue, but an issue of the science of your model. I suggest you try to contact the authors that contributed the pair style(s) or the authors of the publications that were referencing.
Hi @Gianluca the Sivak paper (as well as many others, including e.g. here and here) deals with differences in how the stochastic forces are discretized. This is a different issue from using or not the velocity-Verlet algorithm, which is assumed the best second-order scheme for the special case without stochastic forces.
Thankfully, GJF and all schemes equivalent to it yield the same configurational sampling, which is the same fact already pointed out by @Germain . There may be some differences in the effective diffusivity (time scales), but not in the position ensemble.
Separately from positions, there is the issue of how to define instantaneous velocities in a way that is numerically stable in the high-friction limit. This is discussed in this paper, and implemented with the vhalf option of fix langevin. I and colleagues have experience this issue ourselves when looking at stress tensor components (which include a velocity term).
I have no experience with applying these schemes on granular bodies or iFLD, but hopefully with the help of others you can figure out if the same implementation applies to your setup. See also this previous discussion.
Giacomo
1 Like
For the sake of completeness, one should also mention fix mvv/dpd command — LAMMPS documentation which allows to tweak which velocity is made available during the velocity verlet timestepping at the point of the force computation.
2 Likes
Thanks Giacomo for the clarifications and thanks all for the suggestions !
All of great help !
Cheers,
Gianluca