Fix CG minimization with DD
authorBerk Hess <hess@kth.se>
Sat, 28 Sep 2019 21:22:26 +0000 (23:22 +0200)
committerBerk Hess <hess@kth.se>
Sat, 28 Sep 2019 21:23:32 +0000 (23:23 +0200)
A force for determining the CG direction was not communicated
correctly with domain decomposition to recent refactoring change
5ed943bc.

Fixes #3112

Change-Id: Ic33062f4aad0ccd7c1296c7351bc71d735ebba24

src/gromacs/mdrun/minimize.cpp

index ab29f1e3155af3afff16e6f3f440d9dfac824b41..211377835691c55a6215aa4effb9b5bcfe7fa7e1 100644 (file)
@@ -948,11 +948,12 @@ static double reorder_partsum(const t_commrec *cr, t_grpopts *opts,
     rvec      *fmg;
     snew(fmg, natoms);
 
-    gmx::ArrayRef<const int> indicesMin = s_b->s.cg_gl;
+    gmx::ArrayRef<const int> indicesMin = s_min->s.cg_gl;
     int i = 0;
     for (int a : indicesMin)
     {
         copy_rvec(fm[i], fmg[a]);
+        i++;
     }
     gmx_sum(top_global->natoms*3, fmg[0], cr);