Changeset c6355f

Timestamp:

Feb 15, 2013, 9:51:49 AM (12 years ago)

Author:

Frederik Heber <heber@…>

Branches:

Action_Thermostats, Add_AtomRandomPerturbation, Add_FitFragmentPartialChargesAction, Add_RotateAroundBondAction, Add_SelectAtomByNameAction, Added_ParseSaveFragmentResults, AddingActions_SaveParseParticleParameters, Adding_Graph_to_ChangeBondActions, Adding_MD_integration_tests, Adding_ParticleName_to_Atom, Adding_StructOpt_integration_tests, AtomFragments, Automaking_mpqc_open, AutomationFragmentation_failures, Candidate_v1.5.4, Candidate_v1.6.0, Candidate_v1.6.1, ChangeBugEmailaddress, ChangingTestPorts, ChemicalSpaceEvaluator, CombiningParticlePotentialParsing, Combining_Subpackages, Debian_Package_split, Debian_package_split_molecuildergui_only, Disabling_MemDebug, Docu_Python_wait, EmpiricalPotential_contain_HomologyGraph, EmpiricalPotential_contain_HomologyGraph_documentation, Enable_parallel_make_install, Enhance_userguide, Enhanced_StructuralOptimization, Enhanced_StructuralOptimization_continued, Example_ManyWaysToTranslateAtom, Exclude_Hydrogens_annealWithBondGraph, FitPartialCharges_GlobalError, Fix_BoundInBox_CenterInBox_MoleculeActions, Fix_ChargeSampling_PBC, Fix_ChronosMutex, Fix_FitPartialCharges, Fix_FitPotential_needs_atomicnumbers, Fix_ForceAnnealing, Fix_IndependentFragmentGrids, Fix_ParseParticles, Fix_ParseParticles_split_forward_backward_Actions, Fix_PopActions, Fix_QtFragmentList_sorted_selection, Fix_Restrictedkeyset_FragmentMolecule, Fix_StatusMsg, Fix_StepWorldTime_single_argument, Fix_Verbose_Codepatterns, Fix_fitting_potentials, Fixes, ForceAnnealing_goodresults, ForceAnnealing_oldresults, ForceAnnealing_tocheck, ForceAnnealing_with_BondGraph, ForceAnnealing_with_BondGraph_continued, ForceAnnealing_with_BondGraph_continued_betteresults, ForceAnnealing_with_BondGraph_contraction-expansion, FragmentAction_writes_AtomFragments, FragmentMolecule_checks_bonddegrees, GeometryObjects, Gui_Fixes, Gui_displays_atomic_force_velocity, ImplicitCharges, IndependentFragmentGrids, IndependentFragmentGrids_IndividualZeroInstances, IndependentFragmentGrids_IntegrationTest, IndependentFragmentGrids_Sole_NN_Calculation, JobMarket_RobustOnKillsSegFaults, JobMarket_StableWorkerPool, JobMarket_unresolvable_hostname_fix, MoreRobust_FragmentAutomation, ODR_violation_mpqc_open, PartialCharges_OrthogonalSummation, PdbParser_setsAtomName, PythonUI_with_named_parameters, QtGui_reactivate_TimeChanged_changes, Recreated_GuiChecks, Rewrite_FitPartialCharges, RotateToPrincipalAxisSystem_UndoRedo, SaturateAtoms_findBestMatching, SaturateAtoms_singleDegree, StoppableMakroAction, Subpackage_CodePatterns, Subpackage_JobMarket, Subpackage_LinearAlgebra, Subpackage_levmar, Subpackage_mpqc_open, Subpackage_vmg, Switchable_LogView, ThirdParty_MPQC_rebuilt_buildsystem, TrajectoryDependenant_MaxOrder, TremoloParser_IncreasedPrecision, TremoloParser_MultipleTimesteps, TremoloParser_setsAtomName, Ubuntu_1604_changes, stable

Children:

64bafe0

Parents:

de6dfb

git-author:

Frederik Heber <heber@…> (12/15/12 08:33:06)

git-committer:

Frederik Heber <heber@…> (02/15/13 09:51:49)

Message:

FIX: New non-zero window operations of SamplingGrid are now fully working.

Location:

src/Jobs/Grid

Files:

• SamplingGrid.cpp (modified) (10 diffs)
• SamplingGrid.hpp (modified) (4 diffs)
• unittests/SamplingGridUnitTest.cpp (modified) (12 diffs)
• unittests/SamplingGridUnitTest.hpp (modified) (2 diffs)

src/Jobs/Grid/SamplingGrid.cpp ¶

@@ -60 +60 @@
 SamplingGrid::SamplingGrid(const double _begin[3],
     const double _end[3],
+    const int _level) :
+  SamplingGridProperties(_begin, _end, _level)
+{
+  setWindowSize(zeroOffset, zeroOffset);
+  ASSERT( getWindowGridPoints() == (size_t)0,
+      "SamplingGrid::SamplingGrid() - incorrect number of samples given for the window.");
+}
+
+SamplingGrid::SamplingGrid(const double _begin[3],
+    const double _end[3],
     const int _level,
     const sampledvalues_t &_sampled_grid) :
@@ -74 +84 @@
   sampled_grid(_grid.sampled_grid)
 {
-  setWindowSize(_grid.begin, _grid.end);
+  setWindowSize(_grid.begin_window, _grid.end_window);
+  ASSERT( getWindowGridPoints() == _grid.getWindowGridPoints(),
+      "SamplingGrid::SamplingGrid() - incorrect number of samples given for the window.");
 }
 
@@ -80 +92 @@
   SamplingGridProperties(_props)
 {
-  setWindowSize(_props.begin, _props.end);
+  setWindowSize(zeroOffset, zeroOffset);
+  ASSERT( getWindowGridPoints() == (size_t)0,
+      "SamplingGrid::SamplingGrid() - incorrect number of samples given for the window.");
 }
 
@@ -241 +255 @@
   setWindowSize(_begin_window, _end_window);
   const size_t gridpoints_window = getWindowGridPoints();
+  sampled_grid.clear();
   sampled_grid.resize(gridpoints_window, 0.);
 }
@@ -268 +283 @@
     const double _end_window[3])
 {
-  // check that we truly have to extend the window
 #ifndef NDEBUG
   for(size_t index=0;index < 3; ++index) {
+    // check that we truly have to extend the window
     ASSERT ( begin_window[index] >= _begin_window[index],
         "SamplingGrid::extendWindow() - component "+toString(index)+
@@ -277 +292 @@
         "SamplingGrid::extendWindow() - component "+toString(index)+
         " of window end is less than old value.");
+
+    // check that we are still less than domain
+    ASSERT ( _begin_window[index] >= begin[index],
+        "SamplingGrid::extendWindow() - component "+toString(index)+
+        " of window start is less than domain start.");
+    ASSERT ( _end_window[index] <= end[index],
+        "SamplingGrid::extendWindow() - component "+toString(index)+
+        " of window end is greater than domain end.");
   }
 #endif
@@ -311 +334 @@
 #endif
   // the only issue are indices
-  const size_t gridpoints_axis = pow(2, level);
+  const size_t gridpoints_axis = getGridPointsPerAxis();
   size_t pre_offset[3];
   size_t post_offset[3];
@@ -318 +341 @@
   for(size_t index=0;index<3;++index) {
     pre_offset[index] = (_begin_window[index] - begin_window[index])*gridpoints_axis;
+    length[index] = (_end_window[index] - _begin_window[index])*gridpoints_axis;
     post_offset[index] = (end_window[index] - _end_window[index])*gridpoints_axis;
-    length[index] = (_end_window[index] - _begin_window[index])*gridpoints_axis;
     total[index] = (end_window[index] - begin_window[index])*gridpoints_axis;
     ASSERT( pre_offset[index]+post_offset[index]+length[index] == total[index],
-        "SamplingGrid::addOntoWindow() - pre, past, and length don't sum up to total for "
+        "SamplingGrid::addOntoWindow() - pre, post, and length don't sum up to total for "
         +toString(index)+"th component.");
   }
-  ASSERT( length[0]*length[1]*length[2] ==  _sampled_grid.size(),
-      "SamplingGrid::addOntoWindow() - not enough sampled values given.");
-  ASSERT( length[0]*length[1]*length[2] <=  sampled_grid.size(),
-      "SamplingGrid::addOntoWindow() - not enough sampled values available.");
+#ifndef NDEBUG
+  const size_t calculated_size = length[0]*length[1]*length[2];
+  ASSERT( calculated_size == _sampled_grid.size(),
+      "SamplingGrid::addOntoWindow() - not enough sampled values given: "
+      +toString(calculated_size)+" != "+toString(_sampled_grid.size())+".");
+  ASSERT( calculated_size <=  sampled_grid.size(),
+      "SamplingGrid::addOntoWindow() - not enough sampled values available: "
+      +toString(calculated_size)+" <= "+toString(sampled_grid.size())+".");
+  const size_t total_size = total[0]*total[1]*total[2];
+  ASSERT( total_size == sampled_grid.size(),
+      "SamplingGrid::addOntoWindow() - total size is not equal to number of present points: "
+      +toString(total_size)+" != "+toString(sampled_grid.size())+".");
+#endif
   size_t N[3];
+  size_t counter = 0;
   sampledvalues_t::iterator griditer = sampled_grid.begin();
   std::advance(griditer, pre_offset[0]*total[1]*total[2]);
@@ -338 +371 @@
       std::advance(griditer, pre_offset[2]);
       for(N[2]=0; N[2] < length[2]; ++N[2]) {
+        ASSERT( griditer != sampled_grid.end(),
+            "SamplingGrid::addOntoWindow() - griditer is already at end of window.");
+        ASSERT( copyiter != _sampled_grid.end(),
+            "SamplingGrid::addOntoWindow() - griditer is already at end of window.");
         *griditer++ += prefactor*(*copyiter++);
       }
@@ -344 +381 @@
     std::advance(griditer, post_offset[1]*total[2]);
   }
-  //std::advance(griditer, post_offset[0]*total[1]*total[2]);
+#ifndef NDEBUG
+  std::advance(griditer, post_offset[0]*total[1]*total[2]);
+  ASSERT( griditer == sampled_grid.end(),
+      "SamplingGrid::addOntoWindow() - griditer is not at end of window.");
+  ASSERT( copyiter == _sampled_grid.end(),
+      "SamplingGrid::addOntoWindow() - copyiter is not at end of window.");
+#endif
   LOG(2, "DEBUG: Grid after adding other is " << sampled_grid << ".");
 }

src/Jobs/Grid/SamplingGrid.hpp ¶

@@ -40 +40 @@
   typedef std::vector< double > sampledvalues_t;
 
-  /** Constructor for class SamplingGrid.
+  /** Constructor for class SamplingGrid for full window.
+   *
+   * Here, the window of sampled values spans the given domain.
    *
    * \param _begin offset for grid per axis
@@ -52 +54 @@
       const sampledvalues_t &_sampled_grid);
 
-  /** Copy constructor for class SamplingGrid.
-   *
-   * \param _grid grid to copy
-   */
-  SamplingGrid(const SamplingGrid &_grid);
-
-  /** Copy constructor for class SamplingGrid from SamplingGridProperties.
+  /** Constructor for class SamplingGrid for empty window.
+   *
+   * Here, the window is initially of size zero.
+   *
+   * \param _begin offset for grid per axis
+   * \param _end edge length of grid per axis
+   * \param _level number of grid points in \f$2^{\text{level}}\f$
+   */
+  SamplingGrid(const double _begin[3],
+      const double _end[3],
+      const int _level);
+
+  /** Copy constructor for class SamplingGrid with full window from SamplingGridProperties.
+   *
+   * Here, the window is initially empty.
    *
    * \param _props properties to copy
@@ -64 +74 @@
   SamplingGrid(const SamplingGridProperties &_props);
 
-  /** Copy constructor for class SamplingGrid from SamplingGridProperties.
+  /** Copy constructor for class SamplingGrid with empty window from SamplingGridProperties.
+   *
+   * Here, the window must span the whole domain
    *
    * \param _props properties to copy
@@ -72 +84 @@
       const SamplingGridProperties &_props,
       const sampledvalues_t &_sampled_grid);
+
+  /** Copy constructor for class SamplingGrid.
+   *
+   * The window of sampled values corresponds to the one on \a _grid.
+   *
+   * \param _grid grid to copy
+   */
+  SamplingGrid(const SamplingGrid &_grid);
 
   /** default cstor.

src/Jobs/Grid/unittests/SamplingGridUnitTest.cpp ¶

@@ -43 +43 @@
 #include "CodePatterns/Assert.hpp"
 
+#include <boost/assign.hpp>
 #include <cmath>
-#include <boost/assign.hpp>
+#include <numeric>
 
 #ifdef HAVE_TESTRUNNER
@@ -57 +58 @@
 
 #define NUMBEROFSAMPLES(n) (size_t)(pow(pow(2,n),3))
+#define DOMAINVOLUME(l) (size_t)pow(l,3)
 
 // Registers the fixture into the 'registry'
@@ -70 +72 @@
   const double begin[3] = { 0., 0., 0. };
   const double end[3] = { 1., 1., 1. };
-  for (size_t i=0; i< pow(1,3)*NUMBEROFSAMPLES(2); ++i)
+  for (size_t i=0; i< DOMAINVOLUME(1)*NUMBEROFSAMPLES(2); ++i)
     values += grid_value;
   grid = new SamplingGrid(begin, end, 2, values);
@@ -89 +91 @@
   const double begin[3] = { 0., 0., 0. };
   const double end[3] = { 2., 2., 2. };
-  SamplingGridProperties illegal_props(begin, end, 1);
+  SamplingGridProperties illegal_props(begin, end, 5);
   SamplingGridProperties legal_props(begin, end, 2);
   CPPUNIT_ASSERT( !grid->isCompatible(illegal_props) );
@@ -98 +100 @@
   SamplingGrid::sampledvalues_t legal_values;
   for (size_t i=0; i< NUMBEROFSAMPLES(2); ++i)
-    illegal_values += 1.5;
+    legal_values += 1.5;
 #ifndef NDEBUG
   // throws because props and size of illegal_values don't match
@@ -114 +116 @@
 }
 
-/** UnitTest for getVolume_Test
+/** 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()
@@ -121 +132 @@
 }
 
-/** UnitTest for getWindowSize_Test
+/** UnitTest for getWindowSize()
  */
 void SamplingGridTest::getWindowSize_Test()
@@ -147 +158 @@
 void SamplingGridTest::extendWindow_Test()
 {
-  // we have a grid with size of two, extend to twice the size and check
+  // 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 };
-  grid->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() );
+  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()
+            );
+      }
 }
 
@@ -164 +193 @@
 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. };
-  grid->sampled_grid.clear();
-  SamplingGrid::sampledvalues_t::iterator griditer = grid->sampled_grid.begin();
-  for (size_t i=0; i< NUMBEROFSAMPLES(grid->level); ++i)
-     *griditer++ += grid_value*i;
-  grid->extendWindow(begin, end);
-
-  CPPUNIT_ASSERT_EQUAL( grid_value, grid->integral() );
+  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() );
+      }
 }
 
@@ -182 +226 @@
   CPPUNIT_ASSERT_EQUAL( 1.*grid_value, grid->integral() );
 
-  // create another window from (.5,.5,.5) to (1., 1., 1.)
-  const double begin[3] = { .5, .5, .5 };
-  const double end[3] = { 1., 1., 1. };
+  // create values for half-sized window
   values.clear();
-  for (size_t i=0; i< NUMBEROFSAMPLES(2)/pow(2,3); ++i) // cut in half
+  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
-  // now perform working operation
-  grid->addOntoWindow(begin, end, values, +1.);
-
-  // check integral to be one and half times the old value
-  CPPUNIT_ASSERT_EQUAL( (1.+pow(.5,3))*grid_value, grid->integral() );
+
+  // 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() );
+      }
 }
 
@@ -205 +260 @@
 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(NUMBEROFSAMPLES(2), 0.);
+  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)
@@ -217 +273 @@
 
   const double integral = grid->integral();
-  const double begin[3] = { 0., 0., 0. };
-  for(double size=2.; size<=5.; size += 1.) {
-    const double newend[3] = {size, size, size};
-    grid->extendWindow(begin, newend);
-    CPPUNIT_ASSERT_EQUAL( pow(size,3)*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() );
 }
 

src/Jobs/Grid/unittests/SamplingGridUnitTest.hpp ¶

@@ -27 +27 @@
     CPPUNIT_TEST ( operatorPlusEqual_Test );
     CPPUNIT_TEST ( operatorMinusEqual_Test );
+    CPPUNIT_TEST ( integral_Test );
     CPPUNIT_TEST ( getVolume_Test );
     CPPUNIT_TEST ( getWindowSize_Test );
@@ -39 +40 @@
       void tearDown();
       void compatibleGrids_Test();
+      void integral_Test();
       void getVolume_Test();
       void getWindowSize_Test();