Segementation Fault of any fix using atomstyle template

LAMMPS LAMMPS Beginners
1

Dear all,

I’ve been trying using lammps to simulate a caffeine monohydrate crystal. I use the atomstyle template and define the molecules as in the attached files. However, the lammps collapses immediately with energy minimization or fix viscous. Even if I tried “setforce 0. 0. 0.” it goes into segmentation fault. I don’t think it’s a potential issue because even if I tuned epsilons in lj/cut as 0, it still goes into segmention fault immediately. Nor do I think it is a initial configuration problem because the configuration was converted from a npt simulation result using GROMACS. I tried both lammps/02Aug2023 and lammps/04Feb2025 but neither of them works. Could you please help me find the issue? Thank you in advance!

caff.txt (634 Bytes)
w.txt (62 Bytes)
run.in (1.8 KB)
monohydrate.params (1.6 KB)
cg_slab_start2.lmpdata (114.4 KB)

2

Thank you very much for providing a complete(!) input deck to help tracking down the origin of the segmentation fault.

This looks like a bug in the neighbor list code that does not anticipate molecule files without bonds and thus without special atom information.

Please try the following change to the LAMMPS source (either version), recompile and run again:

  diff --git a/src/npair_bin.cpp b/src/npair_bin.cpp
  index 2c6fbbb49b..0dc9211a15 100644
  --- a/src/npair_bin.cpp
  +++ b/src/npair_bin.cpp
  @@ -199,9 +199,9 @@ void NPairBin<HALF, NEWTON, TRI, SIZE, ATOMONLY>::build(NeighList *list)
                 if (molecular != Atom::ATOMIC) {
                   if (!moltemplate)
                     which = find_special(special[i], nspecial[i], tag[j]);
  -                else if (imol >= 0)
  -                  which = find_special(onemols[imol]->special[iatom], onemols[imol]->nspecial[iatom],
  -                                       tag[j] - tagprev);
  +                else if ((imol >= 0) && (onemols[imol]->special))
  +                  which = find_special(onemols[imol]->special[iatom],
  +                                       onemols[imol]->nspecial[iatom], tag[j] - tagprev);
                   else
                     which = 0;
                   if (which == 0)
  @@ -222,7 +222,7 @@ void NPairBin<HALF, NEWTON, TRI, SIZE, ATOMONLY>::build(NeighList *list)
                 if (molecular != Atom::ATOMIC) {
                   if (!moltemplate)
                     which = find_special(special[i], nspecial[i], tag[j]);
  -                else if (imol >= 0)
  +                else if ((imol >= 0) && (onemols[imol]->special))
                     which = find_special(onemols[imol]->special[iatom], onemols[imol]->nspecial[iatom],
                                          tag[j] - tagprev);
                   else
3

Dear Dr. Kohlmeyer,

Thank you very much! It works perfectly.

After the minimization, I would like to “fix rigid/small” the molecule caffeine, specified in caff.txt, so I tried to use command delete_bonds to remove the bond potentials defined in the template. However, whether I used group all or group caf_atoms (defined as “group caf_atoms type 1:7”), the lammps gave the error of “ERROR: Cannot use delete_bonds with non-molecular system (src/delete_bonds.cpp:47)”. I’m wondering if this is because I used atom_style template rather than molecular, or is it also because of the one-atom molecule. Thank you!

4

You cannot. The template is the template and you cannot modify it and the whole point of atom style template is to avoid storing the information individually. it is not necessary to delete the bonds fix rigid will work with them or without them.

If you want to make your life easier in that respect, just use atom style full.

5

I see. Thank you so much! I really appreciate your help. It solves my confusion!