This fix requires that the pairwise weightings defined by the -special_bonds command be 0,1,1 for 1-2, 1-3, and -1-4 neighbors within the bond topology. This effectively means that -the pairwise interaction between atoms I and J is turned off when a -bond between them exists and will be turned on when the bond is -broken. It also means that the pairwise interaction of I with J's -other bond partners is unaffected by the existence of the bond. +special_bonds command be *,1,1 for 1-3 and 1-4 +neighbors within the bond topology (the 1-2 setting is not +restricted). This means that the pairwise interaction of I with J's +other bond partners is unaffected by the breaking of the bond.
Related commands:
diff --git a/doc/fix_bond_break.txt b/doc/fix_bond_break.txt index d6ee6bec0e..c3b1a502c2 100644 --- a/doc/fix_bond_break.txt +++ b/doc/fix_bond_break.txt @@ -132,12 +132,10 @@ for these angles and dihedrals, even if your simulation defines an "dihedral_style"_dihedral_style.html. This fix requires that the pairwise weightings defined by the -"special_bonds"_special_bonds.html command be 0,1,1 for 1-2, 1-3, and -1-4 neighbors within the bond topology. This effectively means that -the pairwise interaction between atoms I and J is turned off when a -bond between them exists and will be turned on when the bond is -broken. It also means that the pairwise interaction of I with J's -other bond partners is unaffected by the existence of the bond. +"special_bonds"_special_bonds.html command be *,1,1 for 1-3 and 1-4 +neighbors within the bond topology (the 1-2 setting is not +restricted). This means that the pairwise interaction of I with J's +other bond partners is unaffected by the breaking of the bond. [Related commands:] diff --git a/doc/fix_bond_create.html b/doc/fix_bond_create.html index a521d01914..079bc359e8 100644 --- a/doc/fix_bond_create.html +++ b/doc/fix_bond_create.html @@ -200,12 +200,10 @@ an angle_style or dihedral_style.This fix requires that the pairwise weightings defined by the -special_bonds command be 0,1,1 for 1-2, 1-3, and -1-4 neighbors within the bond topology. This effectively means that -the pairwise interaction between atoms I and J will be turned off when -a bond between them is created. It also means that the pairwise -interaction of I with J's other bond partners will be unaffected by -the new bond. +special_bonds command be *,1,1 for 1-3 and 1-4 +neighbors within the bond topology (the 1-2 setting is not +restricted). This means that the pairwise interaction of I with J's +other bond partners is unaffected by the creation of a new bond.
Related commands:
diff --git a/doc/fix_bond_create.txt b/doc/fix_bond_create.txt index 17fbb77ede..89fb3a0c88 100644 --- a/doc/fix_bond_create.txt +++ b/doc/fix_bond_create.txt @@ -188,12 +188,10 @@ an "angle_style"_angle_style.html or "dihedral_style"_dihedral_style.html. This fix requires that the pairwise weightings defined by the -"special_bonds"_special_bonds.html command be 0,1,1 for 1-2, 1-3, and -1-4 neighbors within the bond topology. This effectively means that -the pairwise interaction between atoms I and J will be turned off when -a bond between them is created. It also means that the pairwise -interaction of I with J's other bond partners will be unaffected by -the new bond. +"special_bonds"_special_bonds.html command be *,1,1 for 1-3 and 1-4 +neighbors within the bond topology (the 1-2 setting is not +restricted). This means that the pairwise interaction of I with J's +other bond partners is unaffected by the creation of a new bond. [Related commands:] diff --git a/src/MC/fix_bond_break.cpp b/src/MC/fix_bond_break.cpp index 53b10ddba5..7f00bf37d5 100644 --- a/src/MC/fix_bond_break.cpp +++ b/src/MC/fix_bond_break.cpp @@ -128,14 +128,15 @@ int FixBondBreak::setmask() void FixBondBreak::init() { - // require special bonds = 0,1,1 + // require special bonds = *,1,1 + // [0] can be anything b/c I,J are removed from each other's special list + // [1],[2] must be 1.0 b/c only special lists of I,J are updated when + // bond I-J is broken, not special lists of neighbors of I,J,etc int flag = 0; - if (force->special_lj[1] != 0.0 || force->special_lj[2] != 1.0 || - force->special_lj[3] != 1.0) flag = 1; - if (force->special_coul[1] != 0.0 || force->special_coul[2] != 1.0 || - force->special_coul[3] != 1.0) flag = 1; - if (flag) error->all(FLERR,"Fix bond/break requires special_bonds = 0,1,1"); + if (force->special_lj[2] != 1.0 || force->special_lj[3] != 1.0) flag = 1; + if (force->special_coul[2] != 1.0 || force->special_coul[3] != 1.0) flag = 1; + if (flag) error->all(FLERR,"Fix bond/break requires special_bonds = *,1,1"); // warn if angles, dihedrals, impropers are being used diff --git a/src/MC/fix_bond_create.cpp b/src/MC/fix_bond_create.cpp index e0ee0b58be..1fa49e25a7 100644 --- a/src/MC/fix_bond_create.cpp +++ b/src/MC/fix_bond_create.cpp @@ -180,16 +180,17 @@ void FixBondCreate::init() if (force->pair == NULL || cutsq > force->pair->cutsq[iatomtype][jatomtype]) error->all(FLERR,"Fix bond/create cutoff is longer than pairwise cutoff"); - // require special bonds = 0,1,1 + // require special bonds = *,1,1 + // [0] can be anything b/c duplicate bond is checked for + // [1],[2] must be 1.0 b/c only special lists of I,J are updated when + // bond I-J is created, not special lists of neighbors of I,J,etc - if (force->special_lj[1] != 0.0 || force->special_lj[2] != 1.0 || - force->special_lj[3] != 1.0) - error->all(FLERR,"Fix bond/create requires special_bonds lj = 0,1,1"); + if (force->special_lj[2] != 1.0 || force->special_lj[3] != 1.0) + error->all(FLERR,"Fix bond/create requires special_bonds lj = *,1,1"); if (atom->q_flag) - if (force->special_coul[1] != 0.0 || force->special_coul[2] != 1.0 || - force->special_coul[3] != 1.0) - error->all(FLERR,"Fix bond/create requires special_bonds coul = 0,1,1"); + if (force->special_coul[2] != 1.0 || force->special_coul[3] != 1.0) + error->all(FLERR,"Fix bond/create requires special_bonds coul = *,1,1"); // warn if angles, dihedrals, impropers are being used @@ -267,7 +268,7 @@ void FixBondCreate::setup(int vflag) void FixBondCreate::post_integrate() { - int i,j,m,ii,jj,inum,jnum,itype,jtype,n1,n3,possible; + int i,j,k,m,ii,jj,inum,jnum,itype,jtype,n1,n3,possible; double xtmp,ytmp,ztmp,delx,dely,delz,rsq; int *ilist,*jlist,*numneigh,**firstneigh; tagint *slist; @@ -309,6 +310,8 @@ void FixBondCreate::post_integrate() double **x = atom->x; tagint *tag = atom->tag; + tagint **bond_atom = atom->bond_atom; + int *num_bond = atom->num_bond; int *mask = atom->mask; int *type = atom->type; @@ -316,8 +319,6 @@ void FixBondCreate::post_integrate() ilist = list->ilist; numneigh = list->numneigh; firstneigh = list->firstneigh; - tagint **bond_atom = atom->bond_atom; - int *num_bond = atom->num_bond; for (ii = 0; ii < inum; ii++) { i = ilist[ii]; @@ -345,23 +346,17 @@ void FixBondCreate::post_integrate() (imaxbond == 0 || bondcount[j] < imaxbond)) possible = 1; } - - // avoid to attempt to add bonds that already exist. we have to - // check both, tag[i]-tag[j] and tag[j]-tag[i] since we do not - // account for newton_bond or newton_pair at this point and - // use a half neighbor list. - if (possible) { - int k; - for (k=0; k < num_bond[i]; ++k) - if (bond_atom[i][k] == tag[j]) - possible = 0; - - for (k=0; k < num_bond[j]; ++k) - if (bond_atom[j][k] == tag[i]) - possible = 0; - } if (!possible) continue; + // do not allow a duplicate bond to be created + // check existing bonds of both I and J + + for (k = 0; k < num_bond[i]; k++) + if (bond_atom[i][k] == tag[j]) possible = 0; + if (j < nlocal) + for (k = 0; k < num_bond[j]; k++) + if (bond_atom[j][k] == tag[i]) possible = 0; + if (!possible) continue; delx = xtmp - x[j][0]; dely = ytmp - x[j][1]; @@ -412,7 +407,6 @@ void FixBondCreate::post_integrate() for (i = 0; i < nlocal; i++) { if (partner[i] == 0) continue; j = atom->map(partner[i]); - if (j < 0) continue; if (partner[j] != tag[i]) continue; // apply probability constraint using RN for atom with smallest ID @@ -454,11 +448,10 @@ void FixBondCreate::post_integrate() // increment bondcount, convert atom to new type if limit reached bondcount[i]++; - if (type[i] == iatomtype) { + if (type[i] == iatomtype) if (bondcount[i] == imaxbond) type[i] = inewtype; - } else { + else if (bondcount[i] == jmaxbond) type[i] = jnewtype; - } // count the created bond once diff --git a/src/read_data.cpp b/src/read_data.cpp index cf35f81007..ef69875c29 100644 --- a/src/read_data.cpp +++ b/src/read_data.cpp @@ -184,10 +184,9 @@ void ReadData::command(int narg, char **arg) } else error->all(FLERR,"Illegal read_data command"); } - // perform 1-pass read if no molecular topoogy in file - // perform 2-pass read if molecular topology - // 1st pass calculates max topology/atom + // perform 2-pass read if molecular topology, + // first pass calculates max topology/atom int atomflag,topoflag; int bondflag,angleflag,dihedralflag,improperflag; @@ -226,10 +225,6 @@ void ReadData::command(int narg, char **arg) atom->allocate_type_arrays(); atom->avec->grow(n); - // initialize here, in case there are no bonds in the data file yet. - atom->bond_per_atom = atom->extra_bond_per_atom; - n = atom->nmax; - domain->print_box(" "); domain->set_initial_box(); domain->set_global_box(); @@ -497,16 +492,19 @@ void ReadData::command(int narg, char **arg) if (!topoflag) break; firstpass = 0; - // reallocate bond,angle,diehdral,improper arrays via grow(), - // using new bond,angle,dihedral,improper per-atom values from 1st pass + // reallocate bond,angle,diehdral,improper arrays via grow() + // use new bond,angle,dihedral,improper per-atom values from 1st pass // should leave other atom arrays unchanged, since already nmax in length + // if bonds/etc not in data file, initialize per-atom size + // with extra settings before grow() of these topology arrays if (bondflag) { memory->destroy(atom->bond_type); memory->destroy(atom->bond_atom); atom->bond_type = NULL; atom->bond_atom = NULL; - } + } else atom->bond_per_atom = atom->extra_bond_per_atom; + if (angleflag) { memory->destroy(atom->angle_type); memory->destroy(atom->angle_atom1); @@ -514,7 +512,8 @@ void ReadData::command(int narg, char **arg) memory->destroy(atom->angle_atom3); atom->angle_type = NULL; atom->angle_atom1 = atom->angle_atom2 = atom->angle_atom3 = NULL; - } + } else atom->angle_per_atom = atom->extra_angle_per_atom; + if (dihedralflag) { memory->destroy(atom->dihedral_type); memory->destroy(atom->dihedral_atom1); @@ -524,7 +523,8 @@ void ReadData::command(int narg, char **arg) atom->dihedral_type = NULL; atom->dihedral_atom1 = atom->dihedral_atom2 = atom->dihedral_atom3 = atom->dihedral_atom4 = NULL; - } + } else atom->dihedral_per_atom = atom->extra_dihedral_per_atom; + if (improperflag) { memory->destroy(atom->improper_type); memory->destroy(atom->improper_atom1); @@ -534,7 +534,7 @@ void ReadData::command(int narg, char **arg) atom->improper_type = NULL; atom->improper_atom1 = atom->improper_atom2 = atom->improper_atom3 = atom->improper_atom4 = NULL; - } + } else atom->improper_per_atom = atom->extra_improper_per_atom; atom->avec->grow(atom->nmax); } diff --git a/src/version.h b/src/version.h index 3f8506ce84..7ca72240cd 100644 --- a/src/version.h +++ b/src/version.h @@ -1 +1 @@ -#define LAMMPS_VERSION "10 Feb 2014" +#define LAMMPS_VERSION "11 Feb 2014"