[lammps-users] question about dumping forces

LAMMPS LAMMPS Mailing List Mirror
1

Dear LAMMPS users,

The LAMMPS version I use is 3Mar20

I have a question about how LAMMPS dumps forces.
In principal, I have a system with rigid walls that exert some pressure in the “z-direction” using “fix aveforce command”

Then, I am performing a calculation for the pressure in the z-direction using “compute stress/atom command” which with conversion to the pressure units produces the value that I wanted.

To be sure that I am in the thermodynamic conditions which I applied, I have also

performed the calculations of the pressure with the method of planes methodology with an auxiliary external program. It is essentially the eq. 14 of the paper specified in the LAMMPS documentation (B. D. Todd, Denis J. Evans, and Peter J. Daivis: “Pressure tensor for inhomogeneous fluids”, Phys. Rev. E 52, 1627 (1995)), i.e.

P_zz=1/2A \sum_i F_zi(sgn(zi-z))

while the kinetic part is assumed to be P_kin=kT \rho(zi)

where \rho(zi) is the local density taken from the density profile, exactly for the plane I take into account.

When I print the forces with “dump custom” command, and use the “fix aveforce” command, the resulting pressure profile (virial + kinetic parts) has approx. value of 0 atm in the middle of the slit.

On the other, hand, if I have the same system which was pre-equilibrated with “fix aveforce” command but now the walls are fixed (not integrated) and therefore there is no point of using “fix aveforce” (in other words, I assume that their position should correspond to the desired pressure), the resulting pressure profile (virial + kinetic parts) has the desired value in the middle of the slit.

Therefore, my question is: Does LAMMPS take into account the external forces that are applied with, for instance, “fix aveforce/addforce” commands to be printed by the “dump custom” command by the default, or is there some additional “modify” command to do so?

Thanks in advance,

2

how to include contributions from fixes is documented in the manual.
I would advise to use the latest LAMMPS version for that, since there were some inconsistencies in the past and that was resolved in the code not so long ago (you have to look through the changelog/bugs pages to see exactly when).

BTW: have you seen that LAMMPS has a compute stress/mop ?

3

Dear Axel,

Could you please send me where I can find that? I have found something such as fix store/force but it seemed to me to be the opposite to what I want. Another option is a fix modify but again, it is apparently not consistent with aveforce plus on top of that it can be added to the global and not per atom pressure which makes perfectly sense.

According to the method of planes command that is available in LAMMPS. It is very limited as it is only able to account for pairwise forces and not include kspace, intramolecular forces, etc.

To be honest I do not know why you implemented this in a way of eq. 15 of the paper you mentioned in the documentation rather than just a simple sum of forces in the entire system, i.e. eq. 14. In this sense, one could get rid of the limitations and wouldn’t need to afford a bit of additional computational time. On the other hand, kinetic contribution could also be computed as
P_ii=d/2kT \rho(I),

so essentially from the density profile where “d” is a dimension of the system and that would work for any interparticle potential.

Best wishes,

Łukasz Baran

pt., 8 paź 2021, 16:01 użytkownik Axel Kohlmeyer <[email protected]> napisał:

4

Dear Axel,

Could you please send me where I can find that? I have found something such as fix store/force but it seemed to me to be the opposite to what I want. Another option is a fix modify but again, it is apparently not consistent with aveforce plus on top of that it can be added to the global and not per atom pressure which makes perfectly sense.

that would be a limitation of fix aveforcem then.
perhaps you can emulate the fix aveforce behavior with fix addforce and some variables.

According to the method of planes command that is available in LAMMPS. It is very limited as it is only able to account for pairwise forces and not include kspace, intramolecular forces, etc.

To be honest I do not know why you implemented this in a way of eq. 15 of the paper you mentioned in the documentation rather than just a simple sum of forces in the entire system, i.e. eq. 14. In this sense, one could get rid of the limitations and wouldn’t need to afford a bit of additional computational time. On the other hand, kinetic contribution could also be computed as
P_ii=d/2kT \rho(I),

i didn’t implement this. in fact, most of the “user visible” features in LAMMPS are contributed code, not developed by the core LAMMPS developers. So you would have to ask the individual author such specific questions. Please keep in mind that LAMMPS is designed to run efficiently in parallel with domain decomposition, which makes any operation that requires global or otherwise spatially sampled data require extra communication. that is a limitation that you usually don’t have to worry about in post-processing.

so essentially from the density profile where “d” is a dimension of the system and that would work for any interparticle potential.

feel free to consider contributing such an implementation.