Dear Lammps-community,
I would like to apply a uniform shear with a constant velocity in my simulations (in xz) while controlling the pressure in xx, yy and zz. For this purpose I use the following combination of fix_npt and fix_deform
fix 1 all npt temp 0.9 0.9 1 aniso 0 0 10 fixedpoint 0 0 0 scalexz no scaleyz no scalexy no
fix 2 all deform 1 xz vel 0.1 remap v
What I expect to get is a constant increase in the xz tilt factor. However, this is only the case in the time interval up to the first cell flip. As soon as the cell is flipped for the first time, the increase in xz starts to oscillate slightly. The amplitude of these oscillations increases with an increasing number of cell flips. Setting “flip no” in both fix-npt and fix-deform results in the expected behavior. However, as I would like to perform long simulations, this appears to be no possibility due to the extreme tilt value. As soon as flips are allowed (either in npt or deform) the problem appears. Using no barostat causes no problems.
It seems that the flip operation of one of the fixes is not compatible with the other fix.
I attached a simple example that reproduces this behavior. The .png shows the aforementioned oscillations in xz over time starting from a first flip. The version I used is LAMMPS (16 Mar 2018).
Please apologize in case I overlooked something and/or my reasoning is wrong.
Has anyone any suggestions to resolve the issue?
Many thanks in advance for any help/comment/reply.
Best wishes
Thomas Reichenbach
Thomas Reichenbach
Multiscale Modelling and Tribosimulations
Fraunhofer Institute for Mechanics of Materials IWM | www.iwm.fraunhofer.de
struct.lmp (210 KB)
in.lmp (600 Bytes)
What happens if you just run the fix deform xz, without the fix npt command (maybe just fix nvt)?
Does the xz tilt factor then increase monotonically and proceed thru many flips?
I’m thinking the noise in xz is due to the fact that NPT is continuously adjusting the x and z box
lengths (independently) hence the instantaneous xz tilt is non-monotonic. What you really care about are
the box boundaries for the deforming box. Are those more smooth?
Steve
I took a look at this. It turns out that nothing funny is going on. If you define the total box displacement to be xz+nx*lx, you recover the expected behavior delta_xz = 0.03, no fluctuations. In this expression, nx is the number of box flips and lx is the instantaneous value of the box length.
Aidan
Thanks for your help Steve and Aidan!
I see your point and agree that everything behaves as expected.
The problem I am not sure how to deal with is that in long runs, after many flips, the potential energy becomes unstable.
I suspect that the origin is that xz of the box becomes, with increasing number of flips nx, more and more sensitive towards changes in the instantaneous box length lx caused by the barostat in order to linearly increase " xz + nx * lx ".
A different damping-constant in the barostat delays the occurrence of the oscillations but sooner or later the potential energy becomes instable.
For the example I provided in my first mail the step vs energy curve behaves as shown in the attached plot before the system finally explodes.
In the end, my setup should correspond to Lees-Edwards boundary conditions with constant pressure, but I am not certain whether my approach makes sense at all....
Best wishes and thanks once again for your help
Thomas
Thanks for your help Steve and Aidan!
I see your point and agree that everything behaves as expected.
The problem I am not sure how to deal with is that in long runs, after many flips, the potential energy becomes unstable.
I suspect that the origin is that xz of the box becomes, with increasing number of flips nx, more and more sensitive towards changes in the instantaneous box length lx caused by the barostat in order to linearly increase " xz + nx * lx ".
A different damping-constant in the barostat delays the occurrence of the oscillations but sooner or later the potential energy becomes instable.
For the example I provided in my first mail the step vs energy curve behaves as shown in the attached plot before the system finally explodes.
please have a look at the image flags of your atoms. if you have a long flow simulation, they may be overflowing, since those are by default only 10-bit integers, i.e range from -512 to 511.
you’ll have to compile LAMMPS with -DLAMMPS_BIGBIG to increase those to 21-bits and thus have a range of about -1M to +1M.
there are more details about the impact of -DLAMMPS_SMALLSMALL, -DLAMMPS_SMALLBIG (the default) and -DLAMMPS_BIGBIG in the manual.
axel.
I’m still not quite clear on the problem here.
Is the issue perhaps this: After many flips and large changes
to the x box length (due to NPT applied in that dim), that using the
instantaneous box length to compute a new tilt value of xz
causes a large delta (versus the previous step xz) when a flip happens?
This would be a computation that fix deform makes to reset
the tilt factor every step. I.e. in the formula xz
- nx*lx, the effect of a small change in lx (between 2 flips) is amplified
as nx becomes large.
Maybe fix deform should not update the tilt factor that way when
lx is a variable. Instead it could continually increment it. I think
we avoided doing that, thinking there could be round-off issues
over millions of steps. But if we are computing it every step
as if lx had its current value since time 0, maybe that is the issue here?
Steve
Yes I think that is the issue here is the use of the instantaneous lx. However, a large delta in the tilt occurs not only when a flip happens but every step, provided many flips occurred previously and the barostat changes lx, in order to satisfy the linear increase of xz + nx*lx.
Moreover, I am confident that, as you suggested, continually incrementing the current xz likely solves the issue such that not "xz + nx*lx" but rather something like "xz[timestep_i] - xz[timestep_(i-1)]" is controlled.
Is there any way to test this using the current implementation or would a modification of the source code be necessary?
@Axel: Thanks for your comment. So far I have not paid any attention to image flags, for further tests I will use the the lammps version with -DLAMMPS_BIGBIG. Unfortunately , in my test runs the pe explosions prevail nonetheless.
Thanks!
Thomas
This would require a code change in fix deform.
Perhaps the code should detect when lx is
not static and do the xz update (each step)
differently. I’ll take a look …
For dynamic lx, it seems like some of the
fix deform options for the requested rate or amount
of final deformation will then be ambiguous.
I.e. what the user requests and what amount
of deformation is produced will be different?
Steve