Changeset 709

Timestamp:

09/29/15 14:01:53 (9 years ago)

Author:

mhnguyen

Message:

Correcting a bug in interpolation domain

+) Replace shared send buffer by seperate buffer for each proc
+) Remove some redundant codes

Test
+) On Curie
+) test_client, test_complete and test_remap pass

Location:

XIOS/trunk

Files:

• bld.cfg (modified) (1 diff)
• inputs/REMAP/iodef.xml (modified) (2 diffs)
• src/filter/spatial_transform_filter.cpp (modified) (5 diffs)
• src/io/nc4_data_output.cpp (modified) (1 diff)
• src/test/test_remap.f90 (modified) (1 diff)
• src/transformation/domain_algorithm_interpolate.cpp (modified) (8 diffs)
• src/transformation/domain_algorithm_interpolate.hpp (modified) (1 diff)
• src/transformation/generic_algorithm_transformation.cpp (modified) (1 diff)

XIOS/trunk/bld.cfg

@@ -30 +30 @@
 src::remap $PWD/extern/remap/src
 bld::lib xios
-bld::lib::remap mapper
-bld::target libxios.a libmapper.a
+bld::target libxios.a 
 #bld::target generate_fortran_interface.exe 
 bld::target xios_server.exe test_remap.exe
-#bld::target test_new_features.exe test_unstruct_complete.exe 
+bld::target test_new_features.exe test_unstruct_complete.exe 
 bld::target test_client.exe test_complete.exe
 bld::exe_dep

XIOS/trunk/inputs/REMAP/iodef.xml

@@ -13 +13 @@
 
 
-   <file_definition type="multiple_file" par_access="collective" output_freq="1ts" output_level="10" enabled=".TRUE.">
+   <file_definition type="one_file" par_access="collective" output_freq="1ts" output_level="10" enabled=".TRUE.">
      <file id="output" name="output">
         <field field_ref="src_field" name="field" />
@@ -32 +32 @@
      <domain id="src_domain" />
      <domain id="dst_domain" domain_src="src_domain">
-       <interpolate_domain file="weight.nc"/>
+<!--       <interpolate_domain file="weight.nc"/>-->
+       <interpolate_domain/>
      </domain>
-     <domain id="dst_domain_regular" domain_src="src_domain">
+     <domain id="dst_domain_regular" domain_src="src_domain" ni_glo="180" nj_glo="90">
        <generate_rectilinear_domain />
        <interpolate_domain/>

XIOS/trunk/src/filter/spatial_transform_filter.cpp

@@ -80 +80 @@
     std::map<int, CArray<int,1> >::const_iterator itbSend = localIndexToSend.begin(), itSend,
                                                   iteSend = localIndexToSend.end();
-    int sendBuffSize = 0;
-    for (itSend = itbSend; itSend != iteSend; ++itSend) sendBuffSize = (sendBuffSize < itSend->second.numElements())
-                                                                     ? itSend->second.numElements(): sendBuffSize;
-    double* sendBuff;
-    if (0 != sendBuffSize) sendBuff = new double[sendBuffSize];
+    int idxSendBuff = 0;
+    std::vector<double*> sendBuff(localIndexToSend.size());
+    for (itSend = itbSend; itSend != iteSend; ++itSend, ++idxSendBuff)
+    {
+      if (0 != itSend->second.numElements())
+        sendBuff[idxSendBuff] = new double[itSend->second.numElements()];
+    }
+
+    idxSendBuff = 0;
     std::vector<MPI_Request> sendRequest;
-    for (itSend = itbSend; itSend != iteSend; ++itSend)
+    for (itSend = itbSend; itSend != iteSend; ++itSend, ++idxSendBuff)
     {
       int destRank = itSend->first;
@@ -93 +97 @@
       for (int idx = 0; idx < countSize; ++idx)
       {
-        sendBuff[idx] = dataSrc(localIndex_p(idx));
+        sendBuff[idxSendBuff][idx] = dataSrc(localIndex_p(idx));
       }
       sendRequest.push_back(MPI_Request());
-      MPI_Isend(sendBuff, countSize, MPI_DOUBLE, destRank, 12, client->intraComm, &sendRequest.back());
+      MPI_Isend(sendBuff[idxSendBuff], countSize, MPI_DOUBLE, destRank, 12, client->intraComm, &sendRequest.back());
     }
 
@@ -113 +117 @@
       int countSize = itRecv->second.size();
       MPI_Recv(recvBuff, recvBuffSize, MPI_DOUBLE, srcRank, 12, client->intraComm, &status);
+      int countBuff = 0;
+      MPI_Get_count(&status, MPI_DOUBLE, &countBuff);
+      if (countBuff != countSize)
+        ERROR("CSpatialTransformFilterEngine::apply(const CArray<double, 1>& dataSrc, CArray<double,1>& dataDest)",
+              "Incoherent between the received size and expected size");
       for (int idx = 0; idx < countSize; ++idx)
       {
@@ -119 +128 @@
         for (int i = 0; i < numIndex; ++i)
         {
-//        if (localIndex_p[i].first >= dataDest.numElements())
           dataDest(localIndex_p[i].first) += recvBuff[idx] * localIndex_p[i].second;
         }
@@ -125 +133 @@
     }
 
-    std::vector<MPI_Status> requestStatus(sendRequest.size());
-    if (!sendRequest.empty()) MPI_Wait(&sendRequest[0], &requestStatus[0]);
-    if (0 != sendBuffSize) delete [] sendBuff;
+
+    if (!sendRequest.empty()) MPI_Waitall(sendRequest.size(), &sendRequest[0], MPI_STATUSES_IGNORE);
+    idxSendBuff = 0;
+    for (itSend = itbSend; itSend != iteSend; ++itSend, ++idxSendBuff)
+    {
+      if (0 != itSend->second.numElements())
+        delete [] sendBuff[idxSendBuff];
+    }
     if (0 != recvBuffSize) delete [] recvBuff;
   }

XIOS/trunk/src/io/nc4_data_output.cpp

@@ -1625 +1625 @@
                   }
 /*
-                  int ssize = nZoomBeginGlobal.size();
                   for (int i = numElement - 1; i >= 0; --i)
 
 
-                  start.resize(ssize);
-                  count.resize(ssize);
-
-                  for (int i = 0; i < ssize; ++i)
+                  start.reserve(nZoomBeginGlobal.size());
+                  count.reserve(nZoomBeginGlobal.size());
+
+
+                  for (int i = numElement - 1; i >= 0; --i)
                   {
-                    start[i] = nZoomBeginServer[ssize - i - 1] - nZoomBeginGlobal[ssize - i - 1];
-                    count[i] = nZoomSizeServer[ssize - i - 1];
+                    if (axisDomainOrder(i))
+                    {
+                      CDomain* domain = CDomain::get(domainList[idxDomain]);
+
+                      if (CDomain::type_attr::unstructured != domain->type)
+                      {
+                        start.push_back(nZoomBeginServer[idx] - nZoomBeginGlobal[idx]);
+                        count.push_back(nZoomSizeServer[idx]);
+                      }
+                      --idx;
+                      start.push_back(nZoomBeginServer[idx] - nZoomBeginGlobal[idx]);
+                      count.push_back(nZoomSizeServer[idx]);
+                      --idx;
+
+                      --idxDomain;
+                      //idx -= 2;
+                    }
+                    else
+                    {
+                      CAxis* axis = CAxis::get(axisList[idxAxis]);
+
+                      start.push_back(nZoomBeginServer[idx] - nZoomBeginGlobal[idx]);
+                      count.push_back(nZoomSizeServer[idx]);
+
+                      --idxAxis;
+                      --idx;
+                    }
                   }
-*/
+
+//                  int ssize = nZoomBeginGlobal.size();
+//
+//                  start.resize(ssize);
+//                  count.resize(ssize);
+//
+//                  for (int i = 0; i < ssize; ++i)
+//                  {
+//                    start[i] = nZoomBeginServer[ssize - i - 1] - nZoomBeginGlobal[ssize - i - 1];
+//                    count[i] = nZoomSizeServer[ssize - i - 1];
+//                  }
                 }
 

XIOS/trunk/src/test/test_remap.f90

@@ -119 +119 @@
                             bounds_lon_1D=dst_boundslon, bounds_lat_1D=dst_boundslat, nvertex=dst_nvertex)
 
-  CALL xios_set_domain_attr("dst_domain_regular", ni_glo=180, nj_glo=90, type="rectilinear")
+  CALL xios_set_domain_attr("dst_domain_regular", type="rectilinear")
 
   dtime%second = 3600

XIOS/trunk/src/transformation/domain_algorithm_interpolate.cpp

@@ -76 +76 @@
   {
     nVertexSrc = constNVertex;
+    domainSrc_->fillInRectilinearBoundLonLat(boundsLonSrc, boundsLatSrc);
   }
 
@@ -114 +115 @@
   }
 
+
+
   // Ok, now use mapper to calculate
   int nSrcLocal = domainSrc_->i_index.numElements();
   int nDstLocal = domainDest_->i_index.numElements();
+  long int * globalSrc = new long int [nSrcLocal];
+  long int * globalDst = new long int [nDstLocal];
+
+  long int globalIndex;
+  int i_ind, j_ind;
+  for (int idx = 0; idx < nSrcLocal; ++idx)
+  {
+    i_ind=domainSrc_->i_index(idx) ;
+    j_ind=domainSrc_->j_index(idx) ;
+
+    globalIndex = i_ind + j_ind * domainSrc_->ni_glo;
+    globalSrc[idx] = globalIndex;
+  }
+
+  for (int idx = 0; idx < nDstLocal; ++idx)
+  {
+    i_ind=domainDest_->i_index(idx) ;
+    j_ind=domainDest_->j_index(idx) ;
+
+    globalIndex = i_ind + j_ind * domainDest_->ni_glo;
+    globalDst[idx] = globalIndex;
+  }
+
+
+  // Calculate weight index
   Mapper mapper(client->intraComm);
   mapper.setVerbosity(PROGRESS) ;
-  mapper.setSourceMesh(boundsLonSrc.dataFirst(), boundsLatSrc.dataFirst(), nVertexSrc, nSrcLocal, &srcPole[0]);
-  mapper.setTargetMesh(boundsLonDest.dataFirst(), boundsLatDest.dataFirst(), nVertexDest, nDstLocal, &dstPole[0]);
+  mapper.setSourceMesh(boundsLonSrc.dataFirst(), boundsLatSrc.dataFirst(), nVertexSrc, nSrcLocal, &srcPole[0], globalSrc);
+  mapper.setTargetMesh(boundsLonDest.dataFirst(), boundsLatDest.dataFirst(), nVertexDest, nDstLocal, &dstPole[0], globalDst);
   std::vector<double> timings = mapper.computeWeights(orderInterp);
 
+  std::map<int,std::vector<std::pair<int,double> > > interpMapValue;
   for (int idx = 0;  idx < mapper.nWeights; ++idx)
   {
-    transformationMapping_[mapper.targetWeightId[idx]].push_back(mapper.sourceWeightId[idx]);
-    transformationWeight_[mapper.targetWeightId[idx]].push_back(mapper.remapMatrix[idx]);
-  }
+    interpMapValue[mapper.targetWeightId[idx]].push_back(make_pair(mapper.sourceWeightId[idx],mapper.remapMatrix[idx]));
+  }
+  exchangeRemapInfo(interpMapValue);
+
+  delete [] globalSrc;
+  delete [] globalDst;
 }
 
@@ -141 +173 @@
 }
 
+void CDomainAlgorithmInterpolate::readRemapInfo()
+{
+  CContext* context = CContext::getCurrent();
+  CContextClient* client=context->client;
+  int clientRank = client->clientRank;
+
+  std::string filename = interpDomain_->file.getValue();
+  std::map<int,std::vector<std::pair<int,double> > > interpMapValue;
+  readInterpolationInfo(filename, interpMapValue);
+
+  exchangeRemapInfo(interpMapValue);
+}
+
+
 /*!
   Read remap information from file then distribute it among clients
 */
-void CDomainAlgorithmInterpolate::readRemapInfo()
+void CDomainAlgorithmInterpolate::exchangeRemapInfo(const std::map<int,std::vector<std::pair<int,double> > >& interpMapValue)
 {
   CContext* context = CContext::getCurrent();
   CContextClient* client=context->client;
   int clientRank = client->clientRank;
-
-  std::string filename = interpDomain_->file.getValue();
-  std::map<int,std::vector<std::pair<int,double> > > interpMapValue;
-  readInterpolationInfo(filename, interpMapValue);
 
   boost::unordered_map<size_t,int> globalIndexOfDomainDest;
@@ -189 +231 @@
   int* sendBuff = new int[nbClient];
   int* recvBuff = new int[nbClient];
-  for (int i = 0; i < nbClient; ++i) sendBuff[i] = 0;
+  for (int i = 0; i < nbClient; ++i)
+  {
+    sendBuff[i] = 0;
+    recvBuff[i] = 0;
+  }
   int sendBuffSize = 0;
   std::map<int, std::vector<size_t> >::const_iterator itbMap = globalIndexInterpSendToClient.begin(), itMap,
@@ -195 +241 @@
   for (itMap = itbMap; itMap != iteMap; ++itMap)
   {
+    const std::vector<size_t>& tmp = itMap->second;
     int sizeIndex = 0, mapSize = (itMap->second).size();
     for (int idx = 0; idx < mapSize; ++idx)
     {
-      sizeIndex += interpMapValue[(itMap->second)[idx]].size();
+      sizeIndex += interpMapValue.at((itMap->second)[idx]).size();
     }
     sendBuff[itMap->first] = sizeIndex;
     sendBuffSize += sizeIndex;
   }
+
 
   MPI_Allreduce(sendBuff, recvBuff, nbClient, MPI_INT, MPI_SUM, client->intraComm);
@@ -211 +259 @@
 
   std::vector<MPI_Request> sendRequest;
-  // Now send index and weight
-  int sendOffSet = 0;
+
+  int sendOffSet = 0, l = 0;
   for (itMap = itbMap; itMap != iteMap; ++itMap)
   {
+    const std::vector<size_t>& indexToSend = itMap->second;
+    int mapSize = indexToSend.size();
     int k = 0;
-    int mapSize = (itMap->second).size();
-
     for (int idx = 0; idx < mapSize; ++idx)
     {
-      std::vector<std::pair<int,double> >& interpMap = interpMapValue[(itMap->second)[idx]];
+      const std::vector<std::pair<int,double> >& interpMap = interpMapValue.at(indexToSend[idx]);
       for (int i = 0; i < interpMap.size(); ++i)
       {
-        sendIndexDestBuff[k] = (itMap->second)[idx];
-        sendIndexSrcBuff[k]  = interpMap[i].first;
-        sendWeightBuff[k]    = interpMap[i].second;
+        sendIndexDestBuff[l] = indexToSend[idx];
+        sendIndexSrcBuff[l]  = interpMap[i].first;
+        sendWeightBuff[l]    = interpMap[i].second;
         ++k;
+        ++l;
       }
     }
@@ -303 +352 @@
 
   std::vector<MPI_Status> requestStatus(sendRequest.size());
-  MPI_Wait(&sendRequest[0], &requestStatus[0]);
+  MPI_Waitall(sendRequest.size(), &sendRequest[0], MPI_STATUS_IGNORE);
 
   delete [] sendIndexDestBuff;
@@ -322 +371 @@
 */
 void CDomainAlgorithmInterpolate::readInterpolationInfo(std::string& filename,
-                                                                std::map<int,std::vector<std::pair<int,double> > >& interpMapValue)
+                                                        std::map<int,std::vector<std::pair<int,double> > >& interpMapValue)
 {
   int ncid ;

XIOS/trunk/src/transformation/domain_algorithm_interpolate.hpp

@@ -32 +32 @@
   void computeRemap();
   void readRemapInfo();
-
+  void exchangeRemapInfo(const std::map<int,std::vector<std::pair<int,double> > >& interpMapValue);
 private:
   CInterpolateDomain* interpDomain_;

XIOS/trunk/src/transformation/generic_algorithm_transformation.cpp

@@ -32 +32 @@
                                                              std::map<size_t, std::vector<std::pair<size_t,double> > >& globaIndexWeightFromDestToSource)
 {
-  std::map<int, std::vector<int> >::const_iterator itbTransMap = transformationMapping_.begin(),
-                                                   itTransMap = itbTransMap,
+  std::map<int, std::vector<int> >::const_iterator itbTransMap = transformationMapping_.begin(), itTransMap,
                                                    iteTransMap = transformationMapping_.end();
   std::map<int, std::vector<double> >::const_iterator itTransWeight = transformationWeight_.begin();
-  std::map<size_t, std::vector<std::pair<size_t,double> > >::iterator iteWeight, itWeight;
-  std::vector<int>::const_iterator itbVec, itVec, iteVec;
   std::vector<std::vector<size_t> > globalIndexSrcGrid;
   CArray<size_t,1> globalIndexDestGrid;