pair hybrid/overlay

Dear Lammps Developers,

I am going to use " pair_style overlay" to combine two potentials (potential A and potential B) for the system I want to study. I want to get the values of force and energy from each potential individually at each time step and then do some calculation on them in order to calculate the contribution of each potential to total heat flux separately. Then, I am going to output both total heat flux and contribution of each potential to heat flux separately to a file. Could you please guide me how I can do this?
I know ev_tally() function in each pair style adds up the energies and forces of each potential but If I use two or more potential in the simulation, I have know idea how " pair_style overlay" adds up the energies and forces from the contribution of each pair style.

I looked at the “pair_hybrid_overlay.cpp” but compute(int eflag, int vflag) function in the code. I need to know how pair_hybrid/overlay do the calculation of the forces and energies.

Please note that compute heat/flux cannot do this because at each time step and after getting the forces I need to do some heavy calculations on the forces and energis.

Thank you

Hossein

You can try using the “rerun” command with potentials separately.

Dear Lammps Developers,
I am going to use " pair_style overlay" to combine two potentials (potential
A and potential B) for the system I want to study. I want to get the values
of force and energy from each potential individually at each time step and
then do some calculation on them in order to calculate the contribution of
each potential to total heat flux separately. Then, I am going to output
both total heat flux and contribution of each potential to heat flux
separately to a file. Could you please guide me how I can do this?

you will likely have to do some significant programming. LAMMPS tries
to avoid keeping excessive copies of data around that is not needed,
but rather computes it incrementally, or re-computes it after the fact
(as in compute group/group via Pair::single()). there are different
ways how this could be done. it depends a lot on what your potentials
A and B are and what exactly you are going to do with the results. for
example, you need to write your own overlay pair style that will save
and store intermediate force and energy data into separate storage and
then make it accessible to a custom heat flux compute. i don't quite
see what the benefit of outputting this data to a file would be. i'd
rather process the data directly from inside LAMMPS, e.g. with a
compute or a fix. these will likely be very large files and thus
outputting them will have a significant negative performance and
parallel scaling impact.

I know ev_tally() function in each pair style adds up the energies and
forces of each potential but If I use two or more potential in the
simulation, I have know idea how " pair_style overlay" adds up the energies
and forces from the contribution of each pair style.
I looked at the "pair_hybrid_overlay.cpp" but compute(int eflag, int vflag)
function in the code. I need to know how pair_hybrid/overlay do the
calculation of the forces and energies.

it is a simple addition. each pair style (and bond, angle, dihedral,
improper, kspace style) incrementally adds to per atom forces and
energies or stress if requested. thus you cannot find any particular
code in pair_style hybrid or hybrid/overlay that shows how this is
combined, since there is nothing stored separately. the hybrid pair
styles primarily need to take care of maintaining the multiple
instances of the different pair styles, suitable sub-neighbor lists
for the individual invocations of the sub styles and make sure that
the requirements for specifying the input are met.

Please note that compute heat/flux cannot do this because at each time step
and after getting the forces I need to do some heavy calculations on the
forces and energis.

what kind of "heavy" calculations would that be?

Thank you very much for your responses. I can calculate the heat fluxes directly inside of each potential file and output them at each time step. I mean in cpp file for each potential, after the force calculation, I calculate the heat flux array (Q[0], Q[1], Q[2] where 0,1,2 are x, y, and z component of heat flux) and then output them to a file. The problem is that if I use for example N potential using pair_style hybrid, N output files will be generated. I can do a post processing on these N output files to add the value of heat flux for each potential and store them to a new file. But the major problem is that the size of heat flux files are very large and I don’t want to save them even temporary in my hard drive. I want to add the heat fluxes at each time step and then store them to a file. Is there any way to add them in the during the calculation? How and where can I add them?

For example, let’s say If I use lj potential and tersoff and I have the Qlj and Qtersoff at each time step. where and how I can output Qtotal = Qlj+Qtersoff to a another file?

I looked at the pair.cpp where both lj and tersoff are dervied but I couldn’t figure out where I have to define my Qtotal array and output them.

Thank you very much

Hossein

Thank you very much for your responses. I can calculate the heat fluxes
directly inside of each potential file and output them at each time step. I
mean in cpp file for each potential, after the force calculation, I
calculate the heat flux array (Q[0], Q[1], Q[2] where 0,1,2 are x, y, and z
component of heat flux) and then output them to a file. The problem is that
if I use for example N potential using pair_style hybrid, N output files
will be generated. I can do a post processing on these N output files to add
the value of heat flux for each potential and store them to a new file. But
the major problem is that the size of heat flux files are very large and I
don't want to save them even temporary in my hard drive. I want to add the
heat fluxes at each time step and then store them to a file. Is there any
way to add them in the during the calculation? How and where can I add them?
For example, let's say If I use lj potential and tersoff and I have the Qlj
and Qtersoff at each time step. where and how I can output Qtotal =
Qlj+Qtersoff to a another file?
I looked at the pair.cpp where both lj and tersoff are dervied but I
couldn't figure out where I have to define my Qtotal array and output them.

sorry, but i had already outlined some possible steps how to do this.
explaining things at the level of detail that you are requesting here
would take me more effort than actually implementing this. what you
want to do requires some deep understanding of how LAMMPS work and the
only way to get this is through carefully reading thing several of the
style files, checking which of them do similar things and adapting
this for your personal needs. what you seem to be hoping for is a
rather minimal approach that is based on some simplified assumptions
about how LAMMPS works, which are obviously not true. if you don't
take the time to study the different parts of LAMMPS that are relevant
to your needs, you won't get your task done. not unless you can entice
somebody that already has that insight to collaborate for you and do
the programming for you.

some comments and a more detailed outline (as much as i am willing to
do) follows here:

- i consider it a *very* bad idea to modify the individual pair
styles. what you are computing is a general property and thus it
should be implemented outside of them, so that it can be applied to
*all* pair styles and will be future proof.

- there are different ways in LAMMPS to store per atom properties.
they can be part of the atom style (that is for properties that are
integral and commonly used, migrate with atoms between sub-domains and
need to be accessible at different stages of the time integration).
these properties can be augmented by properties maintained through fix
property/atom, and if you need to keep information from old time
steps, there is fix store/state. otherwise, you can keep per style
storage in the individual pair styles, see for example pair style eam.
this is a very educational example, since it also demonstrates forward
and reverse communication.

- as i already mentioned, all force computing styles add their forces
incrementally to the force arrays. as you can see from studying the
verlet integrator, pair styles are computed first and the force arrays
are zeroed before. also, there can be only one pair style at a time.
so if you want to keep two separate copies of the forces from two
different pair styles as you describe it, the most straightforward way
would be to write a "multi" pair style that allows you to run multiple
pair styles in a way similar to pair style hybrid/overlay, but then,
while looping over the substyles, you would copy the computed force
into a local per atom storage array maintained by your pair style (for
newton pair on you will also have to do a reverse communication before
that) and then zero out the force array again, to get the second copy
(if needed, you may have to store this as a second local copy with an
optional reverse communication as well) and then proceed. please note
that for heat flux calculations, you do not have up-to-date velocities
at this point, as for the second half step of the velocity verlet, you
will need the new forces computed at the updated positions. there is
also the question of how to handle bonded interactions at this point,
if you want to write a tool that is actually general and can be used
in multiple ways. now it is up to you how to process this data, e.g.
you can access it from a compute later on that would then verify that
your pair style is the new hybrid variant, then casts the pointer to
the generic pair class instance to a specific class instance, and
access the stashed away data and - if needed - the per atom force
array. then you can compute your heat flux as desired. so the
programming required would be a variant of hybrid/overlay (a rather
difficult task) and a variant of the heat/flux compute that would also
be a bit more complex as the original. this is quite demanding, but
the benefit would be a very generic solution.

- an alternate approach would be using a facility that i have been
working on with a collaborator that is not yet in LAMMPS but already
publicly available in LAMMPS-ICMS. there we have implemented a way to
add a callback from a compute to every(!) call of Pair::ev_tally().
this way the callback would be invoked for each pair of atoms. in the
case of your scenario of LJ and Tersoff, you could just set up an
input for Tersoff and then have your compute maintain an instance of
second pair style. then in the PairTersoff::compute() function if
Pair_ev_tally() is called for the non-bonded part (mind you at this
step the full force is not available), you would use to compute the
additional force for that pair from calling Pair::single() for the
instance that the compute maintains and store it in a suitable way for
later processing (e.g. as per atom forces or per atom stress tensor
contribution) and then do the post-processing of that in a similar
fashion to what compute heat/flux does also in the compute_vector()
method. this is not trivial either, but might be simpler to do, for as
long as you are comfortable using an infrastructure that is rather new
and has only been validated for similar but somewhat different uses.

- in both cases, you have to pay attention to how data is managed in
general in LAMMPS and particularly how parallel computation and the
newton pair flag affects it as well as compare with how different
styles do similar things in LAMMPS. this *will* require significant
reading of existing code and some careful experimentation and
debugging and testing. everybody that does program complex features
into LAMMPS will have to do this. "there ain't no escape from the
blues", or has to entice/convince/hire somebody to do this work for
him/her. this definitely is not something that can be done easily even
for an experienced LAMMPS programmer. however, unless there are
boundary conditions that you have not mentioned (a lot of your
description is still rather vague), it would be straightforward for an
experienced programming to take my suggestions and implement them. i
would also assume there are other ways to do this, e.g. a using a
modified integrator class, at different levels of generality and
complexity. it all depends on individual creativity, perseverance and
skill.

as far as i am concerned, from here on you are on your own. you have
to make a choice of how much energy and time you want to invest
yourself or how you would want to address implementing something into
LAMMPS that is non-trivial.

axel.

Dear Axel,

I really appreciate for your complete explanation. Thank you very much for your time and consideration.

Hossein