I'm trying out both the fix:array and the find_custom options you
I'm now able to define new variables through the lammps_command
function and find them again through the code snippet below.
int tmp = -1;
int n = lmp->atom->find_custom("d_xDuDt",tmp);
double *dvector = lmp->atom->dvector[n];
How can I now write to individual atoms?
Given that I have nlocal atoms on my MPI rank and I know the index
of the atom of interest "atomI", I'd like to write a new double to
dvector[atomI] or to lmp->atom->dvector[n][atomI]. If tried it in
this manner but I get an error and gdb isn't very informative about
it (only gives me the line).
Long story short, problem in chair, not in computer. Thanks for the
elaborate reply though. I had to search for "xDuDt" rather than
"d_xDuDt", otherwise you get -1 as a dvector index and the code crashed
due to out of bound array indexing. Could be prevented by a simple "if(n
< 0)" check.
it is not clear what you are referring to as "index". atoms in LAMMPS
have an "atom id", which is stored in atom->tag. this "tag" value is
the atom id used in a data file or the "id" field in a dump file. these
allow to address atoms globally.
With index I mean the index of the atom in a for loop over local atoms
on the rank/ local domain.
I keep track of them by using a lmpCpuId variable which is just a place
holder for the MPI_Comm_rank() output. This way I only send stuff to the
local domain atoms. The reason this is needed is that OpenFoam has a
mesh based decomposition and a particle may reside on a different processor.
In the coupling lib function I match the "atom id" between the incoming
info and the local atoms by sorting the tag arrays and looping over the
indexes (I hope the way I put it here make sense to you). That way I
don't have to care about the order are in, I jsut make sure the atom
od's are unique (they should be of course).
This was all pre-implemented in the sediFoam framework by Heng Xiao et
al. ( https://github.com/xiaoh/sediFoam ). I just ported the whole thing
to the latest lammps version and modified some things on the openfoam side.
however, the per-atom storage, is associated with atoms and they are
distributed and in a different order. you have the "local" atoms (i.e.
atoms owned by a subdomain == MPI rank) between index 0 and nlocal and
"ghost" atoms (copies of atoms from neighboring subdomains up to the
communication cutoff) between nlocal and nlocal+nghost. to get the
global atom id from a local or ghost atom, you look at atom->tag[i] you
will get the global atom id of an atom, and the function atom->map(id)
will return the local atom index or -1 from a global atom id. index ==
-1 means, this atom is not present on the local MPI rank.
Question though, will the following command give me back ghost atom tags
or only local atom tags ?
LAMMPS *lmp = (LAMMPS * ) ptr;
int *tag = lmp->atom->tag
if you are looping over atom indices from neighbor lists, be aware that
those are using local indices and in the jlist, atoms indices may be
"decorated" with extra bits, when an atom is connected via a bond or
angle or dihedral interaction and "special_bond" scaling factors are
I just do granular stuff without bonds etc, but I'll look into this as a
if you want gdb to be more helpful, please compile without optimization
and with -g. then gcc won't optimize away loop indices or other
temporary variables and you can look up their values as expected.
I still had some optimizing flags in the makefile, lazy of me sorry.
Removed them and now it's more helpful.