Reading Energy in Custom Pair Function

LAMMPS LAMMPS Mailing List Mirror
1

Hi, here is my update:

I am writing a fix (fix_addboost) that identifies the particle in a group possessing the highest pair potential energy, and then applies a repulsive potential to that particle only (this principle is used in the bond-boost method, especially of Miron-Fichthorn).

My code compiles, but I get a "segmentation fault (core dumped)" error upon running. I have located the problem in the following lines (lines 277-278 in the full code attached) that try to read "eatoms[i]".

fprintf(screen,"eatom[i] = %f \n", eatom[i]);
currentEnergy = eatom[i];

I am guessing that eatoms[i] have not been called properly. I have followed the syntax in "compute_pe_atom.cpp" and used the following in my code:

if (force->pair){
double *eatom = force->pair->eatom;
        ....
}

I have also included pair.h and force.h at the beginning. I have gone through the LAMMPS developer manual also, but I still cannot figure out what is going wrong.

Can anybody help me solve this problem?

Thanks,
S

fix_addboost.cpp (14.7 KB)

fix_addboost.h (2.42 KB)

2

Hi, here is my update:

I am writing a fix (fix_addboost) that identifies the particle in a group
possessing the highest pair potential energy, and then applies a repulsive
potential to that particle only (this principle is used in the bond-boost
method, especially of Miron-Fichthorn).

My code compiles, but I get a "segmentation fault (core dumped)" error upon
running. I have located the problem in the following lines (lines 277-278 in
the full code attached) that try to read "eatoms[i]".

fprintf(screen,"eatom[i] = %f \n", eatom[i]);
currentEnergy = eatom[i];

I am guessing that eatoms[i] have not been called properly. I have followed
the syntax in "compute_pe_atom.cpp" and used the following in my code:

if (force->pair){
double *eatom = force->pair->eatom;
        ....
}

I have also included pair.h and force.h at the beginning. I have gone
through the LAMMPS developer manual also, but I still cannot figure out what
is going wrong.

with your current setup, LAMMPS doesn't know, that you want to use
per-atom energies, so the array is not allocated and the energy not
tallied into that array.
the clean way to deal with this would be:
- create an instance of compute pe/atom with the suitable arguments
from within your fix
- read the data from the compute's vector_atom member after having
called its compute_peratom() method
- call modify->clearstep_compute() /
modify->addstep_compute(nextstep), where nextstep is the next
timestep, where you need access to the data from the compute.

see other fixes that maintain their own compute instances for more details.

axel.

3

Thanks, Axel, I have followed the example in “fix_ave_atom.cpp” and implemented necessary computes in my fix. It works as expected now.

Anybody in the mailing list interested in the BOND-BOOST METHOD is welcome to use this fix

On a less important note:

This fix writes the amount of repulsive energy applied (“boostAmount”) to a text file, which is in turn read by other fixes (e.g. fix_nve.cpp). Instead, I am thinking of whether boostAmount can be stored as a variable and compute, which can be accessed by other fixes without the need of a separate text file.

Can anybody help me out here?

Thanks,
S

fix_addboost.cpp (15.8 KB)

fix_addboost.h (2.46 KB)

4

Any fix can produce output which other fixes, computes, output access.

See fix indent as a simple example.

Steve

5

Hi, last week I wrote a fix to apply a boost force to a selected atom in a group. Now, I am trying to create a neighbor list in the same fix. I have used the following syntax (lines 346-353 in attached code):

class NeighList *list;

neighbor->build_one(list->index);

inum = list->inum;
ilist = list->ilist;
numneigh = list->numneigh;
firstneigh = list->firstneigh;

But I get “segmentation fault (core dumped).” I am guessing that the neighbor array is not allocated properly. I have looked at other fixes but have not found any example to follow. Can anybody help me correct my code?

The .cpp file is attached.

P.S. Also, for some reason, this fix functions about three times slower than a similar fix (e.g. fix addforce, which was modified into this fix). I am not sure what makes it so slow, and I doubt that the absence of a neighbor list is the reason (I tested by bypassing the loop). If anybody can spare the time to go through my code, please help me figure out how to make it more efficient.

Thanks,

S

fix_addboost.cpp (16.3 KB)

fix_addboost.h (2.51 KB)

6

Hi, last week I wrote a fix to apply a boost force to a selected atom in a
group. Now, I am trying to create a neighbor list in the same fix. I have
used the following syntax (lines 346-353 in attached code):

class NeighList *list;

neighbor->build_one(list->index);

this line is _wrong_! there is no precedence in any other LAMMPS
source files for it. the segfault is well deserved.

inum = list->inum;
ilist = list->ilist;
numneigh = list->numneigh;
firstneigh = list->firstneigh;

But I get "segmentation fault (core dumped)." I am guessing that the
neighbor array is not allocated properly. I have looked at other fixes but
have not found any example to follow. Can anybody help me correct my code?

you need to study other LAMMPS sources more thoroughly and more
carefully. then you will find:

- there are several examples to follow, e.g. fix orient/bcc or
orient/fcc, or fix qeq
- to obtain a neighbor list, you need to create a neighbor list
request instance and configure it for your needs
- you also need to add an init_list() method to your class to get hold
of the neighbor list pointer
- you won't need to request an explicit neighbor list build, but can
stick with the default perpetual neighbor list, that will likely be a
copy of or derived from an existing list.

The .cpp file is attached.

P.S. Also, for some reason, this fix functions about three times slower than
a similar fix (e.g. fix addforce, which was modified into this fix). I am
not sure what makes it so slow, and I doubt that the absence of a neighbor
list is the reason (I tested by bypassing the loop). If anybody can spare
the time to go through my code, please help me figure out how to make it
more efficient.

you are asking for A LOT here. after all, it is *your* code and *your*
work. please ask yourself, would you do this kind of work on somebody
else's code?
there are plenty tools to assist in determining performance
bottlenecks and to profile applications. there are tutorials, courses
and more. if you care about the performance of your code, you *will*
need to learn what is making code fast and what is slowing you down.
for example, have you checked, how much time you spend opening,
reading and parsing (text format) files?
if you feel unable to do this, you need to find a suitable
collaborator, that you can share the credit with, or hire an expert.

axel.