# source:XIOS/dev/branch_yushan_merged/src/transformation/grid_transformation.cpp@1196

Last change on this file since 1196 was 1196, checked in by yushan, 7 years ago

add request_check. test client and complete OK

• Property svn:executable set to `*`
File size: 28.9 KB
Line
1/*!
2   \file grid_transformation.cpp
3   \author Ha NGUYEN
4   \since 14 May 2015
5   \date 02 Jul 2015
6
7   \brief Interface for all transformations.
8 */
9#include "grid_transformation.hpp"
10#include "grid_transformation_factory_impl.hpp"
11#include "algo_types.hpp"
12#include "context.hpp"
13#include "context_client.hpp"
14#include "distribution_client.hpp"
15#include "mpi_tag.hpp"
16#include "grid.hpp"
17#include <boost/unordered_map.hpp>
18
19namespace xios {
20CGridTransformation::CGridTransformation(CGrid* destination, CGrid* source)
21 : CGridTransformationSelector(destination, source),
22  tmpGridDestination_(destination), originalGridSource_(source),
23  tempGridSrcs_(), tempGridDests_(),
24  dynamicalTransformation_(false), timeStamp_()
25{
26}
27
28CGridTransformation::~CGridTransformation()
29{
30}
31
32/*!
33  Select algorithm of a scalar corresponding to its transformation type and its position in each element
34  \param [in] elementPositionInGrid position of element in grid. e.g: a grid has 1 domain and 1 axis, then position of domain is 0 and position of axis is 1
35  \param [in] transType transformation type, for now we have
36  \param [in] transformationOrder position of the transformation in an element (an element can have several transformation)
37*/
38void CGridTransformation::selectScalarAlgo(int elementPositionInGrid, ETranformationType transType, int transformationOrder)
39{
40  std::vector<CScalar*> scaListDestP = gridDestination_->getScalars();
41  int scalarDstIndex =  elementPositionInGridDst2ScalarPosition_[elementPositionInGrid];
42  CScalar::TransMapTypes trans = scaListDestP[scalarDstIndex]->getAllTransformations();
43  CScalar::TransMapTypes::const_iterator it = trans.begin();
44
45  for (int i = 0; i < transformationOrder; ++i, ++it) {}  // Find the correct transformation
46  CGenericAlgorithmTransformation* algo = 0;
47  algo = CGridTransformationFactory<CScalar>::createTransformation(transType,
48                                                                  gridDestination_,
49                                                                  gridSource_,
50                                                                  it->second,
51                                                                  elementPositionInGrid,
52                                                                  elementPositionInGridSrc2ScalarPosition_,
53                                                                  elementPositionInGridSrc2AxisPosition_,
54                                                                  elementPositionInGridSrc2DomainPosition_,
55                                                                  elementPositionInGridDst2ScalarPosition_,
56                                                                  elementPositionInGridDst2AxisPosition_,
57                                                                  elementPositionInGridDst2DomainPosition_);
58  algoTransformation_.push_back(algo);
59}
60
61/*!
62  Select algorithm of an axis corresponding to its transformation type and its position in each element
63  \param [in] elementPositionInGrid position of element in grid. e.g: a grid has 1 domain and 1 axis, then position of domain is 0 and position of axis is 1
64  \param [in] transType transformation type, for now we have zoom_axis, inverse_axis, interpolate_axis
65  \param [in] transformationOrder position of the transformation in an element (an element can have several transformation)
66*/
67void CGridTransformation::selectAxisAlgo(int elementPositionInGrid, ETranformationType transType, int transformationOrder)
68{
69  std::vector<CAxis*> axisListDestP = gridDestination_->getAxis();
70  int axisDstIndex =  elementPositionInGridDst2AxisPosition_[elementPositionInGrid];
71  CAxis::TransMapTypes trans = axisListDestP[axisDstIndex]->getAllTransformations();
72  CAxis::TransMapTypes::const_iterator it = trans.begin();
73  for (int i = 0; i < transformationOrder; ++i, ++it) {}  // Find the correct transformation
74
75  CGenericAlgorithmTransformation* algo = 0;
76  algo = CGridTransformationFactory<CAxis>::createTransformation(transType,
77                                                                 gridDestination_,
78                                                                 gridSource_,
79                                                                 it->second,
80                                                                 elementPositionInGrid,
81                                                                 elementPositionInGridSrc2ScalarPosition_,
82                                                                 elementPositionInGridSrc2AxisPosition_,
83                                                                 elementPositionInGridSrc2DomainPosition_,
84                                                                 elementPositionInGridDst2ScalarPosition_,
85                                                                 elementPositionInGridDst2AxisPosition_,
86                                                                 elementPositionInGridDst2DomainPosition_);
87  algoTransformation_.push_back(algo);
88}
89
90/*!
91  Select algorithm of a domain corresponding to its transformation type and its position in each element
92  \param [in] elementPositionInGrid position of element in grid. e.g: a grid has 1 domain and 1 axis, then position of domain is 0 and position of axis is 1
93  \param [in] transType transformation type, for now we have zoom_domain, interpolate_domain
94  \param [in] transformationOrder position of the transformation in an element (an element can have several transformation)
95*/
96void CGridTransformation::selectDomainAlgo(int elementPositionInGrid, ETranformationType transType, int transformationOrder)
97{
98  std::vector<CDomain*> domainListDestP = gridDestination_->getDomains();
99  int domainIndex =  elementPositionInGridDst2DomainPosition_[elementPositionInGrid];
100  CDomain::TransMapTypes trans = domainListDestP[domainIndex]->getAllTransformations();
101  CDomain::TransMapTypes::const_iterator it = trans.begin();
102  for (int i = 0; i < transformationOrder; ++i, ++it) {}  // Find the correct transformation
103
104  CGenericAlgorithmTransformation* algo = 0;
105  algo = CGridTransformationFactory<CDomain>::createTransformation(transType,
106                                                                   gridDestination_,
107                                                                   gridSource_,
108                                                                   it->second,
109                                                                   elementPositionInGrid,
110                                                                   elementPositionInGridSrc2ScalarPosition_,
111                                                                   elementPositionInGridSrc2AxisPosition_,
112                                                                   elementPositionInGridSrc2DomainPosition_,
113                                                                   elementPositionInGridDst2ScalarPosition_,
114                                                                   elementPositionInGridDst2AxisPosition_,
115                                                                   elementPositionInGridDst2DomainPosition_);
116  algoTransformation_.push_back(algo);
117}
118
119/*!
120  Find position of element in a grid as well as its type (domain, axis, scalar) and position in its own element list
121  \return element position: map<int,<int,int> > corresponds to <element position in grid, <element type, element position in element list> >
122*/
123std::map<int,std::pair<int,int> > CGridTransformation::getElementPosition(CGrid* grid)
124{
125  std::vector<CScalar*> scalarListP = grid->getScalars();
126  std::vector<CAxis*> axisListP = grid->getAxis();
127  std::vector<CDomain*> domListP = grid->getDomains();
128  CArray<int,1> axisDomainOrder = grid->axis_domain_order;
129  int scalarIndex = 0, axisIndex = 0, domainIndex = 0;
130  int nbElement = axisDomainOrder.numElements(), elementDim;
131  std::map<int,std::pair<int,int> > elementPosition;
132  for (int idx = 0; idx < nbElement; ++idx)
133  {
134    elementDim = axisDomainOrder(idx);
135    switch (elementDim)
136    {
137      case 2:
138        elementPosition[idx] = std::make_pair(elementDim, domainIndex);
139        ++domainIndex;
140        break;
141      case 1:
142        elementPosition[idx] = std::make_pair(elementDim, axisIndex);
143        ++axisIndex;
144        break;
145      case 0:
146        elementPosition[idx] = std::make_pair(elementDim, scalarIndex);
147        ++scalarIndex;
148        break;
149      default:
150        break;
151    }
152  }
153
154  return elementPosition;
155}
156
157/*!
158  If there are more than one transformation, a new temporary grid will be created and it will play the role of grid destination.
159This new created one keeps a pointer to the real transformed element of grid destination and generate new copies of other elements from grid source.
160  \param [in] elementPositionInGrid position of element in grid
161  \param [in] transType transformation type
162*/
163void CGridTransformation::setUpGridDestination(int elementPositionInGrid, ETranformationType transType)
164{
165  if (isSpecialTransformation(transType)) return;
166
167  if (!tempGridDests_.empty() && (getNbAlgo() == tempGridDests_.size()))
168  {
169    tempGridDests_.resize(0);
170  }
171
172  if (1 == getNbAlgo())
173  {
174    tmpGridDestination_ = gridDestination_;
175    return;
176  }
177
178  std::vector<CScalar*> scalarListDestP = gridDestination_->getScalars();
179  std::vector<CScalar*> scalarListSrcP = gridSource_->getScalars(), scalarDst;
180
181  std::vector<CAxis*> axisListDestP = gridDestination_->getAxis();
182  std::vector<CAxis*> axisListSrcP = gridSource_->getAxis(), axisDst;
183
184  std::vector<CDomain*> domListDestP = gridDestination_->getDomains();
185  std::vector<CDomain*> domListSrcP = gridSource_->getDomains(), domainDst;
186
187  CArray<int,1> axisDomainOrderSrc = gridSource_->axis_domain_order;
188  CArray<int,1> axisDomainOrderDst = gridDestination_->axis_domain_order;
189
190  std::map<int,std::pair<int,int> > elementPositionSrc = getElementPosition(gridSource_);
191  std::map<int,std::pair<int,int> > elementPositionDst = getElementPosition(gridDestination_);
192
193  CArray<int,1> elementOrder(axisDomainOrderDst.numElements());
194  for (int idx = 0; idx < axisDomainOrderDst.numElements(); ++idx)
195  {
196    if (elementPositionInGrid == idx)
197    {
198      int dimElementDst = elementPositionDst[idx].first;
199      int elementIndex  = elementPositionDst[idx].second;
200      switch (dimElementDst)
201      {
202        case 2:
203          domainDst.push_back(domListDestP[elementIndex]);
204          break;
205        case 1:
206          axisDst.push_back(axisListDestP[elementIndex]);
207          break;
208        case 0:
209          scalarDst.push_back(scalarListDestP[elementIndex]);
210          break;
211        default:
212          break;
213      }
214      elementOrder(idx) = dimElementDst;
215    }
216    else
217    {
218      int dimElementSrc = elementPositionSrc[idx].first;
219      int elementIndex  = elementPositionSrc[idx].second;
220      switch (dimElementSrc)
221      {
222        case 2:
223          domainDst.push_back(domListSrcP[elementIndex]);
224          break;
225        case 1:
226          axisDst.push_back(axisListSrcP[elementIndex]);
227          break;
228        case 0:
229          scalarDst.push_back(scalarListSrcP[elementIndex]);
230          break;
231        default:
232          break;
233      }
234      elementOrder(idx) = dimElementSrc;
235    }
236  }
237
238  tmpGridDestination_ = CGrid::createGrid(domainDst, axisDst, scalarDst, elementOrder);
239  tempGridDests_.push_back(tmpGridDestination_);
240}
241
242/*!
243  Assign the current grid destination to the grid source in the new transformation.
244The current grid destination plays the role of grid source in next transformation (if any).
245Only element on which the transformation is performed is modified
246  \param [in] elementPositionInGrid position of element in grid
247  \param [in] transType transformation type
248*/
249void CGridTransformation::setUpGridSource(int elementPositionInGrid)
250{
251  if (!tempGridSrcs_.empty() && (getNbAlgo()-1) == tempGridSrcs_.size())
252  {
253    tempGridSrcs_.resize(0);
254  }
255
256  std::vector<CScalar*> scalarListDestP = tmpGridDestination_->getScalars();
257  std::vector<CScalar*> scalarListSrcP = gridSource_->getScalars(), scalarSrc;
258
259  std::vector<CAxis*> axisListDestP = tmpGridDestination_->getAxis();
260  std::vector<CAxis*> axisListSrcP = gridSource_->getAxis(), axisSrc;
261
262  std::vector<CDomain*> domListDestP = tmpGridDestination_->getDomains();
263  std::vector<CDomain*> domListSrcP = gridSource_->getDomains(), domainSrc;
264
265  CArray<int,1> axisDomainOrderSrc = gridSource_->axis_domain_order;
266  CArray<int,1> axisDomainOrderDst = tmpGridDestination_->axis_domain_order;
267
268  std::map<int,std::pair<int,int> > elementPositionSrc = getElementPosition(gridSource_);
269  std::map<int,std::pair<int,int> > elementPositionDst = getElementPosition(tmpGridDestination_);
270
271  for (int idx = 0; idx < axisDomainOrderDst.numElements(); ++idx)
272  {
273    if (elementPositionInGrid == idx)
274    {
275      int dimElementDst = elementPositionDst[idx].first;
276      int elementIndex  = elementPositionDst[idx].second;
277      if (2 == dimElementDst)
278      {
279        CDomain* domain = CDomain::createDomain();
280        domain->domain_ref.setValue(domListDestP[elementIndex]->getId());
281        domain->solveRefInheritance(true);
282        domain->checkAttributesOnClient();
283        domainSrc.push_back(domain);
284      }
285      else if (1 == dimElementDst)
286      {
287        CAxis* axis = CAxis::createAxis();
288        axis->axis_ref.setValue(axisListDestP[elementIndex]->getId());
289        axis->solveRefInheritance(true);
290        axis->checkAttributesOnClient();
291        axisSrc.push_back(axis);
292      }
293      else
294      {
295        CScalar* scalar = CScalar::createScalar();
296        scalar->scalar_ref.setValue(scalarListDestP[elementIndex]->getId());
297        scalar->solveRefInheritance(true);
298        scalar->checkAttributesOnClient();
299        scalarSrc.push_back(scalar);
300      }
301    }
302    else
303    {
304      int dimElementDst = elementPositionDst[idx].first;
305      int elementIndex  = elementPositionDst[idx].second;
306      switch (dimElementDst)
307      {
308        case 2:
309          domainSrc.push_back(domListDestP[elementIndex]);
310          break;
311        case 1:
312          axisSrc.push_back(axisListDestP[elementIndex]);
313          break;
314        case 0:
315          scalarSrc.push_back(scalarListDestP[elementIndex]);
316          break;
317        default:
318          break;
319      }
320    }
321  }
322
323  gridSource_ = CGrid::createGrid(domainSrc, axisSrc, scalarSrc, tmpGridDestination_->axis_domain_order);
324
325  tempGridSrcs_.push_back(gridSource_);
326}
327
328/*!
329  Perform all transformations
330  For each transformation, there are some things to do:
331  -) Chose the correct algorithm by transformation type and position of element
332  -) Calculate the mapping of global index between the current grid source and grid destination
333  -) Calculate the mapping of global index between current grid DESTINATION and grid SOURCE
334  -) Make current grid destination become grid source in the next transformation
335*/
336void CGridTransformation::computeAll(const std::vector<CArray<double,1>* >& dataAuxInputs, Time timeStamp)
337{
338  if (nbNormalAlgos_ < 1) return;
339  if (!auxInputs_.empty() && !dynamicalTransformation_) { dynamicalTransformation_ = true; return; }
340  if (dynamicalTransformation_)
341  {
342    if (timeStamp_.insert(timeStamp).second)   //Reset map
343    {
344      std::list<SendingIndexGridSourceMap>().swap(localIndexToSendFromGridSource_);
346      std::list<size_t>().swap(nbLocalIndexOnGridDest_);
348    }
349    else
350      return;
351  }
352
353  CContext* context = CContext::getCurrent();
354  CContextClient* client = context->client;
355
356  ListAlgoType::const_iterator itb = listAlgos_.begin(),
357                               ite = listAlgos_.end(), it;
358
359  CGenericAlgorithmTransformation* algo = 0;
360  int nbAgloTransformation = 0; // Only count for executed transformation. Generate domain is a special one, not executed in the list
361  for (it = itb; it != ite; ++it)
362  {
363    int elementPositionInGrid = it->first;
364    ETranformationType transType = (it->second).first;
365    int transformationOrder = (it->second).second.first;
366    int algoType = ((it->second).second.second); //algoTypes_[std::distance(itb, it)];
367    SourceDestinationIndexMap globaIndexWeightFromSrcToDst;
368
369
370    // Create a temporary grid destination which contains transformed element of grid destination and
371    // non-transformed elements to grid source
372    setUpGridDestination(elementPositionInGrid, transType);
373
374    // First of all, select an algorithm
375    if (!dynamicalTransformation_ || (algoTransformation_.size() < listAlgos_.size()))
376    {
377      selectAlgo(elementPositionInGrid, transType, transformationOrder, algoType);
378      algo = algoTransformation_.back();
379    }
380    else
381      algo = algoTransformation_[std::distance(itb, it)];
382
383    if ((0 != algo) &&
384        ((CGenericAlgorithmTransformation::ELEMENT_NO_MODIFICATION_WITH_DATA == algo->type()) ||
385        (CGenericAlgorithmTransformation::ELEMENT_MODIFICATION_WITH_DATA == algo->type()))) // Only registered transformation can be executed
386    {
387      algo->computeIndexSourceMapping(dataAuxInputs);
388
389      // ComputeTransformation of global index of each element
390      int elementPosition = it->first;
391      algo->computeGlobalSourceIndex(elementPosition,
392                                     gridSource_,
393                                     tmpGridDestination_,
394                                     globaIndexWeightFromSrcToDst);
395
396      // Compute transformation of global indexes among grids
397      computeTransformationMapping(globaIndexWeightFromSrcToDst);
398
399      if (1 < nbNormalAlgos_)
400      {
401        // Now grid destination becomes grid source in a new transformation
402        if (nbAgloTransformation != (nbNormalAlgos_-1)) setUpGridSource(elementPositionInGrid);
403      }
404      ++nbAgloTransformation;
405    }
406  }
407}
408
409/*!
410  Compute exchange index between grid source and grid destination
411  \param [in] globalIndexWeightFromDestToSource global index mapping between grid destination and grid source
412*/
413void CGridTransformation::computeTransformationMapping(const SourceDestinationIndexMap& globaIndexWeightFromSrcToDst)
414{
415  CContext* context = CContext::getCurrent();
416  CContextClient* client = context->client;
417  int nbClient = client->clientSize;
418  int clientRank = client->clientRank;
419
420  // Recalculate the distribution of grid destination
421  CDistributionClient distributionClientDest(client->clientRank, tmpGridDestination_);
422  CDistributionClient::GlobalLocalDataMap& globalLocalIndexGridDestSendToServer = distributionClientDest.getGlobalLocalDataSendToServer();
423
424  // Update number of local index on each transformation
425  size_t nbLocalIndex = globalLocalIndexGridDestSendToServer.size();
426  nbLocalIndexOnGridDest_.push_back(nbLocalIndex);
430
431  // Find out number of index sent from grid source and number of index received on grid destination
432  SourceDestinationIndexMap::const_iterator itbIndex = globaIndexWeightFromSrcToDst.begin(),
433                                            iteIndex = globaIndexWeightFromSrcToDst.end(), itIndex;
434  typedef boost::unordered_map<size_t, std::vector<std::pair<size_t,double> > > SendIndexMap;
435  std::map<int,int> sendRankSizeMap,recvRankSizeMap;
436  int connectedClient = globaIndexWeightFromSrcToDst.size();
437  int* recvCount=new int[nbClient];
438  int* displ=new int[nbClient];
439  int* sendRankBuff=new int[connectedClient];
440  int* sendSizeBuff=new int[connectedClient];
441  int n = 0;
442  for (itIndex = itbIndex; itIndex != iteIndex; ++itIndex, ++n)
443  {
444    sendRankBuff[n] = itIndex->first;
445    const SendIndexMap& sendIndexMap = itIndex->second;
446    SendIndexMap::const_iterator itbSend = sendIndexMap.begin(), iteSend = sendIndexMap.end(), itSend;
447    int sendSize = 0;
448    for (itSend = itbSend; itSend != iteSend; ++itSend)
449    {
450      sendSize += itSend->second.size();
451    }
452    sendSizeBuff[n] = sendSize;
453    sendRankSizeMap[itIndex->first] = sendSize;
454  }
455  MPI_Allgather(&connectedClient,1,MPI_INT,recvCount,1,MPI_INT,client->intraComm);
456
457  displ[0]=0 ;
458  for(int n=1;n<nbClient;n++) displ[n]=displ[n-1]+recvCount[n-1];
459  int recvSize=displ[nbClient-1]+recvCount[nbClient-1];
460  int* recvRankBuff=new int[recvSize];
461  int* recvSizeBuff=new int[recvSize];
462  MPI_Allgatherv(sendRankBuff,connectedClient,MPI_INT,recvRankBuff,recvCount,displ,MPI_INT,client->intraComm);
463  MPI_Allgatherv(sendSizeBuff,connectedClient,MPI_INT,recvSizeBuff,recvCount,displ,MPI_INT,client->intraComm);
464  for (int i = 0; i < nbClient; ++i)
465  {
466    int currentPos = displ[i];
467    for (int j = 0; j < recvCount[i]; ++j)
468      if (recvRankBuff[currentPos+j] == clientRank)
469      {
470        recvRankSizeMap[i] = recvSizeBuff[currentPos+j];
471      }
472  }
473
474  // Sending global index of grid source to corresponding process as well as the corresponding mask
475  std::vector<ep_lib::MPI_Request> requests;
476  requests.reserve(2*recvRankSizeMap.size()+2*globaIndexWeightFromSrcToDst.size());
477  std::vector<ep_lib::MPI_Status> status;
478  boost::unordered_map<int, unsigned char* > recvMaskDst;
479  boost::unordered_map<int, unsigned long* > recvGlobalIndexSrc;
480  for (std::map<int,int>::const_iterator itRecv = recvRankSizeMap.begin(); itRecv != recvRankSizeMap.end(); ++itRecv)
481  {
482    int recvRank = itRecv->first;
483    int recvSize = itRecv->second;
484    recvMaskDst[recvRank] = new unsigned char [recvSize];
485    recvGlobalIndexSrc[recvRank] = new unsigned long [recvSize];
486
487    requests.push_back(ep_lib::MPI_Request());
488    MPI_Irecv(recvGlobalIndexSrc[recvRank], recvSize, MPI_UNSIGNED_LONG, recvRank, 46, client->intraComm, &requests.back());
489    requests.push_back(ep_lib::MPI_Request());
490    MPI_Irecv(recvMaskDst[recvRank], recvSize, MPI_UNSIGNED_CHAR, recvRank, 47, client->intraComm, &requests.back());
491  }
492
493  boost::unordered_map<int, CArray<size_t,1> > globalIndexDst;
494  boost::unordered_map<int, CArray<double,1> > weightDst;
495  boost::unordered_map<int, unsigned char* > sendMaskDst;
496  boost::unordered_map<int, unsigned long* > sendGlobalIndexSrc;
497  for (itIndex = itbIndex; itIndex != iteIndex; ++itIndex)
498  {
499    int sendRank = itIndex->first;
500    int sendSize = sendRankSizeMap[sendRank];
501    const SendIndexMap& sendIndexMap = itIndex->second;
502    SendIndexMap::const_iterator itbSend = sendIndexMap.begin(), iteSend = sendIndexMap.end(), itSend;
503    globalIndexDst[sendRank].resize(sendSize);
504    weightDst[sendRank].resize(sendSize);
505    sendMaskDst[sendRank] = new unsigned char [sendSize];
506    sendGlobalIndexSrc[sendRank] = new unsigned long [sendSize];
507    int countIndex = 0;
508    for (itSend = itbSend; itSend != iteSend; ++itSend)
509    {
510      const std::vector<std::pair<size_t,double> >& dstWeight = itSend->second;
511      for (int idx = 0; idx < dstWeight.size(); ++idx)
512      {
513        globalIndexDst[sendRank](countIndex) = dstWeight[idx].first;
514        weightDst[sendRank](countIndex) = dstWeight[idx].second;
515        if (0 < globalLocalIndexGridDestSendToServer.count(dstWeight[idx].first))
517        else
519        sendGlobalIndexSrc[sendRank][countIndex] = itSend->first;
520        ++countIndex;
521      }
522    }
523
524    // Send global index source and mask
525    requests.push_back(ep_lib::MPI_Request());
526    MPI_Isend(sendGlobalIndexSrc[sendRank], sendSize, MPI_UNSIGNED_LONG, sendRank, 46, client->intraComm, &requests.back());
527    requests.push_back(ep_lib::MPI_Request());
528    MPI_Isend(sendMaskDst[sendRank], sendSize, MPI_UNSIGNED_CHAR, sendRank, 47, client->intraComm, &requests.back());
529  }
530
531  status.resize(requests.size());
532  MPI_Waitall(requests.size(), &requests[0], &status[0]);
533
534  // Okie, now use the mask to identify which index source we need to send, then also signal the destination which masked index we will return
535  std::vector<ep_lib::MPI_Request>().swap(requests);
536  std::vector<ep_lib::MPI_Status>().swap(status);
537  requests.reserve(sendRankSizeMap.size()+recvRankSizeMap.size());
538  // Okie, on destination side, we will wait for information of masked index of source
539  for (std::map<int,int>::const_iterator itSend = sendRankSizeMap.begin(); itSend != sendRankSizeMap.end(); ++itSend)
540  {
541    int recvRank = itSend->first;
542    int recvSize = itSend->second;
543
544    requests.push_back(ep_lib::MPI_Request());
545    MPI_Irecv(sendMaskDst[recvRank], recvSize, MPI_UNSIGNED_CHAR, recvRank, 48, client->intraComm, &requests.back());
546  }
547
548  // Ok, now we fill in local index of grid source (we even count for masked index)
549  CDistributionClient distributionClientSrc(client->clientRank, gridSource_);
550  CDistributionClient::GlobalLocalDataMap& globalLocalIndexGridSrcSendToServer = distributionClientSrc.getGlobalLocalDataSendToServer();
551  localIndexToSendFromGridSource_.push_back(SendingIndexGridSourceMap());
552  SendingIndexGridSourceMap& tmpSend = localIndexToSendFromGridSource_.back();
553  for (std::map<int,int>::const_iterator itRecv = recvRankSizeMap.begin(); itRecv != recvRankSizeMap.end(); ++itRecv)
554  {
555    int recvRank = itRecv->first;
556    int recvSize = itRecv->second;
558    unsigned long* recvIndexSrc = recvGlobalIndexSrc[recvRank];
559    int realSendSize = 0;
560    for (int idx = 0; idx < recvSize; ++idx)
561    {
562      if (0 != (*(recvMask+idx))) // OKie, now we have a demand from non-masked index destination
563        if (0 < globalLocalIndexGridSrcSendToServer.count(*(recvIndexSrc+idx))) // check whether index source is masked
564         ++realSendSize;
565        else // inform the destination that this index is masked
567    }
568
569    tmpSend[recvRank].resize(realSendSize);
570    realSendSize = 0;
571    for (int idx = 0; idx < recvSize; ++idx)
572    {
573      if (0 != (*(recvMask+idx))) // OKie, now we have a demand from non-masked index destination
574      {
575        tmpSend[recvRank](realSendSize) = globalLocalIndexGridSrcSendToServer[*(recvIndexSrc+idx)];
576         ++realSendSize;
577      }
578    }
579
580    // Okie, now inform the destination which source index are masked
581    requests.push_back(ep_lib::MPI_Request());
582    MPI_Isend(recvMaskDst[recvRank], recvSize, MPI_UNSIGNED_CHAR, recvRank, 48, client->intraComm, &requests.back());
583  }
584  status.resize(requests.size());
585  MPI_Waitall(requests.size(), &requests[0], &status[0]);
586
587  // Cool, now we can fill in local index of grid destination (counted for masked index)
590  for (std::map<int,int>::const_iterator itSend = sendRankSizeMap.begin(); itSend != sendRankSizeMap.end(); ++itSend)
591  {
592    int recvRank = itSend->first;
593    int recvSize = itSend->second;
595
596    CArray<size_t,1>& recvIndexDst = globalIndexDst[recvRank];
597    CArray<double,1>& recvWeightDst = weightDst[recvRank];
598    int realRecvSize = 0;
599    for (int idx = 0; idx < recvSize; ++idx)
600    {
601      if (0 != *(recvMask+idx)) // OKie, now we have a non-masked index destination
602         ++realRecvSize;
603    }
604
605    int localIndexDst;
606    recvTmp[recvRank].resize(realRecvSize);
607    realRecvSize = 0;
608    for (int idx = 0; idx < recvSize; ++idx)
609    {
610      if (0 != *(recvMask+idx)) // OKie, now we have a demand from non-masked index destination
611      {
612        recvTmp[recvRank][realRecvSize].first = globalLocalIndexGridDestSendToServer[recvIndexDst(idx)];
613        recvTmp[recvRank][realRecvSize].second = recvWeightDst(idx);
615         ++realRecvSize;
616      }
617    }
618  }
619
620  delete [] recvCount;
621  delete [] displ;
622  delete [] sendRankBuff;
623  delete [] recvRankBuff;
624  delete [] sendSizeBuff;
625  delete [] recvSizeBuff;
626
627  boost::unordered_map<int, unsigned char* >::const_iterator itChar;
629    delete [] itChar->second;
631    delete [] itChar->second;
632  boost::unordered_map<int, unsigned long* >::const_iterator itLong;
633  for (itLong = sendGlobalIndexSrc.begin(); itLong != sendGlobalIndexSrc.end(); ++itLong)
634    delete [] itLong->second;
635  for (itLong = recvGlobalIndexSrc.begin(); itLong != recvGlobalIndexSrc.end(); ++itLong)
636    delete [] itLong->second;
637
638}
639
640/*!
641  Local index of data which need sending from the grid source
642  \return local index of data
643*/
644const std::list<CGridTransformation::SendingIndexGridSourceMap>& CGridTransformation::getLocalIndexToSendFromGridSource() const
645{
646  return localIndexToSendFromGridSource_;
647}
648
649/*!
650  Local index of data which will be received on the grid destination
651  \return local index of data
652*/
654{
656}
657
658/*!
659 Number of index will be received on the grid destination
660  \return number of index of data
661*/
663{
664  return nbLocalIndexOnGridDest_;
665}
666
667/*!
668  Local mask of data which will be received on the grid destination
669  \return local mask of data
670*/