I came upon a problem that I think I understand and can solve but while I investigating I found an inconsistency that I would like to understand better. I’m hoping that someone can give me an explanation of what is happening…
I am simulating a bead-spring polymer at very high energy/Temperature scales. This can stretch bonds beyond the neighbor list / ghost cutoff and then cause incorrect bond energies to be computed if bonded atoms cross a periodic boundary. attached input script and data file illustrate the problem. The bond energy (E_bond) at step 3 is ten times higher than step 2 because an incorrect bond length is used to compute the bond energy. This can be avoided/fixed by invoking neigh_modify cluster yes or by increasing the comm_modify cutoff.
The part that I don’t understand is that using “compute bond/local eng” and then using compute reduce sum (or just printing for a single bond) gives the correct bond energy whether or not the neighbor list is correct. I put some rudimentary debugging statements in lammps but don’t see how/where/when this calculation is being performed and why it’s being performed correctly (using minimum image) while the ebond calculation is done incorrectly.
Is this expected behavior? Using lammps-7Sep16 if anyone tries to reproduce.
I came upon a problem that I think I understand and can solve but while I
investigating I found an inconsistency that I would like to understand
better. I'm hoping that someone can give me an explanation of what is
happening...
I am simulating a bead-spring polymer at very high energy/Temperature
scales. This can stretch bonds beyond the neighbor list / ghost cutoff and
then cause incorrect bond energies to be computed if bonded atoms cross a
periodic boundary. attached input script and data file illustrate the
problem. The bond energy (E_bond) at step 3 is ten times higher than step 2
because an incorrect bond length is used to compute the bond energy. This
can be avoided/fixed by invoking neigh_modify cluster yes or by increasing
the comm_modify cutoff.
The part that I don't understand is that using "compute bond/local eng" and
then using compute reduce sum (or just printing for a single bond) gives the
correct bond energy whether or not the neighbor list is correct. I put some
rudimentary debugging statements in lammps but don't see how/where/when this
calculation is being performed and why it's being performed correctly (using
minimum image) while the ebond calculation is done incorrectly.
this is because compute bond/local uses the bond->single() function
and passes the atoms and rsq as an argument, while the regular
force/energy computation uses bond->compute().
you are not running in parallel, right?
> Hi all,
>
> I came upon a problem that I think I understand and can solve but while I
> investigating I found an inconsistency that I would like to understand
> better. I'm hoping that someone can give me an explanation of what is
> happening...
>
> I am simulating a bead-spring polymer at very high energy/Temperature
> scales. This can stretch bonds beyond the neighbor list / ghost cutoff
and
> then cause incorrect bond energies to be computed if bonded atoms cross a
> periodic boundary. attached input script and data file illustrate the
> problem. The bond energy (E_bond) at step 3 is ten times higher than
step 2
> because an incorrect bond length is used to compute the bond energy.
This
> can be avoided/fixed by invoking neigh_modify cluster yes or by
increasing
> the comm_modify cutoff.
>
> The part that I don't understand is that using "compute bond/local eng"
and
> then using compute reduce sum (or just printing for a single bond) gives
the
> correct bond energy whether or not the neighbor list is correct. I put
some
> rudimentary debugging statements in lammps but don't see how/where/when
this
> calculation is being performed and why it's being performed correctly
(using
> minimum image) while the ebond calculation is done incorrectly.
this is because compute bond/local uses the bond->single() function
and passes the atoms and rsq as an argument, while the regular
force/energy computation uses bond->compute().
This explains what I was seeing with my basic print statement debugging.
you are not running in parallel, right?
Does it matter whether I run in serial or parallel? At first I ran in
serial but I just tested in parallel and don't seem to see a difference.