Fixed columnwise sums to remove compile-flag-dependent differences in results.

core/src/superpose.cc

@@ -48,6 +48,17 @@ inline Real GetRMSD(const promod3::core::EMatX& pos_one,
   return std::sqrt(rmsd/pos_one.rows());
 }
 
+// calculate average row of the point set
+inline promod3::core::ERVec3 AverageRow(promod3::core::EMatX& pos) {
+  // return pos.colwise().sum()/pos.rows();
+  // NOTE: the colwise-operator is opt-flag-dependent and we hence do it manually
+  promod3::core::ERVec3 avg_row = promod3::core::ERVec3::Zero();
+  for (uint i = 0; i < pos.rows(); ++i) {
+    avg_row += pos.row(i);
+  }
+  return avg_row / pos.rows();
+}
+
 }
 
 
@@ -61,17 +72,13 @@ Real SuperposedRMSD(EMatX& pos_one, EMatX& pos_two){
     throw promod3::Error(ss.str());
   }
 
-  ERVec3 avg_one = pos_one.colwise().sum()/pos_one.rows();
-  ERVec3 avg_two = pos_two.colwise().sum()/pos_two.rows();
+  ERVec3 avg_one = AverageRow(pos_one);
+  ERVec3 avg_two = AverageRow(pos_two);
 
   // SVD only determines the rotational component of the superposition
-  // we need to manually shift the centers of the two point sets on onto 
-  // origin
+  // we need to shift the centers of the two point sets onto origin
   for (uint i = 0; i < pos_one.rows(); ++i){
     pos_one.row(i) -= avg_one;
-  }
-
-  for (uint i = 0; i < pos_two.rows(); ++i){
     pos_two.row(i) -= avg_two;
   }
 
@@ -88,17 +95,13 @@ geom::Mat4 MinRMSDSuperposition(EMatX& pos_one, EMatX& pos_two){
     throw promod3::Error(ss.str());
   }
 
-  ERVec3 avg_one = pos_one.colwise().sum()/pos_one.rows();
-  ERVec3 avg_two = pos_two.colwise().sum()/pos_two.rows();
+  ERVec3 avg_one = AverageRow(pos_one);
+  ERVec3 avg_two = AverageRow(pos_two);
 
   // SVD only determines the rotational component of the superposition
-  // we need to manually shift the centers of the two point sets on onto 
-  // origin
+  // we need to shift the centers of the two point sets onto origin
   for (uint i = 0; i < pos_one.rows(); ++i){
     pos_one.row(i) -= avg_one;
-  }
-
-  for (uint i = 0; i < pos_two.rows(); ++i){
     pos_two.row(i) -= avg_two;
   }
 
@@ -127,17 +130,13 @@ std::pair<geom::Mat4,Real> Superposition(EMatX& pos_one, EMatX& pos_two){
     throw promod3::Error(ss.str());
   }
 
-  ERVec3 avg_one = pos_one.colwise().sum()/pos_one.rows();
-  ERVec3 avg_two = pos_two.colwise().sum()/pos_two.rows();
+  ERVec3 avg_one = AverageRow(pos_one);
+  ERVec3 avg_two = AverageRow(pos_two);
 
   // SVD only determines the rotational component of the superposition
-  // we need to manually shift the centers of the two point sets on onto 
-  // origin
+  // we need to shift the centers of the two point sets onto origin
   for (uint i = 0; i < pos_one.rows(); ++i){
     pos_one.row(i) -= avg_one;
-  }
-
-  for (uint i = 0; i < pos_two.rows(); ++i){
     pos_two.row(i) -= avg_two;
   }