FeSpaceMapping.h 11.4 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
// ============================================================================
// ==                                                                        ==
// == AMDiS - Adaptive multidimensional simulations                          ==
// ==                                                                        ==
// ==  http://www.amdis-fem.org                                              ==
// ==                                                                        ==
// ============================================================================
//
// 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.



/** \file FeSpaceMapping.h */

23
#include <vector>
24
#include <map>
Thomas Witkowski's avatar
FETI-DP  
Thomas Witkowski committed
25
26
#include <set>
#include "parallel/DofComm.h"
27
28
#include "parallel/MpiHelper.h"
#include "parallel/ParallelTypes.h"
29
#include "parallel/StdMpi.h"
30
31
32
33
34
35
36

#ifndef AMDIS_FE_SPACE_MAPPING_H
#define AMDIS_FE_SPACE_MAPPING_H

namespace AMDiS {
  using namespace std;

Thomas Witkowski's avatar
Thomas Witkowski committed
37
38
39
40
41
  struct MultiIndex
  {
    int local, global;
  };

42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
  /** \brief
   * Defines for each system component a mapping for sets of global DOF indices
   * to global matrix indices. The mapping is defined for all DOFs in rank's 
   * subdomain. When periodic boundary conditions are used, then the mapping
   * stores also information for the periodic associations of rank's DOF on
   * periodic boundaries.
   */
  class DofToMatIndex
  {
  public:
    DofToMatIndex() {}

    /// Reset the data structure.
    inline void clear()
    {
      data.clear();
    }

    /// Add a new mapping.
    inline void add(int component, DegreeOfFreedom dof, int globalMatIndex)
    {
      data[component][dof] = globalMatIndex;
    }

    /// Map a global DOF index to the global matrix index for a specific 
    /// system component number.
    inline int get(int component, DegreeOfFreedom dof)
    {
      FUNCNAME("DofToMatIndex::get()");

      TEST_EXIT_DBG(data.count(component))
	("No mapping data for component %d available!\n", component);

      TEST_EXIT_DBG(data[component].count(dof))
	("Mapping for DOF %d in component %d does not exists!\n",
	 dof, component);

      return data[component][dof];
    }

  private:
    /// The mapping data. For each system component there is a specific map that
    /// maps global DOF indices to global matrix indices.
    map<int, map<DegreeOfFreedom, int> > data;
  };


89
90
91
  class GlobalDofMap
  {
  public:
92
93
94
95
96
97
98
    /// This constructor exists only to create std::map of this class and make
    /// use of the operator [] for read access. Should never be called.
    GlobalDofMap() 
    {
      ERROR_EXIT("Should not be called!\n");
    }
      
99
100
    GlobalDofMap(MPI::Intracomm* m)
      : mpiComm(m),
Thomas Witkowski's avatar
Thomas Witkowski committed
101
102
	sendDofs(NULL),
	recvDofs(NULL),
103
	feSpace(NULL),
Thomas Witkowski's avatar
Thomas Witkowski committed
104
	needGlobalMapping(false),
105
	nRankDofs(0),
Thomas Witkowski's avatar
Thomas Witkowski committed
106
	nLocalDofs(0),
107
	nOverallDofs(0),
Thomas Witkowski's avatar
FETI-DP  
Thomas Witkowski committed
108
109
	rStartDofs(0),
	overlap(false)
110
111
    {}
    
112
    void clear();
113
    
Thomas Witkowski's avatar
Thomas Witkowski committed
114
    MultiIndex& operator[](DegreeOfFreedom d)
115
116
117
118
119
120
    {
      TEST_EXIT_DBG(dofMap.count(d))("Should not happen!\n");

      return dofMap[d];
    }
    
Thomas Witkowski's avatar
Thomas Witkowski committed
121
    void insertRankDof(DegreeOfFreedom dof0, DegreeOfFreedom dof1 = -1)
122
123
124
125
126
    {
      FUNCNAME("GlobalDofMap::insertRankDof()");
      
      TEST_EXIT_DBG(dofMap.count(dof0) == 0)("Should not happen!\n");
      
Thomas Witkowski's avatar
Thomas Witkowski committed
127
      dofMap[dof0].local = dof1;
Thomas Witkowski's avatar
Thomas Witkowski committed
128
      nLocalDofs++;
Thomas Witkowski's avatar
Thomas Witkowski committed
129
130
      if (dof1 != -1)
	nRankDofs++;
131
132
    }
    
Thomas Witkowski's avatar
FETI-DP  
Thomas Witkowski committed
133
    void insert(DegreeOfFreedom dof0, DegreeOfFreedom dof1 = -1)
134
135
136
137
138
    {
      FUNCNAME("GlobalDofMap::insert()");
      
      TEST_EXIT_DBG(dofMap.count(dof0) == 0)("Should not happen!\n");
      
Thomas Witkowski's avatar
Thomas Witkowski committed
139
      dofMap[dof0].local = dof1;
Thomas Witkowski's avatar
Thomas Witkowski committed
140
      nLocalDofs++;
Thomas Witkowski's avatar
FETI-DP  
Thomas Witkowski committed
141
      nonRankDofs.insert(dof0);
142
143
144
145
146
147
    }
    
    bool isSet(DegreeOfFreedom dof)
    {
      return static_cast<bool>(dofMap.count(dof));
    }
148
149
150
151
152

    bool isRankDof(DegreeOfFreedom dof)
    {
      return !(static_cast<bool>(nonRankDofs.count(dof)));
    }
153
154
155
156
157
158
    
    unsigned int size()
    {
      return dofMap.size();
    }
    
Thomas Witkowski's avatar
Thomas Witkowski committed
159
    map<DegreeOfFreedom, MultiIndex>& getMap()
160
161
162
163
    {
      return dofMap;
    }
    
Thomas Witkowski's avatar
Thomas Witkowski committed
164
    void update();
165

166
    void addOffset(int offset);
167

Thomas Witkowski's avatar
FETI-DP  
Thomas Witkowski committed
168
169
170
171
172
173
174
175
176
177
178
179
180
181
    void computeNonLocalIndices();

    void print();

    void setFeSpace(const FiniteElemSpace *fe)
    {
      feSpace = fe;
    }

    void setOverlap(bool b)
    {
      overlap = b;
    }

182
    void setNeedGlobalMapping(bool b)
Thomas Witkowski's avatar
FETI-DP  
Thomas Witkowski committed
183
    {
184
      needGlobalMapping = b;
Thomas Witkowski's avatar
FETI-DP  
Thomas Witkowski committed
185
186
    }

187
    void setDofComm(DofComm &pSend, DofComm &pRecv)
Thomas Witkowski's avatar
Thomas Witkowski committed
188
    {
189
190
      sendDofs = &pSend;
      recvDofs = &pRecv;
Thomas Witkowski's avatar
Thomas Witkowski committed
191
192
    }

193
194
  private:
    MPI::Intracomm* mpiComm;
Thomas Witkowski's avatar
FETI-DP  
Thomas Witkowski committed
195

196
197
198
199
    DofComm *sendDofs;
    
    DofComm *recvDofs;

Thomas Witkowski's avatar
FETI-DP  
Thomas Witkowski committed
200
201
    const FiniteElemSpace *feSpace;

202
    /// 
Thomas Witkowski's avatar
Thomas Witkowski committed
203
    map<DegreeOfFreedom, MultiIndex> dofMap;
204

Thomas Witkowski's avatar
FETI-DP  
Thomas Witkowski committed
205
206
    std::set<DegreeOfFreedom> nonRankDofs;

Thomas Witkowski's avatar
Thomas Witkowski committed
207
    bool needGlobalMapping;
208
209
  public:
    /// 
Thomas Witkowski's avatar
Thomas Witkowski committed
210
    int nRankDofs, nLocalDofs, nOverallDofs, rStartDofs;
Thomas Witkowski's avatar
FETI-DP  
Thomas Witkowski committed
211
212

    bool overlap;
213
214
215
216
217
218
219
  };
  

  template<typename T>
  class FeSpaceData
  {
  public:
220
221
    FeSpaceData() 
      : mpiComm(NULL),
222
223
224
	sendDofs(NULL),
	recvDofs(NULL),
	overlap(false),
225
226
227
228
229
	feSpaces(0),
	nRankDofs(-1),
	nOverallDofs(-1),
	rStartDofs(-1)
    {} 
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247

    T& operator[](const FiniteElemSpace* feSpace)
    {
      FUNCNAME("FeSpaceData::operator[]()");

      TEST_EXIT_DBG(data.count(feSpace))("Should not happen!\n");

      return data.find(feSpace)->second;
    }

    void addFeSpace(const FiniteElemSpace* feSpace)
    {
      FUNCNAME("FeSpaceData::addFeSpace()");
      
      if (data.count(feSpace))
	data.find(feSpace)->second.clear();
      else
	data.insert(make_pair(feSpace, T(mpiComm)));
Thomas Witkowski's avatar
FETI-DP  
Thomas Witkowski committed
248
249

      data.find(feSpace)->second.setFeSpace(feSpace);
250
251
    }    

252
    int getRankDofs(vector<const FiniteElemSpace*> &fe)
253
    {
254
255
      FUNCNAME("FeSpaceData::getRankDofs()");

256
      int result = 0;
257
258
259
      for (unsigned int i = 0; i < fe.size(); i++) {
	TEST_EXIT_DBG(data.count(fe[i]))("Cannot find FE space: %p\n", fe[i]);
	result += data[fe[i]].nRankDofs;
260
261
262
263
264
      }

      return result;
    }

265
266
267
268
269
270
271
    inline int getRankDofs()
    {
      TEST_EXIT_DBG(nRankDofs >= 0)("Should not happen!\n");

      return nRankDofs;
    }

Thomas Witkowski's avatar
Thomas Witkowski committed
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
    int getLocalDofs(vector<const FiniteElemSpace*> &fe)
    {
      FUNCNAME("FeSpaceData::getLocalDofs()");

      int result = 0;
      for (unsigned int i = 0; i < fe.size(); i++) {
	TEST_EXIT_DBG(data.count(fe[i]))("Cannot find FE space: %p\n", fe[i]);
	result += data[fe[i]].nLocalDofs;
      }

      return result;
    }

    inline int getLocalDofs()
    {
      TEST_EXIT_DBG(nLocalDofs >= 0)("Should not happen!\n");

      return nLocalDofs;
    }

292
293
294
295
296
297
298
299
300
301
302
    int getOverallDofs(vector<const FiniteElemSpace*> &feSpaces)
    {
      int result = 0;
      for (unsigned int i = 0; i < feSpaces.size(); i++) {
	TEST_EXIT_DBG(data.count(feSpaces[i]))("Should not happen!\n");
	result += data.find(feSpaces[i])->second.nOverallDofs;
      }

      return result;
    }

303
304
305
306
307
308
309
    inline int getOverallDofs()
    {
      TEST_EXIT_DBG(nOverallDofs >= 0)("Should not happen!\n");

      return nOverallDofs;
    }

310
311
312
313
314
315
316
317
318
319
320
    int getStartDofs(vector<const FiniteElemSpace*> &feSpaces)
    {
      int result = 0;
      for (unsigned int i = 0; i < feSpaces.size(); i++) {
	TEST_EXIT_DBG(data.count(feSpaces[i]))("Should not happen!\n");
	result += data.find(feSpaces[i])->second.rStartDofs;
      }

      return result;
    }

321
322
323
324
325
326
327
    int getStartDofs()
    {
      TEST_EXIT_DBG(rStartDofs >= 0)("Should not happen!\n");

      return rStartDofs;
    }

Thomas Witkowski's avatar
Thomas Witkowski committed
328
329
    void init(MPI::Intracomm *m,
	      vector<const FiniteElemSpace*> &fe,
330
331
	      bool needGlobalMapping,
	      bool ol)
332
    {
Thomas Witkowski's avatar
Thomas Witkowski committed
333
      mpiComm = m;
334
      feSpaces = fe;
335
      overlap = ol;
Thomas Witkowski's avatar
Thomas Witkowski committed
336
      for (unsigned int i = 0; i < feSpaces.size(); i++) {
Thomas Witkowski's avatar
Thomas Witkowski committed
337
338
	feSpacesUnique.insert(feSpaces[i]);
	
339
	addFeSpace(feSpaces[i]);
Thomas Witkowski's avatar
Thomas Witkowski committed
340
	data[feSpaces[i]].setNeedGlobalMapping(needGlobalMapping);
341
	data[feSpaces[i]].setOverlap(overlap);
Thomas Witkowski's avatar
Thomas Witkowski committed
342
      }
343
344
345
346
    }

    void update()
    {
Thomas Witkowski's avatar
Thomas Witkowski committed
347
348
349
      for (std::set<const FiniteElemSpace*>::iterator it = feSpacesUnique.begin();
	   it != feSpacesUnique.end(); ++it)
	data[*it].update();
Thomas Witkowski's avatar
Thomas Witkowski committed
350

351
      nRankDofs = getRankDofs(feSpaces);
Thomas Witkowski's avatar
Thomas Witkowski committed
352
      nLocalDofs = getLocalDofs(feSpaces);
353
354
      nOverallDofs = getOverallDofs(feSpaces);
      rStartDofs = getStartDofs(feSpaces);
355
356

      computeMatIndex();
357
358
    }

359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419

    void computeMatIndex()
    {
      dofToMatIndex.clear();

      int offset = rStartDofs;

      for (unsigned int i = 0; i < feSpaces.size(); i++) {

	map<DegreeOfFreedom, MultiIndex>& dofMap = data[feSpaces[i]].getMap();
	typedef map<DegreeOfFreedom, MultiIndex>::iterator ItType;
	for (ItType it = dofMap.begin(); it != dofMap.end(); ++it) {
	  if (data[feSpaces[i]].isRankDof(it->first)) {
	    int globalMatIndex = 
	      it->second.local - data[feSpaces[i]].rStartDofs + offset;
	    dofToMatIndex.add(i, it->first, globalMatIndex);
	  }
	}
	
	if (!overlap)
	  continue;
	
	TEST_EXIT_DBG(sendDofs != NULL && recvDofs != NULL)
	  ("No communicator given!\n");
	
	StdMpi<vector<DegreeOfFreedom> > stdMpi(*mpiComm);
	for (DofComm::Iterator it(*sendDofs, feSpaces[i]); 
	     !it.end(); it.nextRank()) {
	  vector<DegreeOfFreedom> sendGlobalDofs;
	  
	  for (; !it.endDofIter(); it.nextDof())
	    if (dofMap.count(it.getDofIndex()))
	      sendGlobalDofs.push_back(dofToMatIndex.get(i, it.getDofIndex()));
	  
	  stdMpi.send(it.getRank(), sendGlobalDofs);
	}
	
	for (DofComm::Iterator it(*recvDofs, feSpaces[i]); 
	     !it.end(); it.nextRank())
	  stdMpi.recv(it.getRank());
	
	stdMpi.startCommunication();
	
	{
	  int counter = 0;
	  for (DofComm::Iterator it(*recvDofs, feSpaces[i]); 
	       !it.end(); it.nextRank()) {
	    for (; !it.endDofIter(); it.nextDof()) {
	      if (dofMap.count(it.getDofIndex())) {
		DegreeOfFreedom d = stdMpi.getRecvData(it.getRank())[counter++];
		dofToMatIndex.add(i, it.getDofIndex(), d);
	      }
	    }
	  }
	}

	offset += data[feSpaces[i]].nRankDofs;
      }
    }


420
421
422
423
424
    inline int mapLocal(int index, int ithFeSpace)
    {
      int result = 0;
      for (int i = 0; i < ithFeSpace; i++)
	result += data[feSpaces[i]].nRankDofs;
Thomas Witkowski's avatar
Thomas Witkowski committed
425
      result += data[feSpaces[ithFeSpace]][index].local;
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443

      return result;
    }

    inline int mapLocal(int index, const FiniteElemSpace *feSpace)
    {
      for (unsigned int i = 0; i < feSpaces.size(); i++)
	if (feSpaces[i] == feSpace)
	  return mapLocal(index, feSpace, i);

      return -1;
    }

    inline int mapGlobal(int index, int ithFeSpace)
    {
      int result = rStartDofs;
      for (int i = 0; i < ithFeSpace; i++)
	result += data[feSpaces[i]].nRankDofs;
Thomas Witkowski's avatar
Thomas Witkowski committed
444
      result += data[feSpaces[ithFeSpace]][index].local;
445
446
447
448
449
450
451
452
453
454
455
456
457

      return result;
    }

    inline int mapGlobal(int index, const FiniteElemSpace *feSpace)
    {
      for (unsigned int i = 0; i < feSpaces.size(); i++)
	if (feSpaces[i] == feSpace)
	  return mapGlobal(index, feSpace, i);

      return -1;
    }

458
459
460
461
462
463
464
465
466
467
    void setDofComm(DofComm &pSend, DofComm &pRecv)
    {
      sendDofs = &pSend;
      recvDofs = &pRecv;

      for (std::set<const FiniteElemSpace*>::iterator it = feSpacesUnique.begin();
	   it != feSpacesUnique.end(); ++it)
	data[*it].setDofComm(pSend, pRecv);
    }

468
469
  private:
    MPI::Intracomm* mpiComm;
470
471
472
473
474
475

    DofComm *sendDofs;
    
    DofComm *recvDofs;

    bool overlap;
476
477
    
    map<const FiniteElemSpace*, T> data;
478
479
480

    vector<const FiniteElemSpace*> feSpaces;

Thomas Witkowski's avatar
Thomas Witkowski committed
481
482
483
    std::set<const FiniteElemSpace*> feSpacesUnique;

    int nRankDofs, nLocalDofs, nOverallDofs, rStartDofs;
484
485
486
487

    /// Mapping from global DOF indices to global matrix indices under 
    /// consideration of possibly multiple components.
    DofToMatIndex dofToMatIndex;
488
489
490
491
492
  };

}

#endif