ParallelDomainProblem.h 12.4 KB
Newer Older
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19
// ============================================================================
// ==                                                                        ==
// == AMDiS - Adaptive multidimensional simulations                          ==
// ==                                                                        ==
// ============================================================================
// ==                                                                        ==
// ==  crystal growth group                                                  ==
// ==                                                                        ==
// ==  Stiftung caesar                                                       ==
// ==  Ludwig-Erhard-Allee 2                                                 ==
// ==  53175 Bonn                                                            ==
// ==  germany                                                               ==
// ==                                                                        ==
// ============================================================================
// ==                                                                        ==
// ==  http://www.caesar.de/cg/AMDiS                                         ==
// ==                                                                        ==
// ============================================================================

20
/** \file ParallelDomain.h */
21

22 23
#ifndef AMDIS_PARALLELDOMAIN_H
#define AMDIS_PARALLELDOMAIN_H
24 25

#include <map>
26
#include <set>
27 28 29 30
#include <vector>

#include "ProblemTimeInterface.h"
#include "ProblemIterationInterface.h"
31
#include "FiniteElemSpace.h"
32
#include "AdaptInfo.h"
Thomas Witkowski's avatar
Thomas Witkowski committed
33
#include "InteriorBoundary.h"
34 35
#include "AMDiS_fwd.h"

36 37
#include "petsc.h"
#include "petscsys.h"
38
#include "petscao.h"
39 40 41 42 43 44
#include "mpi.h"

namespace AMDiS {

  class ParMetisPartitioner;

45
  class ParallelDomainBase : public ProblemIterationInterface,
46 47
                                    public ProblemTimeInterface
  {
Thomas Witkowski's avatar
Thomas Witkowski committed
48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66
  private:
    /// Defines type for a vector of DOFs.
    typedef std::vector<const DegreeOfFreedom*> DofContainer;

    /// Defines a mapping type from DOFs to rank numbers.
    typedef std::map<const DegreeOfFreedom*, int> DofToRank;

    /// Defines a mapping type from DOFs to a set of rank numbers.
    typedef std::map<const DegreeOfFreedom*, std::set<int> > DofToPartitions;

    /// Defines a mapping type from rank numbers to sets of DOFs.
    typedef std::map<int, DofContainer> RankToDofContainer;

    /// Defines a mapping type from DOF indices to DOF indices.
    typedef std::map<DegreeOfFreedom, DegreeOfFreedom> DofMapping;

    /// Defines a mapping type from DOF indices to boolean values.
    typedef std::map<DegreeOfFreedom, bool> DofToBool;

Thomas Witkowski's avatar
Thomas Witkowski committed
67 68 69 70 71 72
    /// Defines a mapping type from rank numbers to sets of coordinates.
    typedef std::map<int, std::vector<WorldVector<double> > > RankToCoords;

    /// Forward type (it maps rank numbers to the interior boundary objects).
    typedef InteriorBoundary::RankToBoundMap RankToBoundMap;

Thomas Witkowski's avatar
Thomas Witkowski committed
73 74
    typedef std::map<int, DofContainer> ElementIdxToDofs;

75
  public:
76
    ParallelDomainBase(const std::string& name,
Thomas Witkowski's avatar
Thomas Witkowski committed
77 78 79 80
		       ProblemIterationInterface *iterationIF,
		       ProblemTimeInterface *timeIF,
		       FiniteElemSpace *feSpace,
		       RefinementManager *refineManager);
81

82
    virtual ~ParallelDomainBase() {}
83

84
    virtual void initParallelization(AdaptInfo *adaptInfo);
85

86 87 88 89 90 91 92 93 94
    void exitParallelization(AdaptInfo *adaptInfo);

    /** \brief
     * Test, if the mesh consists of macro elements only. The mesh partitioning of
     * the parallelization works for macro meshes only and would fail, if the mesh
     * is already refined in some way. Therefore, this function will exit the program
     * if it finds a non macro element in the mesh.
     */
    void testForMacroMesh();
95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110

    /// Set for each element on the partitioning level the number of leaf elements.
    double setElemWeights(AdaptInfo *adaptInfo);

    void partitionMesh(AdaptInfo *adaptInfo);

    virtual void setTime(AdaptInfo *adaptInfo) {}

    virtual void initTimestep(AdaptInfo *adaptInfo) {}

    virtual void closeTimestep(AdaptInfo *adaptInfo) {}

    virtual void solveInitialProblem(AdaptInfo *adaptInfo) {}
  
    virtual void transferInitialSolution(AdaptInfo *adaptInfo) {}

111 112
    virtual void beginIteration(AdaptInfo *adaptInfo) 
    {
113 114 115
      iterationIF->beginIteration(adaptInfo);
    }

116 117
    virtual Flag oneIteration(AdaptInfo *adaptInfo, Flag toDo = FULL_ITERATION) 
    {
118 119 120
      ERROR_EXIT("Not implemented!\n");
      return 0;
    }
121

122 123
    virtual void endIteration(AdaptInfo *adaptInfo) 
    {
124 125
      iterationIF->endIteration(adaptInfo);
    }
126

127 128 129 130
    virtual int getNumProblems() 
    {
      return 0;
    }
131

132 133
    inline virtual const std::string& getName() 
    { 
134 135 136
      return name; 
    }

137
    /// Returns \ref nRankDOFs, the number of DOFs in the rank mesh.
138 139
    int getNumberRankDOFs() 
    {
140 141 142
      return nRankDOFs;
    }

143 144
    void fillPetscMatrix(DOFMatrix *mat, DOFVector<double> *vec);

145 146
    void solvePetscMatrix(DOFVector<double> *vec);

147 148 149 150
    virtual ProblemStatBase *getProblem(int number = 0) 
    {
      return NULL;
    }
151 152 153 154

    virtual void serialize(std::ostream&) {}

    virtual void deserialize(std::istream&) {}
155 156 157


  protected:
158 159 160 161
    /** \brief
     * Determine the interior boundaries, i.e. boundaries between ranks, and store
     * all information about them in \ref interiorBoundary.
     */
Thomas Witkowski's avatar
Thomas Witkowski committed
162 163
    void createInteriorBoundaryInfo(DofContainer& rankDOFs,
				    DofToRank& boundaryDOFs);
164 165 166 167

    /// Removes all macro elements from the mesh that are not part of ranks partition.
    void removeMacroElements();

168 169 170 171

    /** \brief
     * Creates from a macro mesh a correct local and global DOF index numbering.
     *
Thomas Witkowski's avatar
Thomas Witkowski committed
172
     * \param[out] rankDOFs      Returns all DOFs from the macro mesh, which are owned
173
     *                           by the rank after partitioning the macro mesh.
Thomas Witkowski's avatar
Thomas Witkowski committed
174
     * \param[out] boundaryDOFs  Returns all DOFs from the macro mesh, which lies at
175 176
     *                           an interior boundary of the rank. This object maps
     *                           each such DOF to the rank that owns this DOF.
Thomas Witkowski's avatar
Thomas Witkowski committed
177 178
     * \param[out] nRankDOFs     Number of DOFs owned by rank.
     * \param[out] nOverallDOFs  Number of all DOFs in macro mesh.
179
     */
Thomas Witkowski's avatar
Thomas Witkowski committed
180 181
    void createLocalGlobalNumbering(DofContainer& rankDOFs,
				    DofToRank& boundaryDOFs,
182 183 184 185
				    int& nRankDOFs, 
				    int& nOverallDOFs);

    void updateLocalGlobalNumbering(int& nRankDOFs, int& nOverallDOFs);
186

Thomas Witkowski's avatar
Thomas Witkowski committed
187 188 189
    void addAllVertexDOFs(Element *el, int ithEdge, DofContainer& dofs);

    void addAllEdgeDOFs(Element *el, int ithEdge, DofContainer& dofs);
Thomas Witkowski's avatar
Thomas Witkowski committed
190

191 192
    /** \brief
     * This function traverses the whole mesh, i.e. before it is really partitioned,
Thomas Witkowski's avatar
Thomas Witkowski committed
193 194 195
     * and collects information about which DOF corresponds to which rank. Can only
     * be used, if \ref partitionVec is set correctly. This is only the case, when
     * the macro mesh is partitioned.
196
     *
Thomas Witkowski's avatar
Thomas Witkowski committed
197
     * \param[out] partionDOFs   Stores to each DOF pointer the set of ranks the DOF is
198
     *                           part of.
Thomas Witkowski's avatar
Thomas Witkowski committed
199 200
     * \param[out] rankDOFs      Stores all rank DOFs.
     * \param[out] boundaryDOFs  Stores all DOFs in ranks partition that are on an 
201 202
     *                           interior boundary but correspond to another rank.
     */
Thomas Witkowski's avatar
n  
Thomas Witkowski committed
203 204 205 206
    void createDOFMemberInfo(DofToPartitions& partitionDofs,
			     DofContainer& rankOwnedDofs,
			     DofContainer& rankAllDofs,
			     DofToRank& boundaryDofs);
Thomas Witkowski's avatar
Thomas Witkowski committed
207

Thomas Witkowski's avatar
Thomas Witkowski committed
208 209 210
    void DbgCreateElementMap(ElementIdxToDofs &elMap);
    
    void DbgTestElementMap(ElementIdxToDofs &elMap);
Thomas Witkowski's avatar
Thomas Witkowski committed
211 212

    void DbgTestInteriorBoundary();
213
     
Thomas Witkowski's avatar
Thomas Witkowski committed
214 215 216 217 218 219
    /** \brief
     * This function is used for debugging only. It traverses all interior boundaries
     * and compares the dof indices on them with the dof indices of the boundarys
     * neighbours. The function fails, when dof indices on an interior boundary does
     * not fit together.
     */
Thomas Witkowski's avatar
Thomas Witkowski committed
220
    void DbgTestCommonDofs();
221

Thomas Witkowski's avatar
Thomas Witkowski committed
222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249 250
    inline void orderDOFs(const DegreeOfFreedom* dof1,
			  const DegreeOfFreedom* dof2,
			  const DegreeOfFreedom* dof3,
			  DofContainer &vec)
    {
      vec.resize(3);

      if (*dof1 < *dof2 && *dof1 < *dof3)
	vec[0] = dof1;
      else if (*dof2 < *dof1 && *dof2 < *dof3)
	vec[0] = dof2;
      else 
	vec[0] = dof3;

      if (*dof1 > *dof2 && *dof1 > *dof3)
	vec[2] = dof1;
      else if (*dof2 > *dof1 && *dof2 > *dof3)
	vec[2] = dof2;
      else 
	vec[2] = dof3;

      if (dof1 != vec[0] && dof1 != vec[2]) 
	vec[1] = dof1;
      else if (dof2 != vec[0] && dof2 != vec[2]) 
	vec[1] = dof2;
      else
	vec[1] = dof3;
    }

251
  protected:
252 253 254 255 256 257
    ///
    ProblemIterationInterface *iterationIF;

    ///
    ProblemTimeInterface *timeIF;

258 259 260 261 262 263 264 265 266 267 268 269 270 271 272 273
    /// The rank of the current process.
    int mpiRank;

    /// Overall number of processes.
    int mpiSize;

    /** \brief
     * MPI communicator collected all processes, which should
     * be used for calculation. The Debug procces is not included
     * in this communicator.
     */
    MPI::Intracomm mpiComm;

    /// Name of the problem (as used in the init files)
    std::string name;

274 275 276
    /// Finite element space of the problem.
    FiniteElemSpace *feSpace;

277 278 279
    /// Mesh of the problem.
    Mesh *mesh;

280 281 282
    /// Info level.
    int info;

283 284 285
    /// Refinement manager for the mesh.
    RefinementManager *refinementManager;

286 287 288 289 290 291 292 293 294 295 296 297 298 299 300 301 302 303 304 305
    /// Pointer to the paritioner which is used to devide a mesh into partitions.
    ParMetisPartitioner *partitioner;

    /// Weights for the elements, i.e., the number of leaf elements within this element.
    std::map<int, double> elemWeights;

    /// Is true, if the mesh was not partitioned before, otherwise it's false.
    bool initialPartitionMesh;

    /** \brief
     * Stores to every coarse element index the number of the partition it 
     * corresponds to.
     */
    std::map<int, int> partitionVec;

    /** \brief
     * Stores an old partitioning of elements. To every element index the number
     * of the parition it corresponds to is stored.
     */
    std::map<int, int> oldPartitionVec;    
306

307 308 309 310
    Mat petscMatrix;

    Vec petscRhsVec;
    
311
    Vec petscSolVec;
312

313 314
    /// Number of DOFs in the rank mesh.
    int nRankDOFs;
Thomas Witkowski's avatar
Thomas Witkowski committed
315

316 317 318 319
    /** \brief
     * Set of all interior boundary DOFs in ranks partition. The object maps to
     * each such DOF to the number of the rank that owns this DOF.
     */
Thomas Witkowski's avatar
Thomas Witkowski committed
320
    DofToRank boundaryDOFs;
321

Thomas Witkowski's avatar
Thomas Witkowski committed
322
    /** \brief 
Thomas Witkowski's avatar
Thomas Witkowski committed
323 324 325 326 327 328 329 330 331 332 333 334
     * Defines the interior boundaries of the domain that result from partitioning
     * the whole mesh. Contains only the boundaries, which are owned by the rank, i.e.,
     * the object gives for every neighbour rank i the boundaries this rank owns and 
     * shares with rank i.
     */
    InteriorBoundary myIntBoundary;
    
    /** \brief
     * Defines the interior boundaries of the domain that result from partitioning
     * the whole mesh. Contains only the boundaries, which are not owned by the rank,
     * i.e., the object gives for every neighbour rank i the boundaries that are
     * owned by rank i and are shared with this rank.
Thomas Witkowski's avatar
Thomas Witkowski committed
335
     */
Thomas Witkowski's avatar
Thomas Witkowski committed
336
    InteriorBoundary otherIntBoundary;
Thomas Witkowski's avatar
Thomas Witkowski committed
337

338 339 340 341
    /** \brief
     * This map contains for each rank the list of dofs the current rank must send
     * to exchange solution dofs at the interior boundaries.
     */
Thomas Witkowski's avatar
Thomas Witkowski committed
342
    RankToDofContainer sendDofs;
343 344 345 346 347

    /** \brief
     * This map contains for each rank the list of dofs from which the current rank 
     * must receive solution values of dofs at the interior boundaries.
     */
Thomas Witkowski's avatar
Thomas Witkowski committed
348
    RankToDofContainer recvDofs;
349 350

    /// Maps local to global dof indices.
Thomas Witkowski's avatar
Thomas Witkowski committed
351
    DofMapping mapLocalGlobalDOFs;
352 353

    /// Maps global to local dof indices.
Thomas Witkowski's avatar
Thomas Witkowski committed
354
    DofMapping mapGlobalLocalDOFs;
355 356 357 358 359 360

    /** \brief
     * Maps all DOFs in ranks partition to a bool value. If it is true, the DOF is
     * owned by the rank. Otherwise, its an interior boundary DOF that is owned by
     * another rank.
     */
Thomas Witkowski's avatar
n  
Thomas Witkowski committed
361 362 363
    DofToBool isRankDof;

    int rstart;
364 365
  };

Thomas Witkowski's avatar
Thomas Witkowski committed
366 367 368 369 370
  bool cmpDofsByValue(const DegreeOfFreedom* dof1, const DegreeOfFreedom* dof2)
  {
    return (*dof1 < *dof2);
  }

371
  class ParallelDomainScal : public ParallelDomainBase
372 373
  {
  public:
374
    ParallelDomainScal(const std::string& name,
375 376 377
			      ProblemScal *problem,
			      ProblemInstatScal *problemInstat);

378 379
    void initParallelization(AdaptInfo *adaptInfo);

380 381 382
    virtual Flag oneIteration(AdaptInfo *adaptInfo, Flag toDo = FULL_ITERATION);

  protected:
383 384
    /// Starts the solution of the linear system using Petsc.
    void solve();
385

386 387 388
  protected:
    /// Pointer to the stationary problem.
    ProblemScal *probScal;
389 390 391
  };
}

392
#endif // AMDIS_PARALLELDOMAIN_H