Context Navigation

source: XIOS/dev/XIOS_DEV_CMIP6/src/transformation/generic_algorithm_transformation.cpp @ 1221

Last change on this file since 1221 was 1221, checked in by mhnguyen, 7 years ago

Fixing a minor bug in transformation: No source index at all for transformation

Test
+) On Curie
+) OK

File size: 27.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
16	namespace xios {
17
18	CGenericAlgorithmTransformation::CGenericAlgorithmTransformation()
19	: transformationMapping_(), transformationWeight_(), transformationPosition_(),
20	idAuxInputs_(), type_(ELEMENT_NO_MODIFICATION_WITH_DATA), indexElementSrc_(),
21	computedProcSrcNonTransformedElement_(false)
22	{
23	}
24
25	void CGenericAlgorithmTransformation::updateData(CArray<double,1>& dataOut)
26	{
27
28	}
29
30	void CGenericAlgorithmTransformation::apply(const std::vector<std::pair<int,double> >& localIndex,
31	const double* dataInput,
32	CArray<double,1>& dataOut,
33	std::vector<bool>& flagInitial,
34	bool ignoreMissingValue)
35	{
36	int nbLocalIndex = localIndex.size();
37	double defaultValue = std::numeric_limits<double>::quiet_NaN();
38	if (ignoreMissingValue)
39	{
40	for (int idx = 0; idx < nbLocalIndex; ++idx)
41	{
42	if (NumTraits<double>::isnan(*(dataInput + idx)))
43	{
44	flagInitial[localIndex[idx].first] = false;
45	}
46	else
47	{
48	dataOut(localIndex[idx].first) += (dataInput + idx) localIndex[idx].second;
49	flagInitial[localIndex[idx].first] = true; // Reset flag to indicate not all data source are nan
50	}
51	}
52
53	// If all data source are nan then data destination must be nan
54	for (int idx = 0; idx < nbLocalIndex; ++idx)
55	{
56	if (!flagInitial[localIndex[idx].first])
57	dataOut(localIndex[idx].first) = defaultValue;
58	}
59	}
60	else
61	{
62	for (int idx = 0; idx < nbLocalIndex; ++idx)
63	{
64	dataOut(localIndex[idx].first) += (dataInput + idx) localIndex[idx].second;
65	}
66	}
67	}
68
69	void CGenericAlgorithmTransformation::computePositionElements(CGrid* dst, CGrid* src)
70	{
71	int idxScalar = 0, idxAxis = 0, idxDomain = 0;
72	CArray<int,1> axisDomainOrderDst = dst->axis_domain_order;
73	for (int i = 0; i < axisDomainOrderDst.numElements(); ++i)
74	{
75	int dimElement = axisDomainOrderDst(i);
76	if (2 == dimElement)
77	{
78	elementPositionInGridDst2DomainPosition_[i] = idxDomain;
79	++idxDomain;
80	}
81	else if (1 == dimElement)
82	{
83	elementPositionInGridDst2AxisPosition_[i] = idxAxis;
84	++idxAxis;
85	}
86	else
87	{
88	elementPositionInGridDst2ScalarPosition_[i] = idxScalar;
89	++idxScalar;
90	}
91	}
92
93	idxScalar = idxAxis = idxDomain = 0;
94	CArray<int,1> axisDomainOrderSrc = src->axis_domain_order;
95	for (int i = 0; i < axisDomainOrderSrc.numElements(); ++i)
96	{
97	int dimElement = axisDomainOrderSrc(i);
98	if (2 == dimElement)
99	{
100	elementPositionInGridSrc2DomainPosition_[i] = idxDomain;
101	++idxDomain;
102	}
103	else if (1 == dimElement)
104	{
105	elementPositionInGridSrc2AxisPosition_[i] = idxAxis;
106	++idxAxis;
107	}
108	else
109	{
110	elementPositionInGridSrc2ScalarPosition_[i] = idxScalar;
111	++idxScalar;
112	}
113	}
114	}
115
116	/*!
117	This function computes the global indexes of grid source, which the grid destination is in demand.
118	\param[in] elementPositionInGrid position of an element in a grid .E.g: if grid is composed of domain and axis (in order),
119	then position of axis in grid is 1 and domain is positioned at 0.
120	\param[in] gridSrc Grid source
121	\param[in] gridDst Grid destination
122	\param[in\out] globaIndexWeightFromSrcToDst mapping of each global index source and weight to index destination
123	*/
124	void CGenericAlgorithmTransformation::computeGlobalSourceIndex(int elementPositionInGrid,
125	CGrid* gridSrc,
126	CGrid* gridDst,
127	SourceDestinationIndexMap& globaIndexWeightFromSrcToDst)
128	{
129	CContext* context = CContext::getCurrent();
130	CContextClient* client = context->client;
131	int nbClient = client->clientSize;
132
133	typedef boost::unordered_map<int, std::vector<std::pair<int,double> > > SrcToDstMap;
134	int idx;
135
136	// compute position of elements on grids
137	computePositionElements(gridDst, gridSrc);
138	std::vector<CScalar*> scalarListDestP = gridDst->getScalars();
139	std::vector<CAxis*> axisListDestP = gridDst->getAxis();
140	std::vector<CDomain*> domainListDestP = gridDst->getDomains();
141	CArray<int,1> axisDomainDstOrder = gridDst->axis_domain_order;
142	std::vector<CScalar*> scalarListSrcP = gridSrc->getScalars();
143	std::vector<CAxis*> axisListSrcP = gridSrc->getAxis();
144	std::vector<CDomain*> domainListSrcP = gridSrc->getDomains();
145	CArray<int,1> axisDomainSrcOrder = gridSrc->axis_domain_order;
146
147	bool isTransPosEmpty = transformationPosition_.empty();
148	CArray<size_t,1> transPos;
149	if (!isTransPosEmpty) transPos.resize(transformationMapping_.size());
150	std::set<size_t> allIndexSrc; // All index of source, which can be scattered among processes, need for doing transformation
151
152	for (size_t idxTrans = 0; idxTrans < transformationMapping_.size(); ++idxTrans)
153	{
154	TransformationIndexMap::const_iterator itbTransMap = transformationMapping_[idxTrans].begin(), itTransMap,
155	iteTransMap = transformationMapping_[idxTrans].end();
156	TransformationWeightMap::const_iterator itbTransWeight = transformationWeight_[idxTrans].begin(), itTransWeight;
157
158	// Build mapping between global source element index and global destination element index.
159	itTransWeight = itbTransWeight;
160	for (itTransMap = itbTransMap; itTransMap != iteTransMap; ++itTransMap, ++itTransWeight)
161	{
162	const std::vector<int>& srcIndex = itTransMap->second;
163	for (idx = 0; idx < srcIndex.size(); ++idx)
164	allIndexSrc.insert(srcIndex[idx]);
165	}
166
167	if (!isTransPosEmpty)
168	{
169	TransformationPositionMap::const_iterator itPosMap = transformationPosition_[idxTrans].begin();
170	transPos(idxTrans) = itPosMap->second[0];
171	}
172	}
173
174	size_t indexSrcSize = 0;
175	CArray<size_t,1> indexSrc(allIndexSrc.size());
176	for (std::set<size_t>::iterator it = allIndexSrc.begin(); it != allIndexSrc.end(); ++it, ++indexSrcSize)
177	indexSrc(indexSrcSize) = *it;
178
179	// Flag to indicate whether we will recompute distribution of source global index on processes
180	bool computeGlobalIndexOnProc = false;
181	if (indexElementSrc_.size() != allIndexSrc.size())
182	computeGlobalIndexOnProc = true;
183	else
184	{
185	for (std::set<size_t>::iterator it = allIndexSrc.begin(); it != allIndexSrc.end(); ++it)
186	if (0 == indexElementSrc_.count(*it))
187	{
188	computeGlobalIndexOnProc = true;
189	break;
190	}
191	}
192
193	if (computeGlobalIndexOnProc)
194	indexElementSrc_.swap(allIndexSrc);
195
196	int sendValue = (computeGlobalIndexOnProc) ? 1 : 0;
197	int recvValue = 0;
198	MPI_Allreduce(&sendValue, &recvValue, 1, MPI_INT, MPI_SUM, client->intraComm);
199	computeGlobalIndexOnProc = (0 < recvValue);
200
201	if (computeGlobalIndexOnProc \|\| !computedProcSrcNonTransformedElement_)
202	{
203	// Find out global index source of transformed element on corresponding process.
204	if (globalElementIndexOnProc_.empty())
205	globalElementIndexOnProc_.resize(axisDomainDstOrder.numElements());
206	CClientClientDHTInt::Index2VectorInfoTypeMap globalIndexOfTransformedElementOnProc;
207	for (idx = 0; idx < axisDomainDstOrder.numElements(); ++idx)
208	{
209	if (idx == elementPositionInGrid)
210	computeExchangeGlobalIndex(indexSrc, axisDomainSrcOrder(idx), globalIndexOfTransformedElementOnProc); //globalElementIndexOnProc[idx]);
211	if (!computedProcSrcNonTransformedElement_)
212	{
213	if (2 == axisDomainDstOrder(idx)) // It's domain
214	{
215	if (idx != elementPositionInGrid)
216	computeExchangeDomainIndex(domainListDestP[elementPositionInGridDst2DomainPosition_[idx]],
217	domainListSrcP[elementPositionInGridSrc2DomainPosition_[idx]],
218	transPos,
219	globalElementIndexOnProc_[idx]);
220
221	}
222	else if (1 == axisDomainDstOrder(idx))//it's an axis
223	{
224	if (idx != elementPositionInGrid)
225	computeExchangeAxisIndex(axisListDestP[elementPositionInGridDst2AxisPosition_[idx]],
226	axisListSrcP[elementPositionInGridSrc2AxisPosition_[idx]],
227	transPos,
228	globalElementIndexOnProc_[idx]);
229	}
230	else //it's a scalar
231	{
232	if (idx != elementPositionInGrid)
233	computeExchangeScalarIndex(scalarListDestP[elementPositionInGridDst2ScalarPosition_[idx]],
234	scalarListSrcP[elementPositionInGridSrc2ScalarPosition_[idx]],
235	transPos,
236	globalElementIndexOnProc_[idx]);
237
238	}
239	}
240	}
241
242	if (!isTransPosEmpty && !computedProcSrcNonTransformedElement_)
243	{
244	for (idx = 0; idx < globalElementIndexOnProc_.size(); ++idx)
245	{
246	if (idx != elementPositionInGrid)
247	{
248	boost::unordered_map<int,std::vector<size_t> >::iterator itb = globalElementIndexOnProc_[idx].begin(), it,
249	ite = globalElementIndexOnProc_[idx].end();
250	for (it = itb; it != ite; ++it) it->second.resize(1);
251	}
252	}
253	}
254
255	if (!computedProcSrcNonTransformedElement_)
256	{
257	for (idx = 0; idx < globalElementIndexOnProc_.size(); ++idx)
258	{
259	if (idx != elementPositionInGrid)
260	{
261	boost::unordered_map<int,std::vector<size_t> >::iterator itb = globalElementIndexOnProc_[idx].begin(), it,
262	ite = globalElementIndexOnProc_[idx].end();
263	for (it = itb; it != ite; ++it) procOfNonTransformedElements_.insert(it->first);
264	if (procOfNonTransformedElements_.size() == nbClient)
265	break;
266	}
267	}
268	}
269
270	// Processes contain the source index of transformed element
271	std::set<int> procOfTransformedElement;
272	CClientClientDHTInt::Index2VectorInfoTypeMap::iterator itIdxb = globalIndexOfTransformedElementOnProc.begin(),
273	itIdxe = globalIndexOfTransformedElementOnProc.end(), itIdx;
274	for (itIdx = itIdxb; itIdx != itIdxe; ++itIdx)
275	{
276	std::vector<int>& tmp = itIdx->second;
277	for (int i = 0; i < tmp.size(); ++i)
278	procOfTransformedElement.insert(tmp[i]);
279	if (tmp.size() == nbClient)
280	break;
281	}
282
283	std::set<int>& commonProc = (procOfTransformedElement.size() < procOfNonTransformedElements_.size()) ? procOfTransformedElement
284	: (!procOfNonTransformedElements_.empty() ? procOfNonTransformedElements_ : procOfTransformedElement);
285
286	std::vector<int> procContainSrcElementIdx(commonProc.size());
287	int count = 0;
288	for (std::set<int>::iterator it = commonProc.begin(); it != commonProc.end(); ++it)
289	procContainSrcElementIdx[count++] = *it;
290
291	procContainSrcElementIdx_.swap(procContainSrcElementIdx);
292
293	// For the first time, surely we calculate proc containing non transformed elements
294	if (!computedProcSrcNonTransformedElement_)
295	computedProcSrcNonTransformedElement_ = true;
296	}
297
298	for (size_t idxTrans = 0; idxTrans < transformationMapping_.size(); ++idxTrans)
299	{
300	TransformationIndexMap::const_iterator itbTransMap = transformationMapping_[idxTrans].begin(), itTransMap,
301	iteTransMap = transformationMapping_[idxTrans].end();
302	TransformationWeightMap::const_iterator itbTransWeight = transformationWeight_[idxTrans].begin(), itTransWeight;
303	SrcToDstMap src2DstMap;
304	src2DstMap.rehash(std::ceil(transformationMapping_[idxTrans].size()/src2DstMap.max_load_factor()));
305
306	// Build mapping between global source element index and global destination element index.
307	boost::unordered_map<int,std::vector<size_t> >().swap(globalElementIndexOnProc_[elementPositionInGrid]);
308	std::set<int> tmpCounter;
309	itTransWeight = itbTransWeight;
310	for (itTransMap = itbTransMap; itTransMap != iteTransMap; ++itTransMap, ++itTransWeight)
311	{
312	const std::vector<int>& srcIndex = itTransMap->second;
313	const std::vector<double>& weight = itTransWeight->second;
314	for (idx = 0; idx < srcIndex.size(); ++idx)
315	{
316	src2DstMap[srcIndex[idx]].push_back(make_pair(itTransMap->first, weight[idx]));
317	if (0 == tmpCounter.count(srcIndex[idx]))
318	{
319	tmpCounter.insert(srcIndex[idx]);
320	for (int j = 0; j < procContainSrcElementIdx_.size(); ++j)
321	globalElementIndexOnProc_[elementPositionInGrid][procContainSrcElementIdx_[j]].push_back(srcIndex[idx]);
322	}
323	}
324	}
325
326	if (!isTransPosEmpty)
327	{
328	for (idx = 0; idx < globalElementIndexOnProc_.size(); ++idx)
329	{
330	if (idx != elementPositionInGrid)
331	{
332	boost::unordered_map<int,std::vector<size_t> >::iterator itb = globalElementIndexOnProc_[idx].begin(), it,
333	ite = globalElementIndexOnProc_[idx].end();
334	for (it = itb; it != ite; ++it) it->second[0] = transPos(idxTrans);
335	}
336	}
337	}
338
339	// Ok, now compute global index of grid source and ones of grid destination
340	computeGlobalGridIndexMapping(elementPositionInGrid,
341	procContainSrcElementIdx_, //srcRank,
342	src2DstMap,
343	gridSrc,
344	gridDst,
345	globalElementIndexOnProc_,
346	globaIndexWeightFromSrcToDst);
347	}
348	}
349
350	/*!
351	Compute mapping of global index of grid source and grid destination
352	\param [in] elementPositionInGrid position of element in grid. E.x: grid composed of domain and axis, domain has position 0 and axis 1.
353	\param [in] srcRank rank of client from which we demand global index of element source
354	\param [in] src2DstMap mapping of global index of element source and global index of element destination
355	\param[in] gridSrc Grid source
356	\param[in] gridDst Grid destination
357	\param[in] globalElementIndexOnProc Global index of element source on different client rank
358	\param[out] globaIndexWeightFromSrcToDst Mapping of global index of grid source and grid destination
359	*/
360	void CGenericAlgorithmTransformation::computeGlobalGridIndexMapping(int elementPositionInGrid,
361	const std::vector<int>& srcRank,
362	boost::unordered_map<int, std::vector<std::pair<int,double> > >& src2DstMap,
363	CGrid* gridSrc,
364	CGrid* gridDst,
365	std::vector<boost::unordered_map<int,std::vector<size_t> > >& globalElementIndexOnProc,
366	SourceDestinationIndexMap& globaIndexWeightFromSrcToDst)
367	{
368	std::vector<CDomain*> domainListSrcP = gridSrc->getDomains();
369	std::vector<CAxis*> axisListSrcP = gridSrc->getAxis();
370	std::vector<CScalar*> scalarListSrcP = gridSrc->getScalars();
371	CArray<int,1> axisDomainSrcOrder = gridSrc->axis_domain_order;
372
373	size_t nbElement = axisDomainSrcOrder.numElements();
374	std::vector<size_t> nGlobSrc(nbElement);
375	size_t globalSrcSize = 1;
376	int domainIndex = 0, axisIndex = 0, scalarIndex = 0;
377	for (int idx = 0; idx < nbElement; ++idx)
378	{
379	nGlobSrc[idx] = globalSrcSize;
380	int elementDimension = axisDomainSrcOrder(idx);
381
382	// If this is a domain
383	if (2 == elementDimension)
384	{
385	globalSrcSize = domainListSrcP[domainIndex]->nj_glo.getValue() domainListSrcP[domainIndex]->ni_glo.getValue();
386	++domainIndex;
387	}
388	else if (1 == elementDimension) // So it's an axis
389	{
390	globalSrcSize *= axisListSrcP[axisIndex]->n_glo.getValue();
391	++axisIndex;
392	}
393	else
394	{
395	globalSrcSize *= 1;
396	++scalarIndex;
397	}
398	}
399
400	std::vector<CDomain*> domainListDestP = gridDst->getDomains();
401	std::vector<CAxis*> axisListDestP = gridDst->getAxis();
402	std::vector<CScalar*> scalarListDestP = gridDst->getScalars();
403	CArray<int,1> axisDomainDstOrder = gridDst->axis_domain_order;
404
405	std::vector<size_t> nGlobDst(nbElement);
406	size_t globalDstSize = 1;
407	domainIndex = axisIndex = scalarIndex = 0;
408	for (int idx = 0; idx < nbElement; ++idx)
409	{
410	nGlobDst[idx] = globalDstSize;
411	int elementDimension = axisDomainDstOrder(idx);
412
413	// If this is a domain
414	if (2 == elementDimension)
415	{
416	globalDstSize = domainListDestP[domainIndex]->nj_glo.getValue() domainListDestP[domainIndex]->ni_glo.getValue();
417	++domainIndex;
418	}
419	else if (1 == elementDimension) // So it's an axis
420	{
421	globalDstSize *= axisListDestP[axisIndex]->n_glo.getValue();
422	++axisIndex;
423	}
424	else
425	{
426	globalDstSize *= 1;
427	++scalarIndex;
428	}
429	}
430
431	for (int i = 0; i < srcRank.size(); ++i)
432	{
433	size_t ssize = 1;
434	int rankSrc = srcRank[i];
435	for (int idx = 0; idx < nbElement; ++idx)
436	{
437	ssize *= (globalElementIndexOnProc[idx][rankSrc]).size();
438	}
439
440	std::vector<int> idxLoop(nbElement,0);
441	std::vector<int> currentIndexSrc(nbElement, 0);
442	std::vector<int> currentIndexDst(nbElement, 0);
443	int innnerLoopSize = (globalElementIndexOnProc[0])[rankSrc].size();
444	size_t idx = 0;
445	while (idx < ssize)
446	{
447	for (int ind = 0; ind < nbElement; ++ind)
448	{
449	if (idxLoop[ind] == (globalElementIndexOnProc[ind])[rankSrc].size())
450	{
451	idxLoop[ind] = 0;
452	++idxLoop[ind+1];
453	}
454
455	currentIndexDst[ind] = currentIndexSrc[ind] = (globalElementIndexOnProc[ind])[rankSrc][idxLoop[ind]];
456	}
457
458	for (int ind = 0; ind < innnerLoopSize; ++ind)
459	{
460	currentIndexDst[0] = currentIndexSrc[0] = (globalElementIndexOnProc[0])[rankSrc][ind];
461	int globalElementDstIndexSize = 0;
462	if (1 == src2DstMap.count(currentIndexSrc[elementPositionInGrid]))
463	{
464	globalElementDstIndexSize = src2DstMap[currentIndexSrc[elementPositionInGrid]].size();
465	}
466
467	std::vector<size_t> globalDstVecIndex(globalElementDstIndexSize,0);
468	size_t globalSrcIndex = 0;
469	for (int idxElement = 0; idxElement < nbElement; ++idxElement)
470	{
471	if (idxElement == elementPositionInGrid)
472	{
473	for (int k = 0; k < globalElementDstIndexSize; ++k)
474	{
475	globalDstVecIndex[k] += src2DstMap[currentIndexSrc[elementPositionInGrid]][k].first * nGlobDst[idxElement];
476	}
477	}
478	else
479	{
480	for (int k = 0; k < globalElementDstIndexSize; ++k)
481	{
482	globalDstVecIndex[k] += currentIndexDst[idxElement] * nGlobDst[idxElement];
483	}
484	}
485	globalSrcIndex += currentIndexSrc[idxElement] * nGlobSrc[idxElement];
486	}
487
488	for (int k = 0; k < globalElementDstIndexSize; ++k)
489	{
490	globaIndexWeightFromSrcToDst[rankSrc][globalSrcIndex].push_back(make_pair(globalDstVecIndex[k],src2DstMap[currentIndexSrc[elementPositionInGrid]][k].second ));
491	}
492	++idxLoop[0];
493	}
494	idx += innnerLoopSize;
495	}
496	}
497	}
498
499	/*!
500	Find out proc and global index of axis source which axis destination is on demande
501	\param[in] scalar Scalar destination
502	\param[in] scalar Scalar source
503	\param[in] destGlobalIndexPositionInGrid Relative position of axis corresponds to other element of grid.
504	\param[out] globalScalarIndexOnProc Global index of axis source on different procs
505	*/
506	void CGenericAlgorithmTransformation::computeExchangeScalarIndex(CScalar* scalarDst,
507	CScalar* scalarSrc,
508	CArray<size_t,1>& destGlobalIndexPositionInGrid,
509	boost::unordered_map<int,std::vector<size_t> >& globalScalarIndexOnProc)
510	{
511	CContext* context = CContext::getCurrent();
512	CContextClient* client=context->client;
513	int clientRank = client->clientRank;
514	int clientSize = client->clientSize;
515
516	globalScalarIndexOnProc.rehash(std::ceil(clientSize/globalScalarIndexOnProc.max_load_factor()));
517	for (int idx = 0; idx < clientSize; ++idx)
518	{
519	globalScalarIndexOnProc[idx].push_back(0);
520	}
521	}
522
523	/*!
524	Find out proc and global index of axis source which axis destination is on demande
525	\param[in] axisDst Axis destination
526	\param[in] axisSrc Axis source
527	\param[in] destGlobalIndexPositionInGrid Relative position of axis corresponds to other element of grid.
528	\param[out] globalAxisIndexOnProc Global index of axis source on different procs
529	*/
530	void CGenericAlgorithmTransformation::computeExchangeAxisIndex(CAxis* axisDst,
531	CAxis* axisSrc,
532	CArray<size_t,1>& destGlobalIndexPositionInGrid,
533	boost::unordered_map<int,std::vector<size_t> >& globalAxisIndexOnProc)
534	{
535	CContext* context = CContext::getCurrent();
536	CContextClient* client=context->client;
537	int clientRank = client->clientRank;
538	int clientSize = client->clientSize;
539
540	size_t globalIndex;
541	int nIndexSize = axisSrc->index.numElements();
542	CClientClientDHTInt::Index2VectorInfoTypeMap globalIndex2ProcRank;
543	globalIndex2ProcRank.rehash(std::ceil(nIndexSize/globalIndex2ProcRank.max_load_factor()));
544	for (int idx = 0; idx < nIndexSize; ++idx)
545	{
546	globalIndex = axisSrc->index(idx);
547	globalIndex2ProcRank[globalIndex].push_back(clientRank);
548	}
549
550	CClientClientDHTInt dhtIndexProcRank(globalIndex2ProcRank, client->intraComm);
551	CArray<size_t,1> globalAxisIndex(axisDst->index.numElements());
552	for (int idx = 0; idx < globalAxisIndex.numElements(); ++idx)
553	{
554	globalAxisIndex(idx) = axisDst->index(idx);
555	}
556	dhtIndexProcRank.computeIndexInfoMapping(globalAxisIndex);
557
558	std::vector<int> countIndex(clientSize,0);
559	const CClientClientDHTInt::Index2VectorInfoTypeMap& computedGlobalIndexOnProc = dhtIndexProcRank.getInfoIndexMap();
560	CClientClientDHTInt::Index2VectorInfoTypeMap::const_iterator itb = computedGlobalIndexOnProc.begin(), it,
561	ite = computedGlobalIndexOnProc.end();
562	for (it = itb; it != ite; ++it)
563	{
564	const std::vector<int>& procList = it->second;
565	for (int idx = 0; idx < procList.size(); ++idx) ++countIndex[procList[idx]];
566	}
567
568	globalAxisIndexOnProc.rehash(std::ceil(clientSize/globalAxisIndexOnProc.max_load_factor()));
569	for (int idx = 0; idx < clientSize; ++idx)
570	{
571	if (0 != countIndex[idx])
572	{
573	globalAxisIndexOnProc[idx].resize(countIndex[idx]);
574	countIndex[idx] = 0;
575	}
576	}
577
578	for (it = itb; it != ite; ++it)
579	{
580	const std::vector<int>& procList = it->second;
581	for (int idx = 0; idx < procList.size(); ++idx)
582	{
583	globalAxisIndexOnProc[procList[idx]][countIndex[procList[idx]]] = it->first;
584	++countIndex[procList[idx]];
585	}
586	}
587	}
588
589	/*!
590	Find out proc and global index of domain source which domain destination is on demande
591	\param[in] domainDst Domain destination
592	\param[in] domainSrc Domain source
593	\param[in] destGlobalIndexPositionInGrid Relative position of domain corresponds to other element of grid.
594	\param[out] globalDomainIndexOnProc Global index of domain source on different procs
595	*/
596	void CGenericAlgorithmTransformation::computeExchangeDomainIndex(CDomain* domainDst,
597	CDomain* domainSrc,
598	CArray<size_t,1>& destGlobalIndexPositionInGrid,
599	boost::unordered_map<int,std::vector<size_t> >& globalDomainIndexOnProc)
600	{
601	CContext* context = CContext::getCurrent();
602	CContextClient* client=context->client;
603	int clientRank = client->clientRank;
604	int clientSize = client->clientSize;
605
606	int niGlobSrc = domainSrc->ni_glo.getValue();
607	size_t globalIndex;
608	int i_ind, j_ind;
609	int nIndexSize = (destGlobalIndexPositionInGrid.isEmpty()) ? domainSrc->i_index.numElements()
610	: destGlobalIndexPositionInGrid.numElements();
611	CClientClientDHTInt::Index2VectorInfoTypeMap globalIndex2ProcRank;
612	globalIndex2ProcRank.rehash(std::ceil(nIndexSize/globalIndex2ProcRank.max_load_factor()));
613	if (destGlobalIndexPositionInGrid.isEmpty())
614	{
615	for (int idx = 0; idx < nIndexSize; ++idx)
616	{
617	i_ind=domainSrc->i_index(idx) ;
618	j_ind=domainSrc->j_index(idx) ;
619
620	globalIndex = i_ind + j_ind * niGlobSrc;
621	globalIndex2ProcRank[globalIndex].resize(1);
622	globalIndex2ProcRank[globalIndex][0] = clientRank;
623	}
624	}
625	else
626	{
627	for (int idx = 0; idx < nIndexSize; ++idx)
628	{
629	globalIndex2ProcRank[destGlobalIndexPositionInGrid(idx)].push_back(clientRank);
630	}
631	}
632
633	CArray<size_t,1> globalDomainIndex;
634	if (destGlobalIndexPositionInGrid.isEmpty())
635	{
636	int niGlobDst = domainDst->ni_glo.getValue();
637	globalDomainIndex.resize(domainDst->i_index.numElements());
638	nIndexSize = domainDst->i_index.numElements();
639
640	for (int idx = 0; idx < nIndexSize; ++idx)
641	{
642	i_ind=domainDst->i_index(idx) ;
643	j_ind=domainDst->j_index(idx) ;
644
645	globalDomainIndex(idx) = i_ind + j_ind * niGlobDst;
646	}
647	}
648	else
649	{
650	globalDomainIndex.reference(destGlobalIndexPositionInGrid);
651	}
652
653	CClientClientDHTInt dhtIndexProcRank(globalIndex2ProcRank, client->intraComm);
654	dhtIndexProcRank.computeIndexInfoMapping(globalDomainIndex);
655
656	std::vector<int> countIndex(clientSize,0);
657	const CClientClientDHTInt::Index2VectorInfoTypeMap& computedGlobalIndexOnProc = dhtIndexProcRank.getInfoIndexMap();
658	CClientClientDHTInt::Index2VectorInfoTypeMap::const_iterator itb = computedGlobalIndexOnProc.begin(), it,
659	ite = computedGlobalIndexOnProc.end();
660	for (it = itb; it != ite; ++it)
661	{
662	const std::vector<int>& procList = it->second;
663	for (int idx = 0; idx < procList.size(); ++idx) ++countIndex[procList[idx]];
664	}
665
666	globalDomainIndexOnProc.rehash(std::ceil(clientSize/globalDomainIndexOnProc.max_load_factor()));
667	for (int idx = 0; idx < clientSize; ++idx)
668	{
669	if (0 != countIndex[idx])
670	{
671	globalDomainIndexOnProc[idx].resize(countIndex[idx]);
672	countIndex[idx] = 0;
673	}
674	}
675
676	for (it = itb; it != ite; ++it)
677	{
678	const std::vector<int>& procList = it->second;
679	for (int idx = 0; idx < procList.size(); ++idx)
680	{
681	globalDomainIndexOnProc[procList[idx]][countIndex[procList[idx]]] = it->first;
682	++countIndex[procList[idx]];
683	}
684	}
685	}
686
687	/*!
688	Compute index mapping between element source and element destination with an auxiliary inputs which determine
689	position of each mapped index in global index of grid destination.
690	\param [in] dataAuxInputs auxiliary inputs
691	*/
692	void CGenericAlgorithmTransformation::computeIndexSourceMapping(const std::vector<CArray<double,1>* >& dataAuxInputs)
693	{
694	computeIndexSourceMapping_(dataAuxInputs);
695	}
696
697	std::vector<StdString> CGenericAlgorithmTransformation::getIdAuxInputs()
698	{
699	return idAuxInputs_;
700	}
701
702	CGenericAlgorithmTransformation::AlgoTransType CGenericAlgorithmTransformation::type()
703	{
704	return type_;
705	}
706
707	}

Note: See TracBrowser for help on using the repository browser.

Download in other formats: