Global atom in fix

All,

Sorry if this is not the appropriate place for this.

I’m trying to write a fix that requires an interaction between a single atom and all the others in the system. I’ve tried to find an example of another src file that does something similar, but cannot. I suspect it will require manually adding the atom to every processor and then using atom_map to identify it in the local lists?

Does anyone know if there’s a similar functionality in another LAMMPS SRC file I may be able to use as a prototype, or an easier solution?

Thanks for your help,

the closest thing to what you describe is probably pair style list. however, there is a conceptual difference. pair style list computes interactions between explicitly listed pairs of atoms that can be computed even when they are not stored on the same sub-domain (as usually required by pair styles). however, in your case, that scenario is only well defined, if there are no periodic boundaries, while pair style list has only the limitation that minimum image conventions need to be met (which is not required for regular pair styles where the interaction cutoff may be larger than half the box diameter, since LAMMPS will simply generate more ghost atoms as needed).

what you can do is to have two double[3] arrays for the “global” atom’s positions and forces. in every step you set them to all zero, then set the positions on the MPI rank where that “global” atom is local (i.e. the result from atom->map(global_atom_id) is >= 0 and < atom->nlocal) to the actual positions, and then call MPI_Allreduce() with MPI_SUM to replicate positions on all MPI ranks. then you can loop over all local atoms on all MPI ranks, compute the interaction, add it to the local atom forces and second storage array. then do another MPI_Allreduce() with MPI_SUM on the global atom forces to get the total that force on the “global” atom and then add those force components where that global atom is local.

again, as mentioned before, with periodic boundaries, your scheme becomes ambiguous, since “all” atoms without a cutoff would include all periodic images and that cannot be done without using something equivalent to ewald summation.

HTH,
axel.

I wonder if, as long as it is a single “global” atom we are talking about, you could “fake it” by considering it to be in the center of the box (or anywhere) and using fix addforce with variables on all the other particles. Something along the lines of
variable x0 equal

variable y0 equal

variable z0 equal

variable r2 atom “(v_x - v_x0)^2 + (v_y - v_y0)^2 + (v_z - v_z0)^2”

variable fx atom
etc.
fix global_particle_force all addforce v_fx v_fy v_fz

As far as performance goes though, this will probably not be very efficient because of all the variable evaluations. If you are willing to write a bit of code, I would follow Axel’s approach.

I was actually thinking about doing something similar, and passing the necessary information about the “global” atom in to the fix as a variable. The boundary conditions need to be periodic, as I’m doing long range forces in my simulation as well.

I was actually thinking about doing something similar, and passing the necessary information about the “global” atom in to the fix as a variable. The boundary conditions need to be periodic, as I’m doing long range forces in my simulation as well.

as pointed out before, with periodic boundaries you need a cutoff (and make certain, that the communication cutoff is set properly, and determine the interaction with the closest periodic image) because otherwise your interactions not consistent and you cannot really do “all” atoms, since you are cutting things off arbitrarily.

axel.

Axel,

Thank you very much for all of your help.

I guess I’m a little confused as to the way periodic images are handled in LAMMPS. Are they treated as normal atoms in the atom array? I would have guessed they’d be treated as ghost particles wrapping from the boundary -X -> +X.

To be specific about what I’m trying to do. I have a code for induced magnetic dipoles in permeable materials. Previously I just used a background field to get the initial magnetization, then used the mutual dipole method to take into account the fields from near neighbors as well. Now instead of a uniform background field, I want to have one strong permanently magnetized atom create the field that initially magnetizes the other grains.

I don’t need to calculate any global forces, I just need to know where the atoms on the local list are in relation to the single seed atom.

Oh, I see what you mean, if the seed atom is closest to an atom across the periodic boundary as opposed to direct. I agree this is potentially a problem, in my current simulation however, the seed atom should remain very close to the center of the domain, and the domain is large compared to the grain distribution in order to minimize periodic effects on the long range solver.

Axel,

Thank you very much for all of your help.

I guess I’m a little confused as to the way periodic images are handled in LAMMPS. Are they treated as normal atoms in the atom array? I would have guessed they’d be treated as ghost particles wrapping from the boundary -X → +X.

in LAMMPS periodic boundaries are implemented as part of the domain decomposition. you have to understand that before looking into how to handle periodic boundaries.
domain decomposition breaks the box down to subdomains. the atoms in those are distributed and indexed 0 <= i < atom->nlocal in the per-atom arrays. the atoms within the communication cutoff are communicated to the neighboring subdomains and copied. those atoms are then referred to as ghost atoms and are indexed atom->nlocal <= i < atom->nlocal + atom->nghost

this also implements periodic boundaries, as you look at subdomains across periodic boundaries for generating/communicating ghost atoms. because ghost atoms are copies, LAMMPS is not subject to minimum image conventions (at the expense of maintaining more ghost atoms with longer communication cutoffs).
for finding the closest image, the Domain class has APIs to figure out which atom to access.

axel.

Oh, I see what you mean, if the seed atom is closest to an atom across the periodic boundary as opposed to direct. I agree this is potentially a problem, in my current simulation however, the seed atom should remain very close to the center of the domain, and the domain is large compared to the grain distribution in order to minimize periodic effects on the long range solver.

you should use the utility functions in the Domain class to find the most relevant atom (local or ghost).

axel.

That is as I understood it, I completely missed what you were getting at in terms of the ambiguity however, because I’m not thinking of the system as periodic, I’m using the periodic boundaries to allow for spectral long range solvers, while using a large domain to ensure low error from the periodic assumption.

I agree that the resultant ambiguity would make it unacceptable for wide application, however in the case I am examining I believe the naive approximation will not introduce unacceptable error. I will try to write a more robust solution into my code in the future, or simply leave that functionality out of any code I release.

Thanks again, you have been extremely helpful.