source: XIOS/trunk/src/transformation/generic_algorithm_transformation.cpp @ 888

Last change on this file since 888 was 888, checked in by mhnguyen, 6 years ago

Adding new transformation for scalar: Reducing an axis to a scalar

+) Add new xml node for new transformation
+) Add new algorithms for axis reduction
+) Make change in some place to make sure everything work fine

Test
+) On Curie
+) Tests pass and are correct

File size: 24.9 KB
Line 
1/*!
2   \file generic_algorithm_transformation.hpp
3   \author Ha NGUYEN
4   \since 14 May 2015
5   \date 21 Mars 2016
6
7   \brief Interface for all transformation algorithms.
8 */
9#include "generic_algorithm_transformation.hpp"
10#include "context.hpp"
11#include "context_client.hpp"
12#include "client_client_dht_template.hpp"
13
14namespace xios {
15
16CGenericAlgorithmTransformation::CGenericAlgorithmTransformation()
17 : transformationMapping_(), transformationWeight_(), transformationPosition_(), idAuxInputs_()
18{
19}
20
21void CGenericAlgorithmTransformation::apply(const std::vector<std::pair<int,double> >& localIndex,
22                                            const double* dataInput,
23                                            CArray<double,1>& dataOut,
24                                            std::vector<bool>& flagInitial)
25{
26  int nbLocalIndex = localIndex.size();
27  for (int idx = 0; idx < nbLocalIndex; ++idx)
28  {
29    dataOut(localIndex[idx].first) = *(dataInput + idx) * localIndex[idx].second;
30  }
31}
32
33void CGenericAlgorithmTransformation::computePositionElements(CGrid* dst, CGrid* src)
34{
35  int idxScalar = 0, idxAxis = 0, idxDomain = 0;
36  CArray<int,1> axisDomainOrderDst = dst->axis_domain_order;
37  for (int i = 0; i < axisDomainOrderDst.numElements(); ++i)
38  {
39    int dimElement = axisDomainOrderDst(i);
40    if (2 == dimElement)
41    {
42      elementPositionInGridDst2DomainPosition_[i] = idxDomain;
43      ++idxDomain;
44    }
45    else if (1 == dimElement)
46    {
47      elementPositionInGridDst2AxisPosition_[i] = idxAxis;
48      ++idxAxis;
49    }
50    else
51    {
52      elementPositionInGridDst2ScalarPosition_[i] = idxScalar;
53      ++idxScalar;
54    }
55  }
56
57  idxScalar = idxAxis = idxDomain = 0;
58  CArray<int,1> axisDomainOrderSrc = src->axis_domain_order;
59  for (int i = 0; i < axisDomainOrderSrc.numElements(); ++i)
60  {
61    int dimElement = axisDomainOrderSrc(i);
62    if (2 == dimElement)
63    {
64      elementPositionInGridSrc2DomainPosition_[i] = idxDomain;
65      ++idxDomain;
66    }
67    else if (1 == dimElement)
68    {
69      elementPositionInGridSrc2AxisPosition_[i] = idxAxis;
70      ++idxAxis;
71    }
72    else
73    {
74      elementPositionInGridSrc2ScalarPosition_[i] = idxScalar;
75      ++idxScalar;
76    }
77  }
78}
79
80/*!
81  This function computes the global indexes of grid source, which the grid destination is in demand.
82  \param[in] elementPositionInGrid position of an element in a grid .E.g: if grid is composed of domain and axis (in order),
83                then position of axis in grid is 1 and domain is positioned at 0.
84  \param[in] gridSrc Grid source
85  \param[in] gridDst Grid destination
86  \param[in\out] globaIndexWeightFromSrcToDst mapping of each global index source and weight to index destination
87*/
88void CGenericAlgorithmTransformation::computeGlobalSourceIndex(int elementPositionInGrid,
89                                                               CGrid* gridSrc,
90                                                               CGrid* gridDst,
91                                                               SourceDestinationIndexMap& globaIndexWeightFromSrcToDst)
92 {
93  CContext* context = CContext::getCurrent();
94  CContextClient* client = context->client;
95  int nbClient = client->clientSize;
96
97  typedef boost::unordered_map<int, std::vector<std::pair<int,double> > > SrcToDstMap;
98
99  size_t indexSrcSize = 0;
100  for (size_t idxTrans = 0; idxTrans < transformationMapping_.size(); ++idxTrans)
101  {
102    TransformationIndexMap::const_iterator itbTransMap = transformationMapping_[idxTrans].begin(), itTransMap,
103                                           iteTransMap = transformationMapping_[idxTrans].end();
104    TransformationWeightMap::const_iterator itbTransWeight = transformationWeight_[idxTrans].begin(), itTransWeight;
105
106    itTransWeight = itbTransWeight;
107    for (itTransMap = itbTransMap; itTransMap != iteTransMap; ++itTransMap, ++itTransWeight)
108    {
109       indexSrcSize += (itTransMap->second).size();
110    }
111  }
112
113  bool isTransPosEmpty = transformationPosition_.empty();
114  CArray<size_t,1> transPos;
115  if (!isTransPosEmpty) transPos.resize(transformationMapping_.size());
116  CArray<size_t,1> indexSrc(indexSrcSize);
117  indexSrcSize = 0;
118  for (size_t idxTrans = 0; idxTrans < transformationMapping_.size(); ++idxTrans)
119  {
120    TransformationIndexMap::const_iterator itbTransMap = transformationMapping_[idxTrans].begin(), itTransMap,
121                                           iteTransMap = transformationMapping_[idxTrans].end();
122    TransformationWeightMap::const_iterator itbTransWeight = transformationWeight_[idxTrans].begin(), itTransWeight;
123
124    // Build mapping between global source element index and global destination element index.
125    itTransWeight = itbTransWeight;
126    for (itTransMap = itbTransMap; itTransMap != iteTransMap; ++itTransMap, ++itTransWeight)
127    {
128      const std::vector<int>& srcIndex = itTransMap->second;
129      for (int idx = 0; idx < srcIndex.size(); ++idx)
130      {
131        indexSrc(indexSrcSize) = srcIndex[idx];
132        ++indexSrcSize;
133      }
134    }
135
136    if (!isTransPosEmpty)
137    {
138      TransformationPositionMap::const_iterator itPosMap = transformationPosition_[idxTrans].begin();
139      transPos(idxTrans) = itPosMap->second[0];
140    }
141  }
142
143  // compute position of elements on grids
144  computePositionElements(gridDst, gridSrc);
145  std::vector<CScalar*> scalarListDestP = gridDst->getScalars();
146  std::vector<CAxis*> axisListDestP = gridDst->getAxis();
147  std::vector<CDomain*> domainListDestP = gridDst->getDomains();
148  CArray<int,1> axisDomainDstOrder = gridDst->axis_domain_order;
149  std::vector<CScalar*> scalarListSrcP  = gridSrc->getScalars();
150  std::vector<CAxis*> axisListSrcP = gridSrc->getAxis();
151  std::vector<CDomain*> domainListSrcP = gridSrc->getDomains();
152  CArray<int,1> axisDomainSrcOrder = gridSrc->axis_domain_order;
153
154  // Find out global index source of transformed element on corresponding process.
155  std::vector<boost::unordered_map<int,std::vector<size_t> > > globalElementIndexOnProc(axisDomainDstOrder.numElements());
156  CClientClientDHTInt::Index2VectorInfoTypeMap globalIndexOfTransformedElementOnProc;
157  int axisIndex = 0, domainIndex = 0;
158  for (int idx = 0; idx < axisDomainDstOrder.numElements(); ++idx)
159  {
160    if (idx == elementPositionInGrid)
161      computeExchangeGlobalIndex(indexSrc, axisDomainSrcOrder(idx), globalIndexOfTransformedElementOnProc); //globalElementIndexOnProc[idx]);
162    if (2 == axisDomainDstOrder(idx)) // It's domain
163    {
164      if (idx != elementPositionInGrid)
165        computeExchangeDomainIndex(domainListDestP[elementPositionInGridDst2DomainPosition_[idx]],
166                                   domainListSrcP[elementPositionInGridSrc2DomainPosition_[idx]],
167                                   transPos,
168                                   globalElementIndexOnProc[idx]);
169      ++domainIndex;
170
171    }
172    else if (1 == axisDomainDstOrder(idx))//it's an axis
173    {
174      if (idx != elementPositionInGrid)
175        computeExchangeAxisIndex(axisListDestP[elementPositionInGridDst2AxisPosition_[idx]],
176                                 axisListSrcP[elementPositionInGridSrc2AxisPosition_[idx]],
177                                 transPos,
178                                 globalElementIndexOnProc[idx]);
179      ++axisIndex;
180
181    }
182    else //it's a scalar
183    {
184      if (idx != elementPositionInGrid)
185        computeExchangeScalarIndex(scalarListDestP[elementPositionInGridDst2ScalarPosition_[idx]],
186                                   scalarListSrcP[elementPositionInGridSrc2ScalarPosition_[idx]],
187                                   transPos,
188                                   globalElementIndexOnProc[idx]);
189
190    }
191  }
192
193  if (!isTransPosEmpty)
194  {
195    for (int idx = 0; idx < globalElementIndexOnProc.size(); ++idx)
196    {
197      if (idx != elementPositionInGrid)
198      {
199        boost::unordered_map<int,std::vector<size_t> >::iterator itb = globalElementIndexOnProc[idx].begin(), it,
200                                                                 ite = globalElementIndexOnProc[idx].end();
201        for (it = itb; it != ite; ++it) it->second.resize(1);
202      }
203    }
204  }
205
206  for (size_t idxTrans = 0; idxTrans < transformationMapping_.size(); ++idxTrans)
207  {
208    TransformationIndexMap::const_iterator itbTransMap = transformationMapping_[idxTrans].begin(), itTransMap,
209                                           iteTransMap = transformationMapping_[idxTrans].end();
210    TransformationWeightMap::const_iterator itbTransWeight = transformationWeight_[idxTrans].begin(), itTransWeight;
211    SrcToDstMap src2DstMap;
212    src2DstMap.rehash(std::ceil(transformationMapping_[idxTrans].size()/src2DstMap.max_load_factor()));
213
214    // Build mapping between global source element index and global destination element index.
215    boost::unordered_map<int,std::vector<size_t> >().swap(globalElementIndexOnProc[elementPositionInGrid]);
216    boost::unordered_map<int,int> tmpCounter;
217    itTransWeight = itbTransWeight;
218    for (itTransMap = itbTransMap; itTransMap != iteTransMap; ++itTransMap, ++itTransWeight)
219    {
220      const std::vector<int>& srcIndex = itTransMap->second;
221      const std::vector<double>& weight = itTransWeight->second;
222      for (int idx = 0; idx < srcIndex.size(); ++idx)
223      {
224        src2DstMap[srcIndex[idx]].push_back(make_pair(itTransMap->first, weight[idx]));
225        if (1 == globalIndexOfTransformedElementOnProc.count(srcIndex[idx]) && (0 == tmpCounter.count(srcIndex[idx])))
226        {
227          tmpCounter[srcIndex[idx]] = 1;
228          std::vector<int>& srcProc = globalIndexOfTransformedElementOnProc[srcIndex[idx]];
229          for (int j = 0; j < srcProc.size(); ++j)
230            globalElementIndexOnProc[elementPositionInGrid][srcProc[j]].push_back(srcIndex[idx]);
231        }
232      }
233    }
234
235    if (!isTransPosEmpty)
236    {
237      for (int idx = 0; idx < globalElementIndexOnProc.size(); ++idx)
238      {
239        if (idx != elementPositionInGrid)
240        {
241          boost::unordered_map<int,std::vector<size_t> >::iterator itb = globalElementIndexOnProc[idx].begin(), it,
242                                                                   ite = globalElementIndexOnProc[idx].end();
243          for (it = itb; it != ite; ++it) it->second[0] = transPos(idxTrans);
244        }
245      }
246    }
247
248    std::vector<std::vector<bool> > elementOnProc(axisDomainDstOrder.numElements(), std::vector<bool>(nbClient, false));
249    boost::unordered_map<int,std::vector<size_t> >::const_iterator it, itb, ite;
250    for (int idx = 0; idx < globalElementIndexOnProc.size(); ++idx)
251    {
252      itb = globalElementIndexOnProc[idx].begin();
253      ite = globalElementIndexOnProc[idx].end();
254      for (it = itb; it != ite; ++it) elementOnProc[idx][it->first] = true;
255    }
256
257    // Determine procs which contain global source index
258    std::vector<bool> intersectedProc(nbClient, true);
259    for (int idx = 0; idx < axisDomainDstOrder.numElements(); ++idx)
260    {
261      std::transform(elementOnProc[idx].begin(), elementOnProc[idx].end(),
262                     intersectedProc.begin(), intersectedProc.begin(),
263                     std::logical_and<bool>());
264    }
265
266    std::vector<int> srcRank;
267    for (int idx = 0; idx < nbClient; ++idx)
268    {
269      if (intersectedProc[idx]) srcRank.push_back(idx);
270    }
271
272    // Ok, now compute global index of grid source and ones of grid destination
273    computeGlobalGridIndexMapping(elementPositionInGrid,
274                                  srcRank,
275                                  src2DstMap,
276                                  gridSrc,
277                                  gridDst,
278                                  globalElementIndexOnProc,
279                                  globaIndexWeightFromSrcToDst);
280  }
281 }
282
283/*!
284  Compute mapping of global index of grid source and grid destination
285  \param [in] elementPositionInGrid position of element in grid. E.x: grid composed of domain and axis, domain has position 0 and axis 1.
286  \param [in] srcRank rank of client from which we demand global index of element source
287  \param [in] src2DstMap mapping of global index of element source and global index of element destination
288  \param[in] gridSrc Grid source
289  \param[in] gridDst Grid destination
290  \param[in] globalElementIndexOnProc Global index of element source on different client rank
291  \param[out] globaIndexWeightFromSrcToDst Mapping of global index of grid source and grid destination
292*/
293void CGenericAlgorithmTransformation::computeGlobalGridIndexMapping(int elementPositionInGrid,
294                                                                   const std::vector<int>& srcRank,
295                                                                   boost::unordered_map<int, std::vector<std::pair<int,double> > >& src2DstMap,
296                                                                   CGrid* gridSrc,
297                                                                   CGrid* gridDst,
298                                                                   std::vector<boost::unordered_map<int,std::vector<size_t> > >& globalElementIndexOnProc,
299                                                                   SourceDestinationIndexMap& globaIndexWeightFromSrcToDst)
300{
301  std::vector<CDomain*> domainListSrcP = gridSrc->getDomains();
302  std::vector<CAxis*> axisListSrcP = gridSrc->getAxis();
303  std::vector<CScalar*> scalarListSrcP = gridSrc->getScalars();
304  CArray<int,1> axisDomainSrcOrder = gridSrc->axis_domain_order;
305
306  size_t nbElement = axisDomainSrcOrder.numElements();
307  std::vector<size_t> nGlobSrc(nbElement);
308  size_t globalSrcSize = 1;
309  int domainIndex = 0, axisIndex = 0, scalarIndex = 0;
310  for (int idx = 0; idx < nbElement; ++idx)
311  {
312    nGlobSrc[idx] = globalSrcSize;
313    int elementDimension = axisDomainSrcOrder(idx);
314
315    // If this is a domain
316    if (2 == elementDimension)
317    {
318      globalSrcSize *= domainListSrcP[domainIndex]->nj_glo.getValue() * domainListSrcP[domainIndex]->ni_glo.getValue();
319      ++domainIndex;
320    }
321    else if (1 == elementDimension) // So it's an axis
322    {
323      globalSrcSize *= axisListSrcP[axisIndex]->n_glo.getValue();
324      ++axisIndex;
325    }
326    else
327    {
328      globalSrcSize *= 1;
329      ++scalarIndex;
330    }
331  }
332
333  std::vector<CDomain*> domainListDestP = gridDst->getDomains();
334  std::vector<CAxis*> axisListDestP = gridDst->getAxis();
335  std::vector<CScalar*> scalarListDestP = gridDst->getScalars();
336  CArray<int,1> axisDomainDstOrder = gridDst->axis_domain_order;
337
338  std::vector<size_t> nGlobDst(nbElement);
339  size_t globalDstSize = 1;
340  domainIndex = axisIndex = scalarIndex = 0;
341  for (int idx = 0; idx < nbElement; ++idx)
342  {
343    nGlobDst[idx] = globalDstSize;
344    int elementDimension = axisDomainDstOrder(idx);
345
346    // If this is a domain
347    if (2 == elementDimension)
348    {
349      globalDstSize *= domainListDestP[domainIndex]->nj_glo.getValue() * domainListDestP[domainIndex]->ni_glo.getValue();
350      ++domainIndex;
351    }
352    else if (1 == elementDimension) // So it's an axis
353    {
354      globalDstSize *= axisListDestP[axisIndex]->n_glo.getValue();
355      ++axisIndex;
356    }
357    else
358    {
359      globalDstSize *= 1;
360      ++scalarIndex;
361    }
362  }
363
364  for (int i = 0; i < srcRank.size(); ++i)
365  {
366    size_t ssize = 1;
367    int rankSrc = srcRank[i];
368    for (int idx = 0; idx < nbElement; ++idx)
369    {
370      ssize *= (globalElementIndexOnProc[idx][rankSrc]).size();
371    }
372
373    std::vector<int> idxLoop(nbElement,0);
374    std::vector<int> currentIndexSrc(nbElement, 0);
375    std::vector<int> currentIndexDst(nbElement, 0);
376    int innnerLoopSize = (globalElementIndexOnProc[0])[rankSrc].size();
377    size_t idx = 0;
378    while (idx < ssize)
379    {
380      for (int ind = 0; ind < nbElement; ++ind)
381      {
382        if (idxLoop[ind] == (globalElementIndexOnProc[ind])[rankSrc].size())
383        {
384          idxLoop[ind] = 0;
385          ++idxLoop[ind+1];
386        }
387
388        currentIndexDst[ind] = currentIndexSrc[ind] = (globalElementIndexOnProc[ind])[rankSrc][idxLoop[ind]];
389      }
390
391      for (int ind = 0; ind < innnerLoopSize; ++ind)
392      {
393        currentIndexDst[0] = currentIndexSrc[0] = (globalElementIndexOnProc[0])[rankSrc][ind];
394        int globalElementDstIndexSize = 0;
395        if (1 == src2DstMap.count(currentIndexSrc[elementPositionInGrid]))
396        {
397          globalElementDstIndexSize = src2DstMap[currentIndexSrc[elementPositionInGrid]].size();
398        }
399
400        std::vector<size_t> globalDstVecIndex(globalElementDstIndexSize,0);
401        size_t globalSrcIndex = 0;
402        for (int idxElement = 0; idxElement < nbElement; ++idxElement)
403        {
404          if (idxElement == elementPositionInGrid)
405          {
406            for (int k = 0; k < globalElementDstIndexSize; ++k)
407            {
408              globalDstVecIndex[k] += src2DstMap[currentIndexSrc[elementPositionInGrid]][k].first * nGlobDst[idxElement];
409            }
410          }
411          else
412          {
413            for (int k = 0; k < globalElementDstIndexSize; ++k)
414            {
415              globalDstVecIndex[k] += currentIndexDst[idxElement] * nGlobDst[idxElement];
416            }
417          }
418          globalSrcIndex += currentIndexSrc[idxElement] * nGlobSrc[idxElement];
419        }
420
421        for (int k = 0; k < globalElementDstIndexSize; ++k)
422        {
423          globaIndexWeightFromSrcToDst[rankSrc][globalSrcIndex].push_back(make_pair(globalDstVecIndex[k],src2DstMap[currentIndexSrc[elementPositionInGrid]][k].second ));
424        }
425        ++idxLoop[0];
426      }
427      idx += innnerLoopSize;
428    }
429  }
430}
431
432/*!
433  Find out proc and global index of axis source which axis destination is on demande
434  \param[in] scalar Scalar destination
435  \param[in] scalar Scalar source
436  \param[in] destGlobalIndexPositionInGrid Relative position of axis corresponds to other element of grid.
437  \param[out] globalScalarIndexOnProc Global index of axis source on different procs
438*/
439void CGenericAlgorithmTransformation::computeExchangeScalarIndex(CScalar* scalarDst,
440                                                                 CScalar* scalarSrc,
441                                                                 CArray<size_t,1>& destGlobalIndexPositionInGrid,
442                                                                 boost::unordered_map<int,std::vector<size_t> >& globalScalarIndexOnProc)
443{
444  CContext* context = CContext::getCurrent();
445  CContextClient* client=context->client;
446  int clientRank = client->clientRank;
447  int clientSize = client->clientSize;
448
449  globalScalarIndexOnProc.rehash(std::ceil(clientSize/globalScalarIndexOnProc.max_load_factor()));
450  for (int idx = 0; idx < clientSize; ++idx)
451  {
452    globalScalarIndexOnProc[idx].push_back(0);
453  }
454}
455
456/*!
457  Find out proc and global index of axis source which axis destination is on demande
458  \param[in] axisDst Axis destination
459  \param[in] axisSrc Axis source
460  \param[in] destGlobalIndexPositionInGrid Relative position of axis corresponds to other element of grid.
461  \param[out] globalAxisIndexOnProc Global index of axis source on different procs
462*/
463void CGenericAlgorithmTransformation::computeExchangeAxisIndex(CAxis* axisDst,
464                                                               CAxis* axisSrc,
465                                                               CArray<size_t,1>& destGlobalIndexPositionInGrid,
466                                                               boost::unordered_map<int,std::vector<size_t> >& globalAxisIndexOnProc)
467{
468  CContext* context = CContext::getCurrent();
469  CContextClient* client=context->client;
470  int clientRank = client->clientRank;
471  int clientSize = client->clientSize;
472
473  size_t globalIndex;
474  int nIndexSize = axisSrc->index.numElements();
475  CClientClientDHTInt::Index2VectorInfoTypeMap globalIndex2ProcRank;
476  globalIndex2ProcRank.rehash(std::ceil(nIndexSize/globalIndex2ProcRank.max_load_factor()));
477  for (int idx = 0; idx < nIndexSize; ++idx)
478  {
479    globalIndex = axisSrc->index(idx);
480    globalIndex2ProcRank[globalIndex].push_back(clientRank);
481  }
482
483  CClientClientDHTInt dhtIndexProcRank(globalIndex2ProcRank, client->intraComm);
484  CArray<size_t,1> globalAxisIndex(axisDst->index.numElements());
485  for (int idx = 0; idx < globalAxisIndex.numElements(); ++idx)
486  {
487    globalAxisIndex(idx) = axisDst->index(idx);
488  }
489  dhtIndexProcRank.computeIndexInfoMapping(globalAxisIndex);
490
491  std::vector<int> countIndex(clientSize,0);
492  const CClientClientDHTInt::Index2VectorInfoTypeMap& computedGlobalIndexOnProc = dhtIndexProcRank.getInfoIndexMap();
493  CClientClientDHTInt::Index2VectorInfoTypeMap::const_iterator itb = computedGlobalIndexOnProc.begin(), it,
494                                                               ite = computedGlobalIndexOnProc.end();
495  for (it = itb; it != ite; ++it)
496  {
497    const std::vector<int>& procList = it->second;
498    for (int idx = 0; idx < procList.size(); ++idx) ++countIndex[procList[idx]];
499  }
500
501  globalAxisIndexOnProc.rehash(std::ceil(clientSize/globalAxisIndexOnProc.max_load_factor()));
502  for (int idx = 0; idx < clientSize; ++idx)
503  {
504    if (0 != countIndex[idx])
505    {
506      globalAxisIndexOnProc[idx].resize(countIndex[idx]);
507      countIndex[idx] = 0;
508    }
509  }
510
511  for (it = itb; it != ite; ++it)
512  {
513    const std::vector<int>& procList = it->second;
514    for (int idx = 0; idx < procList.size(); ++idx)
515    {
516      globalAxisIndexOnProc[procList[idx]][countIndex[procList[idx]]] = it->first;
517      ++countIndex[procList[idx]];
518    }
519  }
520}
521
522/*!
523  Find out proc and global index of domain source which domain destination is on demande
524  \param[in] domainDst Domain destination
525  \param[in] domainSrc Domain source
526  \param[in] destGlobalIndexPositionInGrid Relative position of domain corresponds to other element of grid.
527  \param[out] globalDomainIndexOnProc Global index of domain source on different procs
528*/
529void CGenericAlgorithmTransformation::computeExchangeDomainIndex(CDomain* domainDst,
530                                                                 CDomain* domainSrc,
531                                                                 CArray<size_t,1>& destGlobalIndexPositionInGrid,
532                                                                 boost::unordered_map<int,std::vector<size_t> >& globalDomainIndexOnProc)
533{
534  CContext* context = CContext::getCurrent();
535  CContextClient* client=context->client;
536  int clientRank = client->clientRank;
537  int clientSize = client->clientSize;
538
539  int niGlobSrc = domainSrc->ni_glo.getValue();
540  size_t globalIndex;
541  int i_ind, j_ind;
542  int nIndexSize = (destGlobalIndexPositionInGrid.isEmpty()) ? domainSrc->i_index.numElements()
543                                                             : destGlobalIndexPositionInGrid.numElements();
544  CClientClientDHTInt::Index2VectorInfoTypeMap globalIndex2ProcRank;
545  globalIndex2ProcRank.rehash(std::ceil(nIndexSize/globalIndex2ProcRank.max_load_factor()));
546  if (destGlobalIndexPositionInGrid.isEmpty())
547  {
548    for (int idx = 0; idx < nIndexSize; ++idx)
549    {
550      i_ind=domainSrc->i_index(idx) ;
551      j_ind=domainSrc->j_index(idx) ;
552
553      globalIndex = i_ind + j_ind * niGlobSrc;
554      globalIndex2ProcRank[globalIndex].resize(1);
555      globalIndex2ProcRank[globalIndex][0] = clientRank;
556    }
557  }
558  else
559  {
560    for (int idx = 0; idx < nIndexSize; ++idx)
561    {
562      globalIndex2ProcRank[destGlobalIndexPositionInGrid(idx)].push_back(clientRank);
563    }
564  }
565
566  CArray<size_t,1> globalDomainIndex;
567  if (destGlobalIndexPositionInGrid.isEmpty())
568  {
569    int niGlobDst = domainDst->ni_glo.getValue();
570    globalDomainIndex.resize(domainDst->i_index.numElements());
571    nIndexSize = domainDst->i_index.numElements();
572
573    for (int idx = 0; idx < nIndexSize; ++idx)
574    {
575      i_ind=domainDst->i_index(idx) ;
576      j_ind=domainDst->j_index(idx) ;
577
578      globalDomainIndex(idx) = i_ind + j_ind * niGlobDst;
579    }
580  }
581  else
582  {
583    globalDomainIndex.reference(destGlobalIndexPositionInGrid);
584  }
585
586  CClientClientDHTInt dhtIndexProcRank(globalIndex2ProcRank, client->intraComm);
587  dhtIndexProcRank.computeIndexInfoMapping(globalDomainIndex);
588
589  std::vector<int> countIndex(clientSize,0);
590  const CClientClientDHTInt::Index2VectorInfoTypeMap& computedGlobalIndexOnProc = dhtIndexProcRank.getInfoIndexMap();
591  CClientClientDHTInt::Index2VectorInfoTypeMap::const_iterator itb = computedGlobalIndexOnProc.begin(), it,
592                                                               ite = computedGlobalIndexOnProc.end();
593  for (it = itb; it != ite; ++it)
594  {
595    const std::vector<int>& procList = it->second;
596    for (int idx = 0; idx < procList.size(); ++idx) ++countIndex[procList[idx]];
597  }
598
599  globalDomainIndexOnProc.rehash(std::ceil(clientSize/globalDomainIndexOnProc.max_load_factor()));
600  for (int idx = 0; idx < clientSize; ++idx)
601  {
602    if (0 != countIndex[idx])
603    {
604      globalDomainIndexOnProc[idx].resize(countIndex[idx]);
605      countIndex[idx] = 0;
606    }
607  }
608
609  for (it = itb; it != ite; ++it)
610  {
611    const std::vector<int>& procList = it->second;
612    for (int idx = 0; idx < procList.size(); ++idx)
613    {
614      globalDomainIndexOnProc[procList[idx]][countIndex[procList[idx]]] = it->first;
615      ++countIndex[procList[idx]];
616    }
617  }
618}
619
620/*!
621  Compute index mapping between element source and element destination with an auxiliary inputs which determine
622position of each mapped index in global index of grid destination.
623  \param [in] dataAuxInputs auxiliary inputs
624*/
625void CGenericAlgorithmTransformation::computeIndexSourceMapping(const std::vector<CArray<double,1>* >& dataAuxInputs)
626{
627  computeIndexSourceMapping_(dataAuxInputs);
628}
629
630std::vector<StdString> CGenericAlgorithmTransformation::getIdAuxInputs()
631{
632  return idAuxInputs_;
633}
634
635}
Note: See TracBrowser for help on using the repository browser.