How to dump some energy into the system efficiently

Dear Steve, Users,

I plan to add thermal enregy into my system of interest in only selected region (think of a local photon absorbtion on the surface region). I am running my simulations using lammps as library. One of the ways to do it is to user lammps_gather_atoms to extract the cooridnates and velocities, Update the velocities (using some factor) of the atoms which are in certain region and then use lammps_scatter_atom function to feed it back to lammps.

As the lammps_gather_atoms command extract the information of “all” the atoms, I will need to calculate the distance of each atom from the center of photon absorption site to decide whether to update its velocity. If that region is small and my overall system is large, my concern is that it might not be efficient.

I am thinking if I could write a routine in lammps itself which will do a similar thing locally on requesting in input script. The inputs will be group (based on region) and the energy (to be distributed) and possibly how the energy needs to be distributed (gaussian, uniform), etc. As this is not going to be happen at constant time intervals, write a fix probably won’t work.

Another thing could be to write something along the lines of lammps_gather_atom_group in which we specify the group and extract the information of atoms only for that group. For preallocation, we can get the array size from using variable count(group) feature. But after going through library.cpp, it seem to me that atoms “tagging” (which atom goes where in copy array) could create problems when MPI_Allreduce is called?

Your suggestions are welcome. What is more efficient? What templete routines should I look over. (displace_atoms.cpp?) Please guide.

Best Regards,

There are several ways to approach this problem.

(1) look at fix heat, which may do exactly what you want,
or which could be adapted by you to deposit its energy
according to some other formula

(2) The lib interface has a lammps_extract_atom() function
which returns the ptr on each processor to the data for
x,v,whatever. There is no gather/scatter overhead to
using that. Then each proc can loop over its local atoms,
determine if you want to rescale its velocity, and just do
directly. This is the same way something like fix heat operates.

(3) It is not hard to give a fix (fix heat or one you write) an
arg that is an equal-style variable that returns a value for
the next timestep on which you want to invoke it. Thus
if you can write an equal-style variable that generates
the sequence of timesteps you want to invoke your operation,
you can use it. See the dump_modify every v_foo command
as an example. And there are several variable functions
that allow you to generate irregular sequences of timesteps.
Or you could simply code the calculalation of the next timestep
into the fix, have the fix be invoked every step, and simply
return if it is not the "next" timestep.


Dear Steve,

Thank you for all of your suggestions. All of them are makes sense now. I will start from fix heat as you suggested and will try to modify it according to your suggestions.

Best Regards,

Dear Steve,

As per your suggestion, I looked over fix heat and the code by itself is almost clear to me. After thinking a bit, one issue that I figured in using “fix” based methodology is that my group of atoms over which heat is being dissipated varies over time. In addition, the region over which I distribute the energy will change with time (think of a photon whose absortion point (X,Y,Z) will be detemined in stochastic way). In that case, I cannot even use group “all” with some region id as region itself will change. So, I was thinking to adapt the fix heat code to something which is not invoked repeatedly but only when I asked it to with different group definitions or region definitions. The idea is as follows. Lets say we call it dump_energy command

  1. Run A steps (where A is variable)
  2. Define a region of radius R around stochastically chosen point x,y,z point.
  3. Define the group of atoms.
  4. Use dump energy command as follows
    dump_energy group energy (energy units) distribution (unform,gaussian,others)
  5. Delete the group
  6. Delete the region
  7. Go to step 0

Please guide if this makes more sense.

Your second suggestion regarding using lib interface also has merit. But in that case as well, I would need to calculate/extract KE/mass of the group of atoms, calculate/extract center of mass velocity and then rescale the K.E. energy accordingly. In my opinion, writing a module based on fix heat and calling it from the C++ as a command line would be a more systematic way to go.

Please guide.

Best Regards,

Regions can be dynamic - see the region doc page.