Changeset e6a2ef for src/Potentials/PartialNucleiChargeFitter.cpp

Timestamp:

Mar 24, 2017, 10:12:23 AM (9 years ago)

Author:

Frederik Heber <heber@…>

Branches:

FitPartialCharges_GlobalError

Children:

fa2c89

Parents:

a3d08c

git-author:

Frederik Heber <heber@…> (10/10/16 12:34:32)

git-committer:

Frederik Heber <heber@…> (03/24/17 10:12:23)

Message:

Extracted PartialNucleiChargeFitter::calculateResiduum().

File:

• src/Potentials/PartialNucleiChargeFitter.cpp (modified) (4 diffs)

src/Potentials/PartialNucleiChargeFitter.cpp

@@ -305 +305 @@
 }
 
-double PartialNucleiChargeFitter::operator()()
+VectorContent PartialNucleiChargeFitter::calculateResiduum()
+{
+  constructMatrix();
+
+  // calculate residual vector
+  VectorContent residuum = (*PotentialFromCharges) * (*PartialCharges) - SampledPotential;
+
+  return residuum;
+}
+
+void PartialNucleiChargeFitter::prepareCharges(const size_t _size)
 {
   // prepare PartialCharges
@@ -312 +322 @@
     PartialCharges = NULL;
   }
-  PartialCharges = new VectorContent(positions.size());
+  PartialCharges = new VectorContent(_size);
+}
+
+double PartialNucleiChargeFitter::operator()()
+{
+  prepareCharges(positions.size());
 
   // set up over-determined system's problem matrix A for Ax=b
@@ -327 +342 @@
   LOG(2, "DEBUG: Norm of right-hand side is " << SampledPotential.Norm());
 
-  // need to reconstruct matrix as it was overwritten by SVD
-  constructMatrix();
-
-  // calculate residual vector
-  VectorContent residuum = (*PotentialFromCharges) * (*PartialCharges) - SampledPotential;
+  // calculate residuum (forces matrix reconstruction)
+  VectorContent residuum = calculateResiduum();
 
 #ifndef NDEBUG
@@ -367 +379 @@
 
   return residuum.Norm();
+}
+
+bool PartialNucleiChargeFitter::setCharges(const charges_t &_charges)
+{
+  // check sizes
+  if (positions.size() != _charges.size()) {
+    return false;
+  }
+  // (re-)allocate memory
+  prepareCharges(positions.size());
+  // and place charges in vector
+  for(size_t i=0;i<_charges.size();++i)
+    (*PartialCharges)[i] = _charges[i];
+
+  return true;
 }