[lammps-users] freely jointed chain model

Hello,

for test purposes i would like to run a simulation of a freely jointed,
or of a fene bonded chain.
That means no pairwise interactions whatsoever, only langevin dynamics
restricted by bond length (if I fix the bond length) or restrictions due
to the fene bond.

When I try to set this up, I run into following problems.
Doing a simulation with pair_style none,
my first error message is:
ERROR: Atom sorting has bin size = 0.0

Looking it up in the manual, I find out that
"The neighbor cutoff is being used as the bin size, but it is zero. Thus
you must explicitly list a bin size in
the atom_modify sort command or turn off sorting."

Well, the neighbor cutoff is truly zero, as there is no neighbor
interaction calculation at all.
Following the manual, I can set the bin size using

neighbor bin 1

But this command has no effect: even when this command is in effect, I
get the same error "ERROR: Atom sorting has bin size = 0.0",
so clearly the
"neighbor bin 1" command did not set the bin size to 1 distance units.
Is this a bug, a feature, or am I doing something completely wrong?
Anyways, I try to coninue setting
"neighbor nsq 2", but I still get the same error.

After these trials, I conclude that I am unable to set the bin size.
I therefore have to resort to follow the second advice in the error
message, using the atom_modify command.

I first try to turn off sorting (as there is no neighbor interaction
going on, that should be safe, shouldn't it?):

atom_modify sort 0 1

Now I have a new error:
ERROR on proc 1: Bond atoms 651 652 missing on proc 1 at step 2000

(This happens as soon as I switch from soft potential to pair_style
none; note that no error happens when using lj pair_style with extremely
low energy settings)

Looking up the error in the manual, I find that
"
One or both of 2 atoms needed to compute a particular bond are missing
on this processor. Typically this
is because the pairwise cutoff is set too short or the bond has blown
apart and an atom is too far away"

I do not understand this message. First, setting an extremely weak
pairwise interaction did not result in any error.
There is, however, again the mentioning of the pairwise cutoff. For some
reason, it is too short. I still do not understand what the pairwise
cutoff has to do with the calculation of bonds.
And there is of course no reason to believe that a bond has been blown
apart, as the extremely weak pairwise lj interaction should result in
the same simulation (but slower, because calculation still takes place),
than the pair_style none simulation.

So at this place, I am left again with the hint to the pairwise cutoff
being too short. But as I have shown above, I cannot set this cutoff
manually. Switching off sorting doesn't help either.

There also is the neigh_modify command that does something with the
neigbor lists. I tried to switch off any such calculation:
neigh_modify exclude type 1 1
(there is only one kind of atom in my simulation)
or I try to switch it on as often as possible:
neigh_modify every 1 check yes

still the same result, always loosing a bond.

Now I am lost.
How do I set up a simulation with no pairwise interactions (pair_style
none), fene bonds, and langevin dynamics?
It seems the latter is without importance. It does not work for
non-langevin dynamics either.
So the question is:
How do I set up a simulation with no pairwise interactions, but fene
bonds?

"The neighbor cutoff is being used as the bin size, but it is zero. Thus
you must explicitly list a bin size in
the atom_modify sort command or turn off sorting."

This will be fixed by
atom_modify sort 0 0.0\

Now I have a new error:
ERROR on proc 1: Bond atoms 651 652 missing on proc 1 at step 2000

This will be fixed by
communicate single cutoff 1.5

Even with no neighbor lists and pairwise interactions, you still have
to have to get ghost atoms from a bond-length distance away
in order to compute bonds at processor boundaries. This is normally
done by the pair cutoff, which you don't have. So the comm cutoff
command sets it explicitly. 1.5 is the usual max FENE bond length.
You can adjust if that isn't enough for your model.

Steve