[newfield_revision Patch] Revise nearestNeighbors

Jeffrey Oldham oldham at codesourcery.com
Wed Aug 15 18:10:02 UTC 2001


The statigraphic flow code demonstrated the need to generalize nearest
neighbor computations.  These changes permit computing nearest
neighbors given a centering, i.e., all values in one cell; a
FieldOffset, i.e., one value; and a FieldOffsetList, i.e., a set of
values.

Applied to newfield_revision branch.
 
2001-08-15  Jeffrey D. Oldham  <oldham at codesourcery.com>
 
        * FieldOffset.h (FieldOffsetList::FieldOffsetList): Add
        constructor converting from a std::vector to a FieldOffsetList.
        * NearestNeighbors.h (NearestNeighborClass::NearestNeighborClass):
        Changed to an empty constructor.
        (NearestNeighborClass::operator()()): Removed.
        (NearestNeighborClass::operator()(Centering,Centering)): New function.
        (NearestNeighborClass::operator()(Centering,FieldOffset,Centering)):
        New function.
        (NearestNeighborClass::operator()(Centering,FieldOffsetList,Centering)):
        New function.
        (NearestNeighborClass::nearestNeighbors): New function formed from
        previous operator()().
        (manhattanDistance): Initialize variable.
        (NearestNeighborClass): Remove centering data members.
        (nearestNeighbors(Centering,Centering)): Revise to use
        operator()(Centering,Centering).
        (nearestNeighbors(Centering,Centering,bool)): Revise to use
        operator()(Centering,Centering).
        (nearestNeighbors(Centering,FieldOffsetList,Centering)): New
        function.
        (nearestNeighbors(Centering,FieldOffsetList,Centering,bool)): New
        function.
        (nearestNeighbors(Centering,FieldOffset,Centering)): New function.
        (nearestNeighbors(Centering,FieldOffset,Centering,bool)): New
        function.
        * tests/NearestNeighbors.cpp (main): Revise to test new functions.

Applied to      newfield_revision branch
Approved by     Stephen Smith
Tested on       sequential Linux using gcc 3.0.1 by compiling Pooma library and NewField tests

Thanks,
Jeffrey D. Oldham
oldham at codesourcery.com
-------------- next part --------------
? LINUXgcc
? include.mk
? FieldCentering.h.patch
? 09Aug.16.2.ChangeLog
? 13Aug.11.8.patch
? replicate.patch.save
? 09Aug.16.2.patch
? differences-2001Aug14
? 13Aug.11.8.ChangeLog
? 14Aug.14.9.patch
? 13Aug.14.4.patch
? 13Aug.14.4.ChangeLog
? replicate.patch.14Aug.14.8
? 14Aug.14.9.ChangeLog
? 14Aug.15.0.patch
? tests/LINUXgcc
? tests/BasicTest1.out
? tests/BasicTest2.out
? tests/FieldReductions2.cpp
? tests/NearestNeighbors3.cpp
? tests/Gradient.out
? tests/Centerings.out
? tests/ExpressionTest.out
? tests/FieldOffset.out
? tests/FieldReductions.out
? tests/FieldTour1.out
? tests/FieldTour2.out
? tests/FieldTour3.out
? tests/MeshTest1.out
? tests/NearestNeighbors.out
? tests/ScalarCode.out
? tests/StencilTests.out
? tests/VectorTest.out
? tests/Gradient.cpp
? tests/WhereTest.out
? tests/Gradient.cpp.13Aug.13.2.patch
? tests/Gradient.cpp.patch
? tests/Gradient.cpp.ChangeLog
Index: FieldOffset.h
===================================================================
RCS file: /home/pooma/Repository/r2/src/NewField/Attic/FieldOffset.h,v
retrieving revision 1.1.2.3
diff -c -p -r1.1.2.3 FieldOffset.h
*** FieldOffset.h	2001/08/10 17:40:52	1.1.2.3
--- FieldOffset.h	2001/08/14 22:06:18
*************** public:
*** 285,290 ****
--- 285,297 ----
      v_m.reserve(sz);
    }
  
+   // Construct from a vector.
+ 
+   FieldOffsetList(const std::vector<FieldOffset<Dim> > &v) {
+     v_m.resize(v.size());
+     std::copy(v.begin(), v.end(), v_m.begin());
+   }
+ 
    // Copy a vector's entries to this FieldOffsetList.
  
    FieldOffsetList &operator=(const std::vector<FieldOffset<Dim> > &v)
Index: NearestNeighbors.h
===================================================================
RCS file: /home/pooma/Repository/r2/src/NewField/Attic/NearestNeighbors.h,v
retrieving revision 1.1.2.3
diff -c -p -r1.1.2.3 NearestNeighbors.h
*** NearestNeighbors.h	2001/08/14 20:24:17	1.1.2.3
--- NearestNeighbors.h	2001/08/14 22:06:18
***************
*** 40,46 ****
  // Overview: 
  // 
  // nearestNeighbors
! //   - yields FieldOffsetLists corresponding to the specified centerings
  // NearestNeighborClass
  //   - the class performing the work
  //-----------------------------------------------------------------------------
--- 40,46 ----
  // Overview: 
  // 
  // nearestNeighbors
! //   - yields FieldOffsets corresponding to the given parameters
  // NearestNeighborClass
  //   - the class performing the work
  //-----------------------------------------------------------------------------
*************** public:
*** 107,218 ****
    typedef std::vector<MinimumPair> MinimumSet;
  
  
!   NearestNeighborClass(const Center &inputCentering,
! 		       const Center &outputCentering)
!     : inputCentering_m(inputCentering),
!       outputCentering_m(outputCentering)
    {
!     PInsist(inputCentering_m.positions().size() > 0,
  	    "The input centering must be non-empty.");
-   }
  
-   inline
-   Answer_t operator()()
-   {
      Answer_t answer;
!     answer.resize(outputCentering_m.size());
!     Positions inputPositions = inputCentering_m.positions();
!     Positions outputPositions = outputCentering_m.positions();
!     Position positionDifference;
!     double minimumDistance;
  
      // Determine nearest neighbors for each output value.
  
      for (Answer_t::size_type outputIndex = 0;
! 	 outputIndex < outputCentering_m.size();
! 	 ++outputIndex)
!       {
! 	// Compute all input values in the first shell.
! 
! 	MinimumSet minimumSet;	// all input values in first shell
! 	const Position outputValue = outputPositions[outputIndex];
! 	typename Positions::size_type inputIndex = 0;
! 
! 	// Use the first input value to start computing the minimum.
! 
! 	positionDifference = inputPositions[inputIndex] - outputValue;
! 	minimumDistance = 
! 	  (IntraCellOnly ?
! 	   manhattanDistance<Manhattan>(positionDifference) :
! 	   manhattanDistance<ManhattanGrid>(positionDifference));
!  	minimumSet.push_back(std::make_pair(inputIndex, positionDifference));
! 
! 	// Compute the minimum over the rest of the input values.
! 
! 	for (++inputIndex;
! 	     inputIndex < inputPositions.size();
! 	     ++inputIndex) {
! 	  positionDifference = inputPositions[inputIndex] - outputValue;
! 	  const double distance = 
! 	    (IntraCellOnly ?
! 	     manhattanDistance<Manhattan>(positionDifference) :
! 	     manhattanDistance<ManhattanGrid>(positionDifference));
! 	  if (distance < minimumDistance + epsilon) {
! 	    if (distance < minimumDistance) {
! 	      minimumSet.clear();
! 	      minimumDistance = distance;
! 	    }
! 	    minimumSet.push_back(std::make_pair(inputIndex,
! 						positionDifference));
! 	  }
! 	}
  
  
! 	// Convert the minimum set to a set of FieldOffsets.
! 	// minimumSet has all the minimum distance locations.
  
! 	FieldOffset_vt answerHolder;
! 	if (IntraCellOnly) {
! 	  for (MinimumSet::size_type minIndex = 0;
! 	       minIndex < minimumSet.size();
! 	       ++minIndex)
! 	    answerHolder.push_back(FieldOffset_t(Loc<Dim>(0),
! 						 minimumSet[minIndex].first));
! 	}
! 	else {
! 	  FieldOffset_vt partialAnswer;
! 	  for (MinimumSet::size_type minIndex = 0;
! 	       minIndex < minimumSet.size();
! 	       ++minIndex)
! 	    {
! 	      // Compute the cell offsets, appending to the set of answers.
! 
! 	      partialAnswer = computeCellOffsets(minimumSet[minIndex].first,
! 						 minimumSet[minIndex].second);
! 	      answerHolder.insert(answerHolder.end(),
! 				  partialAnswer.begin(), partialAnswer.end());
! 	    }
! 
! 	  // Remove all duplicates from the answer set.
! 
! 	  std::sort(answerHolder.begin(), answerHolder.end(),
! 		    CompareFieldOffset());
! 	  answerHolder.erase(std::unique(answerHolder.begin(),
! 					 answerHolder.end(),
! 					 EqualFieldOffset()),
! 			     answerHolder.end());
! 	}
  
! 	// Store the answer.
  
! 	answer[outputIndex] = answerHolder;
!       }
!     
      return answer;
    }
  
  private:
  
    // Given a difference between two positions in logical coordinate
    // space, return the Manhattan norm distance taking into account
    // that input values are repeated in every grid cell.
--- 107,277 ----
    typedef std::vector<MinimumPair> MinimumSet;
  
  
!   // The constructor performs no work.  The function operators do all
!   // the work.
! 
!   NearestNeighborClass() {}
! 
! 
!   // Return nearest neighbors for all value in an output centering.
! 
!   inline Answer_t
!   operator()(const Center &inputCentering, const Center &outputCentering)
    {
!     PInsist(inputCentering.size() > 0,
  	    "The input centering must be non-empty.");
  
      Answer_t answer;
!     answer.resize(outputCentering.size());
!     const Positions inputPositions = inputCentering.positions();
!     const Positions outputPositions = outputCentering.positions();
  
      // Determine nearest neighbors for each output value.
  
      for (Answer_t::size_type outputIndex = 0;
! 	 outputIndex < outputCentering.size();
! 	 ++outputIndex) {
!       answer[outputIndex] = nearestNeighbors(inputPositions,
! 					     outputPositions[outputIndex]);
!       std::cout << answer[outputIndex] << std::endl; // TMP
!     }
!     return answer;
!   }
  
  
!   // Return the nearest neighbors for one output position, specified
!   // by a FieldOffset.
  
!   inline FieldOffsetList_t
!   operator()(const Center &inputCentering,
! 	     const FieldOffset_t &fieldOffset,
! 	     const Center &outputCentering)
!   {
!     PInsist(inputCentering.size() > 0,
! 	    "The input centering must be non-empty.");
!     PInsist(fieldOffset.subFieldNumber() < outputCentering.size(),
! 	    "The field offset must correspond to the output centering.");
  
!     return nearestNeighbors(inputCentering.positions(),
! 			    outputCentering.position(fieldOffset.subFieldNumber()));
!   }
  
! 
!   // Return the nearest neighbors for multiple output positions, specified
!   // by a FieldOffsetList.
! 
!   inline FieldOffsetList_t
!   operator()(const Center &inputCentering,
! 	     const FieldOffsetList_t &fieldOffsetList,
! 	     const Center &outputCentering)
!   {
!     PInsist(inputCentering.size() > 0,
! 	    "The input centering must be non-empty.");
! 
!     Answer_t answer;
!     answer.resize(fieldOffsetList.size());
!     const Positions inputPositions = inputCentering.positions();
! 
!     // Determine nearest neighbors for each field offset.
! 
!     for (FieldOffsetList_t::size_type folIndex = 0;
! 	 folIndex < outputCentering.size();
! 	 ++folIndex) {
!       PInsist(fieldOffsetList[folIndex].subFieldNumber() < outputCentering.size(),
! 	    "The field offset must correspond to the output centering.");
!       answer[folIndex] =
! 	nearestNeighbors(inputPositions,
! 			 outputCentering.position(fieldOffsetList[folIndex].subFieldNumber()));
!     }
! 
      return answer;
    }
  
+ 
  private:
  
+   // Given an input centering and a position in logical
+   // coordinate space, return a FieldOffsetList of the nearest
+   // neighbors according to Manhattan distance.
+ 
+   inline FieldOffsetList_t
+   nearestNeighbors(const Positions &inputPositions,
+ 		   const Position outputValue)
+   {
+ 
+     // Compute all input values in the first shell.
+ 
+     MinimumSet minimumSet;	// all input values in first shell
+     typename Positions::size_type inputIndex = 0;
+ 
+     // Use the first input value to start computing the minimum.
+ 
+     Position positionDifference = inputPositions[inputIndex] - outputValue;
+     double minimumDistance = 
+       (IntraCellOnly ?
+        manhattanDistance<Manhattan>(positionDifference) :
+        manhattanDistance<ManhattanGrid>(positionDifference));
+     minimumSet.push_back(std::make_pair(inputIndex, positionDifference));
+ 
+     // Compute the minimum over the rest of the input values.
+ 
+     for (++inputIndex;
+ 	 inputIndex < inputPositions.size();
+ 	 ++inputIndex) {
+       positionDifference = inputPositions[inputIndex] - outputValue;
+       const double distance = 
+ 	(IntraCellOnly ?
+ 	 manhattanDistance<Manhattan>(positionDifference) :
+ 	 manhattanDistance<ManhattanGrid>(positionDifference));
+       if (distance < minimumDistance + epsilon) {
+ 	if (distance < minimumDistance) {
+ 	  minimumSet.clear();
+ 	  minimumDistance = distance;
+ 	}
+ 	minimumSet.push_back(std::make_pair(inputIndex,
+ 					    positionDifference));
+       }
+     }
+ 
+ 
+     // Convert the minimum set to a set of FieldOffsets.
+     // minimumSet has all the minimum distance locations.
+ 
+     FieldOffset_vt answerHolder;
+     if (IntraCellOnly) {
+       for (MinimumSet::size_type minIndex = 0;
+ 	   minIndex < minimumSet.size();
+ 	   ++minIndex)
+ 	answerHolder.push_back(FieldOffset_t(Loc<Dim>(0),
+ 					     minimumSet[minIndex].first));
+     }
+     else {
+       FieldOffset_vt partialAnswer;
+       for (MinimumSet::size_type minIndex = 0;
+ 	   minIndex < minimumSet.size();
+ 	   ++minIndex)
+ 	{
+ 	  // Compute the cell offsets, appending to the set of answers.
+ 
+ 	  partialAnswer = computeCellOffsets(minimumSet[minIndex].first,
+ 					     minimumSet[minIndex].second);
+ 	  answerHolder.insert(answerHolder.end(),
+ 			      partialAnswer.begin(), partialAnswer.end());
+ 	}
+ 
+       // Remove all duplicates from the answer set.
+ 
+       std::sort(answerHolder.begin(), answerHolder.end(),
+ 		CompareFieldOffset());
+       answerHolder.erase(std::unique(answerHolder.begin(),
+ 				     answerHolder.end(),
+ 				     EqualFieldOffset()),
+ 			 answerHolder.end());
+     }
+ 
+     return answerHolder;
+   }
+ 
    // Given a difference between two positions in logical coordinate
    // space, return the Manhattan norm distance taking into account
    // that input values are repeated in every grid cell.
*************** private:
*** 240,246 ****
    inline static
    double manhattanDistance(const Position &difference)
    {
!     double answer;
      for (int coordinate = Dim-1; coordinate >= 0; --coordinate)
        answer = Distance()(answer, difference(coordinate));
      return answer;
--- 299,305 ----
    inline static
    double manhattanDistance(const Position &difference)
    {
!     double answer = 0.0;;
      for (int coordinate = Dim-1; coordinate >= 0; --coordinate)
        answer = Distance()(answer, difference(coordinate));
      return answer;
*************** private:
*** 356,364 ****
      }
    };
    
-   const Center &inputCentering_m;
-   const Center &outputCentering_m;
- 
    // Use epsilon when comparing floating-point numbers, which cannot
    // be represented precisely.
  
--- 415,420 ----
*************** std::vector<FieldOffsetList<Dim> >
*** 381,387 ****
  nearestNeighbors(const Centering<Dim> &inputCentering,
  		 const Centering<Dim> &outputCentering)
  {
!   return NearestNeighborClass<Dim>(inputCentering, outputCentering)();
  }
  
  template <int Dim>
--- 437,443 ----
  nearestNeighbors(const Centering<Dim> &inputCentering,
  		 const Centering<Dim> &outputCentering)
  {
!   return NearestNeighborClass<Dim>()(inputCentering, outputCentering);
  }
  
  template <int Dim>
*************** nearestNeighbors(const Centering<Dim> &i
*** 390,398 ****
  		 const Centering<Dim> &outputCentering,
  		 const bool)
  {
!   return NearestNeighborClass<Dim, true>(inputCentering, outputCentering)();
  }
  
  
  //-----------------------------------------------------------------------------
  // inputPosition(inputCentering, fieldOffset)
--- 446,494 ----
  		 const Centering<Dim> &outputCentering,
  		 const bool)
  {
!   return NearestNeighborClass<Dim, true>()(inputCentering, outputCentering);
  }
  
+ template <int Dim>
+ std::vector<FieldOffsetList<Dim> >
+ nearestNeighbors(const Centering<Dim> &inputCentering,
+ 		 const FieldOffsetList<Dim> &fOL,
+ 		 const Centering<Dim> &outputCentering)
+ {
+   return NearestNeighborClass<Dim>()(inputCentering, fOL, outputCentering);
+ }
+ 
+ template <int Dim>
+ std::vector<FieldOffsetList<Dim> >
+ nearestNeighbors(const Centering<Dim> &inputCentering,
+ 		 const FieldOffsetList<Dim> &fOL,
+ 		 const Centering<Dim> &outputCentering,
+ 		 const bool)
+ {
+   return NearestNeighborClass<Dim, true>()
+     (inputCentering, fOL, outputCentering);
+ }
+ 
+ template <int Dim>
+ FieldOffsetList<Dim>
+ nearestNeighbors(const Centering<Dim> &inputCentering,
+ 		 const FieldOffset<Dim> &fieldOffset,
+ 		 const Centering<Dim> &outputCentering)
+ {
+   return NearestNeighborClass<Dim>()
+     (inputCentering, fieldOffset, outputCentering);
+ }
+ 
+ template <int Dim>
+ FieldOffsetList<Dim>
+ nearestNeighbors(const Centering<Dim> &inputCentering,
+ 		 const FieldOffset<Dim> &fieldOffset,
+ 		 const Centering<Dim> &outputCentering,
+ 		 const bool)
+ {
+   return NearestNeighborClass<Dim, true>()
+     (inputCentering, fieldOffset, outputCentering);
+ }
  
  //-----------------------------------------------------------------------------
  // inputPosition(inputCentering, fieldOffset)
Index: tests/NearestNeighbors.cpp
===================================================================
RCS file: /home/pooma/Repository/r2/src/NewField/tests/Attic/NearestNeighbors.cpp,v
retrieving revision 1.1.2.2
diff -c -p -r1.1.2.2 NearestNeighbors.cpp
*** tests/NearestNeighbors.cpp	2001/08/14 20:24:18	1.1.2.2
--- tests/NearestNeighbors.cpp	2001/08/14 22:06:20
*************** int main(int argc, char *argv[])
*** 156,161 ****
--- 156,162 ----
  
    Centering<2> inputCenteringTwo, outputCenteringTwo;
    Centering<3> inputCenteringThree, outputCenteringThree;
+   FieldOffsetList<2> fieldOffsetListTwo;
  
    // Test 2D Continuous Cell -> Continuous Cell.
  
*************** int main(int argc, char *argv[])
*** 210,215 ****
--- 211,231 ----
  			     inputCenteringTwo,
  			     outputCenteringTwo));
  
+   fieldOffsetListTwo =
+     nearestNeighbors(inputCenteringTwo,
+ 		     FieldOffset<2>(Loc<2>(0,0)), outputCenteringTwo);
+   tester.check("vertex->cell intercell",
+ 	       fieldOffsetListTwo.size() == 4 && 
+ 	       checkForFieldOffset(fieldOffsetListTwo,
+ 				   FieldOffset<2>(Loc<2>(0,0))));
+ 
+   fieldOffsetListTwo =
+     nearestNeighbors(inputCenteringTwo,
+ 		     FieldOffset<2>(Loc<2>(0,0)), outputCenteringTwo, true);
+   tester.check("vertex->cell intracell",
+ 	       fieldOffsetListTwo.size() == 1 && 
+ 	       checkForFieldOffset(fieldOffsetListTwo,
+ 				   FieldOffset<2>(Loc<2>(0,0))));
  
    // Test 2D Discontinuous Vertex -> Continuous Cell.
    


More information about the pooma-dev mailing list