ParallelDofMapping.h 23.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 53 54 55 56 57 58 59 60
   */
  class DofToMatIndex
  {
  public:
    DofToMatIndex() {}

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

61 62 63 64 65 66 67
    /** 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.
     */
68 69 70 71 72
    inline void add(int component, DegreeOfFreedom dof, int globalMatIndex)
    {
      data[component][dof] = globalMatIndex;
    }

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

89 90 91 92 93
    /// 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);

94 95 96
  private:
    /// The mapping data. For each system component there is a specific map that
    /// maps global DOF indices to global matrix indices.
97
    map<int, boost::container::flat_map<DegreeOfFreedom, int> > data;
98 99 100
  };


101 102 103
  /// 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.
104
  class ComponentDofMap
105 106
  {
  public:
107 108
    /// This constructor exists only to create std::map of this class and make
    /// use of the operator [] for read access. Should never be called.
109
    ComponentDofMap() 
110 111 112
    {
      ERROR_EXIT("Should not be called!\n");
    }
113 114
     
    /// This is the only valid constructur to be used. 
115
    ComponentDofMap(MeshLevelData* ld);
116 117

    /// Clears all data of the mapping.
118
    void clear();
119

120 121
    /// 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
122
    MultiIndex& operator[](DegreeOfFreedom d)
123
    {
124
      TEST_EXIT_DBG(dofMap.count(d))("DOF %d is not in map!\n", d);
125

126
      return dofMap[d];
127
    }
128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146

    /** \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;
    }
147
    
148 149
    /// Inserts a new DOF to rank's mapping. The DOF is assumed to be owend by
    /// the rank.
150
    void insertRankDof(DegreeOfFreedom dof0, DegreeOfFreedom dof1 = -1)
151
    {
152
      FUNCNAME("ComponentDofMap::insertRankDof()");
153
      
154
      TEST_EXIT_DBG(dofMap.count(dof0) == 0)("Should not happen!\n");
155
      
156 157
      dofMap[dof0].local = dof1;
      nLocalDofs++;
Thomas Witkowski's avatar
Thomas Witkowski committed
158
      if (dof1 != -1)
159
	nRankDofs++;
160 161
    }
    
162 163
    /// 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.
164
    void insertNonRankDof(DegreeOfFreedom dof0, DegreeOfFreedom dof1 = -1)
165
    {
166
      FUNCNAME("ComponentDofMap::insertNonRankDof()");
167
      
168
      TEST_EXIT_DBG(dofMap.count(dof0) == 0)("Should not happen!\n");
169
      
170 171 172
      dofMap[dof0].local = dof1;
      nLocalDofs++;
      nonRankDofs.insert(dof0);
173 174
    }
    
175
    /// Checks if a given DOF is in the DOF mapping.
176
    bool isSet(DegreeOfFreedom dof)
177
    {
178
      return static_cast<bool>(dofMap.count(dof));
179
    }
180

Thomas Witkowski's avatar
Blbu  
Thomas Witkowski committed
181 182 183
    /// 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.
184
    bool isRankDof(DegreeOfFreedom dof)
185
    {
186
      return !(static_cast<bool>(nonRankDofs.count(dof)));
187
    }
Thomas Witkowski's avatar
Blbu  
Thomas Witkowski committed
188 189 190 191 192

    bool isRankGlobalDof(int dof)
    {
      return (dof >= rStartDofs && dof < rStartDofs + nRankDofs);
    }
193
    
194
    /// Returns number of DOFs in the mapping.
195
    unsigned int size()
196
    {
197
      return dofMap.size();
198 199
    }
    
200
    /// Returns the raw data of the mapping.
201
    DofMap& getMap()
202
    {
203
      return dofMap;
204
    }
Thomas Witkowski's avatar
FETI-DP  
Thomas Witkowski committed
205

206 207 208 209 210 211 212 213 214 215
    DofMap::iterator begin()
    {
      return dofMap.begin();
    }

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

216 217
    /// Recomputes the mapping.
    void update();
Thomas Witkowski's avatar
FETI-DP  
Thomas Witkowski committed
218

219
    /// Sets the FE space this mapping corresponds to.
Thomas Witkowski's avatar
FETI-DP  
Thomas Witkowski committed
220 221 222 223 224
    void setFeSpace(const FiniteElemSpace *fe)
    {
      feSpace = fe;
    }

225
    /// Informs the mapping whether a global index must be computed.
226
    void setGlobalMapping(bool b)
Thomas Witkowski's avatar
FETI-DP  
Thomas Witkowski committed
227
    {
228
      globalMapping = b;
Thomas Witkowski's avatar
FETI-DP  
Thomas Witkowski committed
229 230
    }

231
    /// Sets the DOF communicator.
232
    void setDofComm(DofComm &dc)
Thomas Witkowski's avatar
Thomas Witkowski committed
233
    {
234
      dofComm = &dc;
Thomas Witkowski's avatar
Thomas Witkowski committed
235 236
    }

237 238 239 240 241 242
    void setMpiComm(MPI::Intracomm &m, int l)
    {
      mpiComm = m;
      meshLevel = l;
    }

243
  private:
244
    /// Computes a global mapping from the local one.
245
    void computeGlobalMapping();
246 247 248

    /// Computes the global indices of all DOFs in the mapping that are not owned
    /// by the rank.
249
    void computeNonLocalIndices();
250 251

  private:
252
    MeshLevelData *levelData;
Thomas Witkowski's avatar
FETI-DP  
Thomas Witkowski committed
253

254 255
    /// DOF communicator for all DOFs on interior boundaries.
    DofComm *dofComm;
256

257 258 259 260
    MPI::Intracomm mpiComm;

    int meshLevel;

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

265
    /// Mapping data from DOF indices to local and global indices.
266
    DofMap dofMap;
267

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

271
    /// If true, a global index mapping will be computed for all DOFs.
272
    bool globalMapping;
273

274 275
  public:
    /// 
276
    int nRankDofs, nLocalDofs, nOverallDofs, rStartDofs;
277
  };
278

279

280 281
  /// Abstract iterator interface to access containrs containing values of
  /// type \ref ComponentDofMap.
282 283 284 285
  class ComponentIterator {
  public:
    virtual ComponentDofMap& operator*() = 0;

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

288 289 290 291 292 293 294
    virtual bool end() = 0;
    
    virtual void next() = 0;

    virtual void reset() = 0;
  };

295

296 297 298
  /// 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.
299 300 301
  class ComponentDataInterface
  {
  public:
302
    /// Access via component number
303 304
    virtual ComponentDofMap& operator[](int compNumber) = 0;

305
    /// Acess via FE space pointer
306 307
    virtual ComponentDofMap& operator[](const FiniteElemSpace *feSpace) = 0;

308 309
    /// Checks whether the DOF mapping is defined for a specific 
    /// component number.
310 311
    virtual bool isDefinedFor(int compNumber) const = 0;

312
    /// Returns iterator which iterates over all DOF mappings.
313 314
    virtual ComponentIterator& getIteratorData() = 0;
    
315 316 317 318
    /// 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.
319 320
    virtual ComponentIterator& getIteratorComponent() = 0;

321 322 323 324 325 326 327 328 329 330 331 332
    /** \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,
333
		      MeshLevelData &levelData) = 0;
334 335 336

  protected:
    /// The FE spaces for all components.
337
    vector<const FiniteElemSpace*> componentSpaces;
338 339 340

    /// The set of all FE spaces. It uniquly contains all different FE spaces
    /// from \ref feSpaces.
341
    vector<const FiniteElemSpace*> feSpaces;
342 343 344
  };


345 346 347
  /// This class concretizes the interface class \ref ComponentDataInterface. A
  /// DOF mapping is implemented for each component.
  class ComponentData : public ComponentDataInterface
348 349
  {
  public:
350 351
    ComponentData()
      : iter(this)
352 353
    {}

354
    /// Returns DOF mapping for a given component number.
355 356
    ComponentDofMap& operator[](int compNumber)
    {
357 358
      TEST_EXIT_DBG(componentData.count(compNumber))
	("No data for component %d!\n", compNumber);
359

360
      return componentData.find(compNumber)->second;
361 362
    }

363 364 365 366 367
    /// 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)
368
    {
369
      ERROR_EXIT("FE Space access is not possible for component wise defined DOF mappings\n");
370
    }
371 372
    
    /// Return data iterator.
373 374
    ComponentIterator& getIteratorData()
    {
375 376
      iter.reset();
      return iter;
377 378
    }
    
379
    /// Return component iterator.
380 381
    ComponentIterator& getIteratorComponent()
    {
382 383
      iter.reset();
      return iter;
384 385
    }

386 387 388 389 390 391 392 393
    /// 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
394 395
    void init(vector<const FiniteElemSpace*> &f0,
	      vector<const FiniteElemSpace*> &f1,
396
	      bool globalMapping,
397
	      MeshLevelData &levelData);
398

399
  protected:
400

401 402
    /// Iterator class to iterate over all parallel DOF mappings.
    class Iterator : public ComponentIterator {
403
    public:
404 405 406
      Iterator(ComponentData *d)
	: data(d),
	  componentCounter(-1)
407 408
      {}

409 410
      ComponentDofMap& operator*()
      {
411
	(*data)[componentCounter];
412 413
      }

Thomas Witkowski's avatar
Thomas Witkowski committed
414 415
      ComponentDofMap* operator->()
      {
416
	&((*data)[componentCounter]);
Thomas Witkowski's avatar
Thomas Witkowski committed
417 418
      }

419

420 421
      bool end()
      {
422
	return (it == data->componentSpaces.end());
423 424 425 426
      }

      void next()
      {
427
	++it;
428 429 430 431
	++componentCounter;

	if (it == data->componentSpaces.end())
	  componentCounter = -1;
432 433 434 435
      }

      void reset()
      {
436 437
	it = data->componentSpaces.begin();
	componentCounter = 0;
438
      }
439 440

    protected:
441 442
      /// Pointer to data class over which the iterator must iterate.
      ComponentData *data;
443

444
      /// Internal iterator of the internal data from \ref ComponentData.
445
      vector<const FiniteElemSpace*>::iterator it;
446

447 448 449
      /// Component number of current iteration. 
      int componentCounter;
    };
450 451


452 453
    /// Data mapping from component numbers to DOF mapping objects.
    map<unsigned int, ComponentDofMap> componentData;
454

455 456
    /// Iterator object.
    Iterator iter;
457 458


459

460
    friend class Iterator;
461 462 463
  };


464 465 466 467 468

  /// 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
469 470
  {
  public:
471 472 473
    FeSpaceData()
      : iterData(this),
	iterComponent(this)
474 475
    {}

476
    /// Returns DOF mapping for a given component number.
477 478
    ComponentDofMap& operator[](int compNumber)
    {
479 480
      const FiniteElemSpace *feSpace = componentSpaces[compNumber];
      return componentData.find(feSpace)->second;
481 482
    }

483
    /// Returns DOF mapping for a given FE space.
484 485
    ComponentDofMap& operator[](const FiniteElemSpace *feSpace)
    {
486 487 488
      TEST_EXIT_DBG(componentData.count(feSpace))("No data for FE space!\n");;
      
      return componentData.find(feSpace)->second;
489
    }
490 491 492 493 494 495 496 497 498 499

    /// 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.
500 501
    ComponentIterator& getIteratorData()
    {
502 503
      iterData.reset();
      return iterData;
504 505
    }
    
506
    /// Return component iterator.
507 508
    ComponentIterator& getIteratorComponent()
    {
509 510
      iterComponent.reset();
      return iterComponent;
511 512
    }

513
    /// Initialization
514 515
    void init(vector<const FiniteElemSpace*> &f0,
	      vector<const FiniteElemSpace*> &f1,
516
	      bool globalMapping,
517 518
	      MeshLevelData &levelData);

519

520
  protected:
521

522 523
    /// Iterator class to iterate over all parallel DOF mappings.
    class IteratorData : public ComponentIterator {
524
    public:
525 526
      IteratorData(FeSpaceData *d)
	: data(d)
527 528
      {}

529 530
      ComponentDofMap& operator*()
      {
531
	(*data)[*it];
532 533
      }

Thomas Witkowski's avatar
Thomas Witkowski committed
534 535
      ComponentDofMap* operator->()
      {
536
	&((*data)[*it]);
Thomas Witkowski's avatar
Thomas Witkowski committed
537 538
      }

539 540
      bool end()
      {
541
	return (it == data->feSpaces.end());
542 543 544 545
      }

      void next()
      {
546
	++it;
547
      }
Thomas Witkowski's avatar
Thomas Witkowski committed
548

549
      void reset()
Thomas Witkowski's avatar
Thomas Witkowski committed
550
      {
551
	it = data->feSpaces.begin();
Thomas Witkowski's avatar
Thomas Witkowski committed
552 553
      }

554
    protected:
555
      FeSpaceData *data;
556

557
      vector<const FiniteElemSpace*>::iterator it;
558
    };
559

560 561 562 563 564 565 566 567 568 569 570 571 572 573 574 575 576 577 578 579 580 581 582 583 584 585 586 587 588 589 590 591 592 593 594 595 596

    /// 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;
597
    };
598
    
599

600
    map<const FiniteElemSpace*, ComponentDofMap> componentData;
601

602
    IteratorData iterData;
603

604 605 606 607 608
    IteratorComponent iterComponent;

    friend class IteratorData;
    
    friend class IteratorComponent;
609 610
  };

611 612


613 614 615 616 617 618
  /// 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
  };
619

620 621 622 623 624
  /**
   * 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
625
  class ParallelDofMapping
626 627
  {
  public:
628 629 630 631 632 633
    /** \brief 
     * Constructur for parallel DOF mapping.
     *
     * \param[in]   mode  Defines if DOF mapping is defined either per 
     *                    component or per FE space.
     */
634
    ParallelDofMapping(DofMappingMode mode);
635

636 637
    /** \brief 
     * Initialize the parallel DOF mapping.
638 639 640
     *
     * \param[in]  m                  MPI communicator.
     * \param[in]  fe                 The FE spaces of all components of the 
641 642 643 644
     *                                PDE to be solved.
     * \param[in]  uniqueFe           Unique list of FE spaces. Thus, two
     *                                arbitrary elements of this list are always
     *                                different.
645
     * \param[in]  globalMapping      If true, at least one rank's mapping con-
646 647
     *                                taines DOFs that are not owend by the rank.
     */
648
    void init(MeshLevelData& mld,
649
	      vector<const FiniteElemSpace*> &fe,
650
	      vector<const FiniteElemSpace*> &uniqueFe,
651
	      bool globalMapping = true);
652 653 654 655

    /// In the case of having only one FE space, this init function can be used.
    void init(MeshLevelData& mld,
	      const FiniteElemSpace *feSpace,
656
	      bool globalMapping = true)
657 658 659
    {
      vector<const FiniteElemSpace*> feSpaces;
      feSpaces.push_back(feSpace);
660
      init(mld, feSpaces, feSpaces, globalMapping);
661 662
    }

663
    void setMpiComm(MPI::Intracomm &m, int l);
664 665 666
    
    /// Clear all data.
    void clear();
667 668 669

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

672 673 674
    /// Returns the DOF communicator.
    DofComm& getDofComm()
    {
675 676 677
      FUNCNAME("ParallelDofMapping::getDofComm()");

      TEST_EXIT_DBG(dofComm)("No DOF communicator object defined!\n");
678 679 680 681

      return *dofComm;
    }

682 683 684
    /// Changes the computation of matrix indices based of either local or
    /// global DOF indices, see \ref needMatIndexFromGlobal
    void setComputeMatIndex(bool global)
685 686 687 688
    {
      needMatIndexFromGlobal = global;
    }

Thomas Witkowski's avatar
Thomas Witkowski committed
689 690 691 692 693
    inline bool isMatIndexFromGlobal()
    {
      return needMatIndexFromGlobal;
    }

694

695 696 697 698
    /// Access the DOF mapping for a given component number.
    inline ComponentDofMap& operator[](int compNumber)
    {
      return (*data)[compNumber];
699 700
    }

701 702
    /// Access the DOF mapping for a given FE space    
    inline ComponentDofMap& operator[](const FiniteElemSpace *feSpace) 
Thomas Witkowski's avatar
Thomas Witkowski committed
703
    {
704
      return (*data)[feSpace];
Thomas Witkowski's avatar
Thomas Witkowski committed
705 706
    }

707 708
    /// Checks whether the DOF mapping is defined for a specific 
    /// component number.
709
    inline bool isDefinedFor(int compNumber) const
Thomas Witkowski's avatar
Thomas Witkowski committed
710
    {
711
      return data->isDefinedFor(compNumber);
Thomas Witkowski's avatar
Thomas Witkowski committed
712 713
    }

714
    /// Returns the number of solution components the mapping is defined on.
Thomas Witkowski's avatar
Thomas Witkowski committed
715 716
    inline int getNumberOfComponents() const
    {      
717
      return static_cast<int>(componentSpaces.size());
Thomas Witkowski's avatar
Thomas Witkowski committed
718
    }
719

720
    /// Returns \ref nRankDofs, thus the number of DOFs owned by the rank.
721
    inline int getRankDofs()
722
    {
723
      TEST_EXIT_DBG(nRankDofs >= 0)("Should not happen!\n");
724

725
      return nRankDofs;
726 727
    }

728
    /// Returns \ref nLocalDofs, thus the number of DOFs in ranks subdomain.
729
    inline int getLocalDofs()
Thomas Witkowski's avatar
Thomas Witkowski committed
730
    {
731
      TEST_EXIT_DBG(nLocalDofs >= 0)("Should not happen!\n");
Thomas Witkowski's avatar
Thomas Witkowski committed
732

733
      return nLocalDofs;
Thomas Witkowski's avatar
Thomas Witkowski committed
734 735
    }

736
    /// Returns \ref nOverallDofs, thus the number of all DOFs in this mapping.
737
    inline int getOverallDofs()
738
    {
739
      TEST_EXIT_DBG(nOverallDofs >= 0)("Should not happen!\n");
740

741
      return nOverallDofs;
742 743
    }

744 745
    /// Returns \ref rStartDofs, thus the smallest global index of a DOF that is
    /// owned by the rank.
746
    inline int getStartDofs()
747
    {
748
      TEST_EXIT_DBG(rStartDofs >= 0)("Should not happen!\n");
749

750
      return rStartDofs;
751 752
    }

753
    /// Update the mapping.
Thomas Witkowski's avatar
Thomas Witkowski committed
754
    void update();
755

756 757
    /// Returns the global matrix index of a given DOF for a given 
    /// component number.
758
    inline int getMatIndex(int ithComponent, DegreeOfFreedom d)
Thomas Witkowski's avatar
Thomas Witkowski committed
759
    {
760
      return dofToMatIndex.get(ithComponent, d);
Thomas Witkowski's avatar
Thomas Witkowski committed
761 762
    }

763 764 765 766 767 768 769
    /// 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);
    }

770 771
    /// Returns the local matrix index of a given DOF for a given 
    /// component number.
772
    inline int getLocalMatIndex(int ithComponent, DegreeOfFreedom d)
773
    {
774
      return dofToMatIndex.get(ithComponent, d) - rStartDofs;
775 776
    }

Thomas Witkowski's avatar
Thomas Witkowski committed
777 778 779 780 781 782 783 784 785 786 787 788
    // 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");
    }

789
    /// Compute local and global matrix indices.
790
    void computeMatIndex(bool globalIndex);
791

792 793 794 795
    void createIndexSet(IS &is, 
			int firstComponent, 
			int nComponents);

796 797 798 799 800 801 802 803 804 805 806 807 808 809 810
    /// 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
811 812 813 814 815
    inline void createLocalVec(Vec &vec)
    {
      VecCreateSeq(PETSC_COMM_SELF, getRankDofs(), &vec);
    }

816 817
  protected:
    /// Compute \ref nRankDofs.
818
    int computeRankDofs();
819 820

    /// Compute \ref nLocalDofs.
821
    int computeLocalDofs();
822 823

    /// Compute \ref nOverallDofs.
824
    int computeOverallDofs();
825 826

    /// Compute \ref rStartDofs.
827
    int computeStartDofs();
828

829
  private:
830 831 832 833
    MPI::Intracomm mpiComm;

    int meshLevel;

834
    MeshLevelData *levelData;
835

836 837
    /// DOF communicator for all DOFs on interior boundaries.
    DofComm *dofComm;
838

839 840 841
    /// 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.
842
    bool globalMapping;
843 844 845 846 847

    /// 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