MFEM  v3.3
Finite element discretization library
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
pfespace.hpp
Go to the documentation of this file.
1 // Copyright (c) 2010, Lawrence Livermore National Security, LLC. Produced at
2 // the Lawrence Livermore National Laboratory. LLNL-CODE-443211. All Rights
3 // reserved. See file COPYRIGHT for details.
4 //
5 // This file is part of the MFEM library. For more information and source code
6 // availability see http://mfem.org.
7 //
8 // MFEM is free software; you can redistribute it and/or modify it under the
9 // terms of the GNU Lesser General Public License (as published by the Free
10 // Software Foundation) version 2.1 dated February 1999.
11 
12 #ifndef MFEM_PFESPACE
13 #define MFEM_PFESPACE
14 
15 #include "../config/config.hpp"
16 
17 #ifdef MFEM_USE_MPI
18 
19 #include "../linalg/hypre.hpp"
20 #include "../mesh/pmesh.hpp"
21 #include "../mesh/nurbs.hpp"
22 #include "fespace.hpp"
23 
24 namespace mfem
25 {
26 
27 /// Abstract parallel finite element space.
29 {
30 private:
31  /// MPI data.
32  MPI_Comm MyComm;
33  int NRanks, MyRank;
34 
35  /// Parallel mesh.
36  ParMesh *pmesh;
37 
38  /// GroupCommunicator on the local VDofs
39  GroupCommunicator *gcomm;
40 
41  /// Number of true dofs in this processor (local true dofs).
42  int ltdof_size;
43 
44  /// The group of each local dof.
45  Array<int> ldof_group;
46 
47  /// For a local dof: the local true dof number in the master of its group.
48  Array<int> ldof_ltdof;
49 
50  /// Offsets for the dofs in each processor in global numbering.
51  Array<HYPRE_Int> dof_offsets;
52 
53  /// Offsets for the true dofs in each processor in global numbering.
54  Array<HYPRE_Int> tdof_offsets;
55 
56  /// Offsets for the true dofs in neighbor processor in global numbering.
57  Array<HYPRE_Int> tdof_nb_offsets;
58 
59  /// Previous 'dof_offsets' (before Update()), column partition of T.
60  Array<HYPRE_Int> old_dof_offsets;
61 
62  /// The sign of the basis functions at the scalar local dofs.
63  Array<int> ldof_sign;
64 
65  /// The matrix P (interpolation from true dof to dof).
66  HypreParMatrix *P;
67 
68  /// The (block-diagonal) matrix R (restriction of dof to true dof)
69  SparseMatrix *R;
70 
71  ParNURBSExtension *pNURBSext() const
72  { return dynamic_cast<ParNURBSExtension *>(NURBSext); }
73 
74  GroupTopology &GetGroupTopo() const
75  { return (NURBSext) ? pNURBSext()->gtopo : pmesh->gtopo; }
76 
77  void Construct();
78  void Destroy();
79 
80  // ldof_type = 0 : DOFs communicator, otherwise VDOFs communicator
81  void GetGroupComm(GroupCommunicator &gcomm, int ldof_type,
82  Array<int> *ldof_sign = NULL);
83 
84  /// Construct dof_offsets and tdof_offsets using global communication.
85  void GenerateGlobalOffsets();
86 
87  /// Construct ldof_group and ldof_ltdof.
88  void ConstructTrueDofs();
89  void ConstructTrueNURBSDofs();
90 
91  void ApplyLDofSigns(Array<int> &dofs) const;
92 
93  /// Helper struct to store DOF dependencies in a parallel NC mesh.
94  struct Dependency
95  {
96  int rank, dof; ///< master DOF, may be on another processor
97  double coef;
98  Dependency(int r, int d, double c) : rank(r), dof(d), coef(c) {}
99  };
100 
101  /// Dependency list for a local vdof.
102  struct DepList
103  {
104  Array<Dependency> list;
105  int type; ///< 0 = independent, 1 = one-to-one (conforming), 2 = slave
106 
107  DepList() : type(0) {}
108 
109  bool IsTrueDof(int my_rank) const
110  { return type == 0 || (type == 1 && list[0].rank == my_rank); }
111  };
112 
113  void AddSlaveDependencies(DepList deps[], int master_rank,
114  const Array<int> &master_dofs, int master_ndofs,
115  const Array<int> &slave_dofs, DenseMatrix& I);
116 
117  void Add1To1Dependencies(DepList deps[], int owner_rank,
118  const Array<int> &owner_dofs, int owner_ndofs,
119  const Array<int> &dependent_dofs);
120 
121  void GetDofs(int type, int index, Array<int>& dofs);
122  void ReorderFaceDofs(Array<int> &dofs, int orient);
123 
124  /// Build the P and R matrices.
125  void Build_Dof_TrueDof_Matrix();
126 
127  // Used when the ParMesh is non-conforming, i.e. pmesh->pncmesh != NULL.
128  // Constructs the matrices P and R. Determines ltdof_size. Calls
129  // GenerateGlobalOffsets(). Constructs ldof_ltdof.
130  void GetParallelConformingInterpolation();
131 
132  /** Calculate a GridFunction migration matrix after mesh load balancing.
133  The result is a parallel permutation matrix that can be used to update
134  all grid functions defined on this space. */
135  HypreParMatrix* RebalanceMatrix(int old_ndofs,
136  const Table* old_elem_dof);
137 
138  /** Calculate a GridFunction restriction matrix after mesh derefinement.
139  The matrix is constructed so that the new grid function interpolates
140  the original function, i.e., the original function is evaluated at the
141  nodes of the coarse function. */
142  HypreParMatrix* ParallelDerefinementMatrix(int old_ndofs,
143  const Table *old_elem_dof);
144 
145 public:
146  // Face-neighbor data
147  // Number of face-neighbor dofs
149  // Face-neighbor-element to face-neighbor dof
151  // Face-neighbor to ldof in the face-neighbor numbering
153  // The global ldof indices of the face-neighbor dofs
155  // Local face-neighbor data: face-neighbor to ldof
157 
159  int dim = 1, int ordering = Ordering::byNODES);
160 
161  MPI_Comm GetComm() { return MyComm; }
162  int GetNRanks() { return NRanks; }
163  int GetMyRank() { return MyRank; }
164 
165  inline ParMesh *GetParMesh() { return pmesh; }
166 
167  int GetDofSign(int i)
168  { return NURBSext || Nonconforming() ? 1 : ldof_sign[VDofToDof(i)]; }
169  HYPRE_Int *GetDofOffsets() { return dof_offsets; }
170  HYPRE_Int *GetTrueDofOffsets() { return tdof_offsets; }
171  HYPRE_Int GlobalVSize()
172  { return Dof_TrueDof_Matrix()->GetGlobalNumRows(); }
173  HYPRE_Int GlobalTrueVSize()
174  { return Dof_TrueDof_Matrix()->GetGlobalNumCols(); }
175 
176  /// Return the number of local vector true dofs.
177  virtual int GetTrueVSize() { return ltdof_size; }
178 
179  /// Returns indexes of degrees of freedom in array dofs for i'th element.
180  virtual void GetElementDofs(int i, Array<int> &dofs) const;
181 
182  /// Returns indexes of degrees of freedom for i'th boundary element.
183  virtual void GetBdrElementDofs(int i, Array<int> &dofs) const;
184 
185  /** Returns the indexes of the degrees of freedom for i'th face
186  including the dofs for the edges and the vertices of the face. */
187  virtual void GetFaceDofs(int i, Array<int> &dofs) const;
188 
189  void GetSharedEdgeDofs(int group, int ei, Array<int> &dofs) const;
190  void GetSharedFaceDofs(int group, int fi, Array<int> &dofs) const;
191 
192  /// The true dof-to-dof interpolation matrix
194  { if (!P) { Build_Dof_TrueDof_Matrix(); } return P; }
195 
196  /** @brief For a non-conforming mesh, construct and return the interpolation
197  matrix from the partially conforming true dofs to the local dofs. The
198  returned pointer must be deleted by the caller. */
200 
201  /** Create and return a new HypreParVector on the true dofs, which is
202  owned by (i.e. it must be destroyed by) the calling function. */
204  { return (new HypreParVector(MyComm,GlobalTrueVSize(),GetTrueDofOffsets()));}
205 
206  /// Scale a vector of true dofs
207  void DivideByGroupSize(double *vec);
208 
209  /// Return a reference to the internal GroupCommunicator (on VDofs)
210  GroupCommunicator &GroupComm() { return *gcomm; }
211 
212  /// Return a new GroupCommunicator on Dofs
214 
215  /** Given an integer array on the local degrees of freedom, perform
216  a bitwise OR between the shared dofs. */
217  void Synchronize(Array<int> &ldof_marker) const;
218 
219  /// Determine the boundary degrees of freedom
220  virtual void GetEssentialVDofs(const Array<int> &bdr_attr_is_ess,
221  Array<int> &ess_dofs) const;
222 
223  /** Get a list of essential true dofs, ess_tdof_list, corresponding to the
224  boundary attributes marked in the array bdr_attr_is_ess. */
225  virtual void GetEssentialTrueDofs(const Array<int> &bdr_attr_is_ess,
226  Array<int> &ess_tdof_list);
227 
228  /** If the given ldof is owned by the current processor, return its local
229  tdof number, otherwise return -1 */
230  int GetLocalTDofNumber(int ldof);
231  /// Returns the global tdof number of the given local degree of freedom
232  HYPRE_Int GetGlobalTDofNumber(int ldof);
233  /** Returns the global tdof number of the given local degree of freedom in
234  the scalar version of the current finite element space. The input should
235  be a scalar local dof. */
236  HYPRE_Int GetGlobalScalarTDofNumber(int sldof);
237 
238  HYPRE_Int GetMyDofOffset() const;
239  HYPRE_Int GetMyTDofOffset() const;
240 
242  { return Dof_TrueDof_Matrix(); }
243  /// Get the R matrix which restricts a local dof vector to true dof vector.
245  { Dof_TrueDof_Matrix(); return R; }
246 
247  // Face-neighbor functions
248  void ExchangeFaceNbrData();
249  int GetFaceNbrVSize() const { return num_face_nbr_dofs; }
250  void GetFaceNbrElementVDofs(int i, Array<int> &vdofs) const;
251  void GetFaceNbrFaceVDofs(int i, Array<int> &vdofs) const;
252  const FiniteElement *GetFaceNbrFE(int i) const;
253  const FiniteElement *GetFaceNbrFaceFE(int i) const;
254  const HYPRE_Int *GetFaceNbrGlobalDofMap() { return face_nbr_glob_dof_map; }
255 
257  void LoseDofOffsets() { dof_offsets.LoseData(); }
258  void LoseTrueDofOffsets() { tdof_offsets.LoseData(); }
259 
260  bool Conforming() const { return pmesh->pncmesh == NULL; }
261  bool Nonconforming() const { return pmesh->pncmesh != NULL; }
262 
263  /** Reflect changes in the mesh. Calculate one of the refinement/derefinement
264  /rebalance matrices, unless want_transform is false. */
265  virtual void Update(bool want_transform = true);
266 
267  /// Free ParGridFunction transformation matrix (if any), to save memory.
268  virtual void UpdatesFinished()
269  {
271  old_dof_offsets.DeleteAll();
272  }
273 
274  virtual ~ParFiniteElementSpace() { Destroy(); }
275 
276  // Obsolete, kept for backward compatibility
277  int TrueVSize() { return ltdof_size; }
278 };
279 
280 }
281 
282 #endif // MFEM_USE_MPI
283 
284 #endif
Abstract class for Finite Elements.
Definition: fe.hpp:46
virtual int GetTrueVSize()
Return the number of local vector true dofs.
Definition: pfespace.hpp:177
HypreParVector * NewTrueDofVector()
Definition: pfespace.hpp:203
Ordering::Type ordering
Definition: fespace.hpp:74
virtual void Update(bool want_transform=true)
Definition: pfespace.cpp:2198
int VDofToDof(int vdof) const
Definition: fespace.hpp:252
HYPRE_Int * GetDofOffsets()
Definition: pfespace.hpp:169
Data type dense matrix using column-major storage.
Definition: densemat.hpp:22
virtual const Operator * GetProlongationMatrix()
Definition: pfespace.hpp:241
void DivideByGroupSize(double *vec)
Scale a vector of true dofs.
Definition: pfespace.cpp:455
Abstract parallel finite element space.
Definition: pfespace.hpp:28
void Synchronize(Array< int > &ldof_marker) const
Definition: pfespace.cpp:491
virtual const SparseMatrix * GetRestrictionMatrix()
Get the R matrix which restricts a local dof vector to true dof vector.
Definition: pfespace.hpp:244
HypreParMatrix * GetPartialConformingInterpolation()
For a non-conforming mesh, construct and return the interpolation matrix from the partially conformin...
Definition: pfespace.cpp:1544
const HYPRE_Int * GetFaceNbrGlobalDofMap()
Definition: pfespace.hpp:254
HYPRE_Int GetMyDofOffset() const
Definition: pfespace.cpp:612
void DeleteAll()
Delete whole array.
Definition: array.hpp:481
void GetSharedEdgeDofs(int group, int ei, Array< int > &dofs) const
Definition: pfespace.cpp:298
ParNCMesh * pncmesh
Definition: pmesh.hpp:133
HYPRE_Int * GetTrueDofOffsets()
Definition: pfespace.hpp:170
virtual void GetBdrElementDofs(int i, Array< int > &dofs) const
Returns indexes of degrees of freedom for i&#39;th boundary element.
Definition: pfespace.cpp:275
void GetFaceNbrElementVDofs(int i, Array< int > &vdofs) const
Definition: pfespace.cpp:842
HYPRE_Int GetGlobalNumRows() const
Definition: hypre.hpp:371
const FiniteElement * GetFaceNbrFE(int i) const
Definition: pfespace.cpp:876
void GetSharedFaceDofs(int group, int fi, Array< int > &dofs) const
Definition: pfespace.cpp:321
int dim
Definition: ex3.cpp:47
HYPRE_Int GetGlobalNumCols() const
Definition: hypre.hpp:374
Data type sparse matrix.
Definition: sparsemat.hpp:38
void LoseData()
NULL-ifies the data.
Definition: array.hpp:103
HYPRE_Int GetMyTDofOffset() const
Definition: pfespace.cpp:617
int GetLocalTDofNumber(int ldof)
Definition: pfespace.cpp:531
virtual void GetFaceDofs(int i, Array< int > &dofs) const
Definition: pfespace.cpp:289
Wrapper for hypre&#39;s parallel vector class.
Definition: hypre.hpp:72
virtual void GetEssentialTrueDofs(const Array< int > &bdr_attr_is_ess, Array< int > &ess_tdof_list)
Definition: pfespace.cpp:521
GroupCommunicator & GroupComm()
Return a reference to the internal GroupCommunicator (on VDofs)
Definition: pfespace.hpp:210
HYPRE_Int GetGlobalScalarTDofNumber(int sldof)
Definition: pfespace.cpp:575
bool Nonconforming() const
Definition: pfespace.hpp:261
virtual void GetElementDofs(int i, Array< int > &dofs) const
Returns indexes of degrees of freedom in array dofs for i&#39;th element.
Definition: pfespace.cpp:261
const FiniteElement * GetFaceNbrFaceFE(int i) const
Definition: pfespace.cpp:890
int GetFaceNbrVSize() const
Definition: pfespace.hpp:249
void GetFaceNbrFaceVDofs(int i, Array< int > &vdofs) const
Definition: pfespace.cpp:848
virtual void UpdatesFinished()
Free GridFunction transformation matrix (if any), to save memory.
Definition: fespace.hpp:367
virtual void GetEssentialVDofs(const Array< int > &bdr_attr_is_ess, Array< int > &ess_dofs) const
Determine the boundary degrees of freedom.
Definition: pfespace.cpp:508
NURBSExtension * NURBSext
Definition: fespace.hpp:87
GroupTopology gtopo
Definition: pmesh.hpp:120
Abstract operator.
Definition: operator.hpp:21
Wrapper for hypre&#39;s ParCSR matrix class.
Definition: hypre.hpp:174
virtual void UpdatesFinished()
Free ParGridFunction transformation matrix (if any), to save memory.
Definition: pfespace.hpp:268
HypreParMatrix * Dof_TrueDof_Matrix()
The true dof-to-dof interpolation matrix.
Definition: pfespace.hpp:193
GroupCommunicator * ScalarGroupComm()
Return a new GroupCommunicator on Dofs.
Definition: pfespace.cpp:471
Class for parallel meshes.
Definition: pmesh.hpp:28
GroupTopology gtopo
Definition: nurbs.hpp:350
ParFiniteElementSpace(ParMesh *pm, const FiniteElementCollection *f, int dim=1, int ordering=Ordering::byNODES)
Definition: pfespace.cpp:27
HYPRE_Int GetGlobalTDofNumber(int ldof)
Returns the global tdof number of the given local degree of freedom.
Definition: pfespace.cpp:552
Array< HYPRE_Int > face_nbr_glob_dof_map
Definition: pfespace.hpp:154