Bad FENE bond and Bond atoms missing on proc when using large number of processors

Dear Dr. Axel Kohlmeyer and LAMMPS Users,

Hello All. I am currently simulating 100 polymer chains that crosslink between the same Type 3 beads and between the same Type 4 beads using the Kremer-Grest bead (coarse-grained) spring model. The non-connected bonds are governed by the LJ potential (cutoff of 2.5 using LJ units) whereas the connecting bonds and including the newly formed crosslinked bonds are governed by the FENE potential. (K = 30 R = 1.5)

I have attached my input script which reads the data file that contains the coordinates of equilibrated 100 polymer chains. (The initial configuration of 100 separate chains is well-equilibrated) When I use only 8 processors to run the input script, it runs fine free of any warnings or errors.

However, when I use 27 processors without changing my input script and the data file, I obtain a FENE bond too long error after 300 timesteps.

(ERROR on proc 55: Bad FENE bond)

I inserted the minimize command:

min_style cg
minimize 1e-25 1e-25 500000 1000000

Inserting this minimize command allowed me to run the simulation for a while but the simulation terminated after 2081 timesteps with

(ERROR on proc 55: Bond atoms 16222 16223 missing on proc 55 at step 2081)

I understand that both of these errors are due to bad dynamics or bad force fields.
My force fields are standard using the typical Kremer-Grest bead spring parameters, and moreover, the simulation works perfectly fine and the crosslinking process occurs when using only 8 processors.

However, without changing anything except for using 64 cores, I obtain the above errors.

After searching through the mailing list, my pairwise cutoff distance + skin distance greatly exceeds the maximum extension of the FENE bonds. Moreover, I frequently build my neighbor list since fix bond_create will require the neighbor list to be updated.

I am aware that using less number of processors for parallel simulations is preferable. Nevertheless, I would greatly appreciate for any advice on how I can resolve this issue when using larger number of processors. Thank you very much for your time.


Masato Koizumi (1.81 KB)

1 Like

As you said, this kind of error occurs all the time, usually due to bad physics in your model.
If you did not see the error when you ran it on on only 8 CPUs (or 1 CPU), it does not mean that your simulation files were correct.
You see the error now because errors like this are easier to detect when you run the simulation in parallel with a large number of processes.

Physics erros are frequently due to large forces between particles that are too close together. (To test this, try turning off the pair forces (for example using “pair_coeff * * 0.0 1.0”, or “neigh_modify exclude group all all”, or “pair_style none”). Alternatively, in polymer systems, physics errors also frequently due to singularities in the dihedral angle force which occur whenever any of the bond-angles approach 180 degrees. (To test this, try using “dihedral_style none”, or removing the “Dihedrals” section from your DATA file, and at the beginning of that file, setting the number of dihedrals to 0.)

Either way, in my experience, to solve these kinds of problems, it was necessary to do a better job minimizing the system.

Minimizing tightly packed polymer melts is especially difficult because they tend to be topologically frustrated. It can take a while to get one of these systems working.

Firstly, if you haven’t yet, try experimenting with different minimizers and minimization settings. Check out:

Keep in mind, this might not work. In fact, I have given up using the existing LAMMPS minimizers for polymer systems. Instead, I typically run a series of highly-damped Langevin dynamics simulations. (Ie, using “fix langevin” with a small “damp” parameter, and a small timestep, together with “fix nve” or “fix nve/limit”)

If that still does not work, then try turning off (“softening”) the repulsion between non-bonded atoms. This will allow the atoms (chains) to pass through each other. Run the simulation for a short time, and then try slowly increasing the repulsion between atoms again. This can help because repulsion between atoms can prevent chain relaxation if the polymer is tied in a knot, or jammed between tightly packed molecules.

There are several different ways to do this:

(You can also try using this code, but it does not work on GPUs, and it may be confusing to use. You can download it here.)

Incidentally “fix langevin” can be used together with “fix nve/limit” (instead of “fix nve”).
Together, these two fixes are very useful for minimization.
(Unfortunately, if memory serves correctly, “fix nve/limit” does not yet work with “fix rigid”)

Hope this points you in the right direction.


Dear Mr. Andrew Jewett,

Hello Mr. Jewett. I deeply am thankful for your meticulous advice. Thank you for pointing out that using larger number of parallel processors detects these errors more explicitly.

I have found out that my polymer melt itself does not experience any physics errors. The physics error occurs when crosslinking occurs using the fix bond/create:

fix 4 all bond/create 100 3 3 1.0 2 iparam 1 5 jparam 1 5
fix 5 all bond/create 100 4 4 1.0 2 iparam 1 6 jparam 1 6

Although my neighbor list is being updated more frequently than when the crosslinking occurs, it is evident that some bad dynamics, as you mentioned, occurs when the bonds are being created.

I am looking into OVITO to visually track down these bad dynamics. But I would like to express my gratitude for your detailed and valuable advice on LAMMPS during your busy times.


Masato Koizumi