Finite Element Domain Decomposition Library
FEDDLib
Loading...
Searching...
No Matches
Map_def.hpp
1#ifndef MAP_DEF_hpp
2#define MAP_DEF_hpp
3
10
11namespace FEDD {
12template < class LO, class GO, class NO>
13Map<LO,GO,NO>::Map():
14map_()
15{
16
17}
18
19template <class LO, class GO, class NO>
20Map<LO,GO,NO>::Map( TpetraMapConstPtrConst_Type& tpetraMapPtrIn ):
21map_( tpetraMapPtrIn )
22{
23
24}
25
26template < class LO, class GO, class NO>
27Map<LO,GO,NO>::Map( const Map_Type& mapIn ):
28map_()
29{
30 map_ = Teuchos::RCP(new TpetraMap_Type(mapIn.getGlobalNumElements(), mapIn.getNodeElementList(),mapIn.getIndexBase(), mapIn.getComm()));
31
32}
33
34template < class LO, class GO, class NO>
35Map<LO,GO,NO>::Map(GO numGlobalElements,
36 const Teuchos::ArrayView<const GO> &elementList,
37 GO indexBase,
38 const CommConstPtr_Type &comm):
39map_()
40{
41 map_ = Teuchos::RCP(new TpetraMap_Type(numGlobalElements, elementList, indexBase, comm));
42}
43
44template < class LO, class GO, class NO>
45Map<LO,GO,NO>::Map(GO numGlobalElements,
46 LO numLocalElements,
47 GO indexBase,
48 const CommConstPtr_Type &comm):
49map_()
50{
51 map_ = Teuchos::RCP(new TpetraMap_Type(numGlobalElements, numLocalElements, indexBase, comm));
52
53}
54
55template < class LO, class GO, class NO>
56Map<LO,GO,NO>::~Map(){
57
58}
59
60
61template < class LO, class GO, class NO>
62GO Map<LO,GO,NO>::getGlobalElement(LO id) const{
63 TEUCHOS_TEST_FOR_EXCEPTION(map_.is_null(),std::runtime_error,"map is null.");
64 return map_->getGlobalElement(id);
65}
66
67template < class LO, class GO, class NO>
68LO Map<LO,GO,NO>::getLocalElement(GO id) const{
69 TEUCHOS_TEST_FOR_EXCEPTION(map_.is_null(),std::runtime_error,"map is null.");
70 return map_->getLocalElement(id);
71}
72
73template < class LO, class GO, class NO>
74LO Map<LO,GO,NO>::getNodeNumElements() const{
75 TEUCHOS_TEST_FOR_EXCEPTION(map_.is_null(),std::runtime_error,"map is null.");
76 return map_->getLocalNumElements();
77}
78
79template < class LO, class GO, class NO>
80GO Map<LO,GO,NO>::getGlobalNumElements() const{
81 TEUCHOS_TEST_FOR_EXCEPTION(map_.is_null(),std::runtime_error,"map is null.");
82 return map_->getGlobalNumElements();
83}
84
85template < class LO, class GO, class NO>
86typename Map<LO,GO,NO>::CommConstPtr_Type Map<LO,GO,NO>::getComm() const{
87
88 return map_->getComm();
89}
90
91template <class LO,class GO,class NO>
92GO Map<LO,GO,NO>::getIndexBase() const{
93 TEUCHOS_TEST_FOR_EXCEPTION(map_.is_null(),std::runtime_error,"map is null.");
94 return map_->getIndexBase();
95}
96
97template < class LO, class GO, class NO>
98Teuchos::ArrayView< const GO > Map<LO,GO,NO>::getNodeElementList() const{
99 TEUCHOS_TEST_FOR_EXCEPTION(map_.is_null(),std::runtime_error,"map is null.");
100 return map_->getLocalElementList();
101}
102
103template < class LO, class GO, class NO>
104typename Map<LO,GO,NO>::MapPtr_Type Map<LO,GO,NO>::buildVecFieldMap(UN numDofs, std::string ordering) const{
105
106 TEUCHOS_TEST_FOR_EXCEPTION(map_.is_null(),std::runtime_error,"map is null.");
107 TEUCHOS_TEST_FOR_EXCEPTION(ordering.compare("NodeWise"), std::logic_error,"Select a valid ordering: NodeWise");
108 Teuchos::ArrayView<const GO> elementList = map_->getLocalElementList();
109 Teuchos::Array<GO> elementListField( numDofs * elementList.size() );
110 for (UN i=0; i<elementList.size(); i++) {
111 for (UN dof=0; dof<numDofs; dof++)
112 elementListField[ numDofs * i + dof ] = numDofs * elementList[i] + dof;
113 }
114 typedef Teuchos::OrdinalTraits<GO> GOOT;
115 return Teuchos::rcp(new Map_Type(GOOT::invalid(), elementListField(), map_->getIndexBase(), map_->getComm() ) );
116}
117
118
119template < class LO, class GO, class NO>
120typename Map<LO,GO,NO>::TpetraMapConstPtr_Type Map<LO,GO,NO>::getTpetraMap() const{
121
122 TEUCHOS_TEST_FOR_EXCEPTION(map_.is_null(),std::runtime_error,"getTpetraMap(): map_ is null.");
123
124 return map_;
125}
126
127template < class LO, class GO, class NO>
128GO Map<LO,GO,NO>::getMaxAllGlobalIndex() const{
129 return map_->getMaxAllGlobalIndex();
130}
131
132template < class LO, class GO, class NO>
133LO Map<LO,GO,NO>::getMaxLocalIndex() const{
134 return map_->getMaxLocalIndex();
135}
136
137template < class LO, class GO, class NO>
138void Map<LO,GO,NO>::print(Teuchos::EVerbosityLevel verbLevel) const{
139
140 Teuchos::RCP<Teuchos::FancyOStream> out = Teuchos::VerboseObjectBase::getDefaultOStream();
141 map_->describe(*out,verbLevel);
142}
143
144template < class LO, class GO, class NO>
145typename Map<LO,GO,NO>::ThyraVSBConstPtr_Type Map<LO,GO,NO>::getThyraVectorSpaceBase() const{
146
147 TEUCHOS_TEST_FOR_EXCEPTION(map_.is_null(),std::runtime_error,"getThyraVectorSpaceBase(): map_ is null.");
148 //RCP< const Tpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > tpMap = tpetraMap->getTpetra_Map();
149 Teuchos::RCP<Thyra::TpetraVectorSpace<default_sc,LO,GO,NO> > thyraTpetraMap = Thyra::tpetraVectorSpace<default_sc,LO,GO,NO>(map_);
150 return thyraTpetraMap;
151}
152
153template < class LO, class GO, class NO>
154typename Map<LO,GO,NO>::CommPtr_Type Map<LO,GO,NO>::getCommNonConst() {
155 CommConstPtr_Type commConst = map_->getComm();
156 return Teuchos::rcp_const_cast<Comm_Type>(commConst);
157}
158
159template <class LO,class GO,class NO>
160Teuchos::RCP<Map<LO,GO,NO> > Map<LO,GO,NO>::buildUniqueMap( int numFreeProcs ) const
161{
162 TEUCHOS_TEST_FOR_EXCEPTION(map_.is_null(),std::runtime_error,"map is null.");
163 if (numFreeProcs==0) {
164
165 Teuchos::RCP<Tpetra::Vector<GO,LO,GO,NO> > myIndices = Teuchos::RCP( new Tpetra::Vector<GO,LO,GO,NO>(map_));
166 myIndices->putScalar(map_->getComm()->getRank()+1);
167
168 TpetraMapPtr_Type linearMap = Teuchos::RCP(new TpetraMap_Type( map_->getMaxAllGlobalIndex()+1, 0, map_->getComm()));
169 Teuchos::RCP<Tpetra::Vector<GO,LO,GO,NO> > globalIndices = Teuchos::RCP( new Tpetra::Vector<GO,LO,GO,NO>(linearMap));
170
171 Teuchos::RCP<Tpetra::Import<LO,GO,NO> > importer = Teuchos::RCP(new Tpetra::Import<LO, GO, NO>( map_,linearMap));
172 Teuchos::RCP<Tpetra::Import<LO,GO,NO> > importer2 = Teuchos::RCP(new Tpetra::Import<LO, GO, NO>(linearMap, map_));
173 globalIndices->doImport(*myIndices,*importer,Tpetra::INSERT);
174 myIndices->putScalar(0);
175 myIndices->doImport(*globalIndices,*importer2,Tpetra::ADD);
176
177 Teuchos::Array<GO> uniqueVector;
178 for (unsigned i=0; i<myIndices->getLocalLength(); i++) {
179 if (myIndices->getData(0)[i] == map_->getComm()->getRank()+1) {
180 uniqueVector.push_back(map_->getGlobalElement(i));
181 }
182 }
183 TpetraMapPtr_Type mapTpetra = Teuchos::RCP(new TpetraMap_Type(-1, uniqueVector(),0,map_->getComm()));
184 Teuchos::RCP<Map<LO,GO,NO> > map = Teuchos::rcp( new Map<LO,GO,NO>( mapTpetra ) );
185 return map;
186 }
187 else{
188 Teuchos::RCP<Tpetra::Vector<GO,LO,GO,NO> > myIndices = Teuchos::RCP( new Tpetra::Vector<GO,LO,GO,NO>(map_));
189 myIndices->putScalar(map_->getComm()->getRank()+1);
190 GO maxGID = map_->getMaxAllGlobalIndex();
191 int numAvailableRanks = map_->getComm()->getSize() - numFreeProcs;
192 int numElementsForAvailRank = (int) ( ( (maxGID+1) / numAvailableRanks ) + 100.*std::numeric_limits<double>::epsilon() );
193
194 int remainingElement = maxGID+1 - numAvailableRanks * numElementsForAvailRank;
195 bool hasOneMoreElement = false;
196
197 if ( remainingElement > map_->getComm()->getRank() ) {
198 numElementsForAvailRank++;
199 hasOneMoreElement = true;
200 }
201
202 if ( map_->getComm()->getRank() + 1 > map_->getComm()->getSize() - numFreeProcs){
203 numElementsForAvailRank = 0;
204 }
205
206 Teuchos::Array<GO> myElements( numElementsForAvailRank );
207 GO offset = numElementsForAvailRank * map_->getComm()->getRank();
208 if (!hasOneMoreElement) {
209 offset += remainingElement;
210 }
211 for (int i=0; i<myElements.size(); i++) {
212 myElements[i] = i + offset;
213 }
214
215 TpetraMapPtr_Type linearMapAvailRanks = Teuchos::RCP(new TpetraMap_Type(Teuchos::OrdinalTraits<GO>::invalid(), myElements(), 0, map_->getComm()));
216 Teuchos::RCP<Tpetra::Vector<GO,LO,GO,NO> > globalIndices = Teuchos::RCP( new Tpetra::Vector<GO,LO,GO,NO>(linearMapAvailRanks));
217
218 Teuchos::RCP<Tpetra::Import<LO,GO,NO> > importer = Teuchos::RCP(new Tpetra::Import<LO, GO, NO>(map_, linearMapAvailRanks ));
219 Teuchos::RCP<Tpetra::Import<LO,GO,NO> > importer2 = Teuchos::RCP(new Tpetra::Import<LO, GO, NO>( linearMapAvailRanks, map_));
220
221 globalIndices->doImport(*myIndices,*importer,Tpetra::INSERT);
222
223 myIndices->putScalar(0);
224 myIndices->doImport(*globalIndices,*importer2,Tpetra::ADD);
225
226 Teuchos::Array<GO> uniqueVector;
227 for (unsigned i=0; i<myIndices->getLocalLength(); i++) {
228 if (myIndices->getData(0)[i] == map_->getComm()->getRank()+1) {
229 uniqueVector.push_back(map_->getGlobalElement(i));
230 }
231 }
232 TpetraMapPtr_Type mapTpetra = Teuchos::RCP(new TpetraMap_Type(-1, uniqueVector(),0,map_->getComm()));
233 Teuchos::RCP<Map<LO,GO,NO> > map = Teuchos::rcp( new Map<LO,GO,NO>( mapTpetra ) );
234 return map;
235 }
236
237}
238
239//merge with above function
240template <class LO,class GO,class NO>
241Teuchos::RCP<Map<LO,GO,NO> > Map<LO,GO,NO>::buildUniqueMap( tuple_intint_Type rankRange ) const
242{
243 TEUCHOS_TEST_FOR_EXCEPTION(map_.is_null(),std::runtime_error,"map is null.");
244 int rank = map_->getComm()->getRank();
245 Teuchos::RCP<Tpetra::Vector<GO,LO,GO,NO> > myIndices = Teuchos::RCP( new Tpetra::Vector<GO,LO,GO,NO>(map_));
246
247 myIndices->putScalar(rank + 1);
248 GO maxGID = map_->getMaxAllGlobalIndex();
249
250 int numAvailableRanks = std::get<1>(rankRange) - std::get<0>(rankRange) + 1;
251 int numElementsForAvailRank = (int) ( ( (maxGID+1) / numAvailableRanks ) + 100.*std::numeric_limits<double>::epsilon() );
252
253 int remainingElement = maxGID+1 - numAvailableRanks * numElementsForAvailRank;
254 bool hasOneMoreElement = false;
255
256 if ( remainingElement > map_->getComm()->getRank() - std::get<0>(rankRange) ) {
257 numElementsForAvailRank++;
258 hasOneMoreElement = true;
259 }
260 if ( map_->getComm()->getRank() < std::get<0>(rankRange) || map_->getComm()->getRank() > std::get<1>(rankRange) )
261 numElementsForAvailRank = 0;
262
263
264 Teuchos::Array<GO> myElements( numElementsForAvailRank );
265 GO offset = numElementsForAvailRank * ( map_->getComm()->getRank() - std::get<0>(rankRange) );
266
267 if (!hasOneMoreElement)
268 offset += remainingElement;
269
270 for (int i=0; i<myElements.size(); i++)
271 myElements[i] = i + offset;
272
273 TpetraMapPtr_Type linearMapAvailRanks = Teuchos::RCP(new TpetraMap_Type( Teuchos::OrdinalTraits<GO>::invalid(), myElements(), 0, map_->getComm()));
274 Teuchos::RCP<Tpetra::Vector<GO,LO,GO,NO> > globalIndices = Teuchos::RCP( new Tpetra::Vector<GO,LO,GO,NO>(linearMapAvailRanks));
275
276 Teuchos::RCP<Tpetra::Import<LO,GO,NO> > importer = Teuchos::RCP(new Tpetra::Import<LO, GO, NO>( map_,linearMapAvailRanks));
277 Teuchos::RCP<Tpetra::Import<LO,GO,NO> > importer2 = Teuchos::RCP(new Tpetra::Import<LO, GO, NO>(linearMapAvailRanks, map_));
278
279 globalIndices->doImport(*myIndices,*importer,Tpetra::INSERT);
280
281 myIndices->putScalar(0);
282 myIndices->doImport(*globalIndices,*importer2,Tpetra::ADD);
283
284 Teuchos::Array<GO> uniqueVector;
285 for (unsigned i=0; i<myIndices->getLocalLength(); i++) {
286 if (myIndices->getData(0)[i] == map_->getComm()->getRank()+1) {
287 uniqueVector.push_back(map_->getGlobalElement(i));
288 }
289 }
290
291 TpetraMapPtr_Type mapTpetra = Teuchos::RCP(new TpetraMap_Type(-1, uniqueVector(),0,map_->getComm()));
292 Teuchos::RCP<Map<LO,GO,NO> > map = Teuchos::rcp( new Map<LO,GO,NO>( mapTpetra ) );
293 return map;
294
295}
296
297
298}
299#endif
FEDD::Map::buildUniqueMap
Teuchos::RCP< Map< LO, GO, NO > > buildUniqueMap(int numFreeProcs=0) const
Definition Map_def.hpp:160
FEDD
Adaptive Mesh Refinement.
Definition AdaptiveMeshRefinement_decl.hpp:36