[lammps-users] question on a how to apply a force!

[_nicola_varini]

Dear Axel and Steve, thanks for the discussion about how to apply a force and thanks to Axel for gimme a new fix.

In fix_spring_pull I found a bug due to this fact

xc += xv * update->dt;
yc += yv * update->dt;
zc += zv * update->dt;
domain->minimum_image(xc,yc,zc);

dx = xcm[0] - xc;
dy = xcm[1] - yc;
dz = xcm[2] - xc;
domain->minimum_image(dx,dy,dz);

because, due to minimum_image procedure, it’s possible that xc<xcm[0] and this is wrong. So I found that instead do this:
xc += xv * update->dt;
yc += yv * update->dt;
zc += zv * update->dt;
//domain->minimum_image(xc,yc,zc);
check = xcm[0] - xcm_old[0] + boxlength/2.0;
if(check < 0.0) counter+=1;
check = xcm[0] - xcm_old[0] - boxlength/2.0;
if(check > 0.0) counter-=1;

dx = xcm[0] + counterboxlength - xc;
dy = xcm[1] + counterboxlength - yc;
dz = xcm[2] + counter*boxlength - xc;

work properly. Of course this work only if the pulling direction is along x-axis but it’s straight to extend to other directions.

Secondly I found some difference in the medium force between rigid and rigid/nve. I attach the some files.
-in is the lammps input script
-out_lammps_rigid is the output of lammps using rigid in the input script
-out_lammps_rigid is the output of lammps using rigid/nve in the input script
-out_MD_code is the ouput of a reference serial MD code

in lammps output the force is in the column 7 while in the other output is in 2 column

the only difference between the two code is the fact that in lammps the langevin thermostat is applied by using a uniform distribution while in the other code is applied by using a gaussian distribution.

Do you think that fix spring pull style can be useful also for other users?

Thanks for discussion.

N.

2010/11/30 nicola varini <nicola.varini@…24…>

lammps_files.tar.gz (402 KB)

[akohlmey]

Dear Axel and Steve, thanks for the discussion about how to apply a force
and thanks to Axel for gimme a new fix.

In fix_spring_pull I found a bug due to this fact

xc += xv * update->dt;
yc += yv * update->dt;
zc += zv * update->dt;
domain->minimum_image(xc,yc,zc);

dx = xcm[0] - xc;
dy = xcm[1] - yc;
dz = xcm[2] - xc;
domain->minimum_image(dx,dy,dz);

because, due to minimum_image procedure, it's possible that xc<xcm[0] and
this is wrong. So I found that instead do this:

i don't agree with this assessment. the second call to
minimum image should catch this. the only way this would
not work is when your dx/dy/dy is larger than half the
cell, i.e. you have too small a force constant or too
small a cell. as i mentioned before, pulling only makes
sense with a fairly stiff spring.

your hack is not very general and will break badly in
other circumstances. please provide a better argument
and some input to reproduce the problem if you think
that i am in error on this.

thanks,
     axel.

[_nicola_varini]

Dear Axel, I runned two example of different sizes. If you look out_lammp_AK_10 you can see that the com are always more or less in the same position; this problem disappear using a larger sample. The parameters of the two system is the same but maybe the sample in out_lammps_AK_10 is too small.
However I’ve not clear a point:
suppose I’m pulling along the x-axis with v>0 when the minimum_image criteria is applied to xc then xprd is subtracted:
xc=xc_old-xprd
so
dx=xcm[0]-xc_old-xprd

in this case I think it’s not possible that the minimum_image criteria is applied in the same way for xc (dx = dx-xprd) but in the opposite way (dx = dx + xprd) and this can produce a discontinuity in the force.

The other point I would to stress is the fact that using rigid or rigid/nve produce a quite difference medium force; different kind integration can produce this kind of effects?

Cheers.

N.

2010/12/9 Axel Kohlmeyer <[email protected]>

lammps_file.tar.gz (180 KB)

[akohlmey]

nicola,

Dear Axel, I runned two example of different sizes. If you look
out_lammp_AK_10 you can see that the com are always more or less in the same
position; this problem disappear using a larger sample. The parameters of

ok, thanks. i will check this out. i have thought about this a little
bit, and i maintain, that i would like to have a justification why
you need to do a pulling with a small force constant.

to implement the pulling cleanly in a general way, the
distance between the pulling reference point and the
center of mass has to be less than half the simulation
cell diameter. there may be a way to extend this to one
simulation cell length, but this is as far as i would want
to go. anything else, you are free to hack, but you are
on your own.

the two system is the same but maybe the sample in out_lammps_AK_10 is too
small.
However I've not clear a point:
suppose I'm pulling along the x-axis with v>0 when the minimum_image
criteria is applied to xc then xprd is subtracted:
xc=xc_old-xprd

your statement is still not very clear.

xprd is only subtracted or added if the reference point moves
outside the box. that is intentional, since the center of mass of
the reference group is also wrapped back.

so
dx=xcm[0]-xc_old-xprd

in this case I think it's not possible that the minimum_image criteria is
applied in the same way for xc (dx = dx-xprd) but in the opposite way (dx =
dx + xprd) and this can produce a discontinuity in the force.

please see my statement from above. for the code to work
you need a stiff enough spring. give me a good physical
reason to use a softer spring and i will reconsider.

The other point I would to stress is the fact that using rigid or rigid/nve
produce a quite difference medium force; different kind integration can
produce this kind of effects?

i cannot comment on this. i suggest to contact the author of rigid/nve,
which is following a different formulation of the rigid integrator, or
check out the literature.

axel.

[_nicola_varini]

2010/12/10 Axel Kohlmeyer <[email protected]>

nicola,

Dear Axel, I runned two example of different sizes. If you look
out_lammp_AK_10 you can see that the com are always more or less in the same
position; this problem disappear using a larger sample. The parameters of

ok, thanks. i will check this out. i have thought about this a little
bit, and i maintain, that i would like to have a justification why
you need to do a pulling with a small force constant.

to implement the pulling cleanly in a general way, the
distance between the pulling reference point and the
center of mass has to be less than half the simulation
cell diameter. there may be a way to extend this to one
simulation cell length, but this is as far as i would want
to go. anything else, you are free to hack, but you are
on your own.

the two system is the same but maybe the sample in out_lammps_AK_10 is too
small.
However I’ve not clear a point:
suppose I’m pulling along the x-axis with v>0 when the minimum_image
criteria is applied to xc then xprd is subtracted:
xc=xc_old-xprd

your statement is still not very clear.

xprd is only subtracted or added if the reference point moves
outside the box. that is intentional, since the center of mass of
the reference group is also wrapped back.

Exactly, but when the point moves outside the box it’s possible that xcm is inside of the box and can happen that, after wrapping back xc, xc<xcm; or you think I ma wrong? Of course xc<xcm it’s not possible from a physical point of view.
In other word what I mean is that when xc is wrapping back also xcm must be wrapping back in the same time. Do you agree with this statement?

Using the minimum_image or the other methods is in principle the same, the only difference is that xc and dx are always in the simulation box if you use the minimum_image approach.

so
dx=xcm[0]-xc_old-xprd

in this case I think it’s not possible that the minimum_image criteria is
applied in the same way for xc (dx = dx-xprd) but in the opposite way (dx =
dx + xprd) and this can produce a discontinuity in the force.

please see my statement from above. for the code to work
you need a stiff enough spring. give me a good physical
reason to use a softer spring and i will reconsider.

The physical reason for applying a pulling force to this kind of system is for study the stick-slip behavior. You can find this kind of behavior only if you use a soft enough spring and the velocity must be not too high. If you use a spring too stiff and a velocity too high you found more or less the same behavior than using a constant force. Of course both the parameters are important because using a small velocity require a very long simulation time.