ParallelDofMapping.cc 11.6 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
//
// Software License for AMDiS
//
// Copyright (c) 2010 Dresden University of Technology 
// All rights reserved.
// Authors: Simon Vey, Thomas Witkowski et al.
//
// This file is part of AMDiS
//
// See also license.opensource.txt in the distribution.


13
#include "parallel/ParallelDofMapping.h"
Thomas Witkowski's avatar
FETI-DP  
Thomas Witkowski committed
14
#include "parallel/StdMpi.h"
15
16
17
18
19

namespace AMDiS {

  using namespace std;

20
  void FeSpaceDofMap::clear()
21
  {
22
    int nLevel = levelData->getLevelNumber();
23
    dofMap.clear();
24
25
    dofMap.resize(nLevel);

26
    nonRankDofs.clear();
27
28
    nonRankDofs.resize(nLevel);

29
30
31
32
33
    nRankDofs.clear();
    nLocalDofs.clear();
    nOverallDofs.clear();
    rStartDofs.clear();

34
35
36
37
    nRankDofs.resize(nLevel, 0);
    nLocalDofs.resize(nLevel, 0);
    nOverallDofs.resize(nLevel, 0);
    rStartDofs.resize(nLevel, 0);
38
39
40
  }


41
  void FeSpaceDofMap::update()
42
  {
43
    FUNCNAME("FeSpaceDofMap::update()");
44

45
    int nLevel = levelData->getLevelNumber();
46
    for (int i = 0; i < nLevel; i++) {
Thomas Witkowski's avatar
FETI-DP  
Thomas Witkowski committed
47

48
      // === Compute local indices for all rank owned DOFs. ===
49
      
50
51
52
53
54
55
56
57
      for (DofMap::iterator it = dofMap[i].begin(); it != dofMap[i].end(); ++it)
	if (it->second.local == -1 && nonRankDofs[i].count(it->first) == 0)
	  it->second.local = nRankDofs[i]++;
      
      // === Compute number of local and global DOFs in the mapping. ===
      
      nOverallDofs[i] = 0;
      rStartDofs[i] = 0;
58
      mpi::getDofNumbering(mpiComm,
59
			   nRankDofs[i], rStartDofs[i], nOverallDofs[i]);
60
61
62
63
64
65
66
67
68
      
      // === If required, compute also the global indices. ===
      
      if (needGlobalMapping) {
	computeGlobalMapping(i);
	
	if (hasNonLocalDofs)
	  computeNonLocalIndices(i);
      }
69
    }
70
71
72
  }


73
  void FeSpaceDofMap::computeGlobalMapping(int level)
74
  {
75
76
77
78
79
    FUNCNAME("FeSpaceDofMap::computeGlobalMapping()");

    for (DofMap::iterator it = dofMap[level].begin(); 
	 it != dofMap[level].end(); ++it)
      it->second.global = it->second.local + rStartDofs[level];
80
81
  }

82

83
  void FeSpaceDofMap::computeNonLocalIndices(int level)
Thomas Witkowski's avatar
FETI-DP  
Thomas Witkowski committed
84
  {
85
    FUNCNAME("FeSpaceDofMap::computeNonLocalIndices()");
Thomas Witkowski's avatar
FETI-DP  
Thomas Witkowski committed
86
87
88

    typedef map<int, map<const FiniteElemSpace*, DofContainer> >::iterator it_type;

89
90
91
    // === Send all global indices of DOFs that are owned by the rank to all ===
    // === other ranks that also include this DOF.                           ===

92
    StdMpi<vector<int> > stdMpi(mpiComm);
93

94
    for (DofComm::Iterator it(dofComm->getSendDofs(), level, feSpace); 
95
96
97
98
99
	 !it.end(); it.nextRank()) {
      int rank = it.getRank();
      if (meshLevel > 0)
	rank = levelData->mapRank(rank, 0, meshLevel);

Thomas Witkowski's avatar
FETI-DP  
Thomas Witkowski committed
100
      for (; !it.endDofIter(); it.nextDof())
101
102
103
104
105
	if (dofMap[level].count(it.getDofIndex()) && 
	    !nonRankDofs[level].count(it.getDofIndex()))
	  stdMpi.getSendData(rank).
	    push_back(dofMap[level][it.getDofIndex()].global);
    }
Thomas Witkowski's avatar
FETI-DP  
Thomas Witkowski committed
106
107
108

    stdMpi.updateSendDataSize();

109
110
111

    // === Check from which ranks this rank must receive some data. ===

112
    for (DofComm::Iterator it(dofComm->getRecvDofs(), level, feSpace); 
113
	 !it.end(); it.nextRank()) {
Thomas Witkowski's avatar
FETI-DP  
Thomas Witkowski committed
114
115
      bool recvFromRank = false;
      for (; !it.endDofIter(); it.nextDof()) {
116
	if (nonRankDofs[level].count(it.getDofIndex())) {
Thomas Witkowski's avatar
FETI-DP  
Thomas Witkowski committed
117
118
119
120
121
	  recvFromRank = true;
	  break;
	}
      }

122
123
124
125
126
127
128
      if (recvFromRank) {
	int rank = it.getRank();
	if (meshLevel > 0)
	  rank = levelData->mapRank(rank, 0, meshLevel);

	stdMpi.recv(rank);
      }
Thomas Witkowski's avatar
FETI-DP  
Thomas Witkowski committed
129
130
    }

131
132
133

    // === Start communication to exchange global indices. ===

Thomas Witkowski's avatar
FETI-DP  
Thomas Witkowski committed
134
135
    stdMpi.startCommunication();

136
137
138

    // === And set the global indices for all DOFs that are not owned by rank. ===
    
139
    for (DofComm::Iterator it(dofComm->getRecvDofs(), level, feSpace);
Thomas Witkowski's avatar
FETI-DP  
Thomas Witkowski committed
140
	 !it.end(); it.nextRank()) {
141
142
143
144
      int rank = it.getRank();
      if (meshLevel > 0)
	rank = levelData->mapRank(rank, 0, meshLevel);

Thomas Witkowski's avatar
FETI-DP  
Thomas Witkowski committed
145
146
      int i = 0;
      for (; !it.endDofIter(); it.nextDof())
147
	if (nonRankDofs[level].count(it.getDofIndex()))
148
	  dofMap[level][it.getDofIndex()].global = stdMpi.getRecvData(rank)[i++];
Thomas Witkowski's avatar
FETI-DP  
Thomas Witkowski committed
149
150
151
152
    }
  }


153
  void ParallelDofMapping::init(MeshLevelData &ldata,
154
				vector<const FiniteElemSpace*> &fe,
155
				vector<const FiniteElemSpace*> &uniqueFe,
156
157
				bool needGlobalMapping,
				bool bNonLocalDofs)
Thomas Witkowski's avatar
FETI-DP  
Thomas Witkowski committed
158
  {
159
160
    FUNCNAME("ParallelDofMapping::init()");

161
    levelData = &ldata;
162
    feSpaces = fe;
163
    feSpacesUnique = uniqueFe;
164
    hasNonLocalDofs = bNonLocalDofs;
Thomas Witkowski's avatar
FETI-DP  
Thomas Witkowski committed
165
    
166
167
168

    // === Init the mapping for all different FE spaces. ===

169
    for (vector<const FiniteElemSpace*>::iterator it = feSpacesUnique.begin();
170
171
172
173
174
	 it != feSpacesUnique.end(); ++it) {
      addFeSpace(*it);
      data[*it].setNeedGlobalMapping(needGlobalMapping);
      data[*it].setNonLocalDofs(hasNonLocalDofs);
    }
Thomas Witkowski's avatar
FETI-DP  
Thomas Witkowski committed
175
  }
Thomas Witkowski's avatar
Thomas Witkowski committed
176
177


178
179
180
181
  void ParallelDofMapping::clear()
  {
    FUNCNAME("ParallelDofMapping::clear()");

182
183
    TEST_EXIT_DBG(levelData)("No mesh level data object defined!\n");

184
    for (vector<const FiniteElemSpace*>::iterator it = feSpacesUnique.begin();
185
186
187
	 it != feSpacesUnique.end(); ++it)
      data[*it].clear();

188
189
190
191
192
    int nLevel = levelData->getLevelNumber();
    nRankDofs.resize(nLevel);
    nLocalDofs.resize(nLevel);
    nOverallDofs.resize(nLevel);
    rStartDofs.resize(nLevel);
193
    dofToMatIndex.resize(nLevel);
194

195
196
197
198
199
    for (int i = 0; i < nLevel; i++) {
      nRankDofs[i] = -1;
      nLocalDofs[i] = -1;
      nOverallDofs[i] = -1;
      rStartDofs[i] = -1;
200
      dofToMatIndex[i].clear();
201
    }
202
203
  }

204
205
206
207
208
209
210
211
212
213
214
  
  void ParallelDofMapping::setMpiComm(MPI::Intracomm &m, int l)
  {
    mpiComm = m;
    meshLevel = l;

    for (vector<const FiniteElemSpace*>::iterator it = feSpacesUnique.begin();
	 it != feSpacesUnique.end(); ++it)
      data[*it].setMpiComm(m, l);
  }

215

216
  void ParallelDofMapping::setDofComm(DofComm &dc)
217
218
219
  {
    FUNCNAME("ParallelDofMapping::setDofComm()");

220
    dofComm = &dc;
221
222

    // Add the DOF communicator also to all FE space DOF mappings.
223
    for (vector<const FiniteElemSpace*>::iterator it = feSpacesUnique.begin();
224
	 it != feSpacesUnique.end(); ++it)
225
      data[*it].setDofComm(dc);
226
  }
227

228

Thomas Witkowski's avatar
Thomas Witkowski committed
229
230
231
232
233
234
235
  void ParallelDofMapping::addFeSpace(const FiniteElemSpace* feSpace)
  {
    FUNCNAME("ParallelDofMapping::addFeSpace()");
    
    if (data.count(feSpace))
      data.find(feSpace)->second.clear();
    else
236
      data.insert(make_pair(feSpace, FeSpaceDofMap(levelData)));
Thomas Witkowski's avatar
Thomas Witkowski committed
237
238
239
240
241
    
    data.find(feSpace)->second.setFeSpace(feSpace);
  }    


242
  int ParallelDofMapping::computeRankDofs(int level)
Thomas Witkowski's avatar
Thomas Witkowski committed
243
  {
244
    FUNCNAME("ParallelDofMapping::computeRankDofs()");
Thomas Witkowski's avatar
Thomas Witkowski committed
245
246
    
    int result = 0;
247
248
    for (unsigned int i = 0; i < feSpaces.size(); i++) {
      TEST_EXIT_DBG(data.count(feSpaces[i]))("Should not happen!\n");
249
      result += data[feSpaces[i]].nRankDofs[level];
Thomas Witkowski's avatar
Thomas Witkowski committed
250
251
252
253
254
255
    }
    
    return result;
  }


256
  int ParallelDofMapping::computeLocalDofs(int level)
Thomas Witkowski's avatar
Thomas Witkowski committed
257
  {
258
    FUNCNAME("ParallelDofMapping::computeLocalDofs()");
Thomas Witkowski's avatar
Thomas Witkowski committed
259
260
    
    int result = 0;
261
262
    for (unsigned int i = 0; i < feSpaces.size(); i++) {
      TEST_EXIT_DBG(data.count(feSpaces[i]))("Should not happen!\n");
263
      result += data[feSpaces[i]].nLocalDofs[level];
Thomas Witkowski's avatar
Thomas Witkowski committed
264
265
266
267
268
269
    }
    
    return result;
  }


270
  int ParallelDofMapping::computeOverallDofs(int level)
Thomas Witkowski's avatar
Thomas Witkowski committed
271
  {
272
    FUNCNAME("ParallelDofMapping::computeOverallDofs()");
Thomas Witkowski's avatar
Thomas Witkowski committed
273
274
275
276

    int result = 0;
    for (unsigned int i = 0; i < feSpaces.size(); i++) {
      TEST_EXIT_DBG(data.count(feSpaces[i]))("Should not happen!\n");
277
      result += data.find(feSpaces[i])->second.nOverallDofs[level];
Thomas Witkowski's avatar
Thomas Witkowski committed
278
279
280
281
282
283
    }
    
    return result;
  }


284
  int ParallelDofMapping::computeStartDofs(int level)
Thomas Witkowski's avatar
Thomas Witkowski committed
285
  {
286
    FUNCNAME("ParallelDofMapping::computeStartDofs()");
Thomas Witkowski's avatar
Thomas Witkowski committed
287
288
289
290

    int result = 0;
    for (unsigned int i = 0; i < feSpaces.size(); i++) {
      TEST_EXIT_DBG(data.count(feSpaces[i]))("Should not happen!\n");
291
      result += data.find(feSpaces[i])->second.rStartDofs[level];
Thomas Witkowski's avatar
Thomas Witkowski committed
292
293
294
295
296
297
298
299
300
301
    }
    
    return result;
  }


  void ParallelDofMapping::update()
  {
    FUNCNAME("ParallelDofMapping::update()");

302
    // First, update all FE space DOF mappings.
303
    for (vector<const FiniteElemSpace*>::iterator it = feSpacesUnique.begin();
Thomas Witkowski's avatar
Thomas Witkowski committed
304
305
306
	 it != feSpacesUnique.end(); ++it)
      data[*it].update();
    
307
    int nLevel = levelData->getLevelNumber();
308
309
310
311
312
313
314
315
316
317
318
319
    for (int i = 0; i < nLevel; i++) {
      // Compute all numbers from this mappings.

      nRankDofs[i] = computeRankDofs(i);
      nLocalDofs[i] = computeLocalDofs(i);
      nOverallDofs[i] = computeOverallDofs(i);
      rStartDofs[i] = computeStartDofs(i);
      
      // And finally, compute the matrix indices.
      if (needMatIndex)
	computeMatIndex(needMatIndexFromGlobal, i);
    }
Thomas Witkowski's avatar
Thomas Witkowski committed
320
321
  }
  
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336

  void ParallelDofMapping::update(vector<const FiniteElemSpace*>& fe)
  {
    FUNCNAME("ParallelDofMapping::update()");

    for (vector<const FiniteElemSpace*>::iterator it = fe.begin();
	 it != fe.end(); ++it)
      if (find(feSpacesUnique.begin(), feSpacesUnique.end(), *it) == 
	  feSpacesUnique.end())
	ERROR_EXIT("Wrong FE space!\n");

    feSpaces = fe;
    update();
  }

Thomas Witkowski's avatar
Thomas Witkowski committed
337
  
338
  void ParallelDofMapping::computeMatIndex(bool globalIndex, int level)
Thomas Witkowski's avatar
Thomas Witkowski committed
339
340
  {
    FUNCNAME("ParallelDofMapping::computeMatIndex()");
341

342
    dofToMatIndex[level].clear();
Thomas Witkowski's avatar
Thomas Witkowski committed
343
    
344
345
346
    // The offset is always added to the local matrix index. The offset for the
    // DOFs in the first FE spaces is the smalled global index of a DOF that is
    // owned by the rank.
347
    int offset = rStartDofs[level];
Thomas Witkowski's avatar
Thomas Witkowski committed
348
    
349
350
351

    // === Create the matrix indices for all component FE spaces. ===

Thomas Witkowski's avatar
Thomas Witkowski committed
352
    for (unsigned int i = 0; i < feSpaces.size(); i++) {
353
354
355

      // Traverse all DOFs of the FE space and create for all rank owned DOFs
      // a matrix index.
356
      DofMap& dofMap = data[feSpaces[i]].getMap(level);
357
      for (DofMap::iterator it = dofMap.begin(); it != dofMap.end(); ++it) {
358
	if (data[feSpaces[i]].isRankDof(it->first, level)) {
Thomas Witkowski's avatar
Thomas Witkowski committed
359
	  int globalMatIndex = it->second.local + offset;
360
	  if (globalIndex)
361
	    dofToMatIndex[level].add(i, it->second.global, globalMatIndex);
362
	  else
363
	    dofToMatIndex[level].add(i, it->first, globalMatIndex);
Thomas Witkowski's avatar
Thomas Witkowski committed
364
365
366
	}
      }
      
367
368
      // Increase the offset for the next FE space by the number of DOFs owned 
      // by the rank in the current FE space.
369
      offset += data[feSpaces[i]].nRankDofs[level];
Thomas Witkowski's avatar
Thomas Witkowski committed
370
	
371
      // If there are no non local DOFs, continue with the next FE space.
Thomas Witkowski's avatar
Thomas Witkowski committed
372
373
374
      if (!hasNonLocalDofs)
	continue;
      
375
      TEST_EXIT_DBG(dofComm != NULL)("No communicator given!\n");
Thomas Witkowski's avatar
Thomas Witkowski committed
376
      
377
378
379
      // === Communicate the matrix indices for all DOFs that are on some ===
      // === interior boundaries.                                         ===

380
      StdMpi<vector<DegreeOfFreedom> > stdMpi(mpiComm);
381
      for (DofComm::Iterator it(dofComm->getSendDofs(), level, feSpaces[i]); 
Thomas Witkowski's avatar
Thomas Witkowski committed
382
383
384
385
386
	   !it.end(); it.nextRank()) {
	vector<DegreeOfFreedom> sendGlobalDofs;
	
	for (; !it.endDofIter(); it.nextDof())
	  if (dofMap.count(it.getDofIndex()))
387
	    if (globalIndex)
388
	      sendGlobalDofs.push_back(dofToMatIndex[level].get(i, dofMap[it.getDofIndex()].global));
389
	    else
390
	      sendGlobalDofs.push_back(dofToMatIndex[level].get(i, it.getDofIndex()));
Thomas Witkowski's avatar
Thomas Witkowski committed
391
	
392
393
394
395
396
	int rank = it.getRank();
	if (meshLevel > 0)
	  rank = levelData->mapRank(rank, 0, meshLevel);
	
	stdMpi.send(rank, sendGlobalDofs);
Thomas Witkowski's avatar
Thomas Witkowski committed
397
398
      }
      
399
      for (DofComm::Iterator it(dofComm->getRecvDofs(), level, feSpaces[i]); 
400
401
402
403
404
405
406
	   !it.end(); it.nextRank()) {
	int rank = it.getRank();
	if (meshLevel > 0)
	  rank = levelData->mapRank(rank, 0, meshLevel);

	stdMpi.recv(rank);
      }
Thomas Witkowski's avatar
Thomas Witkowski committed
407
408
409
410
      
      stdMpi.startCommunication();
      
      {
411
	for (DofComm::Iterator it(dofComm->getRecvDofs(), level, feSpaces[i]); 
Thomas Witkowski's avatar
Thomas Witkowski committed
412
	     !it.end(); it.nextRank()) {
413
414
415
416
	  int rank = it.getRank();
	  if (meshLevel > 0)
	    rank = levelData->mapRank(rank, 0, meshLevel);
	  
Thomas Witkowski's avatar
Thomas Witkowski committed
417
418
419
	  int counter = 0;
	  for (; !it.endDofIter(); it.nextDof()) {
	    if (dofMap.count(it.getDofIndex())) {
420
	      DegreeOfFreedom d = stdMpi.getRecvData(rank)[counter++];
421
	      if (globalIndex)
422
		dofToMatIndex[level].add(i, dofMap[it.getDofIndex()].global, d);
423
	      else
424
		dofToMatIndex[level].add(i, it.getDofIndex(), d);
Thomas Witkowski's avatar
Thomas Witkowski committed
425
426
427
428
429
430
431
	    }
	  }
	}
      }
    }
  }
  
432
}