[lammps-users] LAMMPS as a library, access 'fix ave/atom' properties

Dear LAMMPS Developers and users, maybe this question is for Steve but everyone who can answer me is welcome. I am continuing my struggle to use LAMMPS as a library with a c++ driver for the moment. Consider the standard example in the couple folder where in the LJ input file we add a line with a fix to average atom positions as in the following:

Input file**********************************************************

3d Lennard-Jones melt

units lj
atom_style atomic
atom_modify map array

lattice fcc 0.8442
region box block 0 4 0 4 0 4
create_box 1 box
create_atoms 1 box
mass 1 1.0

velocity all create 1.44 87287 loop geom

pair_style lj/cut 2.5
pair_coeff 1 1 1.0 1.0 2.5

neighbor 0.3 bin
neigh_modify delay 0 every 20 check no

fix ave all ave/atom 1 5 10 x

fix 1 all nve

run 10

The simplest way to do this is just invoke the compute_scalar()
or compute_vector(int) method of the appropriate fix.

int ifix = modify->find_fix(char *) will file the index of the fix
with ID.

modify->fix[ifix]->compute_scalar() would return the scalar, etc.

You can add a library interface routine to library.cpp to do this
and name it whatever you want. You could also add an extract()
method to fix ave/time which is what lammps_extract() will try
to call.

Steve

Thanks Steve! After having tried all the day I have managed to access the values of the fix using modify->fix[ifix]->array_atom[int][int] as in the following:

for (int i = 0; i < nlocal; i++) {
id = tag[i];
offset = 3*(id-1);

vector1= &lmp->modify->fix[ifix]->array_atom[i][0];
vector2= &lmp->modify->fix[ifix]->array_atom[i][1];
vector3= &lmp->modify->fix[ifix]->array_atom[i][2];

copy[offset+0] = *vector1;
copy[offset+1] = *vector2;
copy[offset+2] = *vector3;
}

Somehow when using modify->fix[ifix]->compute_vector(int) the values in output were always =0.0 … at least now I am matching the results of the dump file… so it should be correct… I will try other possibilities and command for my c++ interface to LAMMPS. Thanks again.

Fabio Pavia

You can of course, do whatever works for you. But by accessing
the data directly instead of thru the provided interface routine,
your code could break if LAMMPS changes. The interface routine
does essentially what your code snippet does, so you must just
have some small error in how you are using it.

Steve

Sorry to bother further, but what do you mean by saying if ‘LAMMPS changes’?

By Checking how these commands are working for example in variable.cpp in the src folder I see that when a variable is taking a scalar value from a vector / array there is a distinction if the vector / array is global or peratom.

With global vector you should get a value only with modify->fix[ifix]->compute_vector(index-1), instead with peratom vector you should get a value only with modify->fix[ifix]->vector_atom[local_index], with local_index is the atom number as locally stored by the owning MPI process.

It looks that “Fix ave/atom” is working with vector_atom[ ] in fix_ave_atom.cpp (peratom values) and “Fix ave/time” is working with compute_vector() in fix_ave_time.cpp (global values).

Thanks again,

Fabio Pavia

Sorry to bother further, but what do you mean by saying if 'LAMMPS
changes'?

This was just a generic comment that you are better off using the
interface routines in fix.h (compute_scalar(), computer_vector()), rather
than accessing the lo-level data structure directly. The interface is
unlikely to change at this point; the lo-level data structures might (though
unlikely).

The rest of your comments are about per-atom vs global datums stored/returned
by fixes, and yes those are different. Not clear if you had a
specific question.

You mentioned earlier that sometimes you were seeing 0.0 returned by the
fix. That is correct behavoir, typical if you are accessing the values before
they have been computed.

Steve