Merge pull request #2406 from jrgissing/bond/react-rmsd_constraint_bugfix

bond/react: rmsd constraint bugfix
This commit is contained in:
Axel Kohlmeyer 2020-10-05 10:43:21 -04:00 committed by GitHub
commit 47e08f63ac
No known key found for this signature in database
GPG Key ID: 4AEE18F83AFDEB23
1 changed files with 12 additions and 2 deletions

View File

@ -1863,6 +1863,8 @@ int FixBondReact::check_constraints()
if (prrhob < rrhandom[(int) constraints[i][2]]->uniform()) return 0; if (prrhob < rrhandom[(int) constraints[i][2]]->uniform()) return 0;
} else if (constraints[i][1] == RMSD) { } else if (constraints[i][1] == RMSD) {
// call superpose // call superpose
int iatom;
int iref = -1; // choose first atom as reference
int n2superpose = 0; int n2superpose = 0;
double **xfrozen; // coordinates for the "frozen" target molecule double **xfrozen; // coordinates for the "frozen" target molecule
double **xmobile; // coordinates for the "mobile" molecule double **xmobile; // coordinates for the "mobile" molecule
@ -1875,20 +1877,28 @@ int FixBondReact::check_constraints()
int myincr = 0; int myincr = 0;
for (int j = 0; j < onemol->natoms; j++) { for (int j = 0; j < onemol->natoms; j++) {
if (onemol->fragmentmask[ifragment][j]) { if (onemol->fragmentmask[ifragment][j]) {
iatom = atom->map(glove[j][1]);
if (iref == -1) iref = iatom;
iatom = domain->closest_image(iref,iatom);
for (int k = 0; k < 3; k++) { for (int k = 0; k < 3; k++) {
xfrozen[myincr][k] = x[atom->map(glove[j][1])][k]; xfrozen[myincr][k] = x[iatom][k];
xmobile[myincr][k] = onemol->x[j][k]; xmobile[myincr][k] = onemol->x[j][k];
} }
myincr++; myincr++;
} }
} }
} else { } else {
int iatom;
int iref = -1; // choose first atom as reference
n2superpose = onemol->natoms; n2superpose = onemol->natoms;
memory->create(xfrozen,n2superpose,3,"bond/react:xfrozen"); memory->create(xfrozen,n2superpose,3,"bond/react:xfrozen");
memory->create(xmobile,n2superpose,3,"bond/react:xmobile"); memory->create(xmobile,n2superpose,3,"bond/react:xmobile");
for (int j = 0; j < n2superpose; j++) { for (int j = 0; j < n2superpose; j++) {
iatom = atom->map(glove[j][1]);
if (iref == -1) iref = iatom;
iatom = domain->closest_image(iref,iatom);
for (int k = 0; k < 3; k++) { for (int k = 0; k < 3; k++) {
xfrozen[j][k] = x[atom->map(glove[j][1])][k]; xfrozen[j][k] = x[iatom][k];
xmobile[j][k] = onemol->x[j][k]; xmobile[j][k] = onemol->x[j][k];
} }
} }