Context Navigation

-              r8ea8c8
+              r4ec18b
 #include <boost/assign.hpp>
+#include <boost/bind.hpp>
 #include <boost/filesystem.hpp>
+#include <boost/function.hpp>
 #include <boost/program_options.hpp>
 …
   else
     return angle/divisor;
+}
+/** Namespace containing all simple extractor functions.
+ *
+ */
+namespace Extractors {
+  /** Simple extractor of all unique pair distances of a given \a fragment.
+   *
+   * \param fragment fragment with all nuclei positions
+   * \param index index refers to the index within the global set of configurations
+   * \return vector of of argument_t containing all found distances
+   */
+  FunctionModel::arguments_t gatherAllDistances(
+      const Fragment& fragment,
+      const size_t index
+      ) {
+    // get distance out of Fragment
+    const Fragment::charges_t charges = fragment.getCharges();
+    const Fragment::positions_t positions = fragment.getPositions();
+    return gatherAllDistanceArguments(charges, positions, index);
+  }
+  /** Gather first distance for the two matching charges.
+   *
+   * \param fragment fragment with all nuclei positions
+   * \param index index refers to the index within the global set of configurations
+   * \param firstelement first element of pair
+   * \param secondelement second element of pair, order is reflected in indices of return argument_t
+   * \return vector of of argument_t containing all found distances
+   */
+  FunctionModel::arguments_t gatherFirstDistance(
+      const Fragment& fragment,
+      const size_t index,
+      const size_t firstelement,
+      const size_t secondelement
+      ) {
+    const Fragment::charges_t charges = fragment.getCharges();
+    const Fragment::positions_t positions = fragment.getPositions();
+    typedef Fragment::charges_t::const_iterator chargeiter_t;
+    std::vector< chargeiter_t > firstpair;
+    firstpair.reserve(2);
+    firstpair +=
+        std::find(charges.begin(), charges.end(), firstelement),
+        std::find(charges.begin(), charges.end(), secondelement);
+    if ((firstpair[0] == charges.end()) || (firstpair[1] == charges.end())) {
+      // complain if tuple not found
+      ELOG(1, "Could not find pair " << firstelement << "," << secondelement
+          << " in fragment " << fragment);
+      return FunctionModel::arguments_t();
+    }
+    // convert position_t to Vector
+    std::vector< std::pair<Vector, size_t> > DistancePair;
+    for (std::vector<chargeiter_t>::const_iterator firstpairiter = firstpair.begin();
+        firstpairiter != firstpair.end(); ++firstpairiter) {
+      Fragment::positions_t::const_iterator positer = positions.begin();
+      const size_t steps = std::distance(charges.begin(), *firstpairiter);
+      std::advance(positer, steps);
+      DistancePair.push_back(
+          std::make_pair(Vector((*positer)[0], (*positer)[1], (*positer)[2]),
+              steps));
+    }
+    // finally convert Vector pair to distance-like argument
+    argument_t arg;
+    arg.indices.first = DistancePair[0].second;
+    arg.indices.second = DistancePair[1].second;
+    arg.distance = DistancePair[0].first.distance(DistancePair[1].first);
+    arg.globalid = index;
+    return FunctionModel::arguments_t(1, arg);
+  }
+}; /* namespace Extractors */
+/** This class encapsulates the training data for a given potential function
+ * to learn.
+ *
+ * The data is added piece-wise by calling the operator() with a specific
+ * Fragment.
+ */
+class TrainingData
+{
+public:
+  //!> typedef for a range within the HomologyContainer at which fragments to look at
+  typedef std::pair<
+      HomologyContainer::const_iterator,
+      HomologyContainer::const_iterator> range_t;
+  //!> Training tuple input vector pair
+  typedef FunctionApproximation::inputs_t InputVector_t;
+  //!> Training tuple output vector pair
+  typedef FunctionApproximation::outputs_t OutputVector_t;
+  //!> Typedef for a function containing how to extract required information from a Fragment.
+  typedef boost::function< FunctionModel::arguments_t (const Fragment &, const size_t)> extractor_t;
+public:
+  /** Constructor for class TrainingData.
+   *
+   */
+  explicit TrainingData(const extractor_t &_extractor) :
+      extractor(extractor)
+  {}
+  /** Destructor for class TrainingData.
+   *
+   */
+  ~TrainingData()
+  {}
+  /** We go through the given \a range of homologous fragments and call
+   * TrainingData::extractor on them in order to gather the distance and
+   * the energy value, stored internally.
+   *
+   * \param range given range within a HomologyContainer of homologous fragments
+   */
+  void operator()(const range_t &range) {
+    for (HomologyContainer::const_iterator iter = range.first; iter != range.second; ++iter) {
+      // get distance out of Fragment
+      const Fragment &fragment = iter->second.first;
+      FunctionModel::arguments_t args = extractor(
+            fragment,
+            DistanceVector.size()
+          );
+      DistanceVector.push_back( args );
+      const double &energy = iter->second.second;
+      EnergyVector.push_back( FunctionModel::results_t(1, energy) );
+    }
+  }
+  /** Getter for const access to internal training data inputs.
+   *
+   * \return const ref to training tuple of input vector
+   */
+  const InputVector_t& getTrainingInputs() const {
+    return DistanceVector;
+  }
+  /** Getter for const access to internal training data outputs.
+   *
+   * \return const ref to training tuple of output vector
+   */
+  const OutputVector_t& getTrainingOutputs() const {
+    return EnergyVector;
+  }
+private:
+  // prohibit use of default constructor, as we always require extraction functor.
+  TrainingData();
+private:
+  //!> private training data vector
+  InputVector_t DistanceVector;
+  OutputVector_t EnergyVector;
+  //!> function to be used for training input data extraction from a fragment
+  const extractor_t extractor;
+};
+// print training data for debugging
+std::ostream &operator<<(std::ostream &out, const TrainingData &data)
+{
+  const TrainingData::InputVector_t &DistanceVector = data.getTrainingInputs();
+  const TrainingData::OutputVector_t &EnergyVector = data.getTrainingOutputs();
+  out << "(" << DistanceVector.size()
+      << "," << EnergyVector.size() << ") data pairs: ";
+  FunctionApproximation::inputs_t::const_iterator initer = DistanceVector.begin();
+  FunctionApproximation::outputs_t::const_iterator outiter = EnergyVector.begin();
+  for (; initer != DistanceVector.end(); ++initer, ++outiter) {
+    for (size_t index = 0; index < (*initer).size(); ++index)
+       out << "(" << (*initer)[index].indices.first << "," << (*initer)[index].indices.second
+          << ") " << (*initer)[index].distance;
+    out << " with energy " << *outiter;
+  }
+  return out;
+}
 …
     // Afterwards we go through all of this type and gather the distance and the energy value
+    typedef std::pair<
+        FunctionApproximation::inputs_t,
+        FunctionApproximation::outputs_t> InputOutputVector_t;
+    InputOutputVector_t DistanceEnergyVector;
+    std::pair<HomologyContainer::const_iterator, HomologyContainer::const_iterator> range =
+        homologies.getHomologousGraphs(graph);
+    for (HomologyContainer::const_iterator iter = range.first; iter != range.second; ++iter) {
+      // get distance out of Fragment
+      const double &energy = iter->second.second;
+      const Fragment &fragment = iter->second.first;
+      const Fragment::charges_t charges = fragment.getCharges();
+      const Fragment::positions_t positions = fragment.getPositions();
+      std::vector< std::pair<Vector, size_t> > DistanceVectors;
+      for (Fragment::charges_t::const_iterator chargeiter = charges.begin();
+          chargeiter != charges.end(); ++chargeiter) {
+        if (*chargeiter == 6) {
+          Fragment::positions_t::const_iterator positer = positions.begin();
+          const size_t steps = std::distance(charges.begin(), chargeiter);
+          std::advance(positer, steps);
+          DistanceVectors.push_back(
+              std::make_pair(Vector((*positer)[0], (*positer)[1], (*positer)[2]),
+                  steps));
+        }
+      }
+      if (DistanceVectors.size() == (size_t)2) {
+        argument_t arg;
+        arg.indices.first = DistanceVectors[0].second;
+        arg.indices.second = DistanceVectors[1].second;
+        arg.distance = DistanceVectors[0].first.distance(DistanceVectors[1].first);
+        arg.globalid = DistanceEnergyVector.first.size();
+        DistanceEnergyVector.first.push_back( FunctionModel::arguments_t(1,arg) );
+        DistanceEnergyVector.second.push_back( FunctionModel::results_t(1,energy) );
+      } else {
+        ELOG(2, "main() - found not exactly two carbon atoms in fragment "
+            << fragment << ".");
+      }
+    }
+    // print training data for debugging
+    {
+      LOG(1, "INFO: I gathered the following (" << DistanceEnergyVector.first.size()
+          << "," << DistanceEnergyVector.second.size() << ") data pairs: ");
+      FunctionApproximation::inputs_t::const_iterator initer = DistanceEnergyVector.first.begin();
+      FunctionApproximation::outputs_t::const_iterator outiter = DistanceEnergyVector.second.begin();
+      for (; initer != DistanceEnergyVector.first.end(); ++initer, ++outiter) {
+        LOG(1, "INFO: (" << (*initer)[0].indices.first << "," << (*initer)[0].indices.second
+            << ") " << (*initer)[0].distance << " with energy " << *outiter);
+      }
+    }
+    TrainingData MorseData(
+        boost::bind(&Extractors::gatherFirstDistance, _1, _2, 6, 6) // gather first carbon pair
+        );
+    MorseData(homologies.getHomologousGraphs(graph));
+    LOG(1, "INFO: I gathered the following training data: " << MorseData);
     // NOTICE that distance are in bohrradi as they come from MPQC!
 …
     morse.setParameters(params);
     FunctionModel &model = morse;
+    FunctionApproximation approximator(1, 1, model);
+    approximator.setTrainingData(DistanceEnergyVector.first,DistanceEnergyVector.second);
+    FunctionApproximation approximator(
+        MorseData.getTrainingInputs().begin()->size(),
+        MorseData.getTrainingOutputs().begin()->size(),
+        model);
+    approximator.setTrainingData(MorseData.getTrainingInputs(),MorseData.getTrainingOutputs());
     if (model.isBoxConstraint() && approximator.checkParameterDerivatives())
       approximator(FunctionApproximation::ParameterDerivative);
 …
         DistanceEnergyVector.first.begin()->size(),
         DistanceEnergyVector.second.begin()->size(),
         model); // CH4 has 5 atoms, hence 5*4/2 distances
+        model);
     approximator.setTrainingData(DistanceEnergyVector.first,DistanceEnergyVector.second);
     if (model.isBoxConstraint() && approximator.checkParameterDerivatives())

Note: See TracChangeset for help on using the changeset viewer.

Context Navigation

Changeset 4ec18b for src/LevMartester.cpp

Legend:

src/LevMartester.cpp

Download in other formats: