ParallelDofMapping.h 24.8 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
24
25
26
#ifndef AMDIS_FE_SPACE_MAPPING_H
#define AMDIS_FE_SPACE_MAPPING_H

#include <mpi.h>
27
#include <vector>
28
#include <map>
Thomas Witkowski's avatar
FETI-DP    
Thomas Witkowski committed
29
#include <set>
30
31
#include <petsc.h>
#include <petscis.h>
32
33
#include <boost/container/flat_map.hpp>
#include <boost/container/flat_set.hpp>
34
35

#include "AMDiS_fwd.h"
Thomas Witkowski's avatar
FETI-DP    
Thomas Witkowski committed
36
#include "parallel/DofComm.h"
37
38
#include "parallel/MpiHelper.h"
#include "parallel/ParallelTypes.h"
39
#include "parallel/StdMpi.h"
40
41

namespace AMDiS {
Thomas Witkowski's avatar
Thomas Witkowski committed
42

43
44
  using namespace std;

45
  /** \brief
46
47
48
   * Is used to store matrix indices to all DOFs in rank's subdomain. Thus, the
   * class defines a mapping from component number and DOF index to a global
   * matrix index. This class does not calculate the indices by itself!
49
50
51
52
   */
  class DofToMatIndex
  {
  public:
53
54
    typedef boost::container::flat_map<DegreeOfFreedom, int> MapType;

55
56
57
58
59
60
61
62
    DofToMatIndex() {}

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

63
64
65
66
67
68
69
    /** Add a new mapping for a given DOF.
     * 
     * \param[in]   component       Component number for which the mapping 
     *                              is defined.
     * \param[in]   dof             DOF index
     * \param[in]   globalMatIndex  Global matrix index.
     */
70
71
72
73
74
    inline void add(int component, DegreeOfFreedom dof, int globalMatIndex)
    {
      data[component][dof] = globalMatIndex;
    }

75
    /// Maps a global DOF index to the global matrix index for a specific 
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
    /// 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];
    }

91
92
93
94
95
96
97
98
99
100
101
102
    /// Returns the number of DOF mappings in one component
    inline int getSize(int component)
    {
      return data[component].size();
    }

    /// Returns the whole mapping for one component
    inline MapType& getData(int component)
    {
      return data[component];
    }
    
103
104
105
106
107
    /// Returns for a given matrix index the component and (local or global) DOF
    /// index. As the data structure is not made for this kind of reverse
    /// search, this is very slow and should be only used for debugging.
    void getReverse(int rowIndex, int &component, int &dofIndex);

108
109
110
  private:
    /// The mapping data. For each system component there is a specific map that
    /// maps global DOF indices to global matrix indices.
111
    map<int, MapType> data;
112
113
114
  };


115
116
117
  /// This class defines the parallel mapping of DOFs for one FE space. It is
  /// used by the class \ref ParallelDofMapping to specifiy the mapping for a
  /// set of  FE spaces.
118
  class ComponentDofMap
119
120
  {
  public:
121
122
    /// This constructor exists only to create std::map of this class and make
    /// use of the operator [] for read access. Should never be called.
123
    ComponentDofMap() 
124
125
126
    {
      ERROR_EXIT("Should not be called!\n");
    }
127
128
     
    /// This is the only valid constructur to be used. 
129
    ComponentDofMap(MeshLevelData* ld);
130
131

    /// Clears all data of the mapping.
132
    void clear();
133

134
135
    /// Maps a DOF index to both, the local and global index of the mapping. The
    /// global index must not be set.
Thomas Witkowski's avatar
Thomas Witkowski committed
136
    MultiIndex& operator[](DegreeOfFreedom d)
137
    {
138
      TEST_EXIT_DBG(dofMap.count(d))("DOF %d is not in map!\n", d);
139

140
      return dofMap[d];
141
    }
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160

    /** \brief
     * Searchs the map for a given DOF. It does not fail, if the DOF is not 
     * mapped by this mapping. In this case, it returns false. If the DOF is 
     * mapped, the result is stored and the function returns true.
     *
     * \param[in]    dof     DOF index for which a mapping is searched.
     * \param[out]   index   In the case that the DOF is mapped, the result
     *                       is stored here.
     */
    inline bool find(DegreeOfFreedom dof, MultiIndex& index)
    {
      DofMap::iterator it = dofMap.find(dof);
      if (it == dofMap.end())
	return false;

      index = it->second;
      return true;
    }
161
    
162
163
    /// Inserts a new DOF to rank's mapping. The DOF is assumed to be owend by
    /// the rank.
Thomas Witkowski's avatar
Thomas Witkowski committed
164
165
    inline void insertRankDof(DegreeOfFreedom dof0, 
			      DegreeOfFreedom dof1 = -1)
166
    {
167
      FUNCNAME("ComponentDofMap::insertRankDof()");
168
      
169
170
      TEST_EXIT_DBG(dofMap.count(dof0) == 0)
	("DOF %d is already defined in mapping!\n", dof0);
171
      
172
173
      dofMap[dof0].local = dof1;
      nLocalDofs++;
Thomas Witkowski's avatar
Thomas Witkowski committed
174
      if (dof1 != -1)
175
	nRankDofs++;
176
177
    }
    
178
179
    /// Inserts a new DOF to rank's mapping. The DOF exists in rank's subdomain
    /// but is owned by a different rank, thus it is part of an interior boundary.
Thomas Witkowski's avatar
Thomas Witkowski committed
180
181
    inline void insertNonRankDof(DegreeOfFreedom dof0, 
				 DegreeOfFreedom dof1 = -1)
182
    {
183
      FUNCNAME("ComponentDofMap::insertNonRankDof()");
184
      
185
      TEST_EXIT_DBG(dofMap.count(dof0) == 0)("Should not happen!\n");
186
      
187
188
189
      dofMap[dof0].local = dof1;
      nLocalDofs++;
      nonRankDofs.insert(dof0);
190
191
    }
    
192
    /// Checks if a given DOF is in the DOF mapping.
193
    bool isSet(DegreeOfFreedom dof)
194
    {
195
      return static_cast<bool>(dofMap.count(dof));
196
    }
197

Thomas Witkowski's avatar
Blbu    
Thomas Witkowski committed
198
199
200
    /// Checks if a given DOF is a rank owned DOF of the DOF mapping. The DOF
    /// must be a DOF of the mapping (this is not checked here), otherwise the
    /// result is meaningsless.
201
    bool isRankDof(DegreeOfFreedom dof)
202
    {
203
      return !(static_cast<bool>(nonRankDofs.count(dof)));
204
    }
Thomas Witkowski's avatar
Blbu    
Thomas Witkowski committed
205
206
207
208
209

    bool isRankGlobalDof(int dof)
    {
      return (dof >= rStartDofs && dof < rStartDofs + nRankDofs);
    }
210
    
211
    /// Returns number of DOFs in the mapping.
212
    unsigned int size()
213
    {
214
      return dofMap.size();
215
216
    }
    
217
    /// Returns the raw data of the mapping.
218
    DofMap& getMap()
219
    {
220
      return dofMap;
221
    }
Thomas Witkowski's avatar
FETI-DP    
Thomas Witkowski committed
222

223
224
225
226
227
228
229
230
231
232
    DofMap::iterator begin()
    {
      return dofMap.begin();
    }

    DofMap::iterator end()
    {
      return dofMap.end();
    }

233
234
    /// Recomputes the mapping.
    void update();
Thomas Witkowski's avatar
FETI-DP    
Thomas Witkowski committed
235

236
    /// Sets the FE space this mapping corresponds to.
Thomas Witkowski's avatar
FETI-DP    
Thomas Witkowski committed
237
238
239
240
241
    void setFeSpace(const FiniteElemSpace *fe)
    {
      feSpace = fe;
    }

242
    /// Informs the mapping whether a global index must be computed.
243
    void setGlobalMapping(bool b)
Thomas Witkowski's avatar
FETI-DP    
Thomas Witkowski committed
244
    {
245
      globalMapping = b;
Thomas Witkowski's avatar
FETI-DP    
Thomas Witkowski committed
246
247
    }

248
    /// Sets the DOF communicator.
249
    void setDofComm(DofComm &dc)
Thomas Witkowski's avatar
Thomas Witkowski committed
250
    {
251
      dofComm = &dc;
Thomas Witkowski's avatar
Thomas Witkowski committed
252
253
    }

254
255
256
257
258
259
    void setMpiComm(MPI::Intracomm &m, int l)
    {
      mpiComm = m;
      meshLevel = l;
    }

260
  private:
261
    /// Computes a global mapping from the local one.
262
    void computeGlobalMapping();
263
264
265

    /// Computes the global indices of all DOFs in the mapping that are not owned
    /// by the rank.
266
    void computeNonLocalIndices();
267
268

  private:
269
    MeshLevelData *levelData;
Thomas Witkowski's avatar
FETI-DP    
Thomas Witkowski committed
270

271
272
    /// DOF communicator for all DOFs on interior boundaries.
    DofComm *dofComm;
273

274
275
276
277
    MPI::Intracomm mpiComm;

    int meshLevel;

278
    /// The FE space this mapping belongs to. This is used only the get the
279
    /// correct DOF communicator in \ref dofComm.
Thomas Witkowski's avatar
FETI-DP    
Thomas Witkowski committed
280
281
    const FiniteElemSpace *feSpace;

282
    /// Mapping data from DOF indices to local and global indices.
283
    DofMap dofMap;
284

285
    /// Set of all DOFs that are in mapping but are not owned by the rank.
Thomas Witkowski's avatar
Thomas Witkowski committed
286
    boost::container::flat_set<DegreeOfFreedom> nonRankDofs;    
Thomas Witkowski's avatar
FETI-DP    
Thomas Witkowski committed
287

288
    /// If true, a global index mapping will be computed for all DOFs.
289
    bool globalMapping;
290

291
292
  public:
    /// 
293
    int nRankDofs, nLocalDofs, nOverallDofs, rStartDofs;
294
  };
295

296

297
298
  /// Abstract iterator interface to access containrs containing values of
  /// type \ref ComponentDofMap.
299
300
301
302
  class ComponentIterator {
  public:
    virtual ComponentDofMap& operator*() = 0;

Thomas Witkowski's avatar
Thomas Witkowski committed
303
304
    virtual ComponentDofMap* operator->() = 0;

305
306
307
308
309
310
311
    virtual bool end() = 0;
    
    virtual void next() = 0;

    virtual void reset() = 0;
  };

312

313
314
315
  /// Abstract interface to acces DOF mapping data for each component. Allows
  /// to hide specific implementations, which allow, e.g., to efficiently map
  /// all components having the same FE space to the same DOF mapping.
316
317
318
  class ComponentDataInterface
  {
  public:
319
    /// Access via component number
320
321
    virtual ComponentDofMap& operator[](int compNumber) = 0;

322
    /// Acess via FE space pointer
323
324
    virtual ComponentDofMap& operator[](const FiniteElemSpace *feSpace) = 0;

325
326
    /// Checks whether the DOF mapping is defined for a specific 
    /// component number.
327
328
    virtual bool isDefinedFor(int compNumber) const = 0;

329
    /// Returns iterator which iterates over all DOF mappings.
330
331
    virtual ComponentIterator& getIteratorData() = 0;
    
332
333
334
335
    /// Returns iterator which iterates over the DOF mappings of all 
    /// components. If the data is defined for each FE space and more than
    /// one commponent is defined on the same FE space, the iterative will
    /// also iterative multple times over the same DOF mapping object.
336
337
    virtual ComponentIterator& getIteratorComponent() = 0;

338
339
340
341
342
343
344
345
346
347
348
349
    /** \brief
     * Initialization of the object.
     *
     * \param[in]   componentSpaces   Set of the FE spaces for each component.
     * \param[in]   feSpaces          Set of all different FE spaces.
     * \param[in]   globalMapping     Mapping is parallel (true) or only 
     *                                local (false).
     * \param[in]   levelData         Data for multi level method.
     */
    virtual void init(vector<const FiniteElemSpace*> &componentSpaces,
		      vector<const FiniteElemSpace*> &feSpaces,
		      bool globalMapping,
350
		      MeshLevelData &levelData) = 0;
351
352
353

  protected:
    /// The FE spaces for all components.
354
    vector<const FiniteElemSpace*> componentSpaces;
355
356
357

    /// The set of all FE spaces. It uniquly contains all different FE spaces
    /// from \ref feSpaces.
358
    vector<const FiniteElemSpace*> feSpaces;
359
360
361
  };


362
363
364
  /// This class concretizes the interface class \ref ComponentDataInterface. A
  /// DOF mapping is implemented for each component.
  class ComponentData : public ComponentDataInterface
365
366
  {
  public:
367
368
    ComponentData()
      : iter(this)
369
370
    {}

371
    /// Returns DOF mapping for a given component number.
372
373
    ComponentDofMap& operator[](int compNumber)
    {
374
375
      TEST_EXIT_DBG(componentData.count(compNumber))
	("No data for component %d!\n", compNumber);
376

377
      return componentData.find(compNumber)->second;
378
379
    }

380
381
382
383
384
    /// Just to implement the corresponding virtual function in \ref 
    /// ComponentDataInterface. Of course it does not work as we have data for
    /// each component. Thus there may be different mappings for the same
    /// FE space.
    ComponentDofMap& operator[](const FiniteElemSpace *feSpace)
385
    {
386
      ERROR_EXIT("FE Space access is not possible for component wise defined DOF mappings\n");
387
    }
388
389
    
    /// Return data iterator.
390
391
    ComponentIterator& getIteratorData()
    {
392
393
      iter.reset();
      return iter;
394
395
    }
    
396
    /// Return component iterator.
397
398
    ComponentIterator& getIteratorComponent()
    {
399
400
      iter.reset();
      return iter;
401
402
    }

403
404
405
406
407
408
409
410
    /// Checks whether the DOF mapping is defined for a specific 
    /// component number.
    bool isDefinedFor(int compNumber) const
    {
      return (static_cast<unsigned int>(compNumber) < componentData.size());
    }

    /// Initialization
411
412
    void init(vector<const FiniteElemSpace*> &f0,
	      vector<const FiniteElemSpace*> &f1,
413
	      bool globalMapping,
414
	      MeshLevelData &levelData);
415

416
  protected:
417

418
419
    /// Iterator class to iterate over all parallel DOF mappings.
    class Iterator : public ComponentIterator {
420
    public:
421
422
423
      Iterator(ComponentData *d)
	: data(d),
	  componentCounter(-1)
424
425
      {}

426
427
      ComponentDofMap& operator*()
      {
428
	(*data)[componentCounter];
429
430
      }

Thomas Witkowski's avatar
Thomas Witkowski committed
431
432
      ComponentDofMap* operator->()
      {
433
	&((*data)[componentCounter]);
Thomas Witkowski's avatar
Thomas Witkowski committed
434
435
      }

436

437
438
      bool end()
      {
439
	return (it == data->componentSpaces.end());
440
441
442
443
      }

      void next()
      {
444
	++it;
445
446
447
448
	++componentCounter;

	if (it == data->componentSpaces.end())
	  componentCounter = -1;
449
450
451
452
      }

      void reset()
      {
453
454
	it = data->componentSpaces.begin();
	componentCounter = 0;
455
      }
456
457

    protected:
458
459
      /// Pointer to data class over which the iterator must iterate.
      ComponentData *data;
460

461
      /// Internal iterator of the internal data from \ref ComponentData.
462
      vector<const FiniteElemSpace*>::iterator it;
463

464
465
466
      /// Component number of current iteration. 
      int componentCounter;
    };
467
468


469
470
    /// Data mapping from component numbers to DOF mapping objects.
    map<unsigned int, ComponentDofMap> componentData;
471

472
473
    /// Iterator object.
    Iterator iter;
474
475


476

477
    friend class Iterator;
478
479
480
  };


481
482
483
484
485

  /// This class concretizes the interface class \ref ComponentDataInterface. A
  /// DOF mapping is implemented for each finite element space. Thus, different
  /// components sharing the same FE space are handled by the same DOF mapping.
  class FeSpaceData : public ComponentDataInterface
486
487
  {
  public:
488
489
490
    FeSpaceData()
      : iterData(this),
	iterComponent(this)
491
492
    {}

493
    /// Returns DOF mapping for a given component number.
494
495
    ComponentDofMap& operator[](int compNumber)
    {
496
497
      const FiniteElemSpace *feSpace = componentSpaces[compNumber];
      return componentData.find(feSpace)->second;
498
499
    }

500
    /// Returns DOF mapping for a given FE space.
501
502
    ComponentDofMap& operator[](const FiniteElemSpace *feSpace)
    {
503
504
505
      TEST_EXIT_DBG(componentData.count(feSpace))("No data for FE space!\n");;
      
      return componentData.find(feSpace)->second;
506
    }
507
508
509
510
511
512
513
514
515
516

    /// Checks whether the DOF mapping is defined for a specific 
    /// component number.
    bool isDefinedFor(int compNumber) const
    {
      const FiniteElemSpace *feSpace = componentSpaces[compNumber];
      return static_cast<bool>(componentData.count(feSpace));
    }

    /// Return data iterator.
517
518
    ComponentIterator& getIteratorData()
    {
519
520
      iterData.reset();
      return iterData;
521
522
    }
    
523
    /// Return component iterator.
524
525
    ComponentIterator& getIteratorComponent()
    {
526
527
      iterComponent.reset();
      return iterComponent;
528
529
    }

530
    /// Initialization
531
532
    void init(vector<const FiniteElemSpace*> &f0,
	      vector<const FiniteElemSpace*> &f1,
533
	      bool globalMapping,
534
535
	      MeshLevelData &levelData);

536

537
  protected:
538

539
540
    /// Iterator class to iterate over all parallel DOF mappings.
    class IteratorData : public ComponentIterator {
541
    public:
542
543
      IteratorData(FeSpaceData *d)
	: data(d)
544
545
      {}

546
547
      ComponentDofMap& operator*()
      {
548
	(*data)[*it];
549
550
      }

Thomas Witkowski's avatar
Thomas Witkowski committed
551
552
      ComponentDofMap* operator->()
      {
553
	&((*data)[*it]);
Thomas Witkowski's avatar
Thomas Witkowski committed
554
555
      }

556
557
      bool end()
      {
558
	return (it == data->feSpaces.end());
559
560
561
562
      }

      void next()
      {
563
	++it;
564
      }
Thomas Witkowski's avatar
Thomas Witkowski committed
565

566
      void reset()
Thomas Witkowski's avatar
Thomas Witkowski committed
567
      {
568
	it = data->feSpaces.begin();
Thomas Witkowski's avatar
Thomas Witkowski committed
569
570
      }

571
    protected:
572
      FeSpaceData *data;
573

574
      vector<const FiniteElemSpace*>::iterator it;
575
    };
576

577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613

    /// Iterator class to iterate over all component DOF mappings.
    class IteratorComponent : public ComponentIterator {
    public:
      IteratorComponent(FeSpaceData *d)
	: data(d)
      {}

      ComponentDofMap& operator*()
      {
	(*data)[*it];
      }

      ComponentDofMap* operator->()
      {
	&((*data)[*it]);
      }

      bool end()
      {
	return (it == data->componentSpaces.end());
      }

      void next()
      {
	++it;
      }

      void reset()
      {
	it = data->componentSpaces.begin();
      }

    protected:
      FeSpaceData *data;

      vector<const FiniteElemSpace*>::iterator it;
614
    };
615
    
616

617
    map<const FiniteElemSpace*, ComponentDofMap> componentData;
618

619
    IteratorData iterData;
620

621
622
623
624
625
    IteratorComponent iterComponent;

    friend class IteratorData;
    
    friend class IteratorComponent;
626
627
  };

628
629


630
631
632
633
634
635
  /// Used to specify whether a parallel DOF mapping is defined for each
  /// specific component or for each FE space.
  enum DofMappingMode {
    COMPONENT_WISE,
    FESPACE_WISE
  };
636

637
638
639
640
641
  /**
   * Implements the mapping from sets of distributed DOF indices to local and
   * global indices. The mapping works for a given set of FE spaces. Furthermore,
   * this class may compute the matrix indices of the set of DOF indices.
   */
Thomas Witkowski's avatar
Thomas Witkowski committed
642
  class ParallelDofMapping
643
644
  {
  public:
645
646
647
    /** \brief 
     * Constructur for parallel DOF mapping.
     *
648
649
650
651
     * \param[in]  mode            Defines if DOF mapping is defined either per
     *                             component or per FE space.
     * \param[in]  matIndexGlobal  If true, the mat index is defined on global
     *                             DOF indices, otherwise on local ones.
652
     */
653
654
    ParallelDofMapping(DofMappingMode mode, 
		       bool matIndexFromGlobal = false);
655

656
657
    /** \brief 
     * Initialize the parallel DOF mapping.
658
659
     *
     * \param[in]  m                  MPI communicator.
660
     * \param[in]  componentSpaces    The FE spaces of all components of the 
661
     *                                PDE to be solved.
662
     * \param[in]  feSpaces           Unique list of FE spaces. Thus, two
663
664
     *                                arbitrary elements of this list are always
     *                                different.
665
     * \param[in]  globalMapping      If true, at least one rank's mapping con-
666
667
     *                                taines DOFs that are not owend by the rank.
     */
668
    void init(MeshLevelData& mld,
669
670
	      vector<const FiniteElemSpace*> &componentSpaces,
	      vector<const FiniteElemSpace*> &feSpaces,
671
	      bool globalMapping = true);
672
673
674
675

    /// In the case of having only one FE space, this init function can be used.
    void init(MeshLevelData& mld,
	      const FiniteElemSpace *feSpace,
676
	      bool globalMapping = true)
677
678
679
    {
      vector<const FiniteElemSpace*> feSpaces;
      feSpaces.push_back(feSpace);
680
      init(mld, feSpaces, feSpaces, globalMapping);
681
682
    }

683
    void setMpiComm(MPI::Intracomm &m, int l);
684
685
686
    
    /// Clear all data.
    void clear();
687
688
689

    /// Set the DOF communicator objects that are required to exchange information
    /// about DOFs that are on interior boundaries.
690
    void setDofComm(DofComm &dofComm);
691

692
693
694
    /// Returns the DOF communicator.
    DofComm& getDofComm()
    {
695
696
697
      FUNCNAME("ParallelDofMapping::getDofComm()");

      TEST_EXIT_DBG(dofComm)("No DOF communicator object defined!\n");
698
699
700
701

      return *dofComm;
    }

Thomas Witkowski's avatar
Thomas Witkowski committed
702
703
704
705
706
    inline bool isMatIndexFromGlobal()
    {
      return needMatIndexFromGlobal;
    }

707
708
709
710
    /// Access the DOF mapping for a given component number.
    inline ComponentDofMap& operator[](int compNumber)
    {
      return (*data)[compNumber];
711
712
    }

713
714
    /// Access the DOF mapping for a given FE space    
    inline ComponentDofMap& operator[](const FiniteElemSpace *feSpace) 
Thomas Witkowski's avatar
Thomas Witkowski committed
715
    {
716
      return (*data)[feSpace];
Thomas Witkowski's avatar
Thomas Witkowski committed
717
718
    }

719
720
    /// Checks whether the DOF mapping is defined for a specific 
    /// component number.
721
    inline bool isDefinedFor(int compNumber) const
Thomas Witkowski's avatar
Thomas Witkowski committed
722
    {
723
      return data->isDefinedFor(compNumber);
Thomas Witkowski's avatar
Thomas Witkowski committed
724
725
    }

726
    /// Returns the number of solution components the mapping is defined on.
Thomas Witkowski's avatar
Thomas Witkowski committed
727
728
    inline int getNumberOfComponents() const
    {      
729
      return static_cast<int>(componentSpaces.size());
Thomas Witkowski's avatar
Thomas Witkowski committed
730
    }
731

732
    /// Returns \ref nRankDofs, thus the number of DOFs owned by the rank.
733
    inline int getRankDofs()
734
    {
735
      TEST_EXIT_DBG(nRankDofs >= 0)("Should not happen!\n");
736

737
      return nRankDofs;
738
739
    }

740
    /// Returns \ref nLocalDofs, thus the number of DOFs in ranks subdomain.
741
    inline int getLocalDofs()
Thomas Witkowski's avatar
Thomas Witkowski committed
742
    {
743
      TEST_EXIT_DBG(nLocalDofs >= 0)("Should not happen!\n");
Thomas Witkowski's avatar
Thomas Witkowski committed
744

745
      return nLocalDofs;
Thomas Witkowski's avatar
Thomas Witkowski committed
746
747
    }

748
    /// Returns \ref nOverallDofs, thus the number of all DOFs in this mapping.
749
    inline int getOverallDofs()
750
    {
751
      TEST_EXIT_DBG(nOverallDofs >= 0)("Should not happen!\n");
752

753
      return nOverallDofs;
754
755
    }

756
757
    /// Returns \ref rStartDofs, thus the smallest global index of a DOF that is
    /// owned by the rank.
758
    inline int getStartDofs()
759
    {
760
      TEST_EXIT_DBG(rStartDofs >= 0)("Should not happen!\n");
761

762
      return rStartDofs;
763
764
    }

765
    /// Update the mapping.
Thomas Witkowski's avatar
Thomas Witkowski committed
766
    void update();
767

768
769
770
    /// Updates only the DOF to matrix index mapping
    void updateMatIndex();

771
772
    /// Returns the global matrix index of a given DOF for a given 
    /// component number.
773
    inline int getMatIndex(int ithComponent, DegreeOfFreedom d)
Thomas Witkowski's avatar
Thomas Witkowski committed
774
    {
775
      return dofToMatIndex.get(ithComponent, d);
Thomas Witkowski's avatar
Thomas Witkowski committed
776
777
    }

778
779
780
781
782
783
784
    /// Returns the component number and local/global DOF index for a given
    /// matrix row index. Should be used for debugging only!
    inline void getReverseMatIndex(int index, int &component, int &dofIndex)
    {
      dofToMatIndex.getReverse(index, component, dofIndex);
    }

785
786
    /// Returns the local matrix index of a given DOF for a given 
    /// component number.
787
    inline int getLocalMatIndex(int ithComponent, DegreeOfFreedom d)
788
    {
789
      return dofToMatIndex.get(ithComponent, d) - rStartDofs;
790
791
    }

792
793
794
795
796
797
    /// Returns the set of unique FE spaces.
    inline vector<const FiniteElemSpace*>& getFeSpaces()
    {
      return feSpaces;
    }

Thomas Witkowski's avatar
Thomas Witkowski committed
798
799
800
801
802
803
804
805
806
807
808
809
    // Writes all data of this object to an output stream.
    void serialize(ostream &out)
    {
      ERROR_EXIT("MUST BE IMPLEMENTED!\n");
    }

    // Reads the object data from an input stream.
    void deserialize(istream &in)
    {
      ERROR_EXIT("MUST BE IMPLEMENTED!\n");
    }

810
    /// Compute local and global matrix indices.
811
    void computeMatIndex(bool globalIndex);
812

813
814
815
816
    void createIndexSet(IS &is, 
			int firstComponent, 
			int nComponents);

817
818
819
820
821
822
823
824
825
826
827
828
829
830
831
    /// Create a parallel distributed PETSc vector based on this mapping.
    inline void createVec(Vec &vec)
    {
      VecCreateMPI(mpiComm, getRankDofs(), getOverallDofs(), &vec);
    }

    /// Create a parallel distributed PETsc vector based on this mapping but
    /// with a different (larger) global size. This is used in multi-level 
    /// method to embed a local vector into a subdomain spaned by several
    /// ranks.
    inline void createVec(Vec &vec, int nGlobalRows)
    {
      VecCreateMPI(mpiComm, getRankDofs(), nGlobalRows, &vec);
    }

Thomas Witkowski's avatar
Thomas Witkowski committed
832
833
834
835
836
    inline void createLocalVec(Vec &vec)
    {
      VecCreateSeq(PETSC_COMM_SELF, getRankDofs(), &vec);
    }

Thomas Witkowski's avatar
bo eh    
Thomas Witkowski committed
837
838
839
840
    /// Prints out some information about the mapping. May be used during
    /// debugging or parallel solver creation.
    void printInfo();   

841
842
  protected:
    /// Compute \ref nRankDofs.
843
    int computeRankDofs();
844
845

    /// Compute \ref nLocalDofs.
846
    int computeLocalDofs();
847
848

    /// Compute \ref nOverallDofs.
849
    int computeOverallDofs();
850
851

    /// Compute \ref rStartDofs.
852
    int computeStartDofs();
853

854
  private:
855
856
857
858
    MPI::Intracomm mpiComm;

    int meshLevel;

859
    MeshLevelData *levelData;
860

861
862
    /// DOF communicator for all DOFs on interior boundaries.
    DofComm *dofComm;
863

864
865
866
    /// Is true if there are DOFs in at least one subdomain that are not owned
    /// by the rank. If the value is false, each rank contains only DOFs that
    /// are also owned by this rank.
867
    bool globalMapping;
868
869
870
871
872

    /// If matrix indices should be computed, this variable defines if the
    /// mapping from DOF indices to matrix row indices is defined on local
    /// or global DOF indices. If true, the mapping is to specify and to use
    /// on global ones, otherwise on local DOF indices.
Thomas Witkowski's avatar
Thomas Witkowski committed
873
874
875
876
877
    /// In most scenarios the mapping stored on local DOF indices is what we
    /// want to have. Only when periodic boundary conditions are used together
    /// with some global matrix approache, the matrix indices must be stored
    /// also for DOFs that are not part of the local subdomain. Thus, the
    /// mapping will be stored on global DOF indices.
878
879
    bool needMatIndexFromGlobal;

880
881
    /// Maps from components to DOF mappings.
    ComponentDataInterface *data;
Thomas Witkowski's avatar
Thomas Witkowski committed
882

883
    /// Number of DOFs owned by rank.
884
    int nRankDofs;
885
886

    /// Number of DOFs in rank's subdomain.
887
    int nLocalDofs;
888
889

    /// Number of global DOFs (this value is thus the same on all ranks).
890
    int nOverallDofs;
891
892

    /// Smallest global index of a DOF owned by the rank.
893
    int rStartDofs;
894
895
896

    /// Mapping from global DOF indices to global matrix indices under 
    /// consideration of possibly multiple components.
897
    DofToMatIndex dofToMatIndex;
898

899
    /// FE spaces of all components.
900
    vector<const FiniteElemSpace*> componentSpaces;
901
902
903

    /// Set of unique FE spaces.
    vector<const FiniteElemSpace*> feSpaces;
Thomas Witkowski's avatar
bo eh    
Thomas Witkowski committed
904
905
906
907
    
    /// Defines the DOF mapping either. The mapping may be defined either for
    /// FE spaces or for component numbers.
    DofMappingMode mode;
908
909
910
911
  };
}

#endif