source: src/Fragmentation/Summation/OrthogonalSummator.hpp@ 14de9a4

/*
 * OrthogonalSummator.hpp
 *
 *  Created on: 28.07.2012
 *      Author: heber
 */

#ifndef ORTHOGONALSUMMATOR_HPP_
#define ORTHOGONALSUMMATOR_HPP_

// include config.h
#ifdef HAVE_CONFIG_H
#include <config.h>
#endif

#include <boost/fusion/sequence.hpp>

#include <iterator>

#include "CodePatterns/Assert.hpp"

#include "Fragmentation/Summation/OrthogonalSummation.hpp"
#include "Fragmentation/Summation/SetValue.hpp"
#include "Fragmentation/Summation/ZeroInstance.hpp"

/** OrthogonalSummator is a general class for making us of OrthogonalSummation.
 *
 * The idea is that we want to sum up not only one value but a whole bunch of.
 * However, in general the types of the values will all differ. We would like
 * to call the summation in a as simple a fashion as possible and at best in
 * a kind of for_each on each type/value.
 * For this purpose we require a list of all the types, represented  by the
 * MPQCDataMap_t map that uses the types enumerated in namespace MPQCDataFused,
 * see the simple example for the associative sequence in boost::fusion.
 * The MPQCDataMap_t then gives the essential mapping from the type to a specific
 * key. Via this key, which is a type and hence can be used for template
 * specification, we may gather all information for launching the
 * OrthogonalSummation. We only need to convert MPQCData into the above
 * MPQCDataMap_t instance and may then access the member variables via the above
 * key and also obtain the type at \b run-time through the key.
 *
 * Note that we then need a conversion from the type stored in the MPQCData,
 * e.g. a vector<double>, to the type in MPQCDataMap_t, e.g. Histogram.
 *
 * Boost::fusion is very cool! If we have a namespace with just structs as keys
 * and the fusion::map from these keys to the true types, then we may do ...
 * \code
 *   MPQCDataMap_t MPQCDataMap;
 *   using namespace MPQCDataFused;
 *   double MPQCData_energy_total = boost::fusion::at_key<MPQCDataFused::energy_total>(MPQCDataMap);
 * \endcode
 * i.e. at_key knows the correct type!
 *
 * Note that you may skip values in the summation by setting their index in
 * _MatrixNrLookup to something outside the range of [0,_data.size()], e.g.
 * -1.
 *
 */
template <typename MapType, typename MapKey>
struct OrthogonalSummator {
  /** We retrieve the type of the MPQCData member variable from the
   * boost::fusion::map and stored it in this typedef. Note that for
   * internal types (e.g. MapValue::const_iterator we need to place
   * \b typename before).
   */
  typedef typename boost::fusion::result_of::value_at_key<MapType, MapKey>::type MapValue;

  /** Constructor for class OrthogonalSummator.
   *
   * \param _subsetmap map with hierarchy of IndexSet's
   * \param _data MPQCData converted to MPQCDataMap_t type, associated to job id
   * \param _container container of IndexSet's such that each set has correct order
   *        to job id and hence to _data.
   * \param _MatrixNrLookup lookup from job id to ordering in above vectors
   */
  OrthogonalSummator(
      SubsetMap::ptr &_subsetmap,
      const std::map<JobId_t, MapType> &_data,
      const IndexSetContainer::Container_t &_container,
      const std::map< JobId_t, size_t > &_MatrixNrLookup) :
    indices(getSubsets(_data.size(),_container)),
    values(createValues(_data, _container, _MatrixNrLookup)),
    OS(indices, values, _subsetmap)
  {
    ASSERT( _data.size() == _MatrixNrLookup.size(),
        "OrthogonalSummator() - ids and MatrixNrLookup don't have same size.");
  }

  /** Summation operator.
   *
   * Initialises instantiated OrthogonalSummation of the respective type via
   * \a OrthogonalSummator::data, uses OrthogonalSummation::operator() to sum and returns
   * the result.
   *
   * \param level up to which level to sum up
   * \return result of OrthogonalSummation for given type from MPQCDataMap_t.
   */
  MapValue operator()(const size_t level)
  {
    // evaluate
    const MapValue result = OS(level);
    return result;
  }

  /** Getter for the contribution of the SetValue<> to a specific \a index.
   *
   * @param index SetValue to this index
   * @return contribution
   */
  MapValue getContributionForIndexSet(const IndexSet::ptr &ptr) const
  {
    typedef SetValueMap<MapValue> setvalues_t;
    const setvalues_t &values = OS.getSetValues();
    return values.getConstValue(ptr)->getContribution();
  }

  /** Getter for the value of the SetValue<> to a specific \a index.
   *
   * @param index SetValue to this index
   * @return value
   */
  MapValue getValueForIndexSet(const IndexSet::ptr &ptr) const
  {
    typedef SetValueMap<MapValue> setvalues_t;
    const setvalues_t &values = OS.getSetValues();
    return values.getConstValue(ptr)->getValue();
  }

private:
  /** Tiny helper to create the indices from a given IndexSetContainer.
   *
   * Basically, we just have to make sure we miss the last one but only
   * if the \a container has more than one set.
   *
   * @param count how many indices of container to use
   * @param container container with IndexSet
   * @return all subsets contained in \a container
   */
  typename OrthogonalSummation<MapValue>::InputSets_t getSubsets(
      const size_t count,
      const IndexSetContainer::Container_t &container)
  {
    IndexSetContainer::Container_t::const_iterator iter = container.begin();
    std::advance(iter, count); // step till after desired element
    typename OrthogonalSummation<MapValue>::InputSets_t indices(container.begin(), iter);
    return indices;
  }

  /** Tiny helper to create the values for the summation in the correct order.
   *
   * @param data
   * @param container
   * @param MatrixNrLookup
   * @return
   */
  typename OrthogonalSummation<MapValue>::InputValues_t createValues(
      const std::map<JobId_t, MapType> &data,
      const IndexSetContainer::Container_t &container,
      const std::map< JobId_t, size_t > &MatrixNrLookup)
  {
    // if the power set of , we don't need to get rid of the "union index set"
    typename OrthogonalSummation<MapValue>::InputValues_t values(data.size(), ZeroInstance<MapValue>());
    for (typename std::map<JobId_t, MapType>::const_iterator dataiter = data.begin();
        dataiter != data.end(); ++dataiter) {
      const MapType &Data = dataiter->second;
      const JobId_t &jobid = dataiter->first;
      const MapValue &value = boost::fusion::at_key<MapKey>(Data);
      const std::map< JobId_t, size_t >::const_iterator nriter = MatrixNrLookup.find(jobid);
      ASSERT( nriter != MatrixNrLookup.end(),
          "OrthogonalSummation<>::createValues() - MatrixNrLookup does not contain id "
          +toString(jobid)+".");
      if ((nriter->second >= 0 ) && (nriter->second < data.size()))
        values[ nriter->second ] = value;
    }
    return values;
  }
private:
  //!> created indices for OS such that we may hand over refs
  typename OrthogonalSummation<MapValue>::InputSets_t indices;
  //!> created values for OS such that we may hand over refs
  typename OrthogonalSummation<MapValue>::InputValues_t values;
  //!> Summation instance to use in operator()(level)
  OrthogonalSummation<MapValue> OS;
};


#endif /* ORTHOGONALSUMMATOR_HPP_ */