MFEM  v4.5.1
Finite element discretization library
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
doftrans.hpp
Go to the documentation of this file.
1 // Copyright (c) 2010-2022, Lawrence Livermore National Security, LLC. Produced
2 // at the Lawrence Livermore National Laboratory. All Rights reserved. See files
3 // LICENSE and NOTICE for details. LLNL-CODE-806117.
4 //
5 // This file is part of the MFEM library. For more information and source code
6 // availability visit https://mfem.org.
7 //
8 // MFEM is free software; you can redistribute it and/or modify it under the
9 // terms of the BSD-3 license. We welcome feedback and contributions, see file
10 // CONTRIBUTING.md for details.
11 
12 #ifndef MFEM_DOFTRANSFORM
13 #define MFEM_DOFTRANSFORM
14 
15 #include "../config/config.hpp"
16 #include "../linalg/linalg.hpp"
17 #include "intrules.hpp"
18 #include "fe.hpp"
19 
20 namespace mfem
21 {
22 
23 /** The DofTransformation class is an abstract base class for a family of
24  transformations that map local degrees of freedom (DoFs), contained within
25  individual elements, to global degrees of freedom, stored within
26  GridFunction objects. These transformations are necessary to ensure that
27  basis functions in neighboring elements align correctly. Closely related but
28  complementary transformations are required for the entries stored in
29  LinearForm and BilinearForm objects. The DofTransformation class is designed
30  to apply the action of both of these types of DoF transformations.
31 
32  Let the "primal transformation" be given by the operator T. This means that
33  given a local element vector v the data that must be placed into a
34  GridFunction object is v_t = T * v.
35 
36  We also need the inverse of the primal transformation T^{-1} so that we can
37  recover the local element vector from data read out of a GridFunction
38  e.g. v = T^{-1} * v_t.
39 
40  We need to preserve the action of our linear forms applied to primal
41  vectors. In other words, if f is the local vector computed by a linear
42  form then f * v = f_t * v_t (where "*" represents an inner product of
43  vectors). This requires that f_t = T^{-T} * f i.e. the "dual transform" is
44  given by the transpose of the inverse of the primal transformation.
45 
46  For bilinear forms we require that v^T * A * v = v_t^T * A_t * v_t. This
47  implies that A_t = T^{-T} * A * T^{-1}. This can be accomplished by
48  performing dual transformations of the rows and columns of the matrix A.
49 
50  For discrete linear operators the range must be modified with the primal
51  transformation rather than the dual transformation because the result is a
52  primal vector rather than a dual vector. This leads to the transformation
53  D_t = T * D * T^{-1}. This can be accomplished by using a primal
54  transformation on the columns of D and a dual transformation on its rows.
55 */
57 {
58 protected:
59  int size_;
60 
62 
64  : size_(size) {}
65 
66 public:
67 
68  inline int Size() const { return size_; }
69  inline int Height() const { return size_; }
70  inline int NumRows() const { return size_; }
71  inline int Width() const { return size_; }
72  inline int NumCols() const { return size_; }
73 
74  /** @brief Configure the transformation using face orientations for the
75  current element. */
76  /// The face_orientation array can be obtained from Mesh::GetElementFaces.
77  inline void SetFaceOrientations(const Array<int> & face_orientation)
78  { Fo = face_orientation; }
79 
80  inline const Array<int> & GetFaceOrientations() const { return Fo; }
81 
82  /** Transform local DoFs to align with the global DoFs. For example, this
83  transformation can be used to map the local vector computed by
84  FiniteElement::Project() to the transformed vector stored within a
85  GridFunction object. */
86  virtual void TransformPrimal(double *v) const = 0;
87  virtual void TransformPrimal(Vector &v) const;
88 
89  /// Transform groups of DoFs stored as dense matrices
90  virtual void TransformPrimalCols(DenseMatrix &V) const;
91 
92  /** Inverse transform local DoFs. Used to transform DoFs from a global vector
93  back to their element-local form. For example, this must be used to
94  transform the vector obtained using GridFunction::GetSubVector before it
95  can be used to compute a local interpolation.
96  */
97  virtual void InvTransformPrimal(double *v) const = 0;
98  virtual void InvTransformPrimal(Vector &v) const;
99 
100  /** Transform dual DoFs as computed by a LinearFormIntegrator before summing
101  into a LinearForm object. */
102  virtual void TransformDual(double *v) const = 0;
103  virtual void TransformDual(Vector &v) const;
104 
105  /** Inverse Transform dual DoFs */
106  virtual void InvTransformDual(double *v) const = 0;
107  virtual void InvTransformDual(Vector &v) const;
108 
109  /** Transform a matrix of dual DoFs entries as computed by a
110  BilinearFormIntegrator before summing into a BilinearForm object. */
111  virtual void TransformDual(DenseMatrix &V) const;
112 
113  /// Transform groups of dual DoFs stored as dense matrices
114  virtual void TransformDualRows(DenseMatrix &V) const;
115  virtual void TransformDualCols(DenseMatrix &V) const;
116 
117  virtual ~DofTransformation() {}
118 };
119 
120 /** Transform a matrix of DoFs entries from different finite element spaces as
121  computed by a DiscreteInterpolator before copying into a
122  DiscreteLinearOperator.
123 */
124 void TransformPrimal(const DofTransformation *ran_dof_trans,
125  const DofTransformation *dom_dof_trans,
126  DenseMatrix &elmat);
127 
128 /** Transform a matrix of dual DoFs entries from different finite element spaces
129  as computed by a BilinearFormIntegrator before summing into a
130  MixedBilinearForm object.
131 */
132 void TransformDual(const DofTransformation *ran_dof_trans,
133  const DofTransformation *dom_dof_trans,
134  DenseMatrix &elmat);
135 
136 /** The VDofTransformation class implements a nested transformation where an
137  arbitrary DofTransformation is replicated with a vdim >= 1.
138 */
140 {
141 private:
142  int vdim_;
143  int ordering_;
144  DofTransformation * doftrans_;
145 
146 public:
147  /** @brief Default constructor which requires that SetDofTransformation be
148  called before use. */
149  VDofTransformation(int vdim = 1, int ordering = 0)
150  : DofTransformation(0),
151  vdim_(vdim), ordering_(ordering),
152  doftrans_(NULL) {}
153 
154  /// Constructor with a known DofTransformation
155  VDofTransformation(DofTransformation & doftrans, int vdim = 1,
156  int ordering = 0)
157  : DofTransformation(vdim * doftrans.Size()),
158  vdim_(vdim), ordering_(ordering),
159  doftrans_(&doftrans) {}
160 
161  /// Set or change the vdim parameter
162  inline void SetVDim(int vdim)
163  {
164  vdim_ = vdim;
165  if (doftrans_)
166  {
167  size_ = vdim_ * doftrans_->Size();
168  }
169  }
170 
171  /// Return the current vdim value
172  inline int GetVDim() const { return vdim_; }
173 
174  /// Set or change the nested DofTransformation object
175  inline void SetDofTransformation(DofTransformation & doftrans)
176  {
177  size_ = vdim_ * doftrans.Size();
178  doftrans_ = &doftrans;
179  }
180 
181  /// Return the nested DofTransformation object
182  inline DofTransformation * GetDofTransformation() const { return doftrans_; }
183 
184  inline void SetFaceOrientation(const Array<int> & face_orientation)
185  { Fo = face_orientation; doftrans_->SetFaceOrientations(face_orientation); }
186 
191 
192  void TransformPrimal(double *v) const;
193  void InvTransformPrimal(double *v) const;
194  void TransformDual(double *v) const;
195  void InvTransformDual(double *v) const;
196 };
197 
198 /** Abstract base class for high-order Nedelec spaces on elements with
199  triangular faces.
200 
201  The Nedelec DoFs on the interior of triangular faces come in pairs which
202  share an interpolation point but have different vector directions. These
203  directions depend on the orientation of the face and can therefore differ in
204  neighboring elements. The mapping required to transform these DoFs can be
205  implemented as series of 2x2 linear transformations. The raw data for these
206  linear transformations is stored in the T_data and TInv_data arrays and can
207  be accessed as DenseMatrices using the GetFaceTransform() and
208  GetFaceInverseTransform() methods.
209 */
211 {
212 protected:
213  static const double T_data[24];
214  static const double TInv_data[24];
215  static const DenseTensor T, TInv;
216  int order;
217  int nedofs; // number of DoFs per edge
218  int nfdofs; // number of DoFs per face
219 
220  ND_DofTransformation(int size, int order);
221 
222 public:
223  // Return the 2x2 transformation operator for the given face orientation
224  static const DenseMatrix & GetFaceTransform(int ori) { return T(ori); }
225 
226  // Return the 2x2 inverse transformation operator
227  static const DenseMatrix & GetFaceInverseTransform(int ori)
228  { return TInv(ori); }
229 };
230 
231 /// DoF transformation implementation for the Nedelec basis on triangles
233 {
234 public:
236 
240 
241  void TransformPrimal(double *v) const;
242 
243  void InvTransformPrimal(double *v) const;
244 
245  void TransformDual(double *v) const;
246 
247  void InvTransformDual(double *v) const;
249 };
250 
251 /// DoF transformation implementation for the Nedelec basis on tetrahedra
253 {
254 public:
256 
261 
262  void TransformPrimal(double *v) const;
263 
264  void InvTransformPrimal(double *v) const;
265 
266  void TransformDual(double *v) const;
267 
268  void InvTransformDual(double *v) const;
269 };
270 
271 /// DoF transformation implementation for the Nedelec basis on wedge elements
273 {
274 public:
276 
281 
282  void TransformPrimal(double *v) const;
283 
284  void InvTransformPrimal(double *v) const;
285 
286  void TransformDual(double *v) const;
287 
288  void InvTransformDual(double *v) const;
289 
290 };
291 
292 } // namespace mfem
293 
294 #endif // MFEM_DOFTRANSFORM
void InvTransformPrimal(double *v) const
Definition: doftrans.cpp:294
void TransformDual(double *v) const
Definition: doftrans.cpp:413
static const DenseTensor TInv
Definition: doftrans.hpp:215
void TransformPrimal(double *v) const
Definition: doftrans.cpp:463
void InvTransformDual(double *v) const
Definition: doftrans.cpp:532
int Height() const
Definition: doftrans.hpp:69
VDofTransformation(int vdim=1, int ordering=0)
Default constructor which requires that SetDofTransformation be called before use.
Definition: doftrans.hpp:149
virtual ~DofTransformation()
Definition: doftrans.hpp:117
void TransformDual(double *v) const
Definition: doftrans.cpp:174
Data type dense matrix using column-major storage.
Definition: densemat.hpp:23
int NumRows() const
Definition: doftrans.hpp:70
DofTransformation(int size)
Definition: doftrans.hpp:63
ND_DofTransformation(int size, int order)
Definition: doftrans.cpp:258
void InvTransformDual(double *v) const
Definition: doftrans.cpp:436
void InvTransformDual(double *v) const
Definition: doftrans.cpp:203
virtual void TransformPrimalCols(DenseMatrix &V) const
Transform groups of DoFs stored as dense matrices.
Definition: doftrans.cpp:22
int GetVDim() const
Return the current vdim value.
Definition: doftrans.hpp:172
void TransformDual(double *v) const
Definition: doftrans.cpp:509
void TransformPrimal(const DofTransformation *ran_dof_trans, const DofTransformation *dom_dof_trans, DenseMatrix &elmat)
Definition: doftrans.cpp:65
void TransformDual(double *v) const
Definition: doftrans.cpp:317
void SetDofTransformation(DofTransformation &doftrans)
Set or change the nested DofTransformation object.
Definition: doftrans.hpp:175
void SetFaceOrientations(const Array< int > &face_orientation)
Configure the transformation using face orientations for the current element.
Definition: doftrans.hpp:77
virtual void TransformDualCols(DenseMatrix &V) const
Definition: doftrans.cpp:52
static const double TInv_data[24]
Definition: doftrans.hpp:214
void TransformPrimal(double *v) const
Definition: doftrans.cpp:367
static const DenseTensor T
Definition: doftrans.hpp:215
int NumCols() const
Definition: doftrans.hpp:72
void InvTransformPrimal(double *v) const
Definition: doftrans.cpp:145
DoF transformation implementation for the Nedelec basis on triangles.
Definition: doftrans.hpp:232
virtual void InvTransformPrimal(double *v) const =0
static const DenseMatrix & GetFaceInverseTransform(int ori)
Definition: doftrans.hpp:227
virtual void TransformPrimal(double *v) const =0
void SetVDim(int vdim)
Set or change the vdim parameter.
Definition: doftrans.hpp:162
void TransformPrimal(double *v) const
Definition: doftrans.cpp:271
virtual void TransformDualRows(DenseMatrix &V) const
Transform groups of dual DoFs stored as dense matrices.
Definition: doftrans.cpp:41
void SetFaceOrientation(const Array< int > &face_orientation)
Definition: doftrans.hpp:184
const Array< int > & GetFaceOrientations() const
Definition: doftrans.hpp:80
static const double T_data[24]
Definition: doftrans.hpp:213
VDofTransformation(DofTransformation &doftrans, int vdim=1, int ordering=0)
Constructor with a known DofTransformation.
Definition: doftrans.hpp:155
static const DenseMatrix & GetFaceTransform(int ori)
Definition: doftrans.hpp:224
Vector data type.
Definition: vector.hpp:60
void TransformDual(const DofTransformation *ran_dof_trans, const DofTransformation *dom_dof_trans, DenseMatrix &elmat)
Definition: doftrans.cpp:93
void InvTransformPrimal(double *v) const
Definition: doftrans.cpp:486
void InvTransformPrimal(double *v) const
Definition: doftrans.cpp:390
void InvTransformDual(double *v) const
Definition: doftrans.cpp:340
virtual void InvTransformDual(double *v) const =0
DoF transformation implementation for the Nedelec basis on tetrahedra.
Definition: doftrans.hpp:252
Rank 3 tensor (array of matrices)
Definition: densemat.hpp:953
int Width() const
Definition: doftrans.hpp:71
virtual void TransformDual(double *v) const =0
DofTransformation * GetDofTransformation() const
Return the nested DofTransformation object.
Definition: doftrans.hpp:182
DoF transformation implementation for the Nedelec basis on wedge elements.
Definition: doftrans.hpp:272
void TransformPrimal(double *v) const
Definition: doftrans.cpp:116