Changeset 630 for XIOS/trunk/src/transformation/grid_transformation.cpp

Timestamp:

07/07/15 10:46:25 (9 years ago)

Author:

mhnguyen

Message:

Implementing interpolation (polynomial) and correct some bugs

+) Implement interpolation (polynomial)
+) Correct some minor bugs relating to memory allocation
+) Clear some redundant codes

Test
+) On Curie
+) test_client and test_complete pass

File:

• XIOS/trunk/src/transformation/grid_transformation.cpp (modified) (21 diffs)

XIOS/trunk/src/transformation/grid_transformation.cpp

@@ -3 +3 @@
    \author Ha NGUYEN
    \since 14 May 2015
-   \date 18 June 2015
+   \date 02 Jul 2015
 
    \brief Interface for all transformations.
@@ -10 +10 @@
 #include "axis_algorithm_inverse.hpp"
 #include "axis_algorithm_zoom.hpp"
+#include "axis_algorithm_interpolate.hpp"
 #include "context.hpp"
 #include "context_client.hpp"
 #include "transformation_mapping.hpp"
-
 #include "axis_algorithm_transformation.hpp"
 
@@ -19 +19 @@
 CGridTransformation::CGridTransformation(CGrid* destination, CGrid* source)
 : gridSource_(source), gridDestination_(destination), originalGridSource_(source),
-  globalIndexOfCurrentGridSource_(0), globalIndexOfOriginalGridSource_(0)
+  globalIndexOfCurrentGridSource_(0), globalIndexOfOriginalGridSource_(0), weightOfGlobalIndexOfOriginalGridSource_(0)
 {
   //Verify the compatibity between two grids
@@ -78 +78 @@
   globalIndexOfCurrentGridSource_   = new CArray<size_t,1>(globalIndexGridDestSendToServer.numElements());
   globalIndexOfOriginalGridSource_  = new CArray<size_t,1>(globalIndexGridDestSendToServer.numElements());
+  weightOfGlobalIndexOfOriginalGridSource_  = new CArray<double,1>(globalIndexGridDestSendToServer.numElements());
   *globalIndexOfCurrentGridSource_  = globalIndexGridDestSendToServer;
   *globalIndexOfOriginalGridSource_ = globalIndexGridDestSendToServer;
+  *weightOfGlobalIndexOfOriginalGridSource_ = 1.0;
 }
 
@@ -87 +89 @@
                                                               ite = algoTransformation_.end();
   for (it = itb; it != ite; ++it) delete (*it);
-
-  std::map<int, std::vector<CArray<int,1>* > >::const_iterator itMapRecv, iteMapRecv;
-  itMapRecv = localIndexToReceiveOnGridDest_.begin();
-  iteMapRecv = localIndexToReceiveOnGridDest_.end();
-  for (; itMapRecv != iteMapRecv; ++itMapRecv)
-  {
-    int numVec = (itMapRecv->second).size();
-    for (int idx = 0; idx < numVec; ++idx) delete (itMapRecv->second)[idx];
-  }
 
   std::map<int, CArray<int,1>* >::const_iterator itMap, iteMap;
@@ -104 +97 @@
   if (0 != globalIndexOfCurrentGridSource_) delete globalIndexOfCurrentGridSource_;
   if (0 != globalIndexOfOriginalGridSource_) delete globalIndexOfOriginalGridSource_;
+  if (0 != weightOfGlobalIndexOfOriginalGridSource_) delete weightOfGlobalIndexOfOriginalGridSource_;
 }
 
@@ -194 +188 @@
 
   CZoomAxis* zoomAxis = 0;
+  CInterpolateAxis* interpAxis = 0;
   CGenericAlgorithmTransformation* algo = 0;
   switch (transType)
   {
+    case TRANS_INTERPOLATE_AXIS:
+      interpAxis = dynamic_cast<CInterpolateAxis*> (it->second);
+      algo = new CAxisAlgorithmInterpolate(axisListDestP[axisIndex], axisListSrcP[axisIndex], interpAxis);
+      break;
     case TRANS_ZOOM_AXIS:
       zoomAxis = dynamic_cast<CZoomAxis*> (it->second);
@@ -237 +236 @@
   switch (transType)
   {
+    case TRANS_INTERPOLATE_AXIS:
     case TRANS_ZOOM_AXIS:
     case TRANS_INVERSE_AXIS:
@@ -263 +263 @@
                                ite = listAlgos_.end(), it;
   CGenericAlgorithmTransformation* algo = 0;
+
   for (it = itb; it != ite; ++it)
   {
-    std::map<size_t, std::set<size_t> > globaIndexMapFromDestToSource;
     int elementPositionInGrid = it->first;
     ETranformationType transType = (it->second).first;
     int transformationOrder = (it->second).second;
+    std::map<size_t, std::vector<std::pair<size_t,double> > > globaIndexWeightFromDestToSource;
 
     // First of all, select an algorithm
@@ -284 +285 @@
                                    gridDestinationDimensionSize,
                                    globalIndexGridDestSendToServer,
-                                   globaIndexMapFromDestToSource);
+                                   globaIndexWeightFromDestToSource);
 
     // Compute transformation of global indexes among grids
-    computeTransformationFromOriginalGridSource(globaIndexMapFromDestToSource);
+    computeTransformationFromOriginalGridSource(globaIndexWeightFromDestToSource);
 
     // Now grid destination becomes grid source in a new transformation
@@ -357 +358 @@
 the final grid destination
 */
-void CGridTransformation::computeTransformationFromOriginalGridSource(const std::map<size_t, std::set<size_t> >& globaIndexMapFromDestToSource)
+void CGridTransformation::computeTransformationFromOriginalGridSource(const std::map<size_t, std::vector<std::pair<size_t,double> > >& globaIndexMapFromDestToSource)
 {
   CContext* context = CContext::getCurrent();
@@ -367 +368 @@
   transformationMap.computeTransformationMapping(globaIndexMapFromDestToSource);
 
-  const std::map<int,std::vector<std::vector<size_t> > >& globalIndexToReceive = transformationMap.getGlobalIndexReceivedOnGridDestMapping();
+  const std::map<int,std::vector<std::vector<std::pair<size_t,double> > > >& globalIndexToReceive = transformationMap.getGlobalIndexReceivedOnGridDestMapping();
   const std::map<int,std::vector<size_t> >& globalIndexToSend = transformationMap.getGlobalIndexSendToGridDestMapping();
 
@@ -404 +405 @@
 
  // Receiving global index of grid source sending from current grid source
- std::map<int,std::vector<std::vector<size_t> > >::const_iterator itbRecv = globalIndexToReceive.begin(), itRecv,
-                                                                  iteRecv = globalIndexToReceive.end();
+ std::map<int,std::vector<std::vector<std::pair<size_t,double> > > >::const_iterator itbRecv = globalIndexToReceive.begin(), itRecv,
+                                                                                     iteRecv = globalIndexToReceive.end();
  int recvBuffSize = 0;
  for (itRecv = itbRecv; itRecv != iteRecv; ++itRecv) recvBuffSize += (itRecv->second).size();
@@ -438 +439 @@
    delete globalIndexOfCurrentGridSource_;
    globalIndexOfCurrentGridSource_ = new CArray<size_t,1>(nbCurrentGridSource);
- }
-
- if (globalIndexOfOriginalGridSource_->numElements() != nbCurrentGridSource)
- {
    delete globalIndexOfOriginalGridSource_;
    globalIndexOfOriginalGridSource_ = new CArray<size_t,1>(nbCurrentGridSource);
+   delete weightOfGlobalIndexOfOriginalGridSource_;
+   weightOfGlobalIndexOfOriginalGridSource_ = new CArray<double,1>(nbCurrentGridSource);
  }
 
@@ -456 +455 @@
      for (int i = 0; i < ssize; ++i)
      {
-       (*globalIndexOfCurrentGridSource_)(k) = (itRecv->second)[idx][i];
+       (*globalIndexOfCurrentGridSource_)(k) = ((itRecv->second)[idx][i]).first;
+       (*weightOfGlobalIndexOfOriginalGridSource_)(k) = ((itRecv->second)[idx][i]).second;
        (*globalIndexOfOriginalGridSource_)(k) = *currentRecvBuff;
        ++k;
@@ -479 +479 @@
   CTransformationMapping transformationMap(gridDestination_, originalGridSource_);
 
-  std::map<size_t, std::set<size_t> > globaIndexMapFromDestToSource;
-
+  std::map<size_t, std::vector<std::pair<size_t,double> > > globaIndexWeightFromDestToSource;
   int nb = globalIndexOfCurrentGridSource_->numElements();
   const size_t sfmax = NumTraits<unsigned long>::sfmax();
@@ -486 +485 @@
   {
     if (sfmax != (*globalIndexOfOriginalGridSource_)(idx))
-      globaIndexMapFromDestToSource[(*globalIndexOfCurrentGridSource_)(idx)].insert((*globalIndexOfOriginalGridSource_)(idx));
+      globaIndexWeightFromDestToSource[(*globalIndexOfCurrentGridSource_)(idx)].push_back(make_pair((*globalIndexOfOriginalGridSource_)(idx),(*weightOfGlobalIndexOfOriginalGridSource_)(idx))) ;
   }
 
   // Then compute transformation mapping among clients
-  transformationMap.computeTransformationMapping(globaIndexMapFromDestToSource);
-
-  const std::map<int,std::vector<std::vector<size_t> > >& globalIndexToReceive = transformationMap.getGlobalIndexReceivedOnGridDestMapping();
+  transformationMap.computeTransformationMapping(globaIndexWeightFromDestToSource);
+
+  const std::map<int,std::vector<std::vector<std::pair<size_t,double> > > >& globalIndexToReceive = transformationMap.getGlobalIndexReceivedOnGridDestMapping();
   const std::map<int,std::vector<size_t> >& globalIndexToSend = transformationMap.getGlobalIndexSendToGridDestMapping();
 
@@ -498 +497 @@
   CDistributionClient distributionClientSrc(client->clientRank, originalGridSource_);
 
-//  const CArray<int, 1>& localIndexOnClientDest = distributionClientDest.getLocalDataIndexOnClient(); //gridDestination_->getDistributionClient()->getLocalDataIndexOnClient();
   const CArray<int, 1>& localIndexOnClientDest = distributionClientDest.getLocalDataIndexSendToServer();
-  const CArray<size_t,1>& globalIndexOnClientDest = distributionClientDest.getGlobalDataIndexSendToServer(); //gridDestination_->getDistributionClient()->getGlobalDataIndexSendToServer();
-
-  const CArray<int, 1>& localIndexOnClientSrc = distributionClientSrc.getLocalDataIndexOnClient(); //gridSource_->getDistributionClient()->getLocalDataIndexOnClient();
-  const CArray<size_t,1>& globalIndexOnClientSrc = distributionClientSrc.getGlobalDataIndexSendToServer(); //gridSource_->getDistributionClient()->getGlobalDataIndexSendToServer();
+  const CArray<size_t,1>& globalIndexOnClientDest = distributionClientDest.getGlobalDataIndexSendToServer();
+
+  const CArray<int, 1>& localIndexOnClientSrc = distributionClientSrc.getLocalDataIndexOnClient();
+  const CArray<size_t,1>& globalIndexOnClientSrc = distributionClientSrc.getGlobalDataIndexSendToServer();
 
   std::vector<size_t>::const_iterator itbVec, itVec, iteVec;
   CArray<size_t, 1>::const_iterator itbArr, itArr, iteArr;
 
-  std::map<int,std::vector<std::vector<size_t> > >::const_iterator itbMapRecv, itMapRecv, iteMapRecv;
+  std::map<int,std::vector<std::vector<std::pair<size_t,double> > > >::const_iterator itbMapRecv, itMapRecv, iteMapRecv;
 
   // Find out local index on grid destination (received)
@@ -522 +520 @@
     {
       int vecSize = ((itMapRecv->second)[i]).size();
-      CArray<int,1>* ptr = new CArray<int,1>(vecSize);
-      localIndexToReceiveOnGridDest_[sourceRank].push_back(ptr);
+      std::vector<std::pair<int,double> > tmpVec;
       for (int idx = 0; idx < vecSize; ++idx)
       {
-        itArr = std::find(itbArr, iteArr, (itMapRecv->second)[i][idx]);
+        size_t globalIndex = (itMapRecv->second)[i][idx].first;
+        double weight = (itMapRecv->second)[i][idx].second;
+        itArr = std::find(itbArr, iteArr, globalIndex);
         if (iteArr != itArr)
         {
           int localIdx = std::distance(itbArr, itArr);
-//          (*localIndexToReceiveOnGridDest_[sourceRank][i])(idx) = localIndexOnClientDest(localIdx); // Local index of un-extracted data (only domain)
-          (*localIndexToReceiveOnGridDest_[sourceRank][i])(idx) = (localIdx); // Local index of extracted data
+          tmpVec.push_back(make_pair(localIdx, weight));
         }
       }
+      localIndexToReceiveOnGridDest_[sourceRank].push_back(tmpVec);
     }
   }
 
+  // Find out local index on grid source (to send)
   std::map<int,std::vector<size_t> >::const_iterator itbMap, itMap, iteMap;
-  // Find out local index on grid source (to send)
   itbMap = globalIndexToSend.begin();
   iteMap = globalIndexToSend.end();
@@ -555 +554 @@
       {
         int localIdx = std::distance(itbArr, itArr);
-//        (*localIndexToSendFromGridSource_[destRank])(idx) = localIndexOnClientSrc(localIdx);
         (*localIndexToSendFromGridSource_[destRank])(idx) = (localIdx);
       }
@@ -575 +573 @@
   \return local index of data
 */
-const std::map<int, std::vector<CArray<int,1>* > >& CGridTransformation::getLocalIndexToReceiveOnGridDest() const
+const std::map<int,std::vector<std::vector<std::pair<int,double> > > >& CGridTransformation::getLocalIndexToReceiveOnGridDest() const
 {
   return localIndexToReceiveOnGridDest_;