ICODE-ADC/lammps
4
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

Solved the dihedral index problem. It was actually the indices used.

Browse Source
This commit is contained in:
Germain Clavier
2022-04-22 00:20:25 +02:00
parent 05a845fafd
commit 0b0947c167
1 changed files with 11 additions and 8 deletions
Download Patch File Download Diff File Expand all files Collapse all files

19
src/EXTRA-COMPUTE/compute_born_matrix.cpp
View File

@ -294,7 +294,7 @@ void ComputeBornMatrix::init()
// need an occasional half neighbor list
neighbor->add_request(this, NeighConst::REQ_OCCASIONAL);
neighbor->add_request(this, NeighConst::REQ_FULL | NeighConst::REQ_OCCASIONAL);
// int irequest = neighbor->request((void *) this);
// neighbor->requests[irequest]->pair = 0;
// neighbor->requests[irequest]->compute = 1;
@ -353,10 +353,12 @@ void ComputeBornMatrix::compute_vector()
// compute Born contribution
if (pairflag) compute_pairs();
if (bondflag) compute_bonds();
if (angleflag) compute_angles();
if (dihedflag) compute_dihedrals();
if (neighbor->ago > 0) {
if (pairflag) compute_pairs();
if (bondflag) compute_bonds();
if (angleflag) compute_angles();
if (dihedflag) compute_dihedrals();
}
// sum Born contributions over all procs
@ -1046,8 +1048,9 @@ void ComputeBornMatrix::compute_dihedrals()
for (atom2 = 0; atom2 < nlocal; atom2++) {
if (!(mask[atom2] & groupbit)) continue;
if (molecular == 1)
nd = num_dihedral[atom2];
// if (molecular == 1)
// nd = num_dihedral[atom2];
if (molecular == Atom::MOLECULAR) nd = num_dihedral[atom2];
else {
if (molindex[atom2] < 0) continue;
imol = molindex[atom2];
@ -1081,7 +1084,7 @@ void ComputeBornMatrix::compute_dihedrals()
// dU(t)/de = dt/dcos(t)*dU(t)/dt*dcos(t)/de
// with dt/dcos(t) = -1/sin(t)
dihedral->born_matrix(nd, atom1, atom2, atom3, atom4, dudih, du2dih);
dihedral->born_matrix(i, atom1, atom2, atom3, atom4, dudih, du2dih);
vb1x = x[atom2][0] - x[atom1][0];
vb1y = x[atom2][1] - x[atom1][1];