1 | #include "grid_remote_connector.hpp" |
---|
2 | #include "client_client_dht_template.hpp" |
---|
3 | #include "leader_process.hpp" |
---|
4 | #include "mpi.hpp" |
---|
5 | |
---|
6 | |
---|
7 | |
---|
8 | namespace xios |
---|
9 | { |
---|
10 | /** |
---|
11 | * \brief class constructor. |
---|
12 | * \param srcView List of sources views. |
---|
13 | * \param dstView List of remotes views. |
---|
14 | * \param localComm Local MPI communicator |
---|
15 | * \param remoteSize Size of the remote communicator |
---|
16 | */ |
---|
17 | CGridRemoteConnector::CGridRemoteConnector(vector<shared_ptr<CLocalView>>& srcView, vector<shared_ptr<CDistributedView>>& dstView, MPI_Comm localComm, int remoteSize) |
---|
18 | : srcView_(srcView), dstView_(dstView), localComm_(localComm), remoteSize_(remoteSize) |
---|
19 | {} |
---|
20 | |
---|
21 | /** |
---|
22 | * \brief class constructor. |
---|
23 | * \param srcView List of sources views. |
---|
24 | * \param dstView List of remotes views. |
---|
25 | * \param localComm Local MPI communicator |
---|
26 | * \param remoteSize Size of the remote communicator |
---|
27 | */ |
---|
28 | CGridRemoteConnector::CGridRemoteConnector(vector<shared_ptr<CLocalView>>& srcView, vector< shared_ptr<CLocalView> >& dstView, MPI_Comm localComm, int remoteSize) |
---|
29 | : srcView_(srcView), localComm_(localComm), remoteSize_(remoteSize) |
---|
30 | { |
---|
31 | for(auto& it : dstView) dstView_.push_back((shared_ptr<CDistributedView>) it) ; |
---|
32 | } |
---|
33 | |
---|
34 | |
---|
35 | /** |
---|
36 | * \brief Compute if each view composing the source grid and the remote grid is distributed or not. |
---|
37 | * Result is stored on internal attributes \b isSrcViewDistributed_ and \b isDstViewDistributed_. |
---|
38 | * \detail To compute this, a hash is computed for each array on indices. The hash must permutable, i.e. |
---|
39 | * the order of the list of global indices doesn't influence the value of the hash. So simply a sum of |
---|
40 | * hash of each indices is used for the whole array. After, the computed hash are compared with each other |
---|
41 | * ranks of \b localComm_ MPI communicator using an MPI_ALLReduce. If, for each ranks, the hash is the same |
---|
42 | * then the view is not distributed |
---|
43 | */ |
---|
44 | void CGridRemoteConnector::computeViewDistribution(void) |
---|
45 | { |
---|
46 | HashXIOS<size_t> hashGlobalIndex; // hash function-object |
---|
47 | |
---|
48 | int nDst = dstView_.size() ; |
---|
49 | vector<size_t> hashRank(remoteSize_) ; |
---|
50 | isDstViewDistributed_.resize(nDst) ; |
---|
51 | |
---|
52 | for(int i=0; i<nDst; i++) |
---|
53 | { |
---|
54 | map<int,CArray<size_t,1>> globalIndexView ; |
---|
55 | dstView_[i]->getGlobalIndexView(globalIndexView) ; |
---|
56 | hashRank.assign(remoteSize_,0) ; // everybody ranks to 0 except rank of the remote view I have |
---|
57 | // that would be assign to my local hash |
---|
58 | for(auto& it : globalIndexView) |
---|
59 | { |
---|
60 | int rank=it.first ; |
---|
61 | CArray<size_t,1>& globalIndex = it.second ; |
---|
62 | size_t globalIndexSize = globalIndex.numElements(); |
---|
63 | size_t hashValue=0 ; |
---|
64 | for(size_t ind=0;ind<globalIndexSize;ind++) hashValue += hashGlobalIndex(globalIndex(ind)) ; |
---|
65 | hashRank[rank] += hashValue ; |
---|
66 | } |
---|
67 | // sum all the hash for every process of the local comm. The reduce is on the size of remote view (remoteSize_) |
---|
68 | // after that for each rank of the remote view, we get the hash |
---|
69 | MPI_Allreduce(MPI_IN_PLACE, hashRank.data(), remoteSize_, MPI_SIZE_T, MPI_SUM, localComm_) ; |
---|
70 | size_t value = hashRank[0] ; |
---|
71 | isDstViewDistributed_[i]=false ; |
---|
72 | for(int j=0 ; j<remoteSize_ ; j++) |
---|
73 | if (value != hashRank[j]) |
---|
74 | { |
---|
75 | isDstViewDistributed_[i]=true ; |
---|
76 | break ; |
---|
77 | } |
---|
78 | } |
---|
79 | |
---|
80 | int nSrc = srcView_.size() ; |
---|
81 | int commSize,commRank ; |
---|
82 | MPI_Comm_size(localComm_,&commSize) ; |
---|
83 | MPI_Comm_rank(localComm_,&commRank) ; |
---|
84 | hashRank.resize(commSize,0) ; |
---|
85 | isSrcViewDistributed_.resize(nSrc) ; |
---|
86 | |
---|
87 | for(int i=0; i<nSrc; i++) |
---|
88 | { |
---|
89 | CArray<size_t,1> globalIndex ; |
---|
90 | srcView_[i]->getGlobalIndexView(globalIndex) ; |
---|
91 | hashRank.assign(commSize,0) ; // 0 for everybody except my rank |
---|
92 | size_t globalIndexSize = globalIndex.numElements() ; |
---|
93 | size_t hashValue=0 ; |
---|
94 | for(size_t ind=0;ind<globalIndexSize;ind++) hashValue += hashGlobalIndex(globalIndex(ind)) ; |
---|
95 | hashRank[commRank] += hashValue ; |
---|
96 | |
---|
97 | // Same method than for remote view |
---|
98 | MPI_Allreduce(MPI_IN_PLACE, hashRank.data(), commSize, MPI_SIZE_T, MPI_SUM, localComm_) ; |
---|
99 | size_t value = hashRank[0] ; |
---|
100 | isSrcViewDistributed_[i]=false ; |
---|
101 | for(int j=0 ; j<commSize ; j++) |
---|
102 | if (value != hashRank[j]) |
---|
103 | { |
---|
104 | isSrcViewDistributed_[i]=true ; |
---|
105 | break ; |
---|
106 | } |
---|
107 | } |
---|
108 | |
---|
109 | } |
---|
110 | |
---|
111 | /** |
---|
112 | * \brief Compute the connector, i.e. compute the \b elements_ attribute. |
---|
113 | * \detail Depending of the distributions of the view computed in the computeViewDistribution() call, the connector is computed in computeConnectorMethods(), and to achieve better optimisation |
---|
114 | * some redondant ranks can be removed from the elements_ map. |
---|
115 | */ |
---|
116 | void CGridRemoteConnector::computeConnector(bool eliminateRedundant) |
---|
117 | { |
---|
118 | if (eliminateRedundant) |
---|
119 | { |
---|
120 | computeViewDistribution() ; |
---|
121 | computeConnectorMethods() ; |
---|
122 | computeRedondantRanks() ; |
---|
123 | for(auto& rank : rankToRemove_) |
---|
124 | for(auto& element : elements_) element.erase(rank) ; |
---|
125 | } |
---|
126 | else |
---|
127 | { |
---|
128 | computeViewDistribution() ; |
---|
129 | computeConnectorRedundant() ; |
---|
130 | } |
---|
131 | } |
---|
132 | |
---|
133 | /** |
---|
134 | * \brief Compute the connector, i.e. compute the \b elements_ attribute. |
---|
135 | * \detail In this routine we don't eliminate redundant cells as it it performed in |
---|
136 | * computeConnectorMethods. It can be usefull to perform reduce operation over processes. |
---|
137 | In future, some optimisation could be done considering full redondance of the |
---|
138 | source view or the destination view. |
---|
139 | */ |
---|
140 | void CGridRemoteConnector::computeConnectorRedundant(void) |
---|
141 | { |
---|
142 | vector<shared_ptr<CLocalView>> srcView ; |
---|
143 | vector<shared_ptr<CDistributedView>> dstView ; |
---|
144 | vector<int> indElements ; |
---|
145 | elements_.resize(srcView_.size()) ; |
---|
146 | |
---|
147 | bool srcViewsNonDistributed=true ; // not usefull now but later for optimization |
---|
148 | for(int i=0;i<srcView_.size();i++) srcViewsNonDistributed = srcViewsNonDistributed && !isSrcViewDistributed_[i] ; |
---|
149 | |
---|
150 | bool dstViewsNonDistributed=true ; // not usefull now but later for optimization |
---|
151 | for(int i=0;i<dstView_.size();i++) dstViewsNonDistributed = dstViewsNonDistributed && !isDstViewDistributed_[i] ; |
---|
152 | |
---|
153 | for(int i=0;i<srcView_.size();i++) |
---|
154 | { |
---|
155 | srcView.push_back(srcView_[i]) ; |
---|
156 | dstView.push_back(dstView_[i]) ; |
---|
157 | indElements.push_back(i) ; |
---|
158 | } |
---|
159 | |
---|
160 | computeGenericMethod(srcView, dstView, indElements) ; |
---|
161 | |
---|
162 | map<int,bool> ranks ; |
---|
163 | for(auto& it : elements_[indElements[0]]) |
---|
164 | { |
---|
165 | if (it.second.numElements()==0) ranks[it.first] = false ; |
---|
166 | else ranks[it.first] = true ; |
---|
167 | } |
---|
168 | |
---|
169 | } |
---|
170 | |
---|
171 | |
---|
172 | /** |
---|
173 | * \brief Compute the connector, i.e. compute the \b elements_ attribute. |
---|
174 | * \detail In order to achive better optimisation, |
---|
175 | * we distingute the case when the grid is not distributed on source grid (\bcomputeSrcNonDistributed), |
---|
176 | * or the remote grid (\b computeDstNonDistributed), or the both (\b computeSrcDstNonDistributed). |
---|
177 | * Otherwise the generic method is called computeGenericMethod. Note that in the case, if one element view |
---|
178 | * is not distributed on the source and on the remote grid, then we can used the tensorial product |
---|
179 | * property to computing it independently using \b computeSrcDstNonDistributed method. |
---|
180 | * After that, we call the \b removeRedondantRanks method to supress blocks of data that can be sent |
---|
181 | * redondantly the the remote servers |
---|
182 | */ |
---|
183 | void CGridRemoteConnector::computeConnectorMethods(void) |
---|
184 | { |
---|
185 | vector<shared_ptr<CLocalView>> srcView ; |
---|
186 | vector<shared_ptr<CDistributedView>> dstView ; |
---|
187 | vector<int> indElements ; |
---|
188 | elements_.resize(srcView_.size()) ; |
---|
189 | |
---|
190 | bool srcViewsNonDistributed=true ; |
---|
191 | for(int i=0;i<srcView_.size();i++) srcViewsNonDistributed = srcViewsNonDistributed && !isSrcViewDistributed_[i] ; |
---|
192 | |
---|
193 | bool dstViewsNonDistributed=true ; |
---|
194 | for(int i=0;i<dstView_.size();i++) dstViewsNonDistributed = dstViewsNonDistributed && !isDstViewDistributed_[i] ; |
---|
195 | |
---|
196 | if (srcViewsNonDistributed) |
---|
197 | { |
---|
198 | int commRank, commSize ; |
---|
199 | MPI_Comm_rank(localComm_,&commRank) ; |
---|
200 | MPI_Comm_size(localComm_,&commSize) ; |
---|
201 | list<int> remoteRanks; |
---|
202 | list<int> notUsed ; |
---|
203 | map<int,bool> ranks ; |
---|
204 | computeLeaderProcess(commRank, commSize, remoteSize_, remoteRanks, notUsed) ; |
---|
205 | for(int rank : remoteRanks) ranks[rank]=true ; |
---|
206 | |
---|
207 | for(int i=0; i<srcView_.size(); i++) |
---|
208 | { |
---|
209 | if (isDstViewDistributed_[i]) computeSrcNonDistributed(i) ; |
---|
210 | else computeSrcDstNonDistributed(i, ranks) ; |
---|
211 | } |
---|
212 | } |
---|
213 | else if (dstViewsNonDistributed) |
---|
214 | { |
---|
215 | map<int,bool> ranks ; |
---|
216 | for(int i=0;i<remoteSize_;i++) ranks[i]=true ; |
---|
217 | for(int i=0; i<srcView_.size(); i++) |
---|
218 | { |
---|
219 | if (isSrcViewDistributed_[i]) computeDstNonDistributed(i,ranks) ; |
---|
220 | else computeSrcDstNonDistributed(i,ranks) ; |
---|
221 | } |
---|
222 | } |
---|
223 | else |
---|
224 | { |
---|
225 | for(int i=0;i<srcView_.size();i++) |
---|
226 | if (isSrcViewDistributed_[i] || isDstViewDistributed_[i]) |
---|
227 | { |
---|
228 | srcView.push_back(srcView_[i]) ; |
---|
229 | dstView.push_back(dstView_[i]) ; |
---|
230 | indElements.push_back(i) ; |
---|
231 | } |
---|
232 | |
---|
233 | computeGenericMethod(srcView, dstView, indElements) ; |
---|
234 | |
---|
235 | map<int,bool> ranks ; |
---|
236 | for(auto& it : elements_[indElements[0]]) |
---|
237 | { |
---|
238 | if (it.second.numElements()==0) ranks[it.first] = false ; |
---|
239 | else ranks[it.first] = true ; |
---|
240 | } |
---|
241 | |
---|
242 | for(int i=0;i<srcView_.size();i++) |
---|
243 | if (!isSrcViewDistributed_[i] && !isDstViewDistributed_[i]) computeSrcDstNonDistributed(i, ranks) ; |
---|
244 | } |
---|
245 | |
---|
246 | } |
---|
247 | |
---|
248 | |
---|
249 | /** |
---|
250 | * \brief Compute the connector for the element \b i when the source view is not distributed. |
---|
251 | * After the call element_[i] is defined. |
---|
252 | * \param i Indice of the element composing the source grid. |
---|
253 | */ |
---|
254 | |
---|
255 | void CGridRemoteConnector::computeSrcNonDistributed(int i) |
---|
256 | { |
---|
257 | auto& element = elements_[i] ; |
---|
258 | map<int,CArray<size_t,1>> globalIndexView ; |
---|
259 | dstView_[i]->getGlobalIndexView(globalIndexView) ; |
---|
260 | |
---|
261 | CClientClientDHTTemplate<int>::Index2InfoTypeMap dataInfo; |
---|
262 | |
---|
263 | for(auto& it : globalIndexView) |
---|
264 | { |
---|
265 | auto& globalIndex=it.second ; |
---|
266 | for(size_t ind : globalIndex) dataInfo[ind]=it.first ; |
---|
267 | } |
---|
268 | |
---|
269 | // First we feed the distributed hash map with key (remote global index) |
---|
270 | // associated with the value of the remote rank |
---|
271 | CClientClientDHTTemplate<int> DHT(dataInfo, localComm_) ; |
---|
272 | // after we feed the DHT with the local global indices of the source view |
---|
273 | |
---|
274 | int commRank, commSize ; |
---|
275 | MPI_Comm_rank(localComm_,&commRank) ; |
---|
276 | MPI_Comm_size(localComm_,&commSize) ; |
---|
277 | CArray<size_t,1> srcIndex ; |
---|
278 | // like the source view is not distributed, then only the rank 0 need to feed the DHT |
---|
279 | if (commRank==0) srcView_[i]->getGlobalIndexView(srcIndex) ; |
---|
280 | |
---|
281 | // compute the mapping |
---|
282 | DHT.computeIndexInfoMapping(srcIndex) ; |
---|
283 | auto& returnInfo = DHT.getInfoIndexMap() ; |
---|
284 | |
---|
285 | // returnInfo contains now the map for each global indices to send to a list of remote rank |
---|
286 | // only for the rank=0 because it is the one to feed the DHT |
---|
287 | // so it need to send the list to each server leader i.e. the local process that handle specifically one or more |
---|
288 | // servers |
---|
289 | |
---|
290 | // rankIndGlo : rankIndGlo[rank][indGlo] : list of indice to send the the remote server of rank "rank" |
---|
291 | vector<vector<size_t>> rankIndGlo(remoteSize_) ; |
---|
292 | if (commRank==0) |
---|
293 | for(auto& it1 : returnInfo) |
---|
294 | for(auto& it2 : it1.second) rankIndGlo[it2].push_back(it1.first) ; |
---|
295 | |
---|
296 | |
---|
297 | vector<MPI_Request> requests ; |
---|
298 | |
---|
299 | if (commRank==0) |
---|
300 | { |
---|
301 | requests.resize(remoteSize_) ; |
---|
302 | for(int i=0 ; i<remoteSize_;i++) |
---|
303 | { |
---|
304 | // ok send only the global indices for a server to the "server leader" |
---|
305 | int rank = getLeaderRank(commSize, remoteSize_, i) ; |
---|
306 | MPI_Isend(rankIndGlo[i].data(), rankIndGlo[i].size(), MPI_SIZE_T, rank, i ,localComm_, &requests[i]) ; |
---|
307 | } |
---|
308 | } |
---|
309 | |
---|
310 | list<int> remoteRanks; |
---|
311 | list<int> notUsed ; |
---|
312 | // I am a server leader of which remote ranks ? |
---|
313 | computeLeaderProcess(commRank, commSize, remoteSize_, remoteRanks, notUsed) ; |
---|
314 | |
---|
315 | for(auto remoteRank : remoteRanks) |
---|
316 | { |
---|
317 | MPI_Status status ; |
---|
318 | int size ; |
---|
319 | MPI_Probe(0,remoteRank,localComm_, &status); |
---|
320 | MPI_Get_count(&status, MPI_SIZE_T, &size) ; |
---|
321 | elements_[i][remoteRank].resize(size) ; |
---|
322 | // for each remote ranks receive the global indices from proc 0 |
---|
323 | MPI_Recv(elements_[i][remoteRank].dataFirst(),size, MPI_SIZE_T,0,remoteRank, localComm_,&status) ; |
---|
324 | } |
---|
325 | |
---|
326 | if (commRank==0) |
---|
327 | { |
---|
328 | vector<MPI_Status> status(remoteSize_) ; |
---|
329 | // asynchronous for sender, wait for completion |
---|
330 | MPI_Waitall(remoteSize_, requests.data(), status.data()) ; |
---|
331 | } |
---|
332 | } |
---|
333 | |
---|
334 | /** |
---|
335 | * \brief Compute the remote connector for the element \b i when the remote view is not distributed. |
---|
336 | * After the call, element_[i] is defined. |
---|
337 | * \param i Indice of the element composing the remote grid. |
---|
338 | * \param ranks The list of rank for which the local proc is in charge to compute the connector |
---|
339 | * (if leader server for exemple). if ranks[rank] == false the corresponding elements_ |
---|
340 | * is set to void array (no data to sent) just in order to notify corresponding remote server |
---|
341 | * that the call is collective with each other one |
---|
342 | */ |
---|
343 | void CGridRemoteConnector::computeDstNonDistributed(int i, map<int,bool>& ranks) |
---|
344 | { |
---|
345 | auto& element = elements_[i] ; |
---|
346 | map<int,CArray<size_t,1>> globalIndexView ; |
---|
347 | dstView_[i]->getGlobalIndexView(globalIndexView) ; |
---|
348 | |
---|
349 | |
---|
350 | CClientClientDHTTemplate<int>::Index2InfoTypeMap dataInfo; |
---|
351 | |
---|
352 | // First we feed the distributed hash map with key (remote global index) |
---|
353 | // associated with the value of the remote rank |
---|
354 | for(auto& it : globalIndexView) |
---|
355 | if (it.first==0) // since the remote view is not distributed, insert only the remote rank 0 |
---|
356 | { |
---|
357 | auto& globalIndex=it.second ; |
---|
358 | for(size_t ind : globalIndex) dataInfo[ind]=0 ; // associated the the rank 0 |
---|
359 | } |
---|
360 | |
---|
361 | CClientClientDHTTemplate<int> DHT(dataInfo, localComm_) ; |
---|
362 | // after we feed the DHT with the local global indices of the source view |
---|
363 | |
---|
364 | CArray<size_t,1> srcIndex ; |
---|
365 | srcView_[i]->getGlobalIndexView(srcIndex) ; |
---|
366 | DHT.computeIndexInfoMapping(srcIndex) ; |
---|
367 | auto& returnInfo = DHT.getInfoIndexMap() ; |
---|
368 | |
---|
369 | // returnInfo contains now the map for each global indices to send to a list of remote rank |
---|
370 | // now construct the element_ list of global indices for each rank in my list except if the erray must be empty |
---|
371 | for (auto& rank : ranks) |
---|
372 | { |
---|
373 | if (rank.second) // non empty array => for rank that have not any data to be received |
---|
374 | { |
---|
375 | int size=0 ; |
---|
376 | for(auto& it : returnInfo) if (!it.second.empty()) size++ ; |
---|
377 | auto& array = element[rank.first] ; |
---|
378 | array.resize(size) ; |
---|
379 | size=0 ; |
---|
380 | for(auto& it : returnInfo) |
---|
381 | if (!it.second.empty()) |
---|
382 | { |
---|
383 | array(size)=it.first ; |
---|
384 | size++ ; |
---|
385 | } |
---|
386 | } |
---|
387 | else element[rank.first] = CArray<size_t,1>(0) ; // empty array => for rank that have not any data to be received |
---|
388 | } |
---|
389 | } |
---|
390 | |
---|
391 | /** |
---|
392 | * \brief Compute the remote connector for the element \b i when the source and the remote view are not distributed. |
---|
393 | * After the call, element_[i] is defined. |
---|
394 | * \param i Indice of the element composing the remote grid. |
---|
395 | * \param ranks The list of rank for which the local proc is in charge to compute the connector |
---|
396 | * (if leader server for exemple). if ranks[rank] == false the corresponding elements_ |
---|
397 | * is set to void array (no data to sent) just in order to notify corresponding remote server |
---|
398 | * that the call is collective with each other one |
---|
399 | */ |
---|
400 | |
---|
401 | void CGridRemoteConnector::computeSrcDstNonDistributed(int i, map<int,bool>& ranks) |
---|
402 | { |
---|
403 | auto& element = elements_[i] ; |
---|
404 | map<int,CArray<size_t,1>> globalIndexView ; |
---|
405 | dstView_[i]->getGlobalIndexView(globalIndexView) ; |
---|
406 | |
---|
407 | |
---|
408 | CClientClientDHTTemplate<int>::Index2InfoTypeMap dataInfo; |
---|
409 | // First we feed the distributed hash map with key (remote global index) |
---|
410 | // associated with the value of the remote rank |
---|
411 | |
---|
412 | for(auto& it : globalIndexView) |
---|
413 | if (it.first==0) // insert only the remote rank 0 since the remote view is not distributed |
---|
414 | { |
---|
415 | auto& globalIndex=it.second ; |
---|
416 | for(size_t ind : globalIndex) dataInfo[ind]=0 ; // associated the the rank 0 |
---|
417 | } |
---|
418 | |
---|
419 | CClientClientDHTTemplate<int> DHT(dataInfo, localComm_) ; |
---|
420 | // after we feed the DHT with the local global indices of the source view |
---|
421 | |
---|
422 | int commRank, commSize ; |
---|
423 | MPI_Comm_rank(localComm_,&commRank) ; |
---|
424 | MPI_Comm_size(localComm_,&commSize) ; |
---|
425 | CArray<size_t,1> srcIndex ; |
---|
426 | |
---|
427 | // like the source view is not distributed, then only the rank 0 need to feed the DHT |
---|
428 | if (commRank==0) srcView_[i]->getGlobalIndexView(srcIndex) ; |
---|
429 | DHT.computeIndexInfoMapping(srcIndex) ; |
---|
430 | auto& returnInfo = DHT.getInfoIndexMap() ; |
---|
431 | |
---|
432 | vector<size_t> indGlo ; |
---|
433 | if (commRank==0) |
---|
434 | for(auto& it1 : returnInfo) |
---|
435 | for(auto& it2 : it1.second) indGlo.push_back(it1.first) ; |
---|
436 | |
---|
437 | // now local rank 0 know which indices to seed to remote rank 0, but all the server |
---|
438 | // must receive the same information. So only the leader rank will sent this. |
---|
439 | // So local rank 0 must broadcast the information to all leader. |
---|
440 | // for this we create a new communicator composed of local process that must send data |
---|
441 | // to a remote rank, data are broadcasted, and element_[i] is construction for each remote |
---|
442 | // rank in charge |
---|
443 | int color=0 ; |
---|
444 | if (ranks.empty()) color=0 ; |
---|
445 | else color=1 ; |
---|
446 | if (commRank==0) color=1 ; |
---|
447 | MPI_Comm newComm ; |
---|
448 | MPI_Comm_split(localComm_, color, commRank, &newComm) ; |
---|
449 | if (color==1) |
---|
450 | { |
---|
451 | // ok, I am part of the process that must send something to one or more remote server |
---|
452 | // so I get the list of global indices from rank 0 |
---|
453 | int dataSize ; |
---|
454 | if (commRank==0) dataSize=indGlo.size() ; |
---|
455 | MPI_Bcast(&dataSize,1,MPI_INT, 0, newComm) ; |
---|
456 | indGlo.resize(dataSize) ; |
---|
457 | MPI_Bcast(indGlo.data(),dataSize,MPI_SIZE_T,0,newComm) ; |
---|
458 | } |
---|
459 | MPI_Comm_free(&newComm) ; |
---|
460 | |
---|
461 | // construct element_[i] from indGlo |
---|
462 | for(auto& rank : ranks) |
---|
463 | { |
---|
464 | if (rank.second) |
---|
465 | { |
---|
466 | int dataSize=indGlo.size(); |
---|
467 | auto& element = elements_[i][rank.first] ; |
---|
468 | element.resize(dataSize) ; |
---|
469 | for(int i=0;i<dataSize; i++) element(i)=indGlo[i] ; |
---|
470 | } |
---|
471 | else element[rank.first] = CArray<size_t,1>(0) ; |
---|
472 | } |
---|
473 | |
---|
474 | } |
---|
475 | |
---|
476 | |
---|
477 | /** |
---|
478 | * \brief Generic method the compute the grid remote connector. Only distributed elements are specifed in the source view and remote view. |
---|
479 | * Connector for non distributed elements are computed separatly to improve performance and memory consumption. After the call, |
---|
480 | * \b elements_ is defined. |
---|
481 | * \param srcView List of the source views composing the grid, without non distributed views |
---|
482 | * \param dstView List of the remote views composing the grid, without non distributed views |
---|
483 | * \param indElements Index of the view making the correspondance between all views and views distributed (that are in input) |
---|
484 | */ |
---|
485 | void CGridRemoteConnector::computeGenericMethod(vector<shared_ptr<CLocalView>>& srcView, vector<shared_ptr<CDistributedView>>& dstView, vector<int>& indElements) |
---|
486 | { |
---|
487 | // generic method, every element can be distributed |
---|
488 | int nDst = dstView.size() ; |
---|
489 | vector<size_t> dstSliceSize(nDst) ; |
---|
490 | dstSliceSize[0] = 1 ; |
---|
491 | for(int i=1; i<nDst; i++) dstSliceSize[i] = dstView[i-1]->getGlobalSize()*dstSliceSize[i-1] ; |
---|
492 | |
---|
493 | CClientClientDHTTemplate<int>::Index2VectorInfoTypeMap dataInfo ; |
---|
494 | CClientClientDHTTemplate<size_t>::Index2VectorInfoTypeMap info ; // info map |
---|
495 | |
---|
496 | // first, we need to feed the DHT with the global index of the remote server |
---|
497 | // for that : |
---|
498 | // First the first element insert the in a DHT with key as the rank and value the list of global index associated |
---|
499 | // Then get the previously stored index associate with the remote rank I am in charge and reinsert the global index |
---|
500 | // corresponding to the position of the element in the remote view suing tensorial product |
---|
501 | // finaly we get only the list of remote global index I am in charge for the whole remote grid |
---|
502 | |
---|
503 | for(int pos=0; pos<nDst; pos++) |
---|
504 | { |
---|
505 | size_t sliceSize=dstSliceSize[pos] ; |
---|
506 | map<int,CArray<size_t,1>> globalIndexView ; |
---|
507 | dstView[pos]->getGlobalIndexView(globalIndexView) ; |
---|
508 | |
---|
509 | CClientClientDHTTemplate<size_t>::Index2VectorInfoTypeMap lastInfo(info) ; |
---|
510 | |
---|
511 | if (pos>0) |
---|
512 | { |
---|
513 | CArray<size_t,1> ranks(globalIndexView.size()) ; |
---|
514 | auto it=globalIndexView.begin() ; |
---|
515 | for(int i=0 ; i<ranks.numElements();i++,it++) ranks(i)=it->first ; |
---|
516 | CClientClientDHTTemplate<size_t> dataRanks(info, localComm_) ; |
---|
517 | dataRanks.computeIndexInfoMapping(ranks) ; |
---|
518 | lastInfo = dataRanks.getInfoIndexMap() ; |
---|
519 | } |
---|
520 | |
---|
521 | info.clear() ; |
---|
522 | for(auto& it : globalIndexView) |
---|
523 | { |
---|
524 | int rank = it.first ; |
---|
525 | auto& globalIndex = it.second ; |
---|
526 | auto& inf = info[rank] ; |
---|
527 | if (pos==0) for(int i=0;i<globalIndex.numElements();i++) inf.push_back(globalIndex(i)) ; |
---|
528 | else |
---|
529 | { |
---|
530 | auto& lastGlobalIndex = lastInfo[rank] ; |
---|
531 | for(size_t lastGlobalInd : lastGlobalIndex) |
---|
532 | { |
---|
533 | for(int i=0;i<globalIndex.numElements();i++) inf.push_back(globalIndex(i)*sliceSize+lastGlobalInd) ; |
---|
534 | } |
---|
535 | } |
---|
536 | } |
---|
537 | |
---|
538 | if (pos==nDst-1) |
---|
539 | { |
---|
540 | for(auto& it : info) |
---|
541 | { |
---|
542 | int rank=it.first ; |
---|
543 | auto& globalIndex = it.second ; |
---|
544 | for(auto globalInd : globalIndex) dataInfo[globalInd].push_back(rank) ; |
---|
545 | } |
---|
546 | } |
---|
547 | } |
---|
548 | |
---|
549 | // we feed the DHT with the remote global index |
---|
550 | CClientClientDHTTemplate<int> dataRanks(dataInfo, localComm_) ; |
---|
551 | |
---|
552 | // generate list of global index for src view |
---|
553 | int nSrc = srcView.size() ; |
---|
554 | vector<size_t> srcSliceSize(nSrc) ; |
---|
555 | |
---|
556 | srcSliceSize[0] = 1 ; |
---|
557 | for(int i=1; i<nSrc; i++) srcSliceSize[i] = srcView[i-1]->getGlobalSize()*srcSliceSize[i-1] ; |
---|
558 | |
---|
559 | vector<size_t> srcGlobalIndex ; |
---|
560 | size_t sliceIndex=0 ; |
---|
561 | srcView[nSrc-1]->getGlobalIndex(srcGlobalIndex, sliceIndex, srcSliceSize.data(), srcView.data(), nSrc-1) ; |
---|
562 | // now we have the global index of the source grid in srcGlobalIndex |
---|
563 | // we feed the DHT with the src global index (if we have) |
---|
564 | if (srcGlobalIndex.size()>0) |
---|
565 | { |
---|
566 | CArray<size_t,1> srcGlobalIndexArray(srcGlobalIndex.data(), shape(srcGlobalIndex.size()),neverDeleteData) ; |
---|
567 | dataRanks.computeIndexInfoMapping(srcGlobalIndexArray) ; |
---|
568 | } |
---|
569 | else |
---|
570 | { |
---|
571 | CArray<size_t,1> srcGlobalIndexArray ; |
---|
572 | dataRanks.computeIndexInfoMapping(srcGlobalIndexArray) ; |
---|
573 | } |
---|
574 | const auto& returnInfo = dataRanks.getInfoIndexMap() ; |
---|
575 | // returnInfo contains now the map for each global indices to send to a list of remote rank |
---|
576 | // but we want to use the tensorial product property to get the same information using only global |
---|
577 | // index of element view. So the idea is to reverse the information : for a global index of the grid |
---|
578 | // to send to the remote server, what is the global index of each element composing the grid ? |
---|
579 | |
---|
580 | vector<map<int, set<size_t>>> elements(nSrc) ; // internal representation of elements composing the grid |
---|
581 | |
---|
582 | for(auto& indRanks : returnInfo) |
---|
583 | { |
---|
584 | size_t gridIndexGlo=indRanks.first ; |
---|
585 | auto& ranks = indRanks.second ; |
---|
586 | for(int i=nSrc-1; i>=0; i--) |
---|
587 | { |
---|
588 | auto& element = elements[i] ; |
---|
589 | size_t localIndGlo = gridIndexGlo / srcSliceSize[i] ; |
---|
590 | gridIndexGlo = gridIndexGlo % srcSliceSize[i] ; |
---|
591 | for(int rank : ranks) element[rank].insert(localIndGlo) ; |
---|
592 | } |
---|
593 | } |
---|
594 | |
---|
595 | // elements_.resize(nSrc) ; |
---|
596 | for(int i=0 ; i<nSrc; i++) |
---|
597 | { |
---|
598 | auto& element=elements[i] ; |
---|
599 | for(auto& rankInd : element) |
---|
600 | { |
---|
601 | int rank=rankInd.first ; |
---|
602 | set<size_t>& indGlo = rankInd.second ; |
---|
603 | CArray<size_t,1>& indGloArray = elements_[indElements[i]][rank] ; |
---|
604 | indGloArray.resize(indGlo.size()) ; |
---|
605 | int j=0 ; |
---|
606 | for (auto index : indGlo) { indGloArray(j) = index ; j++; } |
---|
607 | } |
---|
608 | } |
---|
609 | |
---|
610 | // So what about when there is some server that have no data to receive |
---|
611 | // they must be inform they receive an event with no data. |
---|
612 | // So find remote servers with no data, and one client will take in charge |
---|
613 | // that it receive global index with no data (0-size) |
---|
614 | vector<int> ranks(remoteSize_,0) ; |
---|
615 | for(auto& it : elements_[indElements[0]]) ranks[it.first] = 1 ; |
---|
616 | MPI_Allreduce(MPI_IN_PLACE, ranks.data(), remoteSize_, MPI_INT, MPI_SUM, localComm_) ; |
---|
617 | int commRank, commSize ; |
---|
618 | MPI_Comm_rank(localComm_,&commRank) ; |
---|
619 | MPI_Comm_size(localComm_,&commSize) ; |
---|
620 | int pos=0 ; |
---|
621 | for(int i=0; i<remoteSize_ ; i++) |
---|
622 | if (ranks[i]==0) |
---|
623 | { |
---|
624 | if (pos%commSize==commRank) |
---|
625 | for(int j=0 ; j<nSrc; j++) elements_[indElements[j]][i] = CArray<size_t,1>(0) ; |
---|
626 | pos++ ; |
---|
627 | } |
---|
628 | } |
---|
629 | |
---|
630 | /** |
---|
631 | * \brief Once the connector is computed (compute \b elements_), redondant data can be avoid to be sent to the server. |
---|
632 | * This call compute the redondant rank and store them in \b rankToRemove_ attribute. |
---|
633 | * The goal of this method is to make a hash of each block of indice that determine wich data to send to a |
---|
634 | * of a specific server rank using a hash method. So data to send to a rank is associated to a hash. |
---|
635 | * After we compare hash between local rank and remove redondant data corresponding to the same hash. |
---|
636 | */ |
---|
637 | void CGridRemoteConnector::computeRedondantRanks(void) |
---|
638 | { |
---|
639 | int commRank ; |
---|
640 | MPI_Comm_rank(localComm_,&commRank) ; |
---|
641 | |
---|
642 | set<int> ranks; |
---|
643 | for(auto& element : elements_) |
---|
644 | for(auto& it : element) ranks.insert(it.first) ; |
---|
645 | |
---|
646 | for(auto& element : elements_) |
---|
647 | for(auto& it : element) |
---|
648 | if (ranks.count(it.first)==0) ERROR("void CGridRemoteConnector::removeRedondantRanks(void)",<<"number of ranks in elements is not coherent between each element") ; |
---|
649 | |
---|
650 | HashXIOS<size_t> hashGlobalIndex; |
---|
651 | |
---|
652 | map<int,size_t> hashRanks ; |
---|
653 | for(auto& element : elements_) |
---|
654 | for(auto& it : element) |
---|
655 | { |
---|
656 | auto& globalIndex=it.second ; |
---|
657 | int rank=it.first ; |
---|
658 | size_t hash ; |
---|
659 | hash=0 ; |
---|
660 | for(int i=0; i<globalIndex.numElements(); i++) hash+=hashGlobalIndex(globalIndex(i)) ; |
---|
661 | if (hashRanks.count(rank)==0) hashRanks[rank]=hash ; |
---|
662 | else hashRanks[rank]=hashGlobalIndex.hashCombine(hashRanks[rank],hash) ; |
---|
663 | } |
---|
664 | // a hash is now computed for data block I will sent to the server. |
---|
665 | |
---|
666 | CClientClientDHTTemplate<int>::Index2InfoTypeMap info ; |
---|
667 | |
---|
668 | map<size_t,int> hashRank ; |
---|
669 | HashXIOS<int> hashGlobalIndexRank; |
---|
670 | for(auto& it : hashRanks) |
---|
671 | { |
---|
672 | it.second = hashGlobalIndexRank.hashCombine(it.first,it.second) ; |
---|
673 | info[it.second]=commRank ; |
---|
674 | hashRank[it.second]=it.first ; |
---|
675 | } |
---|
676 | |
---|
677 | // we feed a DHT map with key : hash, value : myrank |
---|
678 | CClientClientDHTTemplate<int> dataHash(info, localComm_) ; |
---|
679 | CArray<size_t,1> hashList(hashRank.size()) ; |
---|
680 | |
---|
681 | int i=0 ; |
---|
682 | for(auto& it : hashRank) { hashList(i)=it.first ; i++; } |
---|
683 | |
---|
684 | // now who are the ranks that have the same hash : feed the DHT with my list of hash |
---|
685 | dataHash.computeIndexInfoMapping(hashList) ; |
---|
686 | auto& hashRankList = dataHash.getInfoIndexMap() ; |
---|
687 | |
---|
688 | |
---|
689 | for(auto& it : hashRankList) |
---|
690 | { |
---|
691 | size_t hash = it.first ; |
---|
692 | auto& ranks = it.second ; |
---|
693 | |
---|
694 | bool first=true ; |
---|
695 | // only the process with the lowest rank get in charge of sendinf data to remote server |
---|
696 | for(int rank : ranks) if (commRank>rank) first=false ; |
---|
697 | if (!first) rankToRemove_.insert(hashRank[hash]) ; |
---|
698 | } |
---|
699 | } |
---|
700 | |
---|
701 | } |
---|