source: XIOS/trunk/src/transformation/axis_algorithm_inverse.cpp @ 829

Last change on this file since 829 was 829, checked in by mhnguyen, 5 years ago

Refactoring transformation code

+) On exchanging information during transformation, not only global index are sent but also local index
+) Correct a bug in distributed hash table (dht)
+) Add new type for dht
+) Clean up some redundant codes

Test
+) On Curie
+) Every test passes
+) Code runs faster in some cases (up to 30%)

File size: 5.8 KB
Line 
1/*!
2   \file axis_algorithm_inverse.hpp
3   \author Ha NGUYEN
4   \since 14 May 2015
5   \date 29 June 2015
6
7   \brief Algorithm for inversing an axis..
8 */
9#include "axis_algorithm_inverse.hpp"
10#include "transformation_mapping.hpp"
11#include "context.hpp"
12#include "context_client.hpp"
13
14namespace xios {
15
16CAxisAlgorithmInverse::CAxisAlgorithmInverse(CAxis* axisDestination, CAxis* axisSource)
17 : CAxisAlgorithmTransformation(axisDestination, axisSource)
18{
19  if (axisDestination->n_glo.getValue() != axisSource->n_glo.getValue())
20  {
21    ERROR("CAxisAlgorithmInverse::CAxisAlgorithmInverse(CAxis* axisDestination, CAxis* axisSource)",
22           << "Two axis have different global size"
23           << "Size of axis source " <<axisSource->getId() << " is " << axisSource->n_glo.getValue()  << std::endl
24           << "Size of axis destionation " <<axisDestination->getId() << " is " << axisDestination->n_glo.getValue());
25  }
26
27//  this->computeIndexSourceMapping();
28}
29
30void CAxisAlgorithmInverse::computeIndexSourceMapping_(const std::vector<CArray<double,1>* >& dataAuxInputs)
31{
32  this->transformationMapping_.resize(1);
33  this->transformationWeight_.resize(1);
34
35  std::map<int, std::vector<int> >& transMap = this->transformationMapping_[0];
36  std::map<int, std::vector<double> >& transWeight = this->transformationWeight_[0];
37
38  int globalIndexSize = axisDestGlobalIndex_.size();
39  for (int idx = 0; idx < globalIndexSize; ++idx)
40  {
41    transMap[axisDestGlobalIndex_[idx]].push_back(axisDestGlobalSize_-axisDestGlobalIndex_[idx]-1);
42    transWeight[axisDestGlobalIndex_[idx]].push_back(1.0);
43  }
44
45  int niSrc   = axisSrc_->n.getValue();
46  int sizeSrc = axisSrc_->n_glo.getValue();
47  if (niSrc != sizeSrc) updateAxisValue();
48  else
49  {
50    for (int idx = 0; idx < sizeSrc; ++idx)
51    {
52      axisDest_->value(idx) = axisSrc_->value(sizeSrc-idx-1);
53    }
54  }
55}
56
57/*!
58  Update value on axis after inversing
59  After an axis is inversed, not only the data on it must be inversed but also the value
60*/
61void CAxisAlgorithmInverse::updateAxisValue()
62{
63  CContext* context = CContext::getCurrent();
64  CContextClient* client=context->client;
65
66  int niSrc     = axisSrc_->n.getValue();
67  int ibeginSrc = axisSrc_->begin.getValue();
68
69  CTransformationMapping transformationMap(axisDest_, axisSrc_);
70
71  std::map<int, std::vector<int> >& transMap = this->transformationMapping_[0];
72  std::map<int, std::vector<double> >& transWeight = this->transformationWeight_[0];
73
74  CTransformationMapping::DestinationIndexMap globaIndexMapFromDestToSource;
75  std::map<int, std::vector<int> >::const_iterator it = transMap.begin(), ite = transMap.end();
76  int localIndex = 0;
77  for (; it != ite; ++it)
78  {
79    globaIndexMapFromDestToSource[it->first].push_back(make_pair(localIndex,make_pair((it->second)[0], (transWeight[it->first])[0])));
80    ++localIndex;
81  }
82
83  transformationMap.computeTransformationMapping(globaIndexMapFromDestToSource);
84
85  const CTransformationMapping::ReceivedIndexMap& globalIndexToReceive = transformationMap.getGlobalIndexReceivedOnGridDestMapping();
86  const CTransformationMapping::SentIndexMap& globalIndexToSend = transformationMap.getGlobalIndexSendToGridDestMapping();
87
88 // Sending global index of original grid source
89  CTransformationMapping::SentIndexMap::const_iterator itbSend = globalIndexToSend.begin(), itSend,
90                                                       iteSend = globalIndexToSend.end();
91 int sendBuffSize = 0;
92 for (itSend = itbSend; itSend != iteSend; ++itSend) sendBuffSize += (itSend->second).size();
93
94 typedef double Scalar;
95 Scalar* sendBuff, *currentSendBuff;
96 if (0 != sendBuffSize) sendBuff = new Scalar[sendBuffSize];
97 for (StdSize idx = 0; idx < sendBuffSize; ++idx) sendBuff[idx] = NumTraits<Scalar>::sfmax();
98
99 int currentBuffPosition = 0;
100 for (itSend = itbSend; itSend != iteSend; ++itSend)
101 {
102   int destRank = itSend->first;
103   const std::vector<std::pair<int,size_t> >& globalIndexOfCurrentGridSourceToSend = itSend->second;
104   int countSize = globalIndexOfCurrentGridSourceToSend.size();
105   for (int idx = 0; idx < (countSize); ++idx)
106   {
107     int index = globalIndexOfCurrentGridSourceToSend[idx].first;
108     sendBuff[idx+currentBuffPosition] = (axisSrc_->value)(index);
109   }
110   currentSendBuff = sendBuff + currentBuffPosition;
111   MPI_Send(currentSendBuff, countSize, MPI_DOUBLE, destRank, 14, client->intraComm);
112   currentBuffPosition += countSize;
113 }
114
115 // Receiving global index of grid source sending from current grid source
116 CTransformationMapping::ReceivedIndexMap::const_iterator itbRecv = globalIndexToReceive.begin(), itRecv,
117                                                          iteRecv = globalIndexToReceive.end();
118 int recvBuffSize = 0;
119 for (itRecv = itbRecv; itRecv != iteRecv; ++itRecv) recvBuffSize += (itRecv->second).size();
120
121 Scalar* recvBuff, *currentRecvBuff;
122 if (0 != recvBuffSize) recvBuff = new Scalar [recvBuffSize];
123 for (StdSize idx = 0; idx < recvBuffSize; ++idx) recvBuff[idx] = NumTraits<Scalar>::sfmax();
124
125 int currentRecvBuffPosition = 0;
126 for (itRecv = itbRecv; itRecv != iteRecv; ++itRecv)
127 {
128   MPI_Status status;
129   int srcRank = itRecv->first;
130   int countSize = (itRecv->second).size();
131   currentRecvBuff = recvBuff + currentRecvBuffPosition;
132   MPI_Recv(currentRecvBuff, countSize, MPI_DOUBLE, srcRank, 14, client->intraComm, &status);
133   currentRecvBuffPosition += countSize;
134 }
135
136 int ibeginDest = axisDest_->begin.getValue();
137 currentRecvBuff = recvBuff;
138 for (itRecv = itbRecv; itRecv != iteRecv; ++itRecv)
139 {
140   int countSize = (itRecv->second).size();
141   for (int idx = 0; idx < countSize; ++idx, ++currentRecvBuff)
142   {
143     int ssize = (itRecv->second)[idx].size();
144     for (int i = 0; i < ssize; ++i)
145     {
146       int index = ((itRecv->second)[idx][i]).first - ibeginDest;
147       (axisDest_->value)(index) = *currentRecvBuff;
148     }
149   }
150 }
151
152 if (0 != sendBuffSize) delete [] sendBuff;
153 if (0 != recvBuffSize) delete [] recvBuff;
154}
155
156}
Note: See TracBrowser for help on using the repository browser.