Index: src/Fragmentation/Summation/SetValues/unittests/SamplingGridUnitTest.cpp
===================================================================
--- src/Fragmentation/Summation/SetValues/unittests/SamplingGridUnitTest.cpp	(revision fbf143b9fb3d409937862a091c997ab0428c8b3e)
+++ src/Fragmentation/Summation/SetValues/unittests/SamplingGridUnitTest.cpp	(revision fbf143b9fb3d409937862a091c997ab0428c8b3e)
@@ -0,0 +1,422 @@
+/*
+ * Project: MoleCuilder
+ * Description: creates and alters molecular systems
+ * Copyright (C)  2012 University of Bonn. All rights reserved.
+ * 
+ *
+ *   This file is part of MoleCuilder.
+ *
+ *    MoleCuilder is free software: you can redistribute it and/or modify
+ *    it under the terms of the GNU General Public License as published by
+ *    the Free Software Foundation, either version 2 of the License, or
+ *    (at your option) any later version.
+ *
+ *    MoleCuilder is distributed in the hope that it will be useful,
+ *    but WITHOUT ANY WARRANTY; without even the implied warranty of
+ *    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ *    GNU General Public License for more details.
+ *
+ *    You should have received a copy of the GNU General Public License
+ *    along with MoleCuilder.  If not, see <http://www.gnu.org/licenses/>.
+ */
+
+/*
+ * SamplingGridUnitTest.cpp
+ *
+ *  Created on: Jul 29, 2012
+ *      Author: heber
+ */
+
+// include config.h
+#ifdef HAVE_CONFIG_H
+#include <config.h>
+#endif
+
+using namespace std;
+
+#include <cppunit/CompilerOutputter.h>
+#include <cppunit/extensions/TestFactoryRegistry.h>
+#include <cppunit/ui/text/TestRunner.h>
+
+#include "SamplingGridUnitTest.hpp"
+
+#include "CodePatterns/Assert.hpp"
+
+#include <boost/assign.hpp>
+#include <cmath>
+#include <numeric>
+
+#ifdef HAVE_TESTRUNNER
+#include "UnitTestMain.hpp"
+#endif /*HAVE_TESTRUNNER*/
+
+using namespace boost::assign;
+
+/********************************************** Test classes **************************************/
+
+const double grid_value=1.;
+
+#define NUMBEROFSAMPLES(n) (size_t)(pow(pow(2,n),3))
+#define DOMAINVOLUME(l) (size_t)pow(l,3)
+
+// Registers the fixture into the 'registry'
+CPPUNIT_TEST_SUITE_REGISTRATION( SamplingGridTest );
+
+
+void SamplingGridTest::setUp()
+{
+  // failing asserts should be thrown
+  ASSERT_DO(Assert::Throw);
+
+  // create the grid
+  const double begin[3] = { 0., 0., 0. };
+  const double end[3] = { 1., 1., 1. };
+  for (size_t i=0; i< DOMAINVOLUME(1)*NUMBEROFSAMPLES(2); ++i)
+    values += grid_value;
+  grid = new SamplingGrid(begin, end, 2, values);
+  CPPUNIT_ASSERT_EQUAL( grid_value, *(grid->sampled_grid.begin()) );
+}
+
+
+void SamplingGridTest::tearDown()
+{
+  delete grid;
+}
+
+/** UnitTest on compatible combination of props and values
+ */
+void SamplingGridTest::compatibleGrids_Test()
+{
+  // check illegal grid
+  const double begin[3] = { 0., 0., 0. };
+  const double end[3] = { 2., 2., 2. };
+  SamplingGridProperties illegal_props(begin, end, 5);
+  SamplingGridProperties legal_props(begin, end, 2);
+  CPPUNIT_ASSERT( !grid->isCompatible(illegal_props) );
+  CPPUNIT_ASSERT( grid->isCompatible(legal_props) );
+  SamplingGrid::sampledvalues_t illegal_values;
+  for (size_t i=0; i< NUMBEROFSAMPLES(1); ++i)
+    illegal_values += 1.5;
+  SamplingGrid::sampledvalues_t legal_values;
+  for (size_t i=0; i< NUMBEROFSAMPLES(2); ++i)
+    legal_values += 1.5;
+#ifndef NDEBUG
+  // throws because props and size of illegal_values don't match
+  std::cout << "The following assertion is intended and does not indicate a failure of the test." << std::endl;
+  CPPUNIT_ASSERT_THROW( SamplingGrid illegal_grid(illegal_props, illegal_values), Assert::AssertionFailure );
+  std::cout << "The following assertion is intended and does not indicate a failure of the test." << std::endl;
+  CPPUNIT_ASSERT_THROW( SamplingGrid illegal_grid(legal_props, illegal_values), Assert::AssertionFailure );
+  std::cout << "The following assertion is intended and does not indicate a failure of the test." << std::endl;
+  CPPUNIT_ASSERT_THROW( SamplingGrid illegal_grid(illegal_props, legal_values), Assert::AssertionFailure );
+#endif
+  // check that grid is still the same
+  for (SamplingGrid::sampledvalues_t::const_iterator iter = grid->sampled_grid.begin();
+      iter != grid->sampled_grid.end(); ++iter)
+    CPPUNIT_ASSERT_EQUAL( grid_value, *iter );
+}
+
+/** UnitTest for isCongruent()
+ */
+void SamplingGridTest::isCongruent_Test()
+{
+  const double begin[3] = { 0., 0., 0. };
+  const double end[3] = { 2., 2., 2. };
+  const double otherbegin[3] = { 0.1, 0.1, 0.1 };
+  const double otherend[3] = { 1., 1., 1. };
+  SamplingGridProperties illegal_begin_props(otherbegin, end, 5);
+  SamplingGridProperties illegal_end_props(begin, otherend, 5);
+  SamplingGridProperties illegal_level_props(begin, end, 5);
+  SamplingGridProperties legal_props(begin, end, 2);
+
+  // differing windows
+  const double begin_window[3] = { 0.5, 0.5, 0.5 };
+  const double end_window[3] = { 1., 1., 1. };
+  const double otherbegin_window[3] = { 0.45, 0.45, 0.45 };
+  const double otherend_window[3] = { 1.05, 1.05, 1.05 };
+
+  // check that incompatible grid are also incongruent
+  SamplingGrid default_grid(legal_props);
+  // note that we always construct a temporary SamplingGrid from given props
+  CPPUNIT_ASSERT( default_grid.isCompatible(illegal_begin_props) == default_grid.isCongruent(illegal_begin_props));
+  CPPUNIT_ASSERT( default_grid.isCompatible(illegal_end_props) == default_grid.isCongruent(illegal_end_props));
+  CPPUNIT_ASSERT( default_grid.isCompatible(illegal_level_props) == default_grid.isCongruent(illegal_level_props));
+  CPPUNIT_ASSERT( default_grid.isCompatible(legal_props) == default_grid.isCongruent(legal_props) );
+
+  default_grid.setWindowSize(begin, end);
+  // same window
+  {
+    SamplingGrid illegal_begin_grid(illegal_begin_props);
+    SamplingGrid illegal_end_grid(illegal_end_props);
+    SamplingGrid illegal_level_grid(illegal_level_props);
+    SamplingGrid legal_grid(legal_props);
+    illegal_begin_grid.setWindowSize(begin, end);
+    illegal_end_grid.setWindowSize(begin, end);
+    illegal_level_grid.setWindowSize(begin, end);
+    legal_grid.setWindowSize(begin, end);
+    CPPUNIT_ASSERT( !illegal_begin_grid.isCongruent(default_grid) );
+    CPPUNIT_ASSERT( !illegal_end_grid.isCongruent(default_grid) );
+    CPPUNIT_ASSERT( !illegal_level_grid.isCongruent(default_grid) );
+    CPPUNIT_ASSERT( legal_grid.isCongruent(default_grid) );
+  }
+
+  // different begin
+  {
+    SamplingGrid illegal_begin_grid(illegal_begin_props);
+    SamplingGrid illegal_end_grid(illegal_end_props);
+    SamplingGrid illegal_level_grid(illegal_level_props);
+    SamplingGrid legal_grid(legal_props);
+    illegal_begin_grid.setWindowSize(otherbegin_window, end);
+    illegal_end_grid.setWindowSize(otherbegin_window, end);
+    illegal_level_grid.setWindowSize(otherbegin_window, end);
+    legal_grid.setWindowSize(begin, end);
+    CPPUNIT_ASSERT( !illegal_begin_grid.isCongruent(legal_grid) );
+    CPPUNIT_ASSERT( !illegal_end_grid.isCongruent(legal_grid) );
+    CPPUNIT_ASSERT( !illegal_level_grid.isCongruent(legal_grid) );
+    CPPUNIT_ASSERT( legal_grid.isCongruent(default_grid) );
+  }
+
+  // different end
+  {
+    SamplingGrid illegal_begin_grid(illegal_begin_props);
+    SamplingGrid illegal_end_grid(illegal_end_props);
+    SamplingGrid illegal_level_grid(illegal_level_props);
+    SamplingGrid legal_grid(legal_props);
+    illegal_begin_grid.setWindowSize(begin, otherend_window);
+    illegal_end_grid.setWindowSize(begin, otherend_window);
+    illegal_level_grid.setWindowSize(begin, otherend_window);
+    legal_grid.setWindowSize(begin, end);
+    CPPUNIT_ASSERT( !illegal_begin_grid.isCongruent(legal_grid) );
+    CPPUNIT_ASSERT( !illegal_end_grid.isCongruent(legal_grid) );
+    CPPUNIT_ASSERT( !illegal_level_grid.isCongruent(legal_grid) );
+    CPPUNIT_ASSERT( legal_grid.isCongruent(default_grid) );
+  }
+
+  // different begin and end
+  {
+    SamplingGrid illegal_begin_grid(illegal_begin_props);
+    SamplingGrid illegal_end_grid(illegal_end_props);
+    SamplingGrid illegal_level_grid(illegal_level_props);
+    SamplingGrid legal_grid(legal_props);
+    illegal_begin_grid.setWindowSize(otherbegin_window, otherend_window);
+    illegal_end_grid.setWindowSize(otherbegin_window, otherend_window);
+    illegal_level_grid.setWindowSize(otherbegin_window, otherend_window);
+    legal_grid.setWindowSize(begin, end);
+    CPPUNIT_ASSERT( !illegal_begin_grid.isCongruent(legal_grid) );
+    CPPUNIT_ASSERT( !illegal_end_grid.isCongruent(legal_grid) );
+    CPPUNIT_ASSERT( !illegal_level_grid.isCongruent(legal_grid) );
+    CPPUNIT_ASSERT( legal_grid.isCongruent(default_grid) );
+  }
+}
+
+/** UnitTest for integral()
+ */
+void SamplingGridTest::integral_Test()
+{
+  double sum = 0.;
+  sum = std::accumulate( grid->sampled_grid.begin(), grid->sampled_grid.end(), sum );
+  CPPUNIT_ASSERT_EQUAL( sum*grid->getVolume()/grid->getWindowGridPoints(), grid->integral() );
+}
+
+/** UnitTest for getVolume()
+ */
+void SamplingGridTest::getVolume_Test()
+{
+  CPPUNIT_ASSERT_EQUAL( 1., grid->getVolume() );
+}
+
+/** UnitTest for getWindowSize()
+ */
+void SamplingGridTest::getWindowSize_Test()
+{
+  // check size of default grid
+  CPPUNIT_ASSERT_EQUAL( (size_t)(NUMBEROFSAMPLES(grid->level)), grid->getWindowGridPoints() );
+
+  // create another one and check its size, too
+  const double begin[3] = { 0., 0., 0. };
+  const double end[3] = { 1., 1., 1. };
+  for (size_t level = 3; level<=6; ++level) {
+    values.clear();
+    for (size_t i=0; i< NUMBEROFSAMPLES(level); ++i)
+      values += grid_value;
+    delete grid;
+    // use other pointer in case something fails
+    SamplingGrid *tmpgrid = new SamplingGrid(begin, end, level, values);
+    grid = tmpgrid;
+    CPPUNIT_ASSERT_EQUAL( (size_t)NUMBEROFSAMPLES(level), grid->getWindowGridPoints() );
+  }
+}
+
+/** UnitTest for extendWindow()
+ */
+void SamplingGridTest::extendWindow_Test()
+{
+  // we have a grid with size of one, extend to twice the size and check
+  const double begin[3] = { 0., 0., 0. };
+  const double size = 2.;
+  const double end[3] = { size, size, size };
+  double offset[3];
+  for (offset[0] = 0.; offset[0] <= 1.; offset[0] += .5)
+    for (offset[1] = 0.; offset[1] <= 1.; offset[1] += .5)
+      for (offset[2] = 0.; offset[2] <= 1.; offset[2] += .5) {
+        const double window_begin[3] = { 0.+offset[0], 0.+offset[1], 0.+offset[2]};
+        const double window_end[3] = { 1.+offset[0], 1.+offset[1], 1.+offset[2]};
+        SamplingGrid newgrid(begin, end, 2);
+        newgrid.setWindowSize(window_begin, window_end);
+        // resize values by hand to new window size. Otherwise they get zero'd.
+        newgrid.sampled_grid = values;
+        newgrid.sampled_grid.resize(NUMBEROFSAMPLES(1));
+        newgrid.extendWindow(begin, end);
+
+        // check integral
+        CPPUNIT_ASSERT_EQUAL(
+            grid_value/(NUMBEROFSAMPLES(grid->level)/NUMBEROFSAMPLES(grid->level)),
+            grid->integral()
+            );
+
+        // check number of points
+        CPPUNIT_ASSERT_EQUAL(
+            (size_t)NUMBEROFSAMPLES(grid->level),
+            grid->getWindowGridPoints()
+            );
+      }
+}
+
+/** UnitTest for extendWindow() with asymmetric values
+ */
+void SamplingGridTest::extendWindow_asymmetric_Test()
+{
+  std::cout << "SamplingGridTest::extendWindow_asymmetric_Test()" << std::endl;
+  const double begin[3] = { 0., 0., 0. };
+  const double end[3] = { 2., 2., 2. };
+  double offset[3];
+  for (offset[0] = 0.; offset[0] <= 1.; offset[0] += .5)
+    for (offset[1] = 0.; offset[1] <= 1.; offset[1] += .5)
+      for (offset[2] = 0.; offset[2] <= 1.; offset[2] += .5) {
+        const double window_begin[3] = { 0.+offset[0], 0.+offset[1], 0.+offset[2]};
+        const double window_end[3] = { 1.+offset[0], 1.+offset[1], 1.+offset[2]};
+        SamplingGrid newgrid(begin, end, 2);
+        CPPUNIT_ASSERT_EQUAL( (size_t)0, newgrid.getWindowGridPoints() );
+        newgrid.setWindowSize(window_begin, window_end);
+        // window size is only half of domain size
+        const size_t max_samples = NUMBEROFSAMPLES(newgrid.level)*pow(0.5,3);
+        for (size_t i=0; i< max_samples; ++i)
+          newgrid.sampled_grid += grid_value*i;
+        const size_t sum_weight = (max_samples)*(max_samples-1)/2;
+        const double integral = newgrid.integral();
+        newgrid.extendWindow(begin, end);
+
+        // check that integral has remained the same
+        CPPUNIT_ASSERT_EQUAL( integral, newgrid.integral() );
+        CPPUNIT_ASSERT_EQUAL( grid_value*sum_weight/DOMAINVOLUME(2), newgrid.integral() );
+      }
+}
+
+/** UnitTest for addOntoWindow()
+ */
+void SamplingGridTest::addOntoWindow_Test()
+{
+  // first window is from (0,0,0) to (1,1,1)
+  CPPUNIT_ASSERT_EQUAL( 1.*grid_value, grid->integral() );
+
+  // create values for half-sized window
+  values.clear();
+  for (size_t i=0; i< (size_t)pow(.5*pow(2,2),3); ++i)
+    values += grid_value;
+
+  // check that too large a window throws
+#ifndef NDEBUG
+  const double begin[3] = { .5, .5, .5 };
+  const double wrongend[3] = { 1.5, 1.5, 1.5 };
+  std::cout << "The following assertion is intended and does not indicate a failure of the test." << std::endl;
+  CPPUNIT_ASSERT_THROW( grid->addOntoWindow(begin, wrongend, values, +1.), Assert::AssertionFailure );
+#endif
+
+  // create another window from (.5,.5,.5) to (1., 1., 1.)
+  double offset[3];
+  for (offset[0] = 0.; offset[0] <= .5; offset[0] += .5)
+    for (offset[1] = 0.; offset[1] <= .5; offset[1] += .5)
+      for (offset[2] = 0.; offset[2] <= .5; offset[2] += .5) {
+        const double window_begin[3] = { 0.+offset[0], 0.+offset[1], 0.+offset[2]};
+        const double window_end[3] = { .5+offset[0], .5+offset[1], .5+offset[2]};
+
+        SamplingGrid newgrid(*grid);
+        // now perform working operation
+        newgrid.addOntoWindow(window_begin, window_end, values, +1.);
+
+        // check integral to be one and one eighth times the old value
+        CPPUNIT_ASSERT_EQUAL( (1.+pow(.5,3))*grid_value, newgrid.integral() );
+      }
+}
+
+/** UnitTest for addOntoWindow() with asymmetric values
+ */
+void SamplingGridTest::addOntoWindow_asymmetric_Test()
+{
+  const size_t size = grid->end[0]-grid->begin[0];
+  // check with asymmetric values
+  grid->sampled_grid.clear();
+  grid->sampled_grid.resize(DOMAINVOLUME(size)*NUMBEROFSAMPLES(grid->level), 0.);
+
+  for (size_t i=0;i<grid->level*(grid->end[0]-grid->begin[0]);++i)
+    grid->sampled_grid[(i*2+0)*2+0] += .5*grid_value;
+  for (size_t i=0;i<grid->level*(grid->end[1]-grid->begin[1]);++i)
+    grid->sampled_grid[(0*2+i)*2+0] += 1.*grid_value;
+  for (size_t i=0;i<grid->level*(grid->end[2]-grid->begin[2]);++i)
+    grid->sampled_grid[(0*2+0)*2+i] += 1.5*grid_value;
+
+  const double integral = grid->integral();
+
+  // now perform working operation
+  grid->addOntoWindow(grid->begin, grid->end, values, +1.);
+  // values is equal to integral of 1.
+  CPPUNIT_ASSERT_EQUAL( 1.+integral, grid->integral() );
+}
+
+/** UnitTest for operator+=()
+ */
+void SamplingGridTest::operatorPlusEqual_Test()
+{
+  // create other grid
+  const double begin[3] = { 0., 0., 0. };
+  const double end[3] = { 1., 1., 1. };
+  SamplingGrid::sampledvalues_t othervalues;
+  const double othergrid_value = 1.5;
+  for (size_t i=0; i< NUMBEROFSAMPLES(2); ++i)
+    othervalues += othergrid_value;
+  SamplingGrid othergrid(begin, end, 2, othervalues);
+  CPPUNIT_ASSERT_EQUAL( othergrid_value, *(othergrid.sampled_grid.begin()) );
+
+  // perform operation
+  CPPUNIT_ASSERT_NO_THROW( *grid += othergrid );
+
+  // check the contents of the grid
+  const double sum = grid_value+othergrid_value;
+  for (SamplingGrid::sampledvalues_t::const_iterator iter = grid->sampled_grid.begin();
+      iter != grid->sampled_grid.end(); ++iter)
+    CPPUNIT_ASSERT_EQUAL( sum, *iter );
+}
+
+/** UnitTest for operator-=()
+ */
+void SamplingGridTest::operatorMinusEqual_Test()
+{
+  // create other grid
+  const double begin[3] = { 0., 0., 0. };
+  const double end[3] = { 1., 1., 1. };
+  SamplingGrid::sampledvalues_t othervalues;
+  const double othergrid_value = 1.5;
+  for (size_t i=0; i< NUMBEROFSAMPLES(2); ++i)
+    othervalues += othergrid_value;
+  SamplingGrid othergrid(begin, end, 2, othervalues);
+  CPPUNIT_ASSERT_EQUAL( othergrid_value, *(othergrid.sampled_grid.begin()) );
+
+  // perform operation
+  CPPUNIT_ASSERT_NO_THROW( *grid -= othergrid );
+
+  // check the contents of the grid
+  const double difference = grid_value-othergrid_value;
+  for (SamplingGrid::sampledvalues_t::const_iterator iter = grid->sampled_grid.begin();
+      iter != grid->sampled_grid.end(); ++iter)
+    CPPUNIT_ASSERT_EQUAL( difference, *iter );
+}
+
