MFEM  v3.3.2
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 class ConformingProlongationOperator;
28 
29 
30 /// Abstract parallel finite element space.
32 {
33 private:
34  /// MPI data.
35  MPI_Comm MyComm;
36  int NRanks, MyRank;
37 
38  /// Parallel mesh.
39  ParMesh *pmesh;
40 
41  /// GroupCommunicator on the local VDofs
42  GroupCommunicator *gcomm;
43 
44  /// Number of true dofs in this processor (local true dofs).
45  mutable int ltdof_size;
46 
47  /// The group of each local dof.
48  Array<int> ldof_group;
49 
50  /// For a local dof: the local true dof number in the master of its group.
51  mutable Array<int> ldof_ltdof;
52 
53  /// Offsets for the dofs in each processor in global numbering.
54  mutable Array<HYPRE_Int> dof_offsets;
55 
56  /// Offsets for the true dofs in each processor in global numbering.
57  mutable Array<HYPRE_Int> tdof_offsets;
58 
59  /// Offsets for the true dofs in neighbor processor in global numbering.
60  mutable Array<HYPRE_Int> tdof_nb_offsets;
61 
62  /// Previous 'dof_offsets' (before Update()), column partition of T.
63  Array<HYPRE_Int> old_dof_offsets;
64 
65  /// The sign of the basis functions at the scalar local dofs.
66  Array<int> ldof_sign;
67 
68  /// The matrix P (interpolation from true dof to dof).
69  mutable HypreParMatrix *P;
71 
72  /// The (block-diagonal) matrix R (restriction of dof to true dof)
73  mutable SparseMatrix *R;
74 
75  ParNURBSExtension *pNURBSext() const
76  { return dynamic_cast<ParNURBSExtension *>(NURBSext); }
77 
78  GroupTopology &GetGroupTopo() const
79  { return (NURBSext) ? pNURBSext()->gtopo : pmesh->gtopo; }
80 
81  void Construct();
82  void Destroy();
83 
84  // ldof_type = 0 : DOFs communicator, otherwise VDOFs communicator
85  void GetGroupComm(GroupCommunicator &gcomm, int ldof_type,
86  Array<int> *ldof_sign = NULL);
87 
88  /// Construct dof_offsets and tdof_offsets using global communication.
89  void GenerateGlobalOffsets() const;
90 
91  /// Construct ldof_group and ldof_ltdof.
92  void ConstructTrueDofs();
93  void ConstructTrueNURBSDofs();
94 
95  void ApplyLDofSigns(Array<int> &dofs) const;
96 
97  /// Helper struct to store DOF dependencies in a parallel NC mesh.
98  struct Dependency
99  {
100  int rank, dof; ///< master DOF, may be on another processor
101  double coef;
102  Dependency(int r, int d, double c) : rank(r), dof(d), coef(c) {}
103  };
104 
105  /// Dependency list for a local vdof.
106  struct DepList
107  {
108  Array<Dependency> list;
109  int type; ///< 0 = independent, 1 = one-to-one (conforming), 2 = slave
110 
111  DepList() : type(0) {}
112 
113  bool IsTrueDof(int my_rank) const
114  { return type == 0 || (type == 1 && list[0].rank == my_rank); }
115  };
116 
117  void AddSlaveDependencies(DepList deps[], int master_rank,
118  const Array<int> &master_dofs, int master_ndofs,
119  const Array<int> &slave_dofs, DenseMatrix& I)
120  const;
121 
122  void Add1To1Dependencies(DepList deps[], int owner_rank,
123  const Array<int> &owner_dofs, int owner_ndofs,
124  const Array<int> &dependent_dofs) const;
125 
126  void GetDofs(int type, int index, Array<int>& dofs) const;
127  void ReorderFaceDofs(Array<int> &dofs, int orient) const;
128 
129  /// Build the P and R matrices.
130  void Build_Dof_TrueDof_Matrix() const;
131 
132  // Used when the ParMesh is non-conforming, i.e. pmesh->pncmesh != NULL.
133  // Constructs the matrices P and R. Determines ltdof_size. Calls
134  // GenerateGlobalOffsets(). Constructs ldof_ltdof.
135  void GetParallelConformingInterpolation() const;
136 
137  /** Calculate a GridFunction migration matrix after mesh load balancing.
138  The result is a parallel permutation matrix that can be used to update
139  all grid functions defined on this space. */
140  HypreParMatrix* RebalanceMatrix(int old_ndofs,
141  const Table* old_elem_dof);
142 
143  /** Calculate a GridFunction restriction matrix after mesh derefinement.
144  The matrix is constructed so that the new grid function interpolates
145  the original function, i.e., the original function is evaluated at the
146  nodes of the coarse function. */
147  HypreParMatrix* ParallelDerefinementMatrix(int old_ndofs,
148  const Table *old_elem_dof);
149 
150 public:
151  // Face-neighbor data
152  // Number of face-neighbor dofs
154  // Face-neighbor-element to face-neighbor dof
156  // Face-neighbor to ldof in the face-neighbor numbering
158  // The global ldof indices of the face-neighbor dofs
160  // Local face-neighbor data: face-neighbor to ldof
162 
164  int dim = 1, int ordering = Ordering::byNODES);
165 
166  MPI_Comm GetComm() const { return MyComm; }
167  int GetNRanks() const { return NRanks; }
168  int GetMyRank() const { return MyRank; }
169 
170  inline ParMesh *GetParMesh() { return pmesh; }
171 
172  int GetDofSign(int i)
173  { return NURBSext || Nonconforming() ? 1 : ldof_sign[VDofToDof(i)]; }
174  HYPRE_Int *GetDofOffsets() const { return dof_offsets; }
175  HYPRE_Int *GetTrueDofOffsets() const { return tdof_offsets; }
176  HYPRE_Int GlobalVSize() const
177  { return Dof_TrueDof_Matrix()->GetGlobalNumRows(); }
178  HYPRE_Int GlobalTrueVSize() const
179  { return Dof_TrueDof_Matrix()->GetGlobalNumCols(); }
180 
181  /// Return the number of local vector true dofs.
182  virtual int GetTrueVSize() const { return ltdof_size; }
183 
184  /// Returns indexes of degrees of freedom in array dofs for i'th element.
185  virtual void GetElementDofs(int i, Array<int> &dofs) const;
186 
187  /// Returns indexes of degrees of freedom for i'th boundary element.
188  virtual void GetBdrElementDofs(int i, Array<int> &dofs) const;
189 
190  /** Returns the indexes of the degrees of freedom for i'th face
191  including the dofs for the edges and the vertices of the face. */
192  virtual void GetFaceDofs(int i, Array<int> &dofs) const;
193 
194  void GetSharedEdgeDofs(int group, int ei, Array<int> &dofs) const;
195  void GetSharedFaceDofs(int group, int fi, Array<int> &dofs) const;
196 
197  /// The true dof-to-dof interpolation matrix
199  { if (!P) { Build_Dof_TrueDof_Matrix(); } return P; }
200 
201  /** @brief For a non-conforming mesh, construct and return the interpolation
202  matrix from the partially conforming true dofs to the local dofs. The
203  returned pointer must be deleted by the caller. */
205 
206  /** Create and return a new HypreParVector on the true dofs, which is
207  owned by (i.e. it must be destroyed by) the calling function. */
209  { return (new HypreParVector(MyComm,GlobalTrueVSize(),GetTrueDofOffsets()));}
210 
211  /// Scale a vector of true dofs
212  void DivideByGroupSize(double *vec);
213 
214  /// Return a reference to the internal GroupCommunicator (on VDofs)
215  GroupCommunicator &GroupComm() { return *gcomm; }
216 
217  /// Return a new GroupCommunicator on Dofs
219 
220  /** Given an integer array on the local degrees of freedom, perform
221  a bitwise OR between the shared dofs. */
222  void Synchronize(Array<int> &ldof_marker) const;
223 
224  /// Determine the boundary degrees of freedom
225  virtual void GetEssentialVDofs(const Array<int> &bdr_attr_is_ess,
226  Array<int> &ess_dofs,
227  int component = -1) const;
228 
229  /** Get a list of essential true dofs, ess_tdof_list, corresponding to the
230  boundary attributes marked in the array bdr_attr_is_ess. */
231  virtual void GetEssentialTrueDofs(const Array<int> &bdr_attr_is_ess,
232  Array<int> &ess_tdof_list,
233  int component = -1);
234 
235  /** If the given ldof is owned by the current processor, return its local
236  tdof number, otherwise return -1 */
237  int GetLocalTDofNumber(int ldof) const;
238  /// Returns the global tdof number of the given local degree of freedom
239  HYPRE_Int GetGlobalTDofNumber(int ldof) const;
240  /** Returns the global tdof number of the given local degree of freedom in
241  the scalar version of the current finite element space. The input should
242  be a scalar local dof. */
243  HYPRE_Int GetGlobalScalarTDofNumber(int sldof);
244 
245  HYPRE_Int GetMyDofOffset() const;
246  HYPRE_Int GetMyTDofOffset() const;
247 
248  virtual const Operator *GetProlongationMatrix();
249  /// Get the R matrix which restricts a local dof vector to true dof vector.
251  { Dof_TrueDof_Matrix(); return R; }
252 
253  // Face-neighbor functions
254  void ExchangeFaceNbrData();
255  int GetFaceNbrVSize() const { return num_face_nbr_dofs; }
256  void GetFaceNbrElementVDofs(int i, Array<int> &vdofs) const;
257  void GetFaceNbrFaceVDofs(int i, Array<int> &vdofs) const;
258  const FiniteElement *GetFaceNbrFE(int i) const;
259  const FiniteElement *GetFaceNbrFaceFE(int i) const;
260  const HYPRE_Int *GetFaceNbrGlobalDofMap() { return face_nbr_glob_dof_map; }
261 
263  void LoseDofOffsets() { dof_offsets.LoseData(); }
264  void LoseTrueDofOffsets() { tdof_offsets.LoseData(); }
265 
266  bool Conforming() const { return pmesh->pncmesh == NULL; }
267  bool Nonconforming() const { return pmesh->pncmesh != NULL; }
268 
269  /** Reflect changes in the mesh. Calculate one of the refinement/derefinement
270  /rebalance matrices, unless want_transform is false. */
271  virtual void Update(bool want_transform = true);
272 
273  /// Free ParGridFunction transformation matrix (if any), to save memory.
274  virtual void UpdatesFinished()
275  {
277  old_dof_offsets.DeleteAll();
278  }
279 
280  virtual ~ParFiniteElementSpace() { Destroy(); }
281 
282  // Obsolete, kept for backward compatibility
283  int TrueVSize() const { return ltdof_size; }
284 };
285 
286 
287 /// Auxiliary class used by ParFiniteElementSpace.
289 {
290 protected:
293 
294 public:
296 
297  virtual void Mult(const Vector &x, Vector &y) const;
298 
299  virtual void MultTranspose(const Vector &x, Vector &y) const;
300 };
301 
302 }
303 
304 #endif // MFEM_USE_MPI
305 
306 #endif
Abstract class for Finite Elements.
Definition: fe.hpp:47
virtual const Operator * GetProlongationMatrix()
Definition: pfespace.cpp:632
virtual void GetEssentialTrueDofs(const Array< int > &bdr_attr_is_ess, Array< int > &ess_tdof_list, int component=-1)
Definition: pfespace.cpp:529
HYPRE_Int GetGlobalTDofNumber(int ldof) const
Returns the global tdof number of the given local degree of freedom.
Definition: pfespace.cpp:562
HypreParVector * NewTrueDofVector()
Definition: pfespace.hpp:208
ConformingProlongationOperator(ParFiniteElementSpace &pfes)
Definition: pfespace.cpp:2299
virtual void GetEssentialVDofs(const Array< int > &bdr_attr_is_ess, Array< int > &ess_dofs, int component=-1) const
Determine the boundary degrees of freedom.
Definition: pfespace.cpp:515
HYPRE_Int * GetDofOffsets() const
Definition: pfespace.hpp:174
Auxiliary class used by ParFiniteElementSpace.
Definition: pfespace.hpp:288
Ordering::Type ordering
Definition: fespace.hpp:74
virtual void Update(bool want_transform=true)
Definition: pfespace.cpp:2225
int VDofToDof(int vdof) const
Definition: fespace.hpp:252
Data type dense matrix using column-major storage.
Definition: densemat.hpp:23
void DivideByGroupSize(double *vec)
Scale a vector of true dofs.
Definition: pfespace.cpp:460
Abstract parallel finite element space.
Definition: pfespace.hpp:31
void Synchronize(Array< int > &ldof_marker) const
Definition: pfespace.cpp:498
virtual const SparseMatrix * GetRestrictionMatrix()
Get the R matrix which restricts a local dof vector to true dof vector.
Definition: pfespace.hpp:250
int GetLocalTDofNumber(int ldof) const
Definition: pfespace.cpp:541
HypreParMatrix * GetPartialConformingInterpolation()
For a non-conforming mesh, construct and return the interpolation matrix from the partially conformin...
Definition: pfespace.cpp:1570
const HYPRE_Int * GetFaceNbrGlobalDofMap()
Definition: pfespace.hpp:260
HYPRE_Int GetMyDofOffset() const
Definition: pfespace.cpp:622
void DeleteAll()
Delete whole array.
Definition: array.hpp:633
void GetSharedEdgeDofs(int group, int ei, Array< int > &dofs) const
Definition: pfespace.cpp:303
ParNCMesh * pncmesh
Definition: pmesh.hpp:134
virtual int GetTrueVSize() const
Return the number of local vector true dofs.
Definition: pfespace.hpp:182
virtual void GetBdrElementDofs(int i, Array< int > &dofs) const
Returns indexes of degrees of freedom for i&#39;th boundary element.
Definition: pfespace.cpp:280
void GetFaceNbrElementVDofs(int i, Array< int > &vdofs) const
Definition: pfespace.cpp:865
HYPRE_Int GetGlobalNumRows() const
Definition: hypre.hpp:374
Communicator performing operations within groups defined by a GroupTopology with arbitrary-size data ...
const FiniteElement * GetFaceNbrFE(int i) const
Definition: pfespace.cpp:899
void GetSharedFaceDofs(int group, int fi, Array< int > &dofs) const
Definition: pfespace.cpp:326
int dim
Definition: ex3.cpp:47
MPI_Comm GetComm() const
Definition: pfespace.hpp:166
HYPRE_Int GetGlobalNumCols() const
Definition: hypre.hpp:377
Data type sparse matrix.
Definition: sparsemat.hpp:38
void LoseData()
NULL-ifies the data.
Definition: array.hpp:104
HypreParMatrix * Dof_TrueDof_Matrix() const
The true dof-to-dof interpolation matrix.
Definition: pfespace.hpp:198
HYPRE_Int GlobalTrueVSize() const
Definition: pfespace.hpp:178
HYPRE_Int GetMyTDofOffset() const
Definition: pfespace.cpp:627
virtual void GetFaceDofs(int i, Array< int > &dofs) const
Definition: pfespace.cpp:294
HYPRE_Int GlobalVSize() const
Definition: pfespace.hpp:176
Wrapper for hypre&#39;s parallel vector class.
Definition: hypre.hpp:73
virtual void MultTranspose(const Vector &x, Vector &y) const
Action of the transpose operator: y=A^t(x). The default behavior in class Operator is to generate an ...
Definition: pfespace.cpp:2369
GroupCommunicator & GroupComm()
Return a reference to the internal GroupCommunicator (on VDofs)
Definition: pfespace.hpp:215
HYPRE_Int GetGlobalScalarTDofNumber(int sldof)
Definition: pfespace.cpp:585
bool Nonconforming() const
Definition: pfespace.hpp:267
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:266
virtual void Mult(const Vector &x, Vector &y) const
Operator application: y=A(x).
Definition: pfespace.cpp:2344
const FiniteElement * GetFaceNbrFaceFE(int i) const
Definition: pfespace.cpp:913
int GetFaceNbrVSize() const
Definition: pfespace.hpp:255
HYPRE_Int * GetTrueDofOffsets() const
Definition: pfespace.hpp:175
void GetFaceNbrFaceVDofs(int i, Array< int > &vdofs) const
Definition: pfespace.cpp:871
Vector data type.
Definition: vector.hpp:41
virtual void UpdatesFinished()
Free GridFunction transformation matrix (if any), to save memory.
Definition: fespace.hpp:373
NURBSExtension * NURBSext
Definition: fespace.hpp:87
GroupTopology gtopo
Definition: pmesh.hpp:121
Abstract operator.
Definition: operator.hpp:21
Wrapper for hypre&#39;s ParCSR matrix class.
Definition: hypre.hpp:175
virtual void UpdatesFinished()
Free ParGridFunction transformation matrix (if any), to save memory.
Definition: pfespace.hpp:274
GroupCommunicator * ScalarGroupComm()
Return a new GroupCommunicator on Dofs.
Definition: pfespace.cpp:478
Class for parallel meshes.
Definition: pmesh.hpp:29
GroupTopology gtopo
Definition: nurbs.hpp:350
ParFiniteElementSpace(ParMesh *pm, const FiniteElementCollection *f, int dim=1, int ordering=Ordering::byNODES)
Definition: pfespace.cpp:27
Array< HYPRE_Int > face_nbr_glob_dof_map
Definition: pfespace.hpp:159