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

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

Change name of several axis attributes and remove some redundant variable of domain

+) Change name of axis attributes to make them consistent with ones of domain
+) Remove zoom_client_* of domain

Test
+) On Curie
+) All tests pass and are correct

File size: 5.7 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  int niSrc   = axisSrc_->n.getValue();
29  int sizeSrc = axisSrc_->n_glo.getValue();
30  if (niSrc != sizeSrc) updateAxisValue();
31  else
32  {
33    for (int idx = 0; idx < sizeSrc; ++idx)
34    {
35      axisDest_->value(idx) = axisSrc_->value(sizeSrc-idx-1);
36    }
37  }
38}
39
40void CAxisAlgorithmInverse::computeIndexSourceMapping()
41{
42  std::map<int, std::vector<int> >& transMap = this->transformationMapping_;
43  std::map<int, std::vector<double> >& transWeight = this->transformationWeight_;
44
45  int globalIndexSize = axisDestGlobalIndex_.size();
46  for (int idx = 0; idx < globalIndexSize; ++idx)
47  {
48    transMap[axisDestGlobalIndex_[idx]].push_back(axisDestGlobalSize_-axisDestGlobalIndex_[idx]-1);
49    transWeight[axisDestGlobalIndex_[idx]].push_back(1.0);
50  }
51}
52
53/*!
54  Update value on axis after inversing
55  After an axis is inversed, not only the data on it must be inversed but also the value
56*/
57void CAxisAlgorithmInverse::updateAxisValue()
58{
59  CContext* context = CContext::getCurrent();
60  CContextClient* client=context->client;
61
62  int niSrc     = axisSrc_->n.getValue();
63  int ibeginSrc = axisSrc_->begin.getValue();
64
65  CTransformationMapping transformationMap(axisDest_, axisSrc_);
66
67  std::map<int, std::vector<int> >& transMap = this->transformationMapping_;
68  std::map<int, std::vector<double> >& transWeight = this->transformationWeight_;
69
70  std::map<size_t, std::vector<std::pair<size_t,double> > > globaIndexMapFromDestToSource;
71  std::map<int, std::vector<int> >::const_iterator it = transMap.begin(), ite = transMap.end();
72  for (; it != ite; ++it)
73  {
74    globaIndexMapFromDestToSource[it->first].push_back(make_pair((it->second)[0], (transWeight[it->first])[0]));
75  }
76
77  transformationMap.computeTransformationMapping(globaIndexMapFromDestToSource);
78
79  const std::map<int,std::vector<std::vector<std::pair<size_t,double> > > >& globalIndexToReceive = transformationMap.getGlobalIndexReceivedOnGridDestMapping();
80  const std::map<int,std::vector<size_t> >& globalIndexToSend = transformationMap.getGlobalIndexSendToGridDestMapping();
81
82 // Sending global index of original grid source
83  std::map<int,std::vector<size_t> >::const_iterator itbSend = globalIndexToSend.begin(), itSend,
84                                                     iteSend = globalIndexToSend.end();
85 int sendBuffSize = 0;
86 for (itSend = itbSend; itSend != iteSend; ++itSend) sendBuffSize += (itSend->second).size();
87
88 typedef double Scalar;
89 Scalar* sendBuff, *currentSendBuff;
90 if (0 != sendBuffSize) sendBuff = new Scalar[sendBuffSize];
91 for (StdSize idx = 0; idx < sendBuffSize; ++idx) sendBuff[idx] = NumTraits<Scalar>::sfmax();
92
93 int currentBuffPosition = 0;
94 for (itSend = itbSend; itSend != iteSend; ++itSend)
95 {
96   int destRank = itSend->first;
97   const std::vector<size_t>& globalIndexOfCurrentGridSourceToSend = itSend->second;
98   int countSize = globalIndexOfCurrentGridSourceToSend.size();
99   for (int idx = 0; idx < (countSize); ++idx)
100   {
101     int index = globalIndexOfCurrentGridSourceToSend[idx] - ibeginSrc;
102     sendBuff[idx+currentBuffPosition] = (axisSrc_->value)(index);
103   }
104   currentSendBuff = sendBuff + currentBuffPosition;
105   MPI_Send(currentSendBuff, countSize, MPI_DOUBLE, destRank, 14, client->intraComm);
106   currentBuffPosition += countSize;
107 }
108
109 // Receiving global index of grid source sending from current grid source
110 std::map<int,std::vector<std::vector<std::pair<size_t,double> > > >::const_iterator itbRecv = globalIndexToReceive.begin(), itRecv,
111                                                                                     iteRecv = globalIndexToReceive.end();
112 int recvBuffSize = 0;
113 for (itRecv = itbRecv; itRecv != iteRecv; ++itRecv) recvBuffSize += (itRecv->second).size();
114
115 Scalar* recvBuff, *currentRecvBuff;
116 if (0 != recvBuffSize) recvBuff = new Scalar [recvBuffSize];
117 for (StdSize idx = 0; idx < recvBuffSize; ++idx) recvBuff[idx] = NumTraits<Scalar>::sfmax();
118
119 int currentRecvBuffPosition = 0;
120 for (itRecv = itbRecv; itRecv != iteRecv; ++itRecv)
121 {
122   MPI_Status status;
123   int srcRank = itRecv->first;
124   int countSize = (itRecv->second).size();
125   currentRecvBuff = recvBuff + currentRecvBuffPosition;
126   MPI_Recv(currentRecvBuff, countSize, MPI_DOUBLE, srcRank, 14, client->intraComm, &status);
127   currentRecvBuffPosition += countSize;
128 }
129
130 int ibeginDest = axisDest_->begin.getValue();
131 currentRecvBuff = recvBuff;
132 for (itRecv = itbRecv; itRecv != iteRecv; ++itRecv)
133 {
134   int countSize = (itRecv->second).size();
135   for (int idx = 0; idx < countSize; ++idx, ++currentRecvBuff)
136   {
137     int ssize = (itRecv->second)[idx].size();
138     for (int i = 0; i < ssize; ++i)
139     {
140       int index = ((itRecv->second)[idx][i]).first - ibeginDest;
141       (axisDest_->value)(index) = *currentRecvBuff;
142     }
143   }
144 }
145
146 if (0 != sendBuffSize) delete [] sendBuff;
147 if (0 != recvBuffSize) delete [] recvBuff;
148}
149
150}
Note: See TracBrowser for help on using the repository browser.