source: XIOS/dev/dev_ym/XIOS_COUPLING/src/transformation/generic_algorithm_transformation.cpp @ 1984

Last change on this file since 1984 was 1984, checked in by ymipsl, 3 years ago

intermediate commit for new tranformation engine?
YM
File size: 43.6 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#include "utils.hpp"
14#include "timer.hpp"
15#include "mpi.hpp"
16#include "transform_connector.hpp"
17#include "weight_transform_connector.hpp"
18
19namespace xios
20{
21
22CGenericAlgorithmTransformation::CGenericAlgorithmTransformation()
23 : transformationMapping_(), transformationWeight_(), transformationPosition_(),
24   idAuxInputs_(), type_(ELEMENT_NO_MODIFICATION_WITH_DATA), indexElementSrc_(),
25   computedProcSrcNonTransformedElement_(false), eliminateRedondantSrc_(true), isDistributedComputed_(false)
26{
27}
28
29void CGenericAlgorithmTransformation::updateData(CArray<double,1>& dataOut)
30{
31
32}
33
34void CGenericAlgorithmTransformation::apply(const std::vector<std::pair<int,double> >& localIndex,
35                                            const double* dataInput,
36                                            CArray<double,1>& dataOut,
37                                            std::vector<bool>& flagInitial,
38                                            bool ignoreMissingValue, bool firstPass  )
39TRY
40{
41  int nbLocalIndex = localIndex.size();   
42  double defaultValue = std::numeric_limits<double>::quiet_NaN();
43   
44  if (ignoreMissingValue)
45  {
46    if (firstPass) dataOut=defaultValue ;
47   
48    for (int idx = 0; idx < nbLocalIndex; ++idx)
49    {
50      if (! NumTraits<double>::isNan(*(dataInput + idx)))
51      {
52        if (flagInitial[localIndex[idx].first]) dataOut(localIndex[idx].first) = *(dataInput + idx) * localIndex[idx].second;
53        else dataOut(localIndex[idx].first) += *(dataInput + idx) * localIndex[idx].second;
54        flagInitial[localIndex[idx].first] = false; // Reset flag to indicate not all data source are nan
55      }
56    }
57
58  }
59  else
60  {
61    for (int idx = 0; idx < nbLocalIndex; ++idx)
62    {
63      dataOut(localIndex[idx].first) += *(dataInput + idx) * localIndex[idx].second;
64    }
65  }
66}
67CATCH
68
69void CGenericAlgorithmTransformation::computePositionElements(CGrid* dst, CGrid* src)
70TRY
71{
72  int idxScalar = 0, idxAxis = 0, idxDomain = 0;
73  CArray<int,1> axisDomainOrderDst = dst->axis_domain_order;
74  for (int i = 0; i < axisDomainOrderDst.numElements(); ++i)
75  {
76    int dimElement = axisDomainOrderDst(i);
77    if (2 == dimElement)
78    {
79      elementPositionInGridDst2DomainPosition_[i] = idxDomain;
80      ++idxDomain;
81    }
82    else if (1 == dimElement)
83    {
84      elementPositionInGridDst2AxisPosition_[i] = idxAxis;
85      ++idxAxis;
86    }
87    else
88    {
89      elementPositionInGridDst2ScalarPosition_[i] = idxScalar;
90      ++idxScalar;
91    }
92  }
93
94  idxScalar = idxAxis = idxDomain = 0;
95  CArray<int,1> axisDomainOrderSrc = src->axis_domain_order;
96  for (int i = 0; i < axisDomainOrderSrc.numElements(); ++i)
97  {
98    int dimElement = axisDomainOrderSrc(i);
99    if (2 == dimElement)
100    {
101      elementPositionInGridSrc2DomainPosition_[i] = idxDomain;
102      ++idxDomain;
103    }
104    else if (1 == dimElement)
105    {
106      elementPositionInGridSrc2AxisPosition_[i] = idxAxis;
107      ++idxAxis;
108    }
109    else
110    {
111      elementPositionInGridSrc2ScalarPosition_[i] = idxScalar;
112      ++idxScalar;
113    }
114  }
115}
116CATCH
117
118bool CGenericAlgorithmTransformation::isDistributedTransformation(int elementPositionInGrid, CGrid* gridSrc, CGrid* gridDst)
119TRY
120{
121
122  if (!isDistributedComputed_)
123  {
124    isDistributedComputed_=true ;
125    if (!eliminateRedondantSrc_) isDistributed_=true ;
126    else
127    {
128      CContext* context = CContext::getCurrent();
129       
130      computePositionElements(gridSrc, gridDst);
131      std::vector<CScalar*> scalarListSrcP  = gridSrc->getScalars();
132      std::vector<CAxis*> axisListSrcP = gridSrc->getAxis();
133      std::vector<CDomain*> domainListSrcP = gridSrc->getDomains();
134      int distributed, distributed_glo ;
135 
136      CArray<int,1> axisDomainSrcOrder = gridSrc->axis_domain_order;
137      if (2 == axisDomainSrcOrder(elementPositionInGrid)) // It's domain
138      {
139        distributed=domainListSrcP[elementPositionInGridSrc2DomainPosition_[elementPositionInGrid]]->isDistributed() ;
140        MPI_Allreduce(&distributed,&distributed_glo, 1, MPI_INT, MPI_LOR, context->intraComm_) ;
141   
142      }
143      else if (1 == axisDomainSrcOrder(elementPositionInGrid))//it's an axis
144      {
145        distributed=axisListSrcP[elementPositionInGridSrc2AxisPosition_[elementPositionInGrid]]->isDistributed() ;
146        MPI_Allreduce(&distributed,&distributed_glo, 1, MPI_INT, MPI_LOR, context->intraComm_) ;
147      }
148      else //it's a scalar
149      {
150        distributed_glo=false ;
151      } 
152      isDistributed_=distributed_glo ;
153    }
154  }
155  return isDistributed_ ;
156}
157CATCH
158
159/*!
160  This function computes the global indexes of grid source, which the grid destination is in demand.
161  \param[in] elementPositionInGrid position of an element in a grid .E.g: if grid is composed of domain and axis (in order),
162                then position of axis in grid is 1 and domain is positioned at 0.
163  \param[in] gridSrc Grid source
164  \param[in] gridDst Grid destination
165  \param[in\out] globaIndexWeightFromSrcToDst mapping of each global index source and weight to index destination
166*/
167void CGenericAlgorithmTransformation::computeGlobalSourceIndex(int elementPositionInGrid,
168                                                               CGrid* gridSrc,
169                                                               CGrid* gridDst,
170                                                               SourceDestinationIndexMap& globaIndexWeightFromSrcToDst)
171TRY
172 {
173  CContext* context = CContext::getCurrent();
174  int nbClient = context->intraCommSize_;
175
176  typedef std::unordered_map<int, std::vector<std::pair<int,double> > > SrcToDstMap;
177  int idx;
178
179  // compute position of elements on grids
180  computePositionElements(gridDst, gridSrc);
181  std::vector<CScalar*> scalarListDestP = gridDst->getScalars();
182  std::vector<CAxis*> axisListDestP = gridDst->getAxis();
183  std::vector<CDomain*> domainListDestP = gridDst->getDomains();
184  CArray<int,1> axisDomainDstOrder = gridDst->axis_domain_order;
185  std::vector<CScalar*> scalarListSrcP  = gridSrc->getScalars();
186  std::vector<CAxis*> axisListSrcP = gridSrc->getAxis();
187  std::vector<CDomain*> domainListSrcP = gridSrc->getDomains();
188  CArray<int,1> axisDomainSrcOrder = gridSrc->axis_domain_order; 
189 
190  bool isTransPosEmpty = transformationPosition_.empty();
191  CArray<size_t,1> transPos;
192  if (!isTransPosEmpty) transPos.resize(transformationMapping_.size());
193  std::set<size_t> allIndexSrc; // All index of source, which can be scattered among processes, need for doing transformation
194 
195  for (size_t idxTrans = 0; idxTrans < transformationMapping_.size(); ++idxTrans)
196  {
197    TransformationIndexMap::const_iterator itbTransMap = transformationMapping_[idxTrans].begin(), itTransMap,
198                                           iteTransMap = transformationMapping_[idxTrans].end();
199    TransformationWeightMap::const_iterator itbTransWeight = transformationWeight_[idxTrans].begin(), itTransWeight;
200
201    // Build mapping between global source element index and global destination element index.
202    itTransWeight = itbTransWeight;
203    for (itTransMap = itbTransMap; itTransMap != iteTransMap; ++itTransMap, ++itTransWeight)
204    {
205      const std::vector<int>& srcIndex = itTransMap->second;
206      for (idx = 0; idx < srcIndex.size(); ++idx)
207        allIndexSrc.insert(srcIndex[idx]);
208    }
209
210    if (!isTransPosEmpty)
211    {
212      TransformationPositionMap::const_iterator itPosMap = transformationPosition_[idxTrans].begin();
213      transPos(idxTrans) = itPosMap->second[0];
214    }
215  }
216
217  size_t indexSrcSize = 0;
218  CArray<size_t,1> indexSrc(allIndexSrc.size());
219  for (std::set<size_t>::iterator it = allIndexSrc.begin(); it != allIndexSrc.end(); ++it, ++indexSrcSize)
220    indexSrc(indexSrcSize) = *it;
221
222  // Flag to indicate whether we will recompute distribution of source global index  on processes
223  bool computeGlobalIndexOnProc = false; 
224  if (indexElementSrc_.size() != allIndexSrc.size())
225    computeGlobalIndexOnProc = true; 
226  else
227  {
228    for (std::set<size_t>::iterator it = allIndexSrc.begin(); it != allIndexSrc.end(); ++it)
229      if (0 == indexElementSrc_.count(*it))
230      {
231        computeGlobalIndexOnProc = true;
232        break;       
233      }
234  }
235
236  if (computeGlobalIndexOnProc)
237    indexElementSrc_.swap(allIndexSrc);
238     
239  int sendValue = (computeGlobalIndexOnProc) ? 1 : 0;
240  int recvValue = 0;
241  MPI_Allreduce(&sendValue, &recvValue, 1, MPI_INT, MPI_SUM, context->intraComm_);
242  computeGlobalIndexOnProc = (0 < recvValue);
243
244//  CClientClientDHTInt::Index2VectorInfoTypeMap globalIndexOfTransformedElementOnProc;
245
246  if (computeGlobalIndexOnProc || !computedProcSrcNonTransformedElement_)
247  {   
248    {
249      CClientClientDHTInt::Index2VectorInfoTypeMap tmp ;
250      globalIndexOfTransformedElementOnProc_.swap(tmp) ;
251    }
252    // Find out global index source of transformed element on corresponding process.   
253    if (globalElementIndexOnProc_.empty())
254      globalElementIndexOnProc_.resize(axisDomainDstOrder.numElements());
255   
256    for (idx = 0; idx < axisDomainDstOrder.numElements(); ++idx)
257    {
258     
259      if (idx == elementPositionInGrid)
260        computeExchangeGlobalIndex(indexSrc, axisDomainSrcOrder(idx), globalIndexOfTransformedElementOnProc_); //globalElementIndexOnProc[idx]);
261      if (!computedProcSrcNonTransformedElement_)
262      {
263        if (2 == axisDomainDstOrder(idx)) // It's domain
264        {
265          if (idx != elementPositionInGrid)
266            computeExchangeDomainIndex(domainListDestP[elementPositionInGridDst2DomainPosition_[idx]],
267                                       domainListSrcP[elementPositionInGridSrc2DomainPosition_[idx]],
268                                       transPos,
269                                       globalElementIndexOnProc_[idx]);     
270
271        }
272        else if (1 == axisDomainDstOrder(idx))//it's an axis
273        {
274          if (idx != elementPositionInGrid)
275            computeExchangeAxisIndex(axisListDestP[elementPositionInGridDst2AxisPosition_[idx]],
276                                     axisListSrcP[elementPositionInGridSrc2AxisPosition_[idx]],
277                                     transPos,
278                                     globalElementIndexOnProc_[idx]);
279        }
280        else //it's a scalar
281        {
282          if (idx != elementPositionInGrid)
283            computeExchangeScalarIndex(scalarListDestP[elementPositionInGridDst2ScalarPosition_[idx]],
284                                       scalarListSrcP[elementPositionInGridSrc2ScalarPosition_[idx]],
285                                       transPos,
286                                       globalElementIndexOnProc_[idx]);
287
288        }
289      }
290    }
291
292    if (!isTransPosEmpty && !computedProcSrcNonTransformedElement_)
293    {
294      for (idx = 0; idx < globalElementIndexOnProc_.size(); ++idx)
295      {
296        if (idx != elementPositionInGrid)
297        {
298          std::unordered_map<int,std::vector<size_t> >::iterator itb = globalElementIndexOnProc_[idx].begin(), it,
299                                                                   ite = globalElementIndexOnProc_[idx].end();
300          for (it = itb; it != ite; ++it) it->second.resize(1);
301        }
302      }
303    }
304
305/*     
306    if (!computedProcSrcNonTransformedElement_)
307    {
308      for (idx = 0; idx < globalElementIndexOnProc_.size(); ++idx)
309      {
310        if (idx != elementPositionInGrid)
311        {
312          std::unordered_map<int,std::vector<size_t> >::iterator itb = globalElementIndexOnProc_[idx].begin(), it,
313                                                                   ite = globalElementIndexOnProc_[idx].end();
314          for (it = itb; it != ite; ++it) procOfNonTransformedElements_.insert(it->first);
315          if (procOfNonTransformedElements_.size() == nbClient)
316            break;
317        }
318      }
319    }
320   
321    // Processes contain the source index of transformed element
322    std::set<int> procOfTransformedElement;
323    CClientClientDHTInt::Index2VectorInfoTypeMap::iterator itIdxb = globalIndexOfTransformedElementOnProc.begin(),
324                                                           itIdxe = globalIndexOfTransformedElementOnProc.end(), itIdx;
325    for (itIdx = itIdxb; itIdx != itIdxe; ++itIdx)
326    {
327      std::vector<int>& tmp = itIdx->second;
328      for (int i = 0; i < tmp.size(); ++i)
329        procOfTransformedElement.insert(tmp[i]);
330      if (tmp.size() == nbClient)
331        break;
332    }                                                           
333   
334    std::set<int>& commonProc = (procOfTransformedElement.size() < procOfNonTransformedElements_.size()) ? procOfTransformedElement
335                              : (!procOfNonTransformedElements_.empty() ? procOfNonTransformedElements_ : procOfTransformedElement);
336   
337    std::vector<int> procContainSrcElementIdx(commonProc.size());
338    int count = 0;
339    for (std::set<int>::iterator it = commonProc.begin(); it != commonProc.end(); ++it)
340      procContainSrcElementIdx[count++] = *it;
341
342    procContainSrcElementIdx_.swap(procContainSrcElementIdx);   
343*/
344   
345    if (procElementList_.empty()) procElementList_.resize(axisDomainDstOrder.numElements()) ;
346    for (idx = 0; idx < axisDomainDstOrder.numElements(); ++idx)
347    {
348      std::set<int>& procList=procElementList_[idx] ;
349      std::set<int> commonTmp ;
350      if (idx == elementPositionInGrid)
351      {
352          set<int> tmpSet ; 
353          procList.swap(tmpSet) ;
354          CClientClientDHTInt::Index2VectorInfoTypeMap::iterator itIdxb = globalIndexOfTransformedElementOnProc_.begin(),
355                                                                 itIdxe = globalIndexOfTransformedElementOnProc_.end(), itIdx;
356          for (itIdx = itIdxb; itIdx != itIdxe; ++itIdx)
357          {
358             std::vector<int>& tmp = itIdx->second;
359             for (int i = 0; i < tmp.size(); ++i) procList.insert(tmp[i]);
360             if (tmp.size() == nbClient)
361             break;
362          }
363      }
364      else
365      {
366        if (!computedProcSrcNonTransformedElement_)
367        {
368          set<int> tmpSet ; 
369          procList.swap(tmpSet) ;
370          std::unordered_map<int,std::vector<size_t> >::iterator itb = globalElementIndexOnProc_[idx].begin(), it,
371                                                                   ite = globalElementIndexOnProc_[idx].end();
372          for (it = itb; it != ite; ++it)
373          {
374            procList.insert(it->first);
375            if (procList.size() == nbClient)  break;
376          }
377        }
378      }
379
380      if (idx==0) commonProc_= procList ;
381      else
382      {
383        for (std::set<int>::iterator it = commonProc_.begin(); it != commonProc_.end(); ++it)
384          if (procList.count(*it)==1) commonTmp.insert(*it) ;
385        commonProc_.swap(commonTmp) ;
386      }
387    }
388    std::vector<int> procContainSrcElementIdx(commonProc_.size());
389    int count = 0;
390    for (std::set<int>::iterator it = commonProc_.begin(); it != commonProc_.end(); ++it) procContainSrcElementIdx[count++] = *it;
391    procContainSrcElementIdx_.swap(procContainSrcElementIdx);
392   
393        // For the first time, surely we calculate proc containing non transformed elements
394    if (!computedProcSrcNonTransformedElement_) computedProcSrcNonTransformedElement_ = true;
395  }
396 
397  for (size_t idxTrans = 0; idxTrans < transformationMapping_.size(); ++idxTrans)
398  {
399    TransformationIndexMap::const_iterator itbTransMap = transformationMapping_[idxTrans].begin(), itTransMap,
400                                           iteTransMap = transformationMapping_[idxTrans].end();
401    TransformationWeightMap::const_iterator itbTransWeight = transformationWeight_[idxTrans].begin(), itTransWeight;
402    SrcToDstMap src2DstMap;
403    src2DstMap.rehash(std::ceil(transformationMapping_[idxTrans].size()/src2DstMap.max_load_factor()));
404
405    // Build mapping between global source element index and global destination element index.
406    std::unordered_map<int,std::vector<size_t> >().swap(globalElementIndexOnProc_[elementPositionInGrid]);
407    std::set<int> tmpCounter;
408    itTransWeight = itbTransWeight;
409    for (itTransMap = itbTransMap; itTransMap != iteTransMap; ++itTransMap, ++itTransWeight)
410    {
411      const std::vector<int>& srcIndex = itTransMap->second;
412      const std::vector<double>& weight = itTransWeight->second;
413      for (idx = 0; idx < srcIndex.size(); ++idx)
414      {
415        src2DstMap[srcIndex[idx]].push_back(make_pair(itTransMap->first, weight[idx]));
416        if (0 == tmpCounter.count(srcIndex[idx]))
417        {         
418          tmpCounter.insert(srcIndex[idx]);
419       
420          vector<int>& rankSrc = globalIndexOfTransformedElementOnProc_[srcIndex[idx]] ;
421          for (int n=0;n<rankSrc.size();++n)
422          {
423            if (commonProc_.count(rankSrc[n])==1) globalElementIndexOnProc_[elementPositionInGrid][rankSrc[n]].push_back(srcIndex[idx]);
424          }
425//          for (int j = 0; j < procContainSrcElementIdx_.size(); ++j)
426//            globalElementIndexOnProc_[elementPositionInGrid][procContainSrcElementIdx_[j]].push_back(srcIndex[idx]);
427        }
428      }
429    }
430 
431    if (!isTransPosEmpty)
432    {
433      for (idx = 0; idx < globalElementIndexOnProc_.size(); ++idx)
434      {
435        if (idx != elementPositionInGrid)
436        {
437          std::unordered_map<int,std::vector<size_t> >::iterator itb = globalElementIndexOnProc_[idx].begin(), it,
438                                                                   ite = globalElementIndexOnProc_[idx].end();
439          for (it = itb; it != ite; ++it) it->second[0] = transPos(idxTrans);
440        }
441      }
442    }
443
444    // Ok, now compute global index of grid source and ones of grid destination
445    computeGlobalGridIndexMapping(elementPositionInGrid,
446                                  procContainSrcElementIdx_, //srcRank,
447                                  src2DstMap,
448                                  gridSrc,
449                                  gridDst,
450                                  globalElementIndexOnProc_,
451                                  globaIndexWeightFromSrcToDst);
452  } 
453 }
454CATCH
455
456/*!
457  Compute mapping of global index of grid source and grid destination
458  \param [in] elementPositionInGrid position of element in grid. E.x: grid composed of domain and axis, domain has position 0 and axis 1.
459  \param [in] srcRank rank of client from which we demand global index of element source
460  \param [in] src2DstMap mapping of global index of element source and global index of element destination
461  \param [in] gridSrc Grid source
462  \param [in] gridDst Grid destination
463  \param [in] globalElementIndexOnProc Global index of element source on different client rank
464  \param [out] globaIndexWeightFromSrcToDst Mapping of global index of grid source and grid destination
465*/
466void CGenericAlgorithmTransformation::computeGlobalGridIndexMapping(int elementPositionInGrid,
467                                                                   const std::vector<int>& srcRank,
468                                                                   std::unordered_map<int, std::vector<std::pair<int,double> > >& src2DstMap,
469                                                                   CGrid* gridSrc,
470                                                                   CGrid* gridDst,
471                                                                   std::vector<std::unordered_map<int,std::vector<size_t> > >& globalElementIndexOnProc,
472                                                                   SourceDestinationIndexMap& globaIndexWeightFromSrcToDst)
473TRY
474{
475  SourceDestinationIndexMap globaIndexWeightFromSrcToDst_tmp ;
476 
477  CContext* context = CContext::getCurrent();
478  int clientRank = context->intraCommRank_;
479 
480  std::vector<CDomain*> domainListSrcP = gridSrc->getDomains();
481  std::vector<CAxis*> axisListSrcP = gridSrc->getAxis();
482  std::vector<CScalar*> scalarListSrcP = gridSrc->getScalars();
483  CArray<int,1> axisDomainSrcOrder = gridSrc->axis_domain_order;
484
485  size_t nbElement = axisDomainSrcOrder.numElements();
486  std::vector<size_t> nGlobSrc(nbElement);
487  size_t globalSrcSize = 1;
488  int domainIndex = 0, axisIndex = 0, scalarIndex = 0;
489  for (int idx = 0; idx < nbElement; ++idx)
490  {
491    nGlobSrc[idx] = globalSrcSize;
492    int elementDimension = axisDomainSrcOrder(idx);
493
494    // If this is a domain
495    if (2 == elementDimension)
496    {
497      globalSrcSize *= domainListSrcP[domainIndex]->nj_glo.getValue() * domainListSrcP[domainIndex]->ni_glo.getValue();
498      ++domainIndex;
499    }
500    else if (1 == elementDimension) // So it's an axis
501    {
502      globalSrcSize *= axisListSrcP[axisIndex]->n_glo.getValue();
503      ++axisIndex;
504    }
505    else
506    {
507      globalSrcSize *= 1;
508      ++scalarIndex;
509    }
510  }
511
512  std::vector<CDomain*> domainListDestP = gridDst->getDomains();
513  std::vector<CAxis*> axisListDestP = gridDst->getAxis();
514  std::vector<CScalar*> scalarListDestP = gridDst->getScalars();
515  CArray<int,1> axisDomainDstOrder = gridDst->axis_domain_order;
516
517  std::vector<size_t> nGlobDst(nbElement);
518  size_t globalDstSize = 1;
519  domainIndex = axisIndex = scalarIndex = 0;
520  set<size_t>  globalSrcStoreIndex ;
521 
522  for (int idx = 0; idx < nbElement; ++idx)
523  {
524    nGlobDst[idx] = globalDstSize;
525    int elementDimension = axisDomainDstOrder(idx);
526
527    // If this is a domain
528    if (2 == elementDimension)
529    {
530      globalDstSize *= domainListDestP[domainIndex]->nj_glo.getValue() * domainListDestP[domainIndex]->ni_glo.getValue();
531      ++domainIndex;
532    }
533    else if (1 == elementDimension) // So it's an axis
534    {
535      globalDstSize *= axisListDestP[axisIndex]->n_glo.getValue();
536      ++axisIndex;
537    }
538    else
539    {
540      globalDstSize *= 1;
541      ++scalarIndex;
542    }
543  }
544
545  std::map< std::pair<size_t,size_t>, int > rankMap ;
546  std::map< std::pair<size_t,size_t>, int >:: iterator rankMapIt ;
547 
548  for (int i = 0; i < srcRank.size(); ++i)
549  {
550    size_t ssize = 1;
551    int rankSrc = srcRank[i];
552    for (int idx = 0; idx < nbElement; ++idx)
553    {
554      ssize *= (globalElementIndexOnProc[idx][rankSrc]).size();
555    }
556
557    std::vector<int> idxLoop(nbElement,0);
558    std::vector<int> currentIndexSrc(nbElement, 0);
559    std::vector<int> currentIndexDst(nbElement, 0);
560    int innnerLoopSize = (globalElementIndexOnProc[0])[rankSrc].size();
561    size_t idx = 0;
562    while (idx < ssize)
563    {
564      for (int ind = 0; ind < nbElement; ++ind)
565      {
566        if (idxLoop[ind] == (globalElementIndexOnProc[ind])[rankSrc].size())
567        {
568          idxLoop[ind] = 0;
569          ++idxLoop[ind+1];
570        }
571
572        currentIndexDst[ind] = currentIndexSrc[ind] = (globalElementIndexOnProc[ind])[rankSrc][idxLoop[ind]];
573      }
574
575      for (int ind = 0; ind < innnerLoopSize; ++ind)
576      {
577        currentIndexDst[0] = currentIndexSrc[0] = (globalElementIndexOnProc[0])[rankSrc][ind];
578        int globalElementDstIndexSize = 0;
579        if (1 == src2DstMap.count(currentIndexSrc[elementPositionInGrid]))
580        {
581          globalElementDstIndexSize = src2DstMap[currentIndexSrc[elementPositionInGrid]].size();
582        }
583
584        std::vector<size_t> globalDstVecIndex(globalElementDstIndexSize,0);
585        size_t globalSrcIndex = 0;
586        for (int idxElement = 0; idxElement < nbElement; ++idxElement)
587        {
588          if (idxElement == elementPositionInGrid)
589          {
590            for (int k = 0; k < globalElementDstIndexSize; ++k)
591            {
592              globalDstVecIndex[k] += src2DstMap[currentIndexSrc[elementPositionInGrid]][k].first * nGlobDst[idxElement];
593            }
594          }
595          else
596          {
597            for (int k = 0; k < globalElementDstIndexSize; ++k)
598            {
599              globalDstVecIndex[k] += currentIndexDst[idxElement] * nGlobDst[idxElement];
600            }
601          }
602          globalSrcIndex += currentIndexSrc[idxElement] * nGlobSrc[idxElement];
603        }
604
605        for (int k = 0; k < globalElementDstIndexSize; ++k)
606        {
607         
608          globaIndexWeightFromSrcToDst_tmp[rankSrc][globalSrcIndex].push_back(make_pair(globalDstVecIndex[k],src2DstMap[currentIndexSrc[elementPositionInGrid]][k].second ));
609          rankMapIt=rankMap.find(make_pair(globalSrcIndex,globalDstVecIndex[k])) ;
610          if (rankMapIt==rankMap.end()) rankMap[make_pair(globalSrcIndex,globalDstVecIndex[k])] = rankSrc ;
611          else if (rankSrc==clientRank) rankMapIt->second = rankSrc ;
612        }
613        ++idxLoop[0];
614      }
615      idx += innnerLoopSize;
616    }
617  }
618
619  // eliminate redondant global src point owned by differrent processes.
620  // Avoid as possible to tranfer data from an other process if the src point is also owned by current process
621      int rankSrc ;
622      size_t globalSrcIndex ;
623      size_t globalDstIndex ;
624      double weight ;
625 
626      SourceDestinationIndexMap::iterator rankIt,rankIte ;
627      std::unordered_map<size_t, std::vector<std::pair<size_t,double> > >::iterator globalSrcIndexIt, globalSrcIndexIte ;
628      std::vector<std::pair<size_t,double> >::iterator vectIt,vectIte ;
629   
630      rankIt=globaIndexWeightFromSrcToDst_tmp.begin() ; rankIte=globaIndexWeightFromSrcToDst_tmp.end() ;
631      for(;rankIt!=rankIte;++rankIt)
632      {
633        rankSrc = rankIt->first ;
634        globalSrcIndexIt = rankIt->second.begin() ; globalSrcIndexIte = rankIt->second.end() ;
635        for(;globalSrcIndexIt!=globalSrcIndexIte;++globalSrcIndexIt)
636        { 
637          globalSrcIndex = globalSrcIndexIt->first ;
638          vectIt = globalSrcIndexIt->second.begin() ; vectIte = globalSrcIndexIt->second.end() ;
639          for(vectIt; vectIt!=vectIte; vectIt++)
640          {
641            globalDstIndex = vectIt->first ;
642            weight = vectIt->second ;
643            if (eliminateRedondantSrc_)
644            {
645              if (rankMap[make_pair(globalSrcIndex,globalDstIndex)] == rankSrc) 
646                globaIndexWeightFromSrcToDst[rankSrc][globalSrcIndex].push_back(make_pair(globalDstIndex,weight)) ;
647            }
648            else globaIndexWeightFromSrcToDst[rankSrc][globalSrcIndex].push_back(make_pair(globalDstIndex,weight)) ;
649         }
650       }
651     }
652}
653CATCH
654
655/*!
656  Find out proc and global index of axis source which axis destination is on demande
657  \param[in] scalar Scalar destination
658  \param[in] scalar Scalar source
659  \param[in] destGlobalIndexPositionInGrid Relative position of axis corresponds to other element of grid.
660  \param[out] globalScalarIndexOnProc Global index of axis source on different procs
661*/
662void CGenericAlgorithmTransformation::computeExchangeScalarIndex(CScalar* scalarDst,
663                                                                 CScalar* scalarSrc,
664                                                                 CArray<size_t,1>& destGlobalIndexPositionInGrid,
665                                                                 std::unordered_map<int,std::vector<size_t> >& globalScalarIndexOnProc)
666TRY
667{
668  CContext* context = CContext::getCurrent();
669  int clientRank = context->intraCommRank_;
670  int clientSize = context->intraCommSize_;
671
672  globalScalarIndexOnProc.rehash(std::ceil(clientSize/globalScalarIndexOnProc.max_load_factor()));
673  for (int idx = 0; idx < clientSize; ++idx)
674  {
675    globalScalarIndexOnProc[idx].push_back(0);
676  }
677}
678CATCH
679
680/*!
681  Find out proc and global index of axis source which axis destination is on demande
682  \param[in] axisDst Axis destination
683  \param[in] axisSrc Axis source
684  \param[in] destGlobalIndexPositionInGrid Relative position of axis corresponds to other element of grid.
685  \param[out] globalAxisIndexOnProc Global index of axis source on different procs
686*/
687void CGenericAlgorithmTransformation::computeExchangeAxisIndex(CAxis* axisDst,
688                                                               CAxis* axisSrc,
689                                                               CArray<size_t,1>& destGlobalIndexPositionInGrid,
690                                                               std::unordered_map<int,std::vector<size_t> >& globalAxisIndexOnProc)
691TRY
692{
693  CContext* context = CContext::getCurrent();
694  int clientRank = context->intraCommRank_;
695  int clientSize = context->intraCommSize_;
696
697  size_t globalIndex;
698  int nIndexSize = axisSrc->index.numElements();
699  CClientClientDHTInt::Index2VectorInfoTypeMap globalIndex2ProcRank;
700  globalIndex2ProcRank.rehash(std::ceil(nIndexSize/globalIndex2ProcRank.max_load_factor()));
701  for (int idx = 0; idx < nIndexSize; ++idx)
702  {
703    if (axisSrc->mask(idx))
704    {
705      globalIndex = axisSrc->index(idx);
706      globalIndex2ProcRank[globalIndex].push_back(clientRank);
707    }
708  }
709
710  CClientClientDHTInt dhtIndexProcRank(globalIndex2ProcRank, context->intraComm_);
711  CArray<size_t,1> globalAxisIndex(axisDst->index.numElements());
712  for (int idx = 0; idx < globalAxisIndex.numElements(); ++idx)
713  {
714    globalAxisIndex(idx) = axisDst->index(idx);
715  }
716  dhtIndexProcRank.computeIndexInfoMapping(globalAxisIndex);
717
718  std::vector<int> countIndex(clientSize,0);
719  const CClientClientDHTInt::Index2VectorInfoTypeMap& computedGlobalIndexOnProc = dhtIndexProcRank.getInfoIndexMap();
720  CClientClientDHTInt::Index2VectorInfoTypeMap::const_iterator itb = computedGlobalIndexOnProc.begin(), it,
721                                                               ite = computedGlobalIndexOnProc.end();
722  for (it = itb; it != ite; ++it)
723  {
724    const std::vector<int>& procList = it->second;
725    for (int idx = 0; idx < procList.size(); ++idx) ++countIndex[procList[idx]];
726  }
727
728  globalAxisIndexOnProc.rehash(std::ceil(clientSize/globalAxisIndexOnProc.max_load_factor()));
729  for (int idx = 0; idx < clientSize; ++idx)
730  {
731    if (0 != countIndex[idx])
732    {
733      globalAxisIndexOnProc[idx].resize(countIndex[idx]);
734      countIndex[idx] = 0;
735    }
736  }
737
738  for (it = itb; it != ite; ++it)
739  {
740    const std::vector<int>& procList = it->second;
741    for (int idx = 0; idx < procList.size(); ++idx)
742    {
743      globalAxisIndexOnProc[procList[idx]][countIndex[procList[idx]]] = it->first;
744      ++countIndex[procList[idx]];
745    }
746  }
747}
748CATCH
749
750/*!
751  Find out proc and global index of domain source which domain destination is on demande
752  \param[in] domainDst Domain destination
753  \param[in] domainSrc Domain source
754  \param[in] destGlobalIndexPositionInGrid Relative position of domain corresponds to other element of grid.
755  \param[out] globalDomainIndexOnProc Global index of domain source on different procs
756*/
757void CGenericAlgorithmTransformation::computeExchangeDomainIndex(CDomain* domainDst,
758                                                                 CDomain* domainSrc,
759                                                                 CArray<size_t,1>& destGlobalIndexPositionInGrid,
760                                                                 std::unordered_map<int,std::vector<size_t> >& globalDomainIndexOnProc)
761TRY
762{
763  CContext* context = CContext::getCurrent();
764  int clientRank = context->intraCommRank_;
765  int clientSize = context->intraCommSize_;
766
767  int niGlobSrc = domainSrc->ni_glo.getValue();
768  size_t globalIndex;
769  int i_ind, j_ind;
770  int nIndexSize = (destGlobalIndexPositionInGrid.isEmpty()) ? domainSrc->i_index.numElements()
771                                                             : destGlobalIndexPositionInGrid.numElements();
772  CClientClientDHTInt::Index2VectorInfoTypeMap globalIndex2ProcRank;
773  globalIndex2ProcRank.rehash(std::ceil(nIndexSize/globalIndex2ProcRank.max_load_factor()));
774 
775  if (destGlobalIndexPositionInGrid.isEmpty())
776  {
777    for (int idx = 0; idx < nIndexSize; ++idx)
778    {
779      i_ind=domainSrc->i_index(idx) ;
780      j_ind=domainSrc->j_index(idx) ;
781
782      if (domainSrc->localMask(idx))
783      {
784        globalIndex = i_ind + j_ind * niGlobSrc;
785        globalIndex2ProcRank[globalIndex].resize(1);
786        globalIndex2ProcRank[globalIndex][0] = clientRank;
787      }
788    }
789  }
790  else
791  {
792    for (int idx = 0; idx < nIndexSize; ++idx)
793    {
794//      if (domainSrc->localMask(idx)) -> not necessairy, mask seem to be included in  destGlobalIndexPositionInGrid(idx)    (ym)
795        globalIndex2ProcRank[destGlobalIndexPositionInGrid(idx)].push_back(clientRank);
796    }
797  }
798
799  CArray<size_t,1> globalDomainIndex;
800  if (destGlobalIndexPositionInGrid.isEmpty())
801  {
802    int niGlobDst = domainDst->ni_glo.getValue();
803    globalDomainIndex.resize(domainDst->i_index.numElements());
804    nIndexSize = domainDst->i_index.numElements();
805
806    for (int idx = 0; idx < nIndexSize; ++idx)
807    {
808      i_ind=domainDst->i_index(idx) ;
809      j_ind=domainDst->j_index(idx) ;
810
811      globalDomainIndex(idx) = i_ind + j_ind * niGlobDst;
812    }
813  }
814  else
815  {
816    globalDomainIndex.reference(destGlobalIndexPositionInGrid);
817  }
818
819  CClientClientDHTInt dhtIndexProcRank(globalIndex2ProcRank, context->intraComm_);
820  dhtIndexProcRank.computeIndexInfoMapping(globalDomainIndex);
821
822  std::vector<int> countIndex(clientSize,0);
823  const CClientClientDHTInt::Index2VectorInfoTypeMap& computedGlobalIndexOnProc = dhtIndexProcRank.getInfoIndexMap();
824  CClientClientDHTInt::Index2VectorInfoTypeMap::const_iterator itb = computedGlobalIndexOnProc.begin(), it,
825                                                               ite = computedGlobalIndexOnProc.end();
826  for (it = itb; it != ite; ++it)
827  {
828    const std::vector<int>& procList = it->second;
829    for (int idx = 0; idx < procList.size(); ++idx) ++countIndex[procList[idx]];
830  }
831
832  globalDomainIndexOnProc.rehash(std::ceil(clientSize/globalDomainIndexOnProc.max_load_factor()));
833  for (int idx = 0; idx < clientSize; ++idx)
834  {
835    if (0 != countIndex[idx])
836    {
837      globalDomainIndexOnProc[idx].resize(countIndex[idx]);
838      countIndex[idx] = 0;
839    }
840  }
841
842  for (it = itb; it != ite; ++it)
843  {
844    const std::vector<int>& procList = it->second;
845    for (int idx = 0; idx < procList.size(); ++idx)
846    {
847      globalDomainIndexOnProc[procList[idx]][countIndex[procList[idx]]] = it->first;
848      ++countIndex[procList[idx]];
849    }
850  }
851}
852CATCH
853
854void CGenericAlgorithmTransformation::computeTransformationMappingNonDistributed(int elementPositionInGrid, CGrid* gridSrc, CGrid* gridDst,
855                                                                                 vector<int>& localSrc, vector<int>& localDst, vector<double>& weight,
856                                                                                 int& nlocalIndexDest)
857TRY
858{
859
860  CContext* context = CContext::getCurrent();
861  int nbClient = context->intraCommSize_;
862
863  computePositionElements(gridDst, gridSrc);
864  std::vector<CScalar*> scalarListDstP = gridDst->getScalars();
865  std::vector<CAxis*> axisListDstP = gridDst->getAxis();
866  std::vector<CDomain*> domainListDstP = gridDst->getDomains();
867  CArray<int,1> axisDomainDstOrder = gridDst->axis_domain_order;
868  std::vector<CScalar*> scalarListSrcP  = gridSrc->getScalars();
869  std::vector<CAxis*> axisListSrcP = gridSrc->getAxis();
870  std::vector<CDomain*> domainListSrcP = gridSrc->getDomains();
871  CArray<int,1> axisDomainSrcOrder = gridSrc->axis_domain_order; 
872
873  int nElement=axisDomainSrcOrder.numElements() ;
874  int indSrc=1 ;
875  int indDst=1 ;
876  vector<int> nIndexSrc(nElement) ;
877  vector<int> nIndexDst(nElement) ;
878  vector< CArray<bool,1>* > maskSrc(nElement) ;
879  vector< CArray<bool,1>* > maskDst(nElement) ;
880   
881  int nlocalIndexSrc=1 ;
882//  int nlocalIndexDest=1 ;
883  nlocalIndexDest=1 ;
884  CArray<bool,1> maskScalar(1) ;
885  maskScalar  = true ;
886 
887 
888  for(int i=0 ; i<nElement; i++)
889  {
890    int dimElement = axisDomainSrcOrder(i);
891    if (2 == dimElement) //domain
892    {
893      CDomain* domain=domainListSrcP[elementPositionInGridSrc2DomainPosition_[i]] ;
894      nIndexSrc[i] = domain->i_index.numElements() ;
895      maskSrc[i]=&domain->localMask ;
896    }
897    else if (1 == dimElement) //axis
898    {
899      CAxis* axis=axisListSrcP[elementPositionInGridSrc2AxisPosition_[i]] ;
900      nIndexSrc[i] = axis->index.numElements() ;
901      maskSrc[i]=&axis->mask ;
902    }
903    else  //scalar
904    {
905      nIndexSrc[i]=1 ;
906      maskSrc[i]=&maskScalar ;
907    }
908    nlocalIndexSrc=nlocalIndexSrc*nIndexSrc[i] ;
909  }
910
911
912
913  int offset=1 ;
914  for(int i=0 ; i<nElement; i++)
915  {
916    int dimElement = axisDomainDstOrder(i);
917    if (2 == dimElement) //domain
918    {
919      CDomain* domain=domainListDstP[elementPositionInGridDst2DomainPosition_[i]] ;
920      int nIndex=domain->i_index.numElements() ;
921      CArray<bool,1>& localMask=domain->localMask ;
922      int nbInd=0 ;
923      for(int j=0;j<nIndex;j++) if (localMask(j)) nbInd++ ;
924      nIndexDst[i] = nbInd ;
925      maskDst[i]=&domain->localMask ;
926    }
927    else if (1 == dimElement) //axis
928    {
929      CAxis* axis = axisListDstP[elementPositionInGridDst2AxisPosition_[i]] ;
930      int nIndex=axis->index.numElements() ;
931      CArray<bool,1>& localMask=axis->mask ;
932      int nbInd=0 ;
933      for(int j=0;j<nIndex;j++) if (localMask(j)) nbInd++ ;
934      nIndexDst[i] = nbInd ;
935      maskDst[i]=&axis->mask ;
936    }
937    else  //scalar
938    {
939      nIndexDst[i]=1 ;
940      maskDst[i]=&maskScalar ;
941    }
942    if (i<elementPositionInGrid) offset=offset*nIndexDst[i] ;
943    nlocalIndexDest=nlocalIndexDest*nIndexDst[i] ;
944  }
945
946  vector<int> dstLocalInd ;
947  int dimElement = axisDomainDstOrder(elementPositionInGrid);
948  if (2 == dimElement) //domain
949  {
950    CDomain* domain = domainListDstP[elementPositionInGridDst2DomainPosition_[elementPositionInGrid]] ;
951    int ni_glo=domain->ni_glo ;
952    int nj_glo=domain->nj_glo ;
953    int nindex_glo=ni_glo*nj_glo ;
954    dstLocalInd.resize(nindex_glo,-1) ;
955    int nIndex=domain->i_index.numElements() ;
956    CArray<bool,1>& localMask=domain->localMask ;
957    int unmaskedInd=0 ;
958    int globIndex ;
959    for(int i=0;i<nIndex;i++)
960    {
961      if (localMask(i))
962      {
963        globIndex=domain->j_index(i)*ni_glo+domain->i_index(i) ;
964        dstLocalInd[globIndex]=unmaskedInd ;
965        unmaskedInd++ ;
966      }
967    }
968  }
969  else if (1 == dimElement) //axis
970  {
971    CAxis* axis = axisListDstP[elementPositionInGridDst2AxisPosition_[elementPositionInGrid]] ;
972    int nindex_glo=axis->n_glo ;
973    dstLocalInd.resize(nindex_glo,-1) ;
974    int nIndex=axis->index.numElements() ;
975    CArray<bool,1>& localMask=axis->mask ; // axis mask must include later data_index
976    int unmaskedInd=0 ;
977    for(int i=0;i<nIndex;i++)
978    {
979      if (localMask(i))
980      {
981        dstLocalInd[axis->index(i)]=unmaskedInd ;
982        unmaskedInd++ ;
983      }
984    }
985  }
986  else  //scalar
987  {
988    dstLocalInd.resize(1) ; 
989    dstLocalInd[0]=0 ; 
990  }
991
992  vector<vector<vector<pair<int,double> > > > dstIndWeight(transformationMapping_.size()) ;
993   
994  for(int t=0;t<transformationMapping_.size();++t)
995  {
996    TransformationIndexMap::const_iterator   itTransMap = transformationMapping_[t].begin(),
997                                             iteTransMap = transformationMapping_[t].end();
998    TransformationWeightMap::const_iterator itTransWeight = transformationWeight_[t].begin();
999    dstIndWeight[t].resize(nIndexSrc[elementPositionInGrid]) ;
1000   
1001    for(;itTransMap!=iteTransMap;++itTransMap,++itTransWeight)
1002    {
1003      int dst=dstLocalInd[itTransMap->first] ;
1004      if (dst!=-1)
1005      {
1006        const vector<int>& srcs=itTransMap->second;
1007        const vector<double>& weights=itTransWeight->second;
1008        for(int i=0;i<srcs.size() ;i++) dstIndWeight[t][srcs[i]].push_back(pair<int,double>(dst*offset+t,weights[i])) ;
1009      }
1010    }
1011  }
1012  int srcInd=0 ; 
1013  int currentInd ;
1014  int t=0 ; 
1015  int srcIndCompressed=0 ;
1016 
1017  nonDistributedrecursiveFunct(nElement-1,true,elementPositionInGrid,maskSrc,maskDst, srcInd, srcIndCompressed, nIndexSrc, t, dstIndWeight, 
1018                               currentInd,localSrc,localDst,weight);
1019               
1020}
1021CATCH
1022
1023
1024void CGenericAlgorithmTransformation::nonDistributedrecursiveFunct(int currentPos, bool masked, int elementPositionInGrid,
1025                                                                   vector< CArray<bool,1>* >& maskSrc, vector< CArray<bool,1>* >& maskDst,
1026                                                                   int& srcInd, int& srcIndCompressed, vector<int>& nIndexSrc,
1027                                                                   int& t, vector<vector<vector<pair<int,double> > > >& dstIndWeight, int currentInd,
1028                                                                   vector<int>& localSrc, vector<int>& localDst, vector<double>& weight)
1029TRY
1030{
1031  int masked_ ;
1032  if (currentPos!=elementPositionInGrid)
1033  {
1034    if (currentPos!=0)
1035    {
1036      CArray<bool,1>& mask = *maskSrc[currentPos] ;
1037     
1038      for(int i=0;i<nIndexSrc[currentPos];i++)
1039      {
1040        masked_=masked ;
1041        if (!mask(i)) masked_=false ;
1042        nonDistributedrecursiveFunct(currentPos-1, masked_, elementPositionInGrid, maskSrc, maskDst, srcInd, srcIndCompressed, nIndexSrc, t,
1043                                     dstIndWeight, currentInd, localSrc, localDst, weight);
1044      }
1045    }
1046    else
1047    {
1048      CArray<bool,1>& mask = *maskSrc[currentPos] ;
1049      for(int i=0;i<nIndexSrc[currentPos];i++)
1050      {
1051        if (masked && mask(i))
1052        {
1053          if (dstIndWeight[t][currentInd].size()>0)
1054          {
1055            for(vector<pair<int,double> >::iterator it = dstIndWeight[t][currentInd].begin(); it!=dstIndWeight[t][currentInd].end(); ++it)
1056            {
1057              localSrc.push_back(srcIndCompressed) ;
1058              localDst.push_back(it->first) ;
1059              weight.push_back(it->second) ;
1060              (it->first)++ ;
1061            }
1062          }
1063          if (t < dstIndWeight.size()-1) t++ ;
1064            srcIndCompressed ++ ;
1065        }
1066        srcInd++ ;
1067      }
1068    }
1069  }
1070  else
1071  {
1072 
1073    if (currentPos!=0)
1074    {
1075
1076      CArray<bool,1>& mask = *maskSrc[currentPos] ;
1077      for(int i=0;i<nIndexSrc[currentPos];i++)
1078      {
1079        t=0 ;
1080        masked_=masked ;
1081        if (!mask(i)) masked_=false ; 
1082        nonDistributedrecursiveFunct(currentPos-1, masked_, elementPositionInGrid, maskSrc, maskDst, srcInd,
1083                                     srcIndCompressed, nIndexSrc, t, dstIndWeight , i,  localSrc, localDst, weight);
1084      }
1085    }
1086    else
1087    {
1088      for(int i=0;i<nIndexSrc[currentPos];i++)
1089      {
1090        if (masked)
1091        {
1092          t=0 ;       
1093          if (dstIndWeight[t][i].size()>0)
1094          {
1095            for(vector<pair<int,double> >::iterator it = dstIndWeight[t][i].begin(); it!=dstIndWeight[t][i].end(); ++it)
1096            {
1097              localSrc.push_back(srcIndCompressed) ;
1098              localDst.push_back(it->first) ;
1099              weight.push_back(it->second) ;
1100              (it->first)++ ;
1101            }
1102           }
1103          if (t < dstIndWeight.size()-1) t++ ;
1104          srcIndCompressed ++ ;
1105        }
1106        srcInd++ ;
1107      }
1108    }
1109  }
1110
1111}
1112CATCH
1113
1114/*!
1115  Compute index mapping between element source and element destination with an auxiliary inputs which determine
1116position of each mapped index in global index of grid destination.
1117  \param [in] dataAuxInputs auxiliary inputs
1118*/
1119void CGenericAlgorithmTransformation::computeIndexSourceMapping(const std::vector<CArray<double,1>* >& dataAuxInputs)
1120TRY
1121{
1122  computeIndexSourceMapping_(dataAuxInputs);
1123}
1124CATCH
1125
1126std::vector<StdString> CGenericAlgorithmTransformation::getIdAuxInputs()
1127TRY
1128{
1129  return idAuxInputs_;
1130}
1131CATCH
1132
1133CGenericAlgorithmTransformation::AlgoTransType CGenericAlgorithmTransformation::type()
1134TRY
1135{
1136  return type_;
1137}
1138CATCH
1139
1140
1141///////////////////////////////////////////////////////////////
1142////////// new algorithm for new method               /////////
1143///////////////////////////////////////////////////////////////
1144
1145
1146
1147void CGenericAlgorithmTransformation::computeAlgorithm(CLocalView* srcView, CLocalView* dstView)
1148{
1149  auto& srcMap = transformationMapping_[0] ;
1150  set<size_t> srcIndex ;
1151  for(auto& it : srcMap)
1152    for(size_t index : it.second) srcIndex.insert(index) ;
1153
1154  CArray<size_t,1> srcArrayIndex(srcIndex.size()) ;
1155  int i=0 ;
1156  for(size_t index : srcIndex) { srcArrayIndex(i) = index ; i++ ;}
1157  CLocalElement recvElement(CContext::getCurrent()->getIntraCommRank(), srcView->getGlobalSize(), srcArrayIndex) ;
1158  recvElement.addFullView() ;
1159
1160  transformConnector_ = new CTransformConnector(srcView, recvElement.getView(CElementView::FULL), CContext::getCurrent()->getIntraComm())  ;
1161  transformConnector_->computeConnector() ;
1162  weightTransformConnector_ = new  CWeightTransformConnector( recvElement.getView(CElementView::FULL), dstView, transformationMapping_[0], transformationWeight_[0]) ; 
1163}
1164 
1165
1166void CGenericAlgorithmTransformation::apply(int dimBefore, int dimAfter, const CArray<double,1>& dataIn, CArray<double,1>& dataOut)
1167{
1168  CArray<double,1> dataOutTmp ;
1169  transformConnector_->transfer(dimBefore, dimAfter, dataIn, dataOutTmp) ;
1170  weightTransformConnector_ -> transfer(dimBefore, dimAfter, dataOutTmp, dataOut) ;
1171}
1172
1173
1174
1175
1176
1177
1178}
Note: See TracBrowser for help on using the repository browser.