source: XIOS/trunk/src/filter/grid_transformation.cpp @ 620

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

Implementing generic transformation algorithm (local commit)

+) Implement 3 important classes:

-gridTransformation to read transformation info from grid and interface with the rest of XIOS
-transformationMapping to be in charge of sending/receiving transformation info among clients
-transformationAlgorithm to represent various algorithms

+) Make some change on field to use the new classes

Test
+) Only test_new_features with inversed axis

File size: 8.1 KB
Line 
1#include "grid_transformation.hpp"
2#include "axis_inverse.hpp"
3#include "transformation_mapping.hpp"
4
5namespace xios {
6CGridTransformation::CGridTransformation(CGrid* destination, CGrid* source)
7: gridSource_(source), gridDestination_(destination)
8{
9  //Verify the compatibity between two grids
10  int numElement = gridDestination_->axis_domain_order.numElements();
11  if (numElement != gridSource_->axis_domain_order.numElements())
12    ERROR("CGridTransformation::CGridTransformation(CGrid* destination, CGrid* source)",
13       << "Two grids have different number of elements"
14       << "Number of elements of grid source " <<gridSource_->getId() << " is " << gridSource_->axis_domain_order.numElements()  << std::endl
15       << "Number of elements of grid destination " <<gridDestination_->getId() << " is " << numElement);
16
17  for (int i = 0; i < numElement; ++i)
18  {
19    if (gridDestination_->axis_domain_order(i) != gridSource_->axis_domain_order(i))
20      ERROR("CGridTransformation::CGridTransformation(CGrid* destination, CGrid* source)",
21         << "Transformed grid and its grid source have incompatible elements"
22         << "Grid source " <<gridSource_->getId() << std::endl
23         << "Grid destination " <<gridDestination_->getId());
24  }
25
26  gridSourceDimensionSize_ = gridSource_->getGlobalDimension();
27  gridDestinationDimensionSize_ = gridDestination_->getGlobalDimension();
28  initializeAlgorithms();
29}
30
31CGridTransformation::~CGridTransformation()
32{
33  std::map<int, std::vector<CGenericAlgorithmTransformation*> >::const_iterator itb = algoTransformation_.begin(), it,
34                                                                                ite = algoTransformation_.end();
35  for (it = itb; it != ite; ++it)
36    for (int idx = 0; idx < (it->second).size(); ++idx)
37      delete (it->second)[idx];
38
39  std::map<int, std::vector<CArray<int,1>* > >::const_iterator itMapRecv, iteMapRecv;
40  itMapRecv = localIndexToReceiveOnGridDest_.begin();
41  iteMapRecv = localIndexToReceiveOnGridDest_.end();
42  for (; itMapRecv != iteMapRecv; ++itMapRecv)
43  {
44    int numVec = (itMapRecv->second).size();
45    for (int idx = 0; idx < numVec; ++idx) delete (itMapRecv->second)[idx];
46  }
47
48  std::map<int, CArray<int,1>* >::const_iterator itMap, iteMap;
49  itMap = localIndexToSendFromGridSource_.begin();
50  iteMap = localIndexToSendFromGridSource_.end();
51  for (; itMap != iteMap; ++itMap) delete (itMap->second);
52}
53
54void CGridTransformation::initializeAlgorithms()
55{
56  initializeAxisAlgorithms();
57  initializeDomainAlgorithms();
58}
59
60/*!
61  Initialize the algorithms corresponding to transformation info contained in each axis.
62If an axis has transformations, these transformations will be represented in form of vector of CTransformation pointers
63In general, each axis can have several transformations performed on itself. However, should they be done seperately or combinely (of course in order)?
64For now, one approach is to do these combinely but maybe this needs changing.
65*/
66void CGridTransformation::initializeAxisAlgorithms()
67{
68  std::vector<int> axisPositionInGrid;
69  std::vector<CAxis*> axisListDestP = gridDestination_->getAxis();
70  std::vector<CAxis*> axisListSrcP = gridSource_->getAxis();
71  if (!axisListDestP.empty())
72  {
73    int idx = 0;
74    for (int i = 0; i < gridDestination_->axis_domain_order.numElements(); ++i)
75    {
76      if (false == (gridDestination_->axis_domain_order)(i))
77      {
78        axisPositionInGrid.push_back(idx);
79        ++idx;
80      }
81      else idx += 2;
82    }
83
84    for (int i = 0; i < axisListDestP.size(); ++i)
85    {
86      if (axisListDestP[i]->hasTransformation())
87      {
88        CGenericAlgorithmTransformation* algo = new CAxisInverse(axisListDestP[i], axisListSrcP[i]);
89        algoTransformation_[axisPositionInGrid[i]].push_back(algo);
90      }
91    }
92  }
93}
94
95void CGridTransformation::initializeDomainAlgorithms()
96{
97
98}
99
100void CGridTransformation::computeTransformation()
101{
102  const CArray<size_t,1>& globalIndexGridDestSendToServer = gridDestination_->getDistributionClient()->getGlobalDataIndexSendToServer();
103  std::map<int, std::vector<CGenericAlgorithmTransformation*> >::const_iterator itbMap, itMap, iteMap;
104  itbMap = algoTransformation_.begin();
105  iteMap = algoTransformation_.end();
106
107  std::vector<CGenericAlgorithmTransformation*>::const_iterator itbVec, itVec, iteVec;
108
109  for (itMap = itbMap; itMap != iteMap; ++itMap)
110  {
111    int elementPosition = itMap->first;
112    itbVec = (itMap->second).begin();
113    iteVec = (itMap->second).end();
114    for (itVec = itbVec; itVec != iteVec; ++itVec)
115    {
116      (*itVec)->computeGlobalSourceIndex(elementPosition,
117                                         gridDestinationDimensionSize_,
118                                         globalIndexGridDestSendToServer,
119                                         globaIndexMapFromDestToSource_);
120    }
121  }
122}
123
124void CGridTransformation::computeTransformationMapping()
125{
126  CTransformationMapping transformationMap(gridDestination_, gridSource_);
127
128  // First of all, need to compute global index mapping representing transformation algorithms
129  computeTransformation();
130
131  // Then compute transformation mapping among clients
132  transformationMap.computeTransformationMapping(globaIndexMapFromDestToSource_);
133
134  const std::map<int,std::vector<std::vector<size_t> > >& globalIndexToReceive = transformationMap.getGlobalIndexReceivedOnGridDestMapping();
135  const std::map<int,std::vector<size_t> >& globalIndexToSend = transformationMap.getGlobalIndexSendToGridDestMapping();
136
137  CArray<int, 1> localIndexOnClientDest = gridDestination_->getDistributionClient()->getLocalDataIndexOnClient();
138  CArray<size_t,1> globalIndexOnClientDest = gridDestination_->getDistributionClient()->getGlobalDataIndexSendToServer();
139
140  CArray<int, 1> localIndexOnClientSrc = gridSource_->getDistributionClient()->getLocalDataIndexOnClient();
141  CArray<size_t,1> globalIndexOnClientSrc = gridSource_->getDistributionClient()->getGlobalDataIndexSendToServer();
142
143
144  std::vector<size_t>::const_iterator itbVec, itVec, iteVec;
145  CArray<size_t, 1>::iterator itbArr, itArr, iteArr;
146
147  std::map<int,std::vector<std::vector<size_t> > >::const_iterator itbMapRecv, itMapRecv, iteMapRecv;
148  // Find out local index on grid destination (received)
149  itbMapRecv = globalIndexToReceive.begin();
150  iteMapRecv = globalIndexToReceive.end();
151  itbArr = globalIndexOnClientDest.begin();
152  iteArr = globalIndexOnClientDest.end();
153  for (itMapRecv = itbMapRecv; itMapRecv != iteMapRecv; ++itMapRecv)
154  {
155    int sourceRank = itMapRecv->first;
156    int numGlobalIndex = (itMapRecv->second).size();
157    for (int i = 0; i < numGlobalIndex; ++i)
158    {
159      int vecSize = ((itMapRecv->second)[i]).size();
160      CArray<int,1>* ptr = new CArray<int,1>(vecSize);
161      localIndexToReceiveOnGridDest_[sourceRank].push_back(ptr);
162      for (int idx = 0; idx < vecSize; ++idx)
163      {
164        itArr = std::find(itbArr, iteArr, (itMapRecv->second)[i][idx]);
165        if (iteArr != itArr)
166        {
167          int localIdx = std::distance(itbArr, itArr);
168          (*localIndexToReceiveOnGridDest_[sourceRank][i])(idx) = localIndexOnClientDest(localIdx);
169        }
170      }
171    }
172
173  }
174
175  std::map<int,std::vector<size_t> >::const_iterator itbMap, itMap, iteMap;
176  // Find out local index on grid source (to send)
177  itbMap = globalIndexToSend.begin();
178  iteMap = globalIndexToSend.end();
179  itbArr = globalIndexOnClientSrc.begin();
180  iteArr = globalIndexOnClientSrc.end();
181  for (itMap = itbMap; itMap != iteMap; ++itMap)
182  {
183    CArray<int,1>* ptr = new CArray<int,1>((itMap->second).size());
184    localIndexToSendFromGridSource_[itMap->first] = ptr;
185    int vecSize = (itMap->second).size();
186    for (int idx = 0; idx < vecSize; ++idx)
187    {
188      itArr = std::find(itbArr, iteArr, (itMap->second)[idx]);
189      if (iteArr != itArr)
190      {
191        int localIdx = std::distance(itbArr, itArr);
192        (*localIndexToSendFromGridSource_[itMap->first])(idx) = localIndexOnClientSrc(localIdx);
193      }
194    }
195  }
196}
197
198std::map<int, CArray<int,1>* > CGridTransformation::getLocalIndexToSendFromGridSource()
199{
200  return localIndexToSendFromGridSource_;
201}
202
203std::map<int, std::vector<CArray<int,1>* > > CGridTransformation::getLocalIndexToReceiveOnGridDest()
204{
205  return localIndexToReceiveOnGridDest_;
206}
207
208}
Note: See TracBrowser for help on using the repository browser.