Molecule command - special bonds generation error

Hello,

I have been playing around with the molecule command. I created a molecule file for SPC/E water but when I try to use it, I get the error.

“Molecule auto special bond generation overflow”

I believe I have accounted for any extra special bonds when creating my simulation box. (I have varied the number for testing.)

“create_box 2 myRegion bond/types 1 angle/types 1 extra/special/per/atom 1”

Am I missing something obvious here?

Tom

Hello,

I have been playing around with the molecule command. I created a molecule
file for SPC/E water but when I try to use it, I get the error.

"Molecule auto special bond generation overflow"

I believe I have accounted for any extra special bonds when creating my
simulation box. (I have varied the number for testing.)

"create_box 2 myRegion bond/types 1 angle/types 1 extra/special/per/atom 1"

1 special atom is too few. you should need at least 2.

Am I missing something obvious here?

the most obvious omission is to tell us which version of LAMMPS you
are using and/or whether you tested this with the latest patchlevel
release. you may be running into an issue, that has already been
resolved.

axel.

Hello Axel and Tom,

the most obvious omission is to tell us which version of LAMMPS you
are using and/or whether you tested this with the latest patchlevel
release. you may be running into an issue, that has already been
resolved.

axel.

The problem is still here with the 20Oct version.

1) Looking into the sources, I've discovered that molecule.cpp looks the
number of special neighbors in atom->maxspecial,
and create_box command sets it in force->special_extra, so atom->maxspecial
seems to stay unchanged.
Changing create_box line 161 to
      atom->maxspecial = force->inumeric(FLERR,arg[iarg+1]);
as it was in older versions removes the error. Not sure if it's not
breaking something else.

2) Molecule command has to go after create_box in order for auto special
bonds generation to work, as I've figured out from the code.

Vasily

Hello Axel and Tom,

the most obvious omission is to tell us which version of LAMMPS you
are using and/or whether you tested this with the latest patchlevel
release. you may be running into an issue, that has already been
resolved.

axel.

The problem is still here with the 20Oct version.

1) Looking into the sources, I've discovered that molecule.cpp looks the
number of special neighbors in atom->maxspecial,
and create_box command sets it in force->special_extra, so atom->maxspecial
seems to stay unchanged.
Changing create_box line 161 to
      atom->maxspecial = force->inumeric(FLERR,arg[iarg+1]);
as it was in older versions removes the error. Not sure if it's not breaking
something else.

yes, it does break something else. that is why it was changed. please
provide a test input to demonstrates where the current code is
failing.

2) Molecule command has to go after create_box in order for auto special
bonds generation to work, as I've figured out from the code.

please explain and - if possible - provide an example demonstrating
that as well.

thanks,
     axel.

>yes, it does break something else. that is why it was changed. please
>provide a test input to demonstrates where the current code is
>failing.

Any molecule bigger than a dimer seems to cause a failure. An example is
attached.
Error comes from molecule.cpp:1142 where the number of derived special
bonds is checked against atom->maxspecial.
atom->maxspecial is 1,

>> 2) Molecule command has to go after create_box in order for auto special
>> bonds generation to work, as I've figured out from the code.

>please explain and - if possible - provide an example demonstrating
>that as well.

The error comes from the same line molecule.cpp:1142

As far as I can see, both atom->maxspecial and force->special_extra are not
set
before create_box command, hence automatic special bonds for more than a
dimer won't work.

Vasily

trimer.txt (175 Bytes)

moltest.in (438 Bytes)

thanks for the input examples.

the following changes should take care of both issues. i.e. disallow
auto-generated exclusions before the box is defined and preserve both,
atom->maxspecial and force->special_extra

axel.

[[email protected] src]$ git diff master create_box.cpp molecule.cpp
diff --git a/src/create_box.cpp b/src/create_box.cpp
index c7fd98b..3a13109 100644
--- a/src/create_box.cpp
+++ b/src/create_box.cpp
@@ -159,6 +159,7 @@ void CreateBox::command(int narg, char **arg)
     } else if (strcmp(arg[iarg],"extra/special/per/atom") == 0) {
       if (iarg+2 > narg) error->all(FLERR,"Illegal create_box command");
       force->special_extra = force->inumeric(FLERR,arg[iarg+1]);
+ atom->maxspecial += force->special_extra;
       iarg += 2;
     } else error->all(FLERR,"Illegal create_box command");
   }
diff --git a/src/molecule.cpp b/src/molecule.cpp
index 5e9aa29..0febb78 100644
--- a/src/molecule.cpp
+++ b/src/molecule.cpp
@@ -595,6 +595,10 @@ void Molecule::read(int flag)
   // set maxspecial on first pass, so allocate() has a size

   if (bondflag && specialflag == 0) {
+ if (domain->box_exist == 0)
+ error->all(FLERR,"Cannot auto-generate special bonds before "
+ "simulation box is defined");

Axel, thanks a lot!

Works for my scripts.

Vasily