Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Add single/improper to create_bonds #1706

Merged
merged 3 commits into from Oct 8, 2019
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Jump to
Jump to file
Failed to load files.
Diff view
Diff view
73 changes: 43 additions & 30 deletions doc/src/create_bonds.txt
Expand Up @@ -23,11 +23,14 @@ style = {many} or {single/bond} or {single/angle} or {single/dihedral} :ule,l
btype = bond type of new bond
batom1,batom2 = atom IDs for two atoms in bond
{single/angle} args = atype aatom1 aatom2 aatom3
atype = bond type of new angle
atype = angle type of new angle
aatom1,aatom2,aatom3 = atom IDs for three atoms in angle
{single/dihedral} args = dtype datom1 datom2 datom3 datom4
dtype = bond type of new dihedral
datom1,datom2,datom3,datom4 = atom IDs for four atoms in dihedral :pre
dtype = dihedral type of new dihedral
datom1,datom2,datom3,datom4 = atom IDs for four atoms in dihedral
{single/improper} args = itype iatom1 iatom2 iatom3 iatom4
itype = improper type of new improper
iatom1,iatom2,iatom3,iatom4 = atom IDs for four atoms in improper :pre
zero or more keyword/value pairs may be appended :l
keyword = {special} :l
{special} value = {yes} or {no} :pre
Expand All @@ -38,51 +41,54 @@ keyword = {special} :l
create_bonds many all all 1 1.0 1.2
create_bonds many surf solvent 3 2.0 2.4
create_bonds single/bond 1 1 2
create_bonds single/angle 5 52 98 107 special no :pre
create_bonds single/angle 5 52 98 107 special no
create_bonds single/dihedral 2 4 19 27 101
create_bonds single/improper 3 23 26 31 57 :pre

[Description:]

Create bonds between pairs of atoms that meet a specified distance
criteria. Or create a single bond, angle, or dihedral between 2, 3,
criteria. Or create a single bond, angle, dihedral or improper between 2, 3,
or 4 specified atoms.

The new bond (angle, dihedral) interactions will then be computed
during a simulation by the bond (angle, dihedral) potential defined by
The new bond (angle, dihedral, improper) interactions will then be computed
during a simulation by the bond (angle, dihedral, improper) potential defined by
the "bond_style"_bond_style.html, "bond_coeff"_bond_coeff.html,
"angle_style"_angle_style.html, "angle_coeff"_angle_coeff.html,
"dihedral_style"_dihedral_style.html,
"dihedral_coeff"_dihedral_coeff.html commands.
"dihedral_coeff"_dihedral_coeff.html, "improper_style"_improper_style.html,
"improper_coeff"_improper_coeff.html commands.

The {many} style is useful for adding bonds to a system, e.g. between
nearest neighbors in a lattice of atoms, without having to enumerate
all the bonds in the data file read by the "read_data"_read_data.html
command.

The {single} styles are useful for adding bonds, angles, dihedrals
The {single} styles are useful for adding bonds, angles, dihedrals, impropers
to a system incrementally, then continuing a simulation.

Note that this command does not auto-create any angle or dihedral
Note that this command does not auto-create any angle, dihedral or improper
interactions when a bond is added. Nor does it auto-create any bonds
when an angle or dihedral is added. Or auto-create any angles when a
dihedral is added. Thus the flexibility of this command is limited.
when an angle, dihedral or improper is added. Or auto-create any angles when a
dihedral or improper is added. Thus the flexibility of this command is limited.
It can be used several times to create different types of bond at
different distances. But it cannot typically auto-create all the
bonds or angles or dihedral that would normally be defined in a data
file for a complex system of molecules.
bonds or angles or dihedrals or impropers that would normally be defined in a
data file for a complex system of molecules.

NOTE: If the system has no bonds (angles, dihedrals) to begin with, or
if more bonds per atom are being added than currently exist, then you
NOTE: If the system has no bonds (angles, dihedrals, impropers) to begin with,
or if more bonds per atom are being added than currently exist, then you
must insure that the number of bond types and the maximum number of
bonds per atom are set to large enough values. And similarly for
angles and dihedrals. Otherwise an error may occur when too many
bonds (angles, dihedrals) are added to an atom. If the
angles, dihedrals and impropers. Otherwise an error may occur when too many
bonds (angles, dihedrals, impropers) are added to an atom. If the
"read_data"_read_data.html command is used to define the system, these
parameters can be set via the "bond types" and "extra bond per atom"
fields in the header section of the data file. If the
"create_box"_create_box.html command is used to define the system,
these 2 parameters can be set via its optional "bond/types" and
"extra/bond/per/atom" arguments. And similarly for angles and
dihedrals. See the doc pages for these 2 commands for details.
"extra/bond/per/atom" arguments. And similarly for angles, dihedrals and
impropers. See the doc pages for these 2 commands for details.

:line

Expand Down Expand Up @@ -137,18 +143,25 @@ ordered linearly within the angle; the central atom is {aatom2}.
{Atype} must be a value between 1 and the number of angle types
defined.

The {single/dihedral} style creates a single dihedral of type {btype}
between two atoms with IDs {batom1} and {batom2}. The ordering of the
atoms is the same as in the {Dihedrals} section of a data file read by
the "read_data"_read_data.html command. I.e. the 4 atoms are ordered
linearly within the dihedral. {Dtype} must be a value between 1 and
The {single/dihedral} style creates a single dihedral of type {dtype}
between four atoms with IDs {datom1}, {datom2}, {datom3}, and {datom4}. The
ordering of the atoms is the same as in the {Dihedrals} section of a data file
read by the "read_data"_read_data.html command. I.e. the 4 atoms are ordered
linearly within the dihedral. {dtype} must be a value between 1 and
the number of dihedral types defined.

The {single/improper} style creates a single improper of type {itype}
between four atoms with IDs {iatom1}, {iatom2}, {iatom3}, and {iatom4}. The
ordering of the atoms is the same as in the {Impropers} section of a data file
read by the "read_data"_read_data.html command. I.e. the 4 atoms are ordered
linearly within the improper. {itype} must be a value between 1 and
the number of improper types defined.

:line

The keyword {special} controls whether an internal list of special
bonds is created after one or more bonds, or a single angle or
dihedral is added to the system.
bonds is created after one or more bonds, or a single angle, dihedral or
improper is added to the system.

The default value is {yes}. A value of {no} cannot be used
with the {many} style.
Expand All @@ -161,16 +174,16 @@ see the "special_bonds"_special_bonds.html command for details.
Thus if you are adding a few bonds or a large list of angles all at
the same time, by using this command repeatedly, it is more efficient
to only trigger the internal list to be created once, after the last
bond (or angle, or dihedral) is added:
bond (or angle, or dihedral, or improper) is added:

create_bonds single/bond 5 52 98 special no
create_bonds single/bond 5 73 74 special no
create_bonds single/bond 5 73 74 special no
...
create_bonds single/bond 5 17 386 special no
create_bonds single/bond 4 112 183 special yes :pre

Note that you MUST insure the internal list is re-built after the last
bond (angle, dihedral) is added, before performing a simulation.
bond (angle, dihedral, improper) is added, before performing a simulation.
Otherwise pairwise interactions will not be properly excluded or
weighted. LAMMPS does NOT check that you have done this correctly.

Expand Down
108 changes: 106 additions & 2 deletions src/create_bonds.cpp
Expand Up @@ -12,7 +12,9 @@
------------------------------------------------------------------------- */

/* ----------------------------------------------------------------------
Contributing authors: Mike Salerno (NRL) added single methods
Contributing authors:
Mike Salerno (NRL) added single methods
Thomas Farmer (ISIS) added single/improper
------------------------------------------------------------------------- */

#include "create_bonds.h"
Expand All @@ -31,7 +33,7 @@

using namespace LAMMPS_NS;

enum{MANY,SBOND,SANGLE,SDIHEDRAL};
enum{MANY,SBOND,SANGLE,SDIHEDRAL,SIMPROPER};

/* ---------------------------------------------------------------------- */

Expand Down Expand Up @@ -100,6 +102,18 @@ void CreateBonds::command(int narg, char **arg)
(datom2 == datom3) || (datom2 == datom4) || (datom3 == datom4))
error->all(FLERR,"Illegal create_bonds command");
iarg = 6;
} else if (strcmp(arg[0],"single/improper") == 0) {
style = SIMPROPER;
if (narg < 6) error->all(FLERR,"Illegal create_bonds command");
dtype = force->inumeric(FLERR,arg[1]);
datom1 = force->tnumeric(FLERR,arg[2]);
datom2 = force->tnumeric(FLERR,arg[3]);
datom3 = force->tnumeric(FLERR,arg[4]);
datom4 = force->tnumeric(FLERR,arg[5]);
if ((datom1 == datom2) || (datom1 == datom3) || (datom1 == datom4) ||
(datom2 == datom3) || (datom2 == datom4) || (datom3 == datom4))
error->all(FLERR,"Illegal create_bonds command");
iarg = 6;
} else error->all(FLERR,"Illegal create_bonds command");

// optional args
Expand Down Expand Up @@ -132,6 +146,9 @@ void CreateBonds::command(int narg, char **arg)
} else if (style == SDIHEDRAL) {
if (dtype <= 0 || dtype > atom->ndihedraltypes)
error->all(FLERR,"Invalid dihedral type in create_bonds command");
} else if (style == SIMPROPER) {
if (dtype <= 0 || dtype > atom->nimpropertypes)
error->all(FLERR,"Invalid improper type in create_bonds command");
}

// invoke creation method
Expand All @@ -140,6 +157,7 @@ void CreateBonds::command(int narg, char **arg)
else if (style == SBOND) single_bond();
else if (style == SANGLE) single_angle();
else if (style == SDIHEDRAL) single_dihedral();
else if (style == SIMPROPER) single_improper();

// trigger special list build

Expand Down Expand Up @@ -512,3 +530,89 @@ void CreateBonds::single_dihedral()
num_dihedral[m]++;
}
}

/* ---------------------------------------------------------------------- */

void CreateBonds::single_improper()
{
int m;

// check that 4 atoms exist

const int nlocal = atom->nlocal;
const int idx1 = atom->map(datom1);
const int idx2 = atom->map(datom2);
const int idx3 = atom->map(datom3);
const int idx4 = atom->map(datom4);

int count = 0;
if ((idx1 >= 0) && (idx1 < nlocal)) count++;
if ((idx2 >= 0) && (idx2 < nlocal)) count++;
if ((idx3 >= 0) && (idx3 < nlocal)) count++;
if ((idx4 >= 0) && (idx4 < nlocal)) count++;

int allcount;
MPI_Allreduce(&count,&allcount,1,MPI_INT,MPI_SUM,world);
if (allcount != 4)
error->all(FLERR,"Create_bonds single/improper atoms do not exist");

// create bond once or 4x if newton_bond set

int *num_improper = atom->num_improper;
int **improper_type = atom->improper_type;
tagint **improper_atom1 = atom->improper_atom1;
tagint **improper_atom2 = atom->improper_atom2;
tagint **improper_atom3 = atom->improper_atom3;
tagint **improper_atom4 = atom->improper_atom4;

if ((m = idx2) >= 0) {
if (num_improper[m] == atom->improper_per_atom)
error->one(FLERR,
"New improper exceeded impropers per atom in create_bonds");
improper_type[m][num_improper[m]] = dtype;
improper_atom1[m][num_improper[m]] = datom1;
improper_atom2[m][num_improper[m]] = datom2;
improper_atom3[m][num_improper[m]] = datom3;
improper_atom4[m][num_improper[m]] = datom4;
num_improper[m]++;
}
atom->nimpropers++;

if (force->newton_bond) return;

if ((m = idx1) >= 0) {
if (num_improper[m] == atom->improper_per_atom)
error->one(FLERR,
"New improper exceeded impropers per atom in create_bonds");
improper_type[m][num_improper[m]] = dtype;
improper_atom1[m][num_improper[m]] = datom1;
improper_atom2[m][num_improper[m]] = datom2;
improper_atom3[m][num_improper[m]] = datom3;
improper_atom4[m][num_improper[m]] = datom4;
num_improper[m]++;
}

if ((m = idx3) >= 0) {
if (num_improper[m] == atom->improper_per_atom)
error->one(FLERR,
"New improper exceeded impropers per atom in create_bonds");
improper_type[m][num_improper[m]] = dtype;
improper_atom1[m][num_improper[m]] = datom1;
improper_atom2[m][num_improper[m]] = datom2;
improper_atom3[m][num_improper[m]] = datom3;
improper_atom4[m][num_improper[m]] = datom4;
num_improper[m]++;
}

if ((m = idx4) >= 0) {
if (num_improper[m] == atom->improper_per_atom)
error->one(FLERR,
"New improper exceeded impropers per atom in create_bonds");
improper_type[m][num_improper[m]] = dtype;
improper_atom1[m][num_improper[m]] = datom1;
improper_atom2[m][num_improper[m]] = datom2;
improper_atom3[m][num_improper[m]] = datom3;
improper_atom4[m][num_improper[m]] = datom4;
num_improper[m]++;
}
}
13 changes: 13 additions & 0 deletions src/create_bonds.h
Expand Up @@ -39,6 +39,7 @@ class CreateBonds : protected Pointers {
void single_bond();
void single_angle();
void single_dihedral();
void single_improper();
};

}
Expand Down Expand Up @@ -87,6 +88,10 @@ E: Invalid dihedral type in create_bonds command

UNDOCUMENTED

E: Invalid improper type in create_bonds command

UNDOCUMENTED

E: Create_bonds requires a pair style be defined

Self-explanatory.
Expand Down Expand Up @@ -135,4 +140,12 @@ E: New dihedral exceeded dihedrals per atom in create_bonds

UNDOCUMENTED

E: Create_bonds single/improper atoms do not exist

UNDOCUMENTED

E: New improper exceeded impropers per atom in create_bonds

UNDOCUMENTED

*/