Changeset a8961d

Timestamp:

Jun 21, 2017, 8:01:27 AM (8 years ago)

Author:

Frederik Heber <frederik.heber@…>

Branches:

IndependentFragmentGrids_Sole_NN_Calculation

Parents:

3fb9ab

git-author:

Frederik Heber <heber@…> (08/19/16 05:27:39)

git-committer:

Frederik Heber <frederik.heber@…> (06/21/17 08:01:27)

Message:

Fixing long-range energy calculations and contribution summation.

• we now calculate e-e+e-n and n-n separately. From the sampled_potential we recover e-e, and vmg delivers e-n. The second calculation yields n-n.
• we need to sum forces from electronic contribution and from second run with nuclei contribution.
• we need to sum the potentials for later fitting partial charges.

Location:

src

Files:

• Actions/FragmentationAction/AnalyseFragmentationResultsAction.cpp (modified) (1 diff)
• Actions/FragmentationAction/FragmentationAutomationAction.cpp (modified) (5 diffs)
• Fragmentation/Summation/Containers/FragmentationLongRangeResults.cpp (modified) (3 diffs)
• Jobs/InterfaceVMGJob.cpp (modified) (2 diffs)

src/Actions/FragmentationAction/AnalyseFragmentationResultsAction.cpp

@@ -392 +392 @@
       value.potential_distribution =
           boost::fusion::at_key<VMGDataFused::both_sampled_potential>(potentialiter->second.second); // contributions
-//      // and re-average to zero (integral is times volume_element which we don't need here)
-//      const double sum =
-//          std::accumulate(
-//              value.potential_distribution.sampled_grid.begin(),
-//              value.potential_distribution.sampled_grid.end(),
-//              0.);
-//      const double offset = sum/(double)value.potential_distribution.sampled_grid.size();
-//      for (SamplingGrid::sampledvalues_t::iterator iter = value.potential_distribution.sampled_grid.begin();
-//          iter != value.potential_distribution.sampled_grid.end();
-//          ++iter)
-//        *iter -= offset;
+      value.potential_distribution +=
+          boost::fusion::at_key<VMGDataFused::sampled_potential>(potentialiter->second.second); // contributions
+      // and re-average to zero (integral is times volume_element which we don't need here)
+      const double sum =
+          std::accumulate(
+              value.potential_distribution.sampled_grid.begin(),
+              value.potential_distribution.sampled_grid.end(),
+              0.);
+      const double offset = sum/(double)value.potential_distribution.sampled_grid.size();
+      for (SamplingGrid::sampledvalues_t::iterator iter = value.potential_distribution.sampled_grid.begin();
+          iter != value.potential_distribution.sampled_grid.end();
+          ++iter)
+        *iter -= offset;
 #else
       ELOG(2, "Long-range information in homology desired but long-range analysis capability not compiled in.");

src/Actions/FragmentationAction/FragmentationAutomationAction.cpp

@@ -170 +170 @@
     )
 {
-  ASSERT( _shortrangedata.size() == _longrangedata.size(),
-      "calculateNucleiNucleiLongRangeContribution() - shortrange and longrange data have inequal size.");
+  // we assume here that the fragmentjobs in short- and longrange case are in the
+  // same order and that the full jobs per level are appended to the end of
+  // longrangedata
   std::map<JobId_t, VMGData>::iterator longrangeiter = _longrangedata.begin();
   std::map<JobId_t, MPQCData>::const_iterator iter = _shortrangedata.begin();
@@ -178 +179 @@
     VMGData &longrange_data = longrangeiter->second;
     // set long-range contributions to zero
+    const double old_energy = longrange_data.nuclei_long;
     longrange_data.nuclei_long = 0.;
     longrange_data.forces.clear();
-    longrange_data.forces.resize(data.positions.size(), FragmentForces::force_t(3,0.));
+    longrange_data.forces.resize(data.positions.size(), FragmentForces::force_t(NDIM,0.));
     longrange_data.hasForces = true;
     // go through positions and evaluate sum naively
@@ -206 +208 @@
         const Vector otherposition((*otherpositer)[0], (*otherpositer)[1], (*otherpositer)[2]);
         const Vector distance = position - otherposition;
-        const double invsqrdist = 1./distance.Norm();
-        const double factor = (*otherchargeiter)*invsqrdist*invsqrdist*invsqrdist;
+        const double invdist = 1./distance.Norm();
+        const double factor = (*otherchargeiter)*invdist*invdist*invdist;
         (*forceiter)[0] = factor*distance[0];
         (*forceiter)[1] = factor*distance[1];
         (*forceiter)[2] = factor*distance[2];
-        longrange_data.nuclei_long += 0.5*(*chargeiter)*(*otherchargeiter)*invsqrdist;
-      }
-    }
+        longrange_data.nuclei_long += 0.5*(*chargeiter)*(*otherchargeiter)*invdist;
+      }
+    }
+    LOG(2, "DEBUG: fragment #" << longrangeiter->first << ": old is " << old_energy
+        << ", new is " << longrange_data.nuclei_long);
   }
 }
@@ -423 +427 @@
           interpolation_degree,
           VMGFragmentController::DoSampleParticles,
-          VMGFragmentController::DoTreatGrid,
+          VMGFragmentController::DontTreatGrid,
           params.DoValenceOnly.get() ? MPQCData::DoSampleValenceOnly : MPQCData::DontSampleValenceOnly,
           params.DoPrintDebug.get(),
@@ -468 +472 @@
         destiter->second.both_sampled_potential = srciter->second.sampled_potential;
         destiter->second.nuclei_long = srciter->second.nuclei_long;
-        destiter->second.forces = srciter->second.forces;
-        destiter->second.hasForces = srciter->second.hasForces;
+        destiter->second.forces += srciter->second.forces;
+        destiter->second.hasForces &= srciter->second.hasForces;
       }
     }

src/Fragmentation/Summation/Containers/FragmentationLongRangeResults.cpp

@@ -121 +121 @@
         Result_perIndexSet_Grid);
 
-    // multiply each short-range potential with the respective charge
+    // multiply each short-range e-e potential with the respective charge
     std::map<JobId_t,MPQCData>::const_iterator mpqciter = fragmentData.begin();
     std::map<JobId_t,VMGData>::iterator vmgiter = longrangeData.begin();
@@ -214 +214 @@
 #endif
 
-      // then, we obtain the e-n+n-n full solution in the same way
+      // then, we obtain the n-n full solution directly
       double nuclei_solution_energy = fullsolutionData[level-1].nuclei_long;
       double nuclei_short_range_energy =
@@ -229 +229 @@
       electron_solution_energy *= .5;
       electron_short_range_energy *= .5;
-
-      // At last, we subtract e-n from n-n+e-n for full solution and short-range
-      // correction.
-      nuclei_solution_energy -= both_solution_energy;
-      nuclei_short_range_energy -= both_short_range_energy;
 
       VMGDataLongRangeMap_t instance;

src/Jobs/InterfaceVMGJob.cpp

@@ -116 +116 @@
 
   Grid& grid = multigrid(multigrid.MaxLevel());
-  grid.Clear();
-  //grid.ClearBoundary(); // we don't have a boundary under periodic boundary conditions
+//  grid.Clear();
+//  grid.ClearBoundary(); // we don't have a boundary under periodic boundary conditions
 
   // print debugging info on grid size
@@ -393 +393 @@
   comm.PrintOnce(Debug, "E_total*:       %e", e_long + e_short_peak + e_short_spline - e_self);
 
-  returndata.nuclei_long = e_long;
-  returndata.electron_long = e_long;
+  returndata.nuclei_long = e;   // nuclei needs self-energy subtracted, ...
+  returndata.electron_long = e_long;  // charge grid gets only long-range interaction
 
   // calculate residual