Merge pull request #1645 from jrgissing/bond/react-max_rxn-bug

Bond/react: max_rxn bugfix + restart reaction counts
This commit is contained in:
Axel Kohlmeyer
2019-09-19 15:21:40 -04:00
committed by GitHub
3 changed files with 89 additions and 26 deletions

View File

@ -186,20 +186,25 @@ reacting atoms.
Some atoms in the pre-reacted template that are not reacting may have
missing topology with respect to the simulation. For example, the
pre-reacted template may contain an atom that would connect to the
rest of a long polymer chain. These are referred to as edge atoms, and
are also specified in the map file. When the pre-reaction template
contains edge atoms, not all atoms, bonds, charges, etc. specified in
the reaction templates will be updated. Specifically, topology that
involves only atoms that are 'too near' to template edges will not be
updated. The definition of 'too near the edge' depends on which
interactions are defined in the simulation. If the simulation has
defined dihedrals, atoms within two bonds of edge atoms are considered
'too near the edge.' If the simulation defines angles, but not
dihedrals, atoms within one bond of edge atoms are considered 'too
near the edge.' If just bonds are defined, only edge atoms are
pre-reacted template may contain an atom that, in the simulation, is
currently connected to the rest of a long polymer chain. These are
referred to as edge atoms, and are also specified in the map file. All
pre-reaction template atoms should be linked to a bonding atom, via at
least one path that does not involve edge atoms. When the pre-reaction
template contains edge atoms, not all atoms, bonds, charges, etc.
specified in the reaction templates will be updated. Specifically,
topology that involves only atoms that are 'too near' to template
edges will not be updated. The definition of 'too near the edge'
depends on which interactions are defined in the simulation. If the
simulation has defined dihedrals, atoms within two bonds of edge atoms
are considered 'too near the edge.' If the simulation defines angles,
but not dihedrals, atoms within one bond of edge atoms are considered
'too near the edge.' If just bonds are defined, only edge atoms are
considered 'too near the edge.'
NOTE: Small molecules, i.e. ones that have all their atoms contained
within the reaction templates, never have edge atoms.
Note that some care must be taken when a building a molecule template
for a given simulation. All atom types in the pre-reacted template
must be the same as those of a potential reaction site in the
@ -392,10 +397,11 @@ local command.
[Restart, fix_modify, output, run start/stop, minimize info:]
No information about this fix is written to "binary restart
files"_restart.html, aside from internally-created per-atom
properties. None of the "fix_modify"_fix_modify.html options are
relevant to this fix.
Cumulative reaction counts for each reaction are written to "binary
restart files"_restart.html. These values are associated with the
reaction name (react-ID). Additionally, internally-created per-atom
properties are stored to allow for smooth restarts. None of the
"fix_modify"_fix_modify.html options are relevant to this fix.
This fix computes one statistic for each {react} argument that it
stores in a global vector, of length 'number of react arguments', that
@ -406,8 +412,8 @@ These is 1 quantity for each react argument:
(1) cumulative # of reactions occurred :ul
No parameter of this fix can be used with the {start/stop} keywords of
the "run"_run.html command. This fix is not invoked during "energy
No parameter of this fix can be used with the {start/stop} keywords
of the "run"_run.html command. This fix is not invoked during "energy
minimization"_minimize.html.
When fix bond/react is 'unfixed,' all internally-created groups are
@ -417,18 +423,20 @@ all other fixes that use any group created by fix bond/react.
[Restrictions:]
This fix is part of the USER-MISC package. It is only enabled if
LAMMPS was built with that package. See the "Build
package"_Build_package.html doc page for more info.
LAMMPS was built with that package. See the
"Build package"_Build_package.html doc page for more info.
[Related commands:]
"fix bond/create"_fix_bond_create.html, "fix
bond/break"_fix_bond_break.html, "fix bond/swap"_fix_bond_swap.html,
"fix bond/create"_fix_bond_create.html,
"fix bond/break"_fix_bond_break.html,
"fix bond/swap"_fix_bond_swap.html,
"dump local"_dump.html, "special_bonds"_special_bonds.html
[Default:]
The option defaults are stabilization = no, prob = 1.0, stabilize_steps = 60, update_edges = none
The option defaults are stabilization = no, prob = 1.0, stabilize_steps = 60,
update_edges = none
:line

View File

@ -56,7 +56,6 @@ static const char cite_fix_bond_react[] =
#define BIG 1.0e20
#define DELTA 16
#define MAXLINE 256
#define MAXGUESS 20 // max # of guesses allowed by superimpose algorithm
#define MAXCONARGS 5 // max # of arguments for any type of constraint
@ -87,6 +86,7 @@ FixBondReact::FixBondReact(LAMMPS *lmp, int narg, char **arg) :
MPI_Comm_size(world,&nprocs);
newton_bond = force->newton_bond;
restart_global = 1;
attempted_rxn = 0;
force_reneighbor = 1;
next_reneighbor = -1;
@ -207,7 +207,9 @@ FixBondReact::FixBondReact(LAMMPS *lmp, int narg, char **arg) :
iarg++;
rxn_name[rxn] = arg[iarg++];
int n = strlen(arg[iarg]) + 1;
if (n > MAXLINE) error->all(FLERR,"Reaction name (react-ID) is too long (limit: 256 characters)");
strncpy(rxn_name[rxn],arg[iarg++],n);
int igroup = group->find(arg[iarg++]);
if (igroup == -1) error->all(FLERR,"Could not find fix group ID");
@ -386,6 +388,10 @@ FixBondReact::FixBondReact(LAMMPS *lmp, int narg, char **arg) :
id_fix3 = NULL;
statted_id = NULL;
custom_exclude_flag = 0;
// used to store restart info
set = new Set[nreacts];
memset(set,0,nreacts*sizeof(Set));
}
/* ---------------------------------------------------------------------- */
@ -469,6 +475,7 @@ FixBondReact::~FixBondReact()
delete [] statted_id;
delete [] guess_branch;
delete [] pioneer_count;
delete [] set;
if (group) {
char **newarg;
@ -1209,7 +1216,7 @@ void FixBondReact::superimpose_algorithm()
rxn_by_proc[j] = -1; // corresponds to ghostly
int itemp = 0;
for (int j = 0; j < nprocs; j++)
for (int k = 0; k < local_rxn_count[j]; k++)
for (int k = 0; k < local_rxncounts[j]; k++)
rxn_by_proc[itemp++] = j;
std::random_shuffle(&rxn_by_proc[0],&rxn_by_proc[delta_rxn]);
for (int j = 0; j < nprocs; j++)
@ -3098,6 +3105,42 @@ void FixBondReact::unpack_reverse_comm(int n, int *list, double *buf)
}
}
/* ----------------------------------------------------------------------
write Set data to restart file
------------------------------------------------------------------------- */
void FixBondReact::write_restart(FILE *fp)
{
set[0].nreacts = nreacts;
for (int i = 0; i < nreacts; i++) {
set[i].reaction_count_total = reaction_count_total[i];
int n = strlen(rxn_name[i]) + 1;
strncpy(set[i].rxn_name,rxn_name[i],n);
}
if (me == 0) {
int size = nreacts*sizeof(Set);
fwrite(&size,sizeof(int),1,fp);
fwrite(set,sizeof(Set),nreacts,fp);
}
}
/* ----------------------------------------------------------------------
use selected state info from restart file to restart the Fix
------------------------------------------------------------------------- */
void FixBondReact::restart(char *buf)
{
Set *set_restart = (Set *) buf;
for (int i = 0; i < set_restart[0].nreacts; i++) {
for (int j = 0; j < nreacts; j++) {
if (strcmp(set_restart[i].rxn_name,rxn_name[j]) == 0) {
reaction_count_total[j] = set_restart[i].reaction_count_total;
}
}
}
}
/* ----------------------------------------------------------------------
memory usage of local atom-based arrays
------------------------------------------------------------------------- */

View File

@ -30,6 +30,9 @@ namespace LAMMPS_NS {
class FixBondReact : public Fix {
public:
enum {MAXLINE=256};
FixBondReact(class LAMMPS *, int, char **);
~FixBondReact();
int setmask();
@ -170,6 +173,15 @@ class FixBondReact : public Fix {
void unlimit_bond();
void limit_bond(int);
void dedup_mega_gloves(int); //dedup global mega_glove
virtual void write_restart(FILE *);
virtual void restart(char *buf);
struct Set {
int nreacts;
char rxn_name[MAXLINE];
int reaction_count_total;
};
Set *set;
// DEBUG