MFEM  v4.6.0
Finite element discretization library
restriction.hpp
Go to the documentation of this file.
1 // Copyright (c) 2010-2023, 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_RESTRICTION
13 #define MFEM_RESTRICTION
14 
15 #include "../linalg/operator.hpp"
16 #include "../mesh/mesh.hpp"
17 
18 namespace mfem
19 {
20 
21 class FiniteElementSpace;
22 enum class ElementDofOrdering;
23 
24 /// Abstract base class that defines an interface for element restrictions.
26 {
27 public:
28  /// @brief Add the E-vector degrees of freedom @a x to the L-vector degrees
29  /// of freedom @a y.
30  void AddMultTranspose(const Vector &x, Vector &y,
31  const double a = 1.0) const override = 0;
32 };
33 
34 /// Operator that converts FiniteElementSpace L-vectors to E-vectors.
35 /** Objects of this type are typically created and owned by FiniteElementSpace
36  objects, see FiniteElementSpace::GetElementRestriction(). */
38 {
39 private:
40  /** This number defines the maximum number of elements any dof can belong to
41  for the FillSparseMatrix method. */
42  static const int MaxNbNbr = 16;
43 
44 protected:
46  const int ne;
47  const int vdim;
48  const bool byvdim;
49  const int ndofs;
50  const int dof;
51  const int nedofs;
55 
56 public:
58  void Mult(const Vector &x, Vector &y) const override;
59  void MultTranspose(const Vector &x, Vector &y) const override;
60  void AddMultTranspose(const Vector &x, Vector &y,
61  const double a = 1.0) const override;
62 
63  /// Compute Mult without applying signs based on DOF orientations.
64  void MultUnsigned(const Vector &x, Vector &y) const;
65  /// Compute MultTranspose without applying signs based on DOF orientations.
66  void MultTransposeUnsigned(const Vector &x, Vector &y) const;
67 
68  /// Compute MultTranspose by setting (rather than adding) element
69  /// contributions; this is a left inverse of the Mult() operation
70  void MultLeftInverse(const Vector &x, Vector &y) const;
71 
72  /// @brief Fills the E-vector y with `boolean` values 0.0 and 1.0 such that each
73  /// each entry of the L-vector is uniquely represented in `y`.
74  /** This means, the sum of the E-vector `y` is equal to the sum of the
75  corresponding L-vector filled with ones. The boolean mask is required to
76  emulate SetSubVector and its transpose on GPUs. This method is running on
77  the host, since the `processed` array requires a large shared memory. */
78  void BooleanMask(Vector& y) const;
79 
80  /// Fill a Sparse Matrix with Element Matrices.
81  void FillSparseMatrix(const Vector &mat_ea, SparseMatrix &mat) const;
82 
83  /** Fill the I array of SparseMatrix corresponding to the sparsity pattern
84  given by this ElementRestriction. */
85  int FillI(SparseMatrix &mat) const;
86  /** Fill the J and Data arrays of SparseMatrix corresponding to the sparsity
87  pattern given by this ElementRestriction, and the values of ea_data. */
88  void FillJAndData(const Vector &ea_data, SparseMatrix &mat) const;
89  /// @private Not part of the public interface (device kernel limitation).
90  ///
91  /// Performs either MultTranspose or AddMultTranspose depending on the
92  /// boolean template parameter @a ADD.
93  template <bool ADD> void TAddMultTranspose(const Vector &x, Vector &y) const;
94 
95  /// @name Low-level access to the underlying element-dof mappings
96  ///@{
97  const Array<int> &GatherMap() const { return gather_map; }
98  const Array<int> &Indices() const { return indices; }
99  const Array<int> &Offsets() const { return offsets; }
100  ///@}
101 };
102 
103 /// Operator that converts L2 FiniteElementSpace L-vectors to E-vectors.
104 /** Objects of this type are typically created and owned by FiniteElementSpace
105  objects, see FiniteElementSpace::GetElementRestriction(). L-vectors
106  corresponding to grid functions in L2 finite element spaces differ from
107  E-vectors only in the ordering of the degrees of freedom. */
109 {
110  const int ne;
111  const int vdim;
112  const bool byvdim;
113  const int ndof;
114  const int ndofs;
115 public:
117  void Mult(const Vector &x, Vector &y) const override;
118  void MultTranspose(const Vector &x, Vector &y) const override;
119  void AddMultTranspose(const Vector &x, Vector &y,
120  const double a = 1.0) const override;
121  /** Fill the I array of SparseMatrix corresponding to the sparsity pattern
122  given by this ElementRestriction. */
123  void FillI(SparseMatrix &mat) const;
124  /** Fill the J and Data arrays of SparseMatrix corresponding to the sparsity
125  pattern given by this L2FaceRestriction, and the values of ea_data. */
126  void FillJAndData(const Vector &ea_data, SparseMatrix &mat) const;
127  /// @private Not part of the public interface (device kernel limitation).
128  ///
129  /// Performs either MultTranspose or AddMultTranspose depending on the
130  /// boolean template parameter @a ADD.
131  template <bool ADD> void TAddMultTranspose(const Vector &x, Vector &y) const;
132 };
133 
134 /** An enum type to specify if only e1 value is requested (SingleValued) or both
135  e1 and e2 (DoubleValued). */
137 
138 /** @brief Base class for operators that extracts Face degrees of freedom.
139 
140  In order to compute quantities on the faces of a mesh, it is often useful to
141  extract the degrees of freedom on the faces of the elements. This class
142  provides an interface for such operations.
143 
144  If the FiniteElementSpace is ordered by Ordering::byVDIM, then the expected
145  format for the L-vector is (vdim x ndofs), otherwise if Ordering::byNODES
146  the expected format is (ndofs x vdim), where ndofs is the total number of
147  degrees of freedom.
148  Since FiniteElementSpace can either be continuous or discontinuous, the
149  degrees of freedom on a face can either be single valued or double valued,
150  this is what we refer to as the multiplicity and is represented by the
151  L2FaceValues enum type.
152  The format of the output face E-vector of degrees of freedom is
153  (face_dofs x vdim x multiplicity x nfaces), where face_dofs is the number of
154  degrees of freedom on each face, and nfaces the number of faces of the
155  requested FaceType (see FiniteElementSpace::GetNFbyType).
156 
157  @note Objects of this type are typically created and owned by
158  FiniteElementSpace objects, see FiniteElementSpace::GetFaceRestriction(). */
159 class FaceRestriction : public Operator
160 {
161 public:
163 
164  FaceRestriction(int h, int w): Operator(h, w) { }
165 
166  virtual ~FaceRestriction() { }
167 
168  /** @brief Extract the face degrees of freedom from @a x into @a y.
169 
170  @param[in] x The L-vector of degrees of freedom.
171  @param[out] y The degrees of freedom on the face, corresponding to a face
172  E-vector.
173  */
174  void Mult(const Vector &x, Vector &y) const override = 0;
175 
176  /** @brief Add the face degrees of freedom @a x to the element degrees of
177  freedom @a y.
178 
179  @param[in] x The face degrees of freedom on the face.
180  @param[in,out] y The L-vector of degrees of freedom to which we add the
181  face degrees of freedom.
182  @param[in] a Scalar coefficient for addition.
183  */
184  virtual void AddMultTranspose(const Vector &x, Vector &y,
185  const double a = 1.0) const override = 0;
186 
187  /** @brief Add the face degrees of freedom @a x to the element degrees of
188  freedom @a y ignoring the signs from DOF orientation. */
189  virtual void AddMultTransposeUnsigned(const Vector &x, Vector &y,
190  const double a = 1.0) const
191  {
192  AddMultTranspose(x, y, a);
193  }
194 
195  /** @brief Add the face degrees of freedom @a x to the element degrees of
196  freedom @a y. Perform the same computation as AddMultTranspose, but
197  @a x is invalid after calling this method.
198 
199  @param[in,out] x The face degrees of freedom on the face.
200  @param[in,out] y The L-vector of degrees of freedom to which we add the
201  face degrees of freedom.
202 
203  @note This method is an optimization of AddMultTranspose where the @a x
204  Vector is used and modified to avoid memory allocation and memcpy.
205  */
206  virtual void AddMultTransposeInPlace(Vector &x, Vector &y) const
207  {
208  AddMultTranspose(x, y);
209  }
210 
211  /** @brief Set the face degrees of freedom in the element degrees of freedom
212  @a y to the values given in @a x.
213 
214  @param[in] x The face degrees of freedom on the face.
215  @param[in,out] y The L-vector of degrees of freedom to which we add the
216  face degrees of freedom.
217  */
218  void MultTranspose(const Vector &x, Vector &y) const override
219  {
220  y = 0.0;
221  AddMultTranspose(x, y);
222  }
223 };
224 
225 /// @brief Operator that extracts face degrees of freedom for H1, ND, or RT
226 /// FiniteElementSpaces.
227 ///
228 /// Objects of this type are typically created and owned by FiniteElementSpace
229 /// objects, see FiniteElementSpace::GetFaceRestriction().
231 {
232 protected:
234  const int nf; // Number of faces of the requested type
235  const int vdim;
236  const bool byvdim;
237  const int face_dofs; // Number of dofs on each face
238  const int elem_dofs; // Number of dofs in each element
239  const int nfdofs; // Total number of face E-vector dofs
240  const int ndofs; // Total number of dofs
241  Array<int> scatter_indices; // Scattering indices for element 1 on each face
242  Array<int> gather_offsets; // offsets for the gathering indices of each dof
243  Array<int> gather_indices; // gathering indices for each dof
244  Array<int> vol_dof_map; // mapping from lexicographic to native ordering
245 
246  /** @brief Construct a ConformingFaceRestriction.
247 
248  @param[in] fes The FiniteElementSpace on which this operates
249  @param[in] f_ordering Request a specific face dof ordering
250  @param[in] type Request internal or boundary faces dofs
251  @param[in] build Request the NCL2FaceRestriction to compute the
252  scatter/gather indices. False should only be used
253  when inheriting from ConformingFaceRestriction.
254  */
256  const ElementDofOrdering f_ordering,
257  const FaceType type,
258  bool build);
259 public:
260  /** @brief Construct a ConformingFaceRestriction.
261 
262  @param[in] fes The FiniteElementSpace on which this operates
263  @param[in] f_ordering Request a specific face dof ordering
264  @param[in] type Request internal or boundary faces dofs */
266  const ElementDofOrdering f_ordering,
267  const FaceType type);
268 
269  /** @brief Scatter the degrees of freedom, i.e. goes from L-Vector to
270  face E-Vector.
271 
272  @param[in] x The L-vector degrees of freedom.
273  @param[out] y The face E-Vector degrees of freedom with the given format:
274  face_dofs x vdim x nf
275  where nf is the number of interior or boundary faces
276  requested by @a type in the constructor.
277  The face_dofs are ordered according to the given
278  ElementDofOrdering. */
279  void Mult(const Vector &x, Vector &y) const override;
280 
282 
283  /** @brief Gather the degrees of freedom, i.e. goes from face E-Vector to
284  L-Vector.
285 
286  @param[in] x The face E-Vector degrees of freedom with the given format:
287  face_dofs x vdim x nf
288  where nf is the number of interior or boundary faces
289  requested by @a type in the constructor.
290  The face_dofs should be ordered according to the given
291  ElementDofOrdering
292  @param[in,out] y The L-vector degrees of freedom.
293  @param[in] a Scalar coefficient for addition. */
294  void AddMultTranspose(const Vector &x, Vector &y,
295  const double a = 1.0) const override;
296 
297  /** @brief Gather the degrees of freedom, i.e. goes from face E-Vector to
298  L-Vector @b not taking into account signs from DOF orientations.
299 
300  @sa AddMultTranspose(). */
301  void AddMultTransposeUnsigned(const Vector &x, Vector &y,
302  const double a = 1.0) const override;
303 
304 private:
305  /** @brief Compute the scatter indices: L-vector to E-vector, and the offsets
306  for the gathering: E-vector to L-vector.
307 
308  @param[in] f_ordering Request a specific face dof ordering.
309  @param[in] type Request internal or boundary faces dofs.
310  */
311  void ComputeScatterIndicesAndOffsets(const ElementDofOrdering f_ordering,
312  const FaceType type);
313 
314  /** @brief Compute the gather indices: E-vector to L-vector.
315 
316  Note: Requires the gather offsets to be computed.
317 
318  @param[in] f_ordering Request a specific face dof ordering.
319  @param[in] type Request internal or boundary faces dofs.
320  */
321  void ComputeGatherIndices(const ElementDofOrdering f_ordering,
322  const FaceType type);
323 
324 protected:
325  mutable Array<int> face_map; // Used in the computation of GetFaceDofs
326 
327  /** @brief Verify that ConformingFaceRestriction is built from a supported
328  finite element space.
329 
330  @param[in] f_ordering The requested face dof ordering.
331  */
332  void CheckFESpace(const ElementDofOrdering f_ordering);
333 
334  /** @brief Set the scattering indices of elem1, and increment the offsets for
335  the face described by the @a face.
336 
337  @param[in] face The face information of the current face.
338  @param[in] face_index The interior/boundary face index.
339  @param[in] f_ordering Request a specific face dof ordering.
340  */
342  const int face_index,
343  const ElementDofOrdering f_ordering);
344 
345  /** @brief Set the gathering indices of elem1 for the interior face described
346  by the @a face.
347 
348  @param[in] face The face information of the current face.
349  @param[in] face_index The interior/boundary face index.
350  @param[in] f_ordering Request a specific face dof ordering.
351  */
353  const int face_index,
354  const ElementDofOrdering f_ordering);
355 };
356 
357 /// @brief Alias for ConformingFaceRestriction, for backwards compatibility and
358 /// as base class for ParNCH1FaceRestriction.
360 
361 /// Operator that extracts Face degrees of freedom for L2 spaces.
362 /** Objects of this type are typically created and owned by FiniteElementSpace
363  objects, see FiniteElementSpace::GetFaceRestriction(). */
365 {
366 protected:
368  const int nf; // Number of faces of the requested type
369  const int ne; // Number of elements
370  const int vdim; // vdim
371  const bool byvdim;
372  const int face_dofs; // Number of dofs on each face
373  const int elem_dofs; // Number of dofs in each element
374  const int nfdofs; // Total number of dofs on the faces
375  const int ndofs; // Total number of dofs
376  const FaceType type;
378  Array<int> scatter_indices1; // Scattering indices for element 1 on each face
379  Array<int> scatter_indices2; // Scattering indices for element 2 on each face
380  Array<int> gather_offsets; // offsets for the gathering indices of each dof
381  Array<int> gather_indices; // gathering indices for each dof
382 
383  /** @brief Constructs an L2FaceRestriction.
384 
385  @param[in] fes The FiniteElementSpace on which this operates
386  @param[in] f_ordering Request a specific face dof ordering
387  @param[in] type Request internal or boundary faces dofs
388  @param[in] m Request the face dofs for elem1, or both elem1 and
389  elem2
390  @param[in] build Request the NCL2FaceRestriction to compute the
391  scatter/gather indices. False should only be used
392  when inheriting from L2FaceRestriction.
393  */
395  const ElementDofOrdering f_ordering,
396  const FaceType type,
397  const L2FaceValues m,
398  bool build);
399 
400 public:
401  /** @brief Constructs an L2FaceRestriction.
402 
403  @param[in] fes The FiniteElementSpace on which this operates
404  @param[in] f_ordering Request a specific face dof ordering
405  @param[in] type Request internal or boundary faces dofs
406  @param[in] m Request the face dofs for elem1, or both elem1 and
407  elem2 */
409  const ElementDofOrdering f_ordering,
410  const FaceType type,
412 
413  /** @brief Scatter the degrees of freedom, i.e. goes from L-Vector to
414  face E-Vector.
415 
416  @param[in] x The L-vector degrees of freedom.
417  @param[out] y The face E-Vector degrees of freedom with the given format:
418  if L2FacesValues::DoubleValued (face_dofs x vdim x 2 x nf)
419  if L2FacesValues::SingleValued (face_dofs x vdim x nf)
420  where nf is the number of interior or boundary faces
421  requested by @a type in the constructor.
422  The face_dofs are ordered according to the given
423  ElementDofOrdering. */
424  void Mult(const Vector &x, Vector &y) const override;
425 
427 
428  /** @brief Gather the degrees of freedom, i.e. goes from face E-Vector to
429  L-Vector.
430 
431  @param[in] x The face E-Vector degrees of freedom with the given format:
432  if L2FacesValues::DoubleValued (face_dofs x vdim x 2 x nf)
433  if L2FacesValues::SingleValued (face_dofs x vdim x nf)
434  where nf is the number of interior or boundary faces
435  requested by @a type in the constructor.
436  The face_dofs should be ordered according to the given
437  ElementDofOrdering
438  @param[in,out] y The L-vector degrees of freedom.
439  @param[in] a Scalar coefficient for addition. */
440  void AddMultTranspose(const Vector &x, Vector &y,
441  const double a = 1.0) const override;
442 
443  /** @brief Fill the I array of SparseMatrix corresponding to the sparsity
444  pattern given by this L2FaceRestriction.
445 
446  @param[in,out] mat The sparse matrix for which we want to initialize the
447  row offsets.
448  @param[in] keep_nbr_block When set to true the SparseMatrix will
449  include the rows (in addition to the columns)
450  corresponding to face-neighbor dofs. The
451  default behavior is to disregard those rows. */
452  virtual void FillI(SparseMatrix &mat, const bool keep_nbr_block = false) const;
453 
454  /** @brief Fill the J and Data arrays of the SparseMatrix corresponding to
455  the sparsity pattern given by this L2FaceRestriction, and the values of
456  fea_data.
457 
458  @param[in] fea_data The dense matrices representing the local operators
459  on each face. The format is:
460  face_dofs x face_dofs x 2 x nf
461  On each face the first local matrix corresponds to
462  the contribution of elem1 on elem2, and the second to
463  the contribution of elem2 on elem1.
464  @param[in,out] mat The sparse matrix that is getting filled.
465  @param[in] keep_nbr_block When set to true the SparseMatrix will
466  include the rows (in addition to the columns)
467  corresponding to face-neighbor dofs. The
468  default behavior is to disregard those rows. */
469  virtual void FillJAndData(const Vector &fea_data,
470  SparseMatrix &mat,
471  const bool keep_nbr_block = false) const;
472 
473  /** @brief This methods adds the DG face matrices to the element matrices.
474 
475  @param[in] fea_data The dense matrices representing the local operators
476  on each face. The format is:
477  face_dofs x face_dofs x 2 x nf
478  On each face the first and second local matrices
479  correspond to the contributions of elem1 and elem2 on
480  themselves respectively.
481  @param[in,out] ea_data The dense matrices representing the element local
482  contributions for each element to which will be
483  added the face contributions.
484  The format is: dofs x dofs x ne, where dofs is the
485  number of dofs per element and ne the number of
486  elements. */
487  virtual void AddFaceMatricesToElementMatrices(const Vector &fea_data,
488  Vector &ea_data) const;
489 
490 private:
491  /** @brief Compute the scatter indices: L-vector to E-vector, and the offsets
492  for the gathering: E-vector to L-vector.
493 
494  @param[in] f_ordering Request a specific face dof ordering.
495  @param[in] type Request internal or boundary faces dofs.
496  */
497  void ComputeScatterIndicesAndOffsets(const ElementDofOrdering f_ordering,
498  const FaceType type);
499 
500  /** @brief Compute the gather indices: E-vector to L-vector.
501 
502  Note: Requires the gather offsets to be computed.
503 
504  @param[in] f_ordering Request a specific face dof ordering.
505  @param[in] type Request internal or boundary faces dofs.
506  */
507  void ComputeGatherIndices(const ElementDofOrdering f_ordering,
508  const FaceType type);
509 
510 protected:
511  mutable Array<int> face_map; // Used in the computation of GetFaceDofs
512 
513  /** @brief Verify that L2FaceRestriction is built from an L2 FESpace.
514 
515  @param[in] f_ordering The requested face dof ordering.
516  */
517  void CheckFESpace(const ElementDofOrdering f_ordering);
518 
519  /** @brief Set the scattering indices of elem1, and increment the offsets for
520  the face described by the @a face. The ordering of the face dofs of elem1
521  is lexicographic relative to elem1.
522 
523  @param[in] face The face information of the current face.
524  @param[in] face_index The interior/boundary face index.
525  */
527  const int face_index);
528 
529  /** @brief Permute and set the scattering indices of elem2, and increment the
530  offsets for the face described by the @a face. The permutation orders the
531  dofs of elem2 lexicographically as the ones of elem1.
532 
533  @param[in] face The face information of the current face.
534  @param[in] face_index The interior/boundary face index.
535  */
537  const int face_index);
538 
539  /** @brief Permute and set the scattering indices of elem2 for the shared
540  face described by the @a face. The permutation orders the dofs of elem2 as
541  the ones of elem1.
542 
543  @param[in] face The face information of the current face.
544  @param[in] face_index The interior/boundary face index.
545  */
547  const Mesh::FaceInformation &face,
548  const int face_index);
549 
550  /** @brief Set the scattering indices of elem2 for the boundary face
551  described by the @a face.
552 
553  @param[in] face The face information of the current face.
554  @param[in] face_index The interior/boundary face index.
555  */
557  const int face_index);
558 
559  /** @brief Set the gathering indices of elem1 for the interior face described
560  by the @a face.
561 
562  Note: This function modifies the offsets.
563 
564  @param[in] face The face information of the current face.
565  @param[in] face_index The interior/boundary face index.
566  */
568  const int face_index);
569 
570  /** @brief Permute and set the gathering indices of elem2 for the interior
571  face described by the @a face. The permutation orders the dofs of elem2 as
572  the ones of elem1.
573 
574  Note: This function modifies the offsets.
575 
576  @param[in] face The face information of the current face.
577  @param[in] face_index The interior/boundary face index.
578  */
580  const int face_index);
581 
582 public:
583  /** @brief Scatter the degrees of freedom, i.e. goes from L-Vector to
584  face E-Vector. Should only be used with conforming faces and when:
585  m == L2FacesValues::SingleValued
586 
587  @param[in] x The L-vector degrees of freedom.
588  @param[out] y The face E-Vector degrees of freedom with the given format:
589  face_dofs x vdim x nf
590  where nf is the number of interior or boundary faces
591  requested by @a type in the constructor.
592  The face_dofs are ordered according to the given
593  ElementDofOrdering. */
594  void SingleValuedConformingMult(const Vector& x, Vector& y) const;
595 
596  /** @brief Scatter the degrees of freedom, i.e. goes from L-Vector to
597  face E-Vector. Should only be used with conforming faces and when:
598  m == L2FacesValues::DoubleValued
599 
600  @param[in] x The L-vector degrees of freedom.
601  @param[out] y The face E-Vector degrees of freedom with the given format:
602  face_dofs x vdim x 2 x nf
603  where nf is the number of interior or boundary faces
604  requested by @a type in the constructor.
605  The face_dofs are ordered according to the given
606  ElementDofOrdering. */
607  virtual void DoubleValuedConformingMult(const Vector& x, Vector& y) const;
608 
609  /** @brief Gather the degrees of freedom, i.e. goes from face E-Vector to
610  L-Vector. Should only be used with conforming faces and when:
611  m == L2FacesValues::SingleValued
612 
613  @param[in] x The face E-Vector degrees of freedom with the given format:
614  face_dofs x vdim x nf
615  where nf is the number of interior or boundary faces
616  requested by @a type in the constructor.
617  The face_dofs should be ordered according to the given
618  ElementDofOrdering
619  @param[in,out] y The L-vector degrees of freedom. */
620  void SingleValuedConformingAddMultTranspose(const Vector& x, Vector& y) const;
621 
622  /** @brief Gather the degrees of freedom, i.e. goes from face E-Vector to
623  L-Vector. Should only be used with conforming faces and when:
624  m == L2FacesValues::DoubleValued
625 
626  @param[in] x The face E-Vector degrees of freedom with the given format:
627  face_dofs x vdim x 2 x nf
628  where nf is the number of interior or boundary faces
629  requested by @a type in the constructor.
630  The face_dofs should be ordered according to the given
631  ElementDofOrdering
632  @param[in,out] y The L-vector degrees of freedom. */
633  void DoubleValuedConformingAddMultTranspose(const Vector& x, Vector& y) const;
634 };
635 
636 /** This struct stores which side is the master nonconforming side and the
637  index of the interpolator, see InterpolationManager class below. */
639 {
640  uint32_t is_non_conforming : 1;
641  uint32_t master_side : 1;
642  uint32_t index : 30;
643 
644  // default constructor, create a conforming face with index 0.
645  InterpConfig() = default;
646 
647  // Non-conforming face
648  InterpConfig(int master_side, int nc_index)
650  { }
651 
652  InterpConfig(const InterpConfig&) = default;
653 
654  InterpConfig &operator=(const InterpConfig &rhs) = default;
655 };
656 
657 /** This struct stores which side is the master nonconforming side and the
658  index of the interpolator, see InterpolationManager class below. */
660 {
662  uint32_t is_non_conforming : 1;
663  uint32_t master_side : 1;
664  uint32_t index : 30;
665 
666  // default constructor.
667  NCInterpConfig() = default;
668 
669  // Non-conforming face
670  NCInterpConfig(int face_index, int master_side, int nc_index)
674  index(nc_index)
675  { }
676 
677  // Non-conforming face
681  master_side(config.master_side),
682  index(config.index)
683  { }
684 
685  NCInterpConfig(const NCInterpConfig&) = default;
686 
687  NCInterpConfig &operator=(const NCInterpConfig &rhs) = default;
688 };
689 
690 /** @brief This class manages the storage and computation of the interpolations
691  from master (coarse) face to slave (fine) face.
692 */
694 {
695 protected:
698  Array<InterpConfig> interp_config; // interpolator index for each face
699  Array<NCInterpConfig> nc_interp_config; // interpolator index for each ncface
700  Vector interpolators; // face_dofs x face_dofs x num_interpolators
701  int nc_cpt; // Counter for interpolators, and used as index.
702 
703  /** The interpolators are associated to a key of containing the address of
704  PointMatrix and a local face identifier. */
705  using Key = std::pair<const DenseMatrix*,int>;
706  /// The temporary map used to store the different interpolators.
707  using Map = std::map<Key, std::pair<int,const DenseMatrix*>>;
708  Map interp_map; // The temporary map that stores the interpolators.
709 
710 public:
711  InterpolationManager() = delete;
712 
713  /** @brief main constructor.
714 
715  @param[in] fes The FiniteElementSpace on which this operates
716  @param[in] ordering Request a specific element ordering.
717  @param[in] type Request internal or boundary faces dofs
718  */
721  FaceType type);
722 
723  /** @brief Register the face with @a face and index @a face_index as a
724  conforming face for the interpolation of the degrees of freedom.
725 
726  @param[in] face The face information of the current face.
727  @param[in] face_index The interior/boundary face index.
728  */
730  int face_index);
731 
732  /** @brief Register the face with @a face and index @a face_index as a
733  conforming face for the interpolation of the degrees of freedom.
734 
735  @param[in] face The face information of the current face.
736  @param[in] face_index The interior/boundary face index.
737  */
739  int face_index);
740 
741  /** @brief Transform the interpolation matrix map into a contiguous memory
742  structure. */
744 
746 
747  /// @brief Return the total number of interpolators.
749  {
750  return nc_cpt;
751  }
752 
753  /** @brief Return an mfem::Vector containing the interpolators in the
754  following format: face_dofs x face_dofs x num_interpolators. */
755  const Vector& GetInterpolators() const
756  {
757  return interpolators;
758  }
759 
760  /** @brief Return an array containing the interpolation configuration for
761  each face registered with RegisterFaceConformingInterpolation and
762  RegisterFaceCoarseToFineInterpolation. */
764  {
765  return interp_config;
766  }
767 
768  /** @brief Return an array containing the interpolation configuration for
769  each face registered with RegisterFaceConformingInterpolation and
770  RegisterFaceCoarseToFineInterpolation. */
772  {
773  return nc_interp_config;
774  }
775 
776 private:
777  /** @brief Returns the interpolation operator from a master (coarse) face to
778  a slave (fine) face.
779 
780  @param[in] face The face information of the current face.
781  @param[in] ptMat The PointMatrix describing the position and orientation
782  of the fine face in the coarse face. This PointMatrix is
783  usually obtained from the mesh through the method
784  GetNCFacesPtMat.
785  @param[in] ordering Request a specific element ordering.
786  @return The dense matrix corresponding to the interpolation of the face
787  degrees of freedom of the master (coarse) face to the slave
788  (fine) face. */
789  const DenseMatrix* GetCoarseToFineInterpolation(
790  const Mesh::FaceInformation &face,
791  const DenseMatrix* ptMat);
792 };
793 
794 /** @brief Operator that extracts face degrees of freedom for L2 nonconforming
795  spaces.
796 
797  In order to support face restrictions on nonconforming meshes, this
798  operator interpolates master (coarse) face degrees of freedom onto the
799  slave (fine) face. This allows face integrators to treat nonconforming
800  faces just as regular conforming faces. */
802 {
803 protected:
805  mutable Vector x_interp;
806 
807  /** @brief Constructs an NCL2FaceRestriction, this is a specialization of a
808  L2FaceRestriction for nonconforming meshes.
809 
810  @param[in] fes The FiniteElementSpace on which this operates
811  @param[in] f_ordering Request a specific face dof ordering
812  @param[in] type Request internal or boundary faces dofs
813  @param[in] m Request the face dofs for elem1, or both elem1 and
814  elem2
815  @param[in] build Request the NCL2FaceRestriction to compute the
816  scatter/gather indices. False should only be used
817  when inheriting from NCL2FaceRestriction.
818  */
820  const ElementDofOrdering f_ordering,
821  const FaceType type,
822  const L2FaceValues m,
823  bool build);
824 public:
825  /** @brief Constructs an NCL2FaceRestriction, this is a specialization of a
826  L2FaceRestriction for nonconforming meshes.
827 
828  @param[in] fes The FiniteElementSpace on which this operates
829  @param[in] f_ordering Request a specific face dof ordering
830  @param[in] type Request internal or boundary faces dofs
831  @param[in] m Request the face dofs for elem1, or both elem1 and
832  elem2
833  */
835  const ElementDofOrdering f_ordering,
836  const FaceType type,
838 
839  /** @brief Scatter the degrees of freedom, i.e. goes from L-Vector to
840  face E-Vector.
841 
842  @param[in] x The L-vector degrees of freedom.
843  @param[out] y The face E-Vector degrees of freedom with the given format:
844  if L2FacesValues::DoubleValued (face_dofs x vdim x 2 x nf),
845  if L2FacesValues::SingleValued (face_dofs x vdim x nf),
846  where nf is the number of interior or boundary faces
847  requested by @a type in the constructor.
848  The face_dofs are ordered according to the given
849  ElementDofOrdering. */
850  void Mult(const Vector &x, Vector &y) const override;
851 
852  /** @brief Gather the degrees of freedom, i.e. goes from face E-Vector to
853  L-Vector.
854 
855  @param[in] x The face E-Vector degrees of freedom with the given format:
856  if L2FacesValues::DoubleValued (face_dofs x vdim x 2 x nf),
857  if L2FacesValues::SingleValued (face_dofs x vdim x nf),
858  where nf is the number of interior or boundary faces
859  requested by @a type in the constructor.
860  The face_dofs should be ordered according to the given
861  ElementDofOrdering
862  @param[in,out] y The L-vector degrees of freedom.
863  @param[in] a Scalar coefficient for addition. */
864  void AddMultTranspose(const Vector &x, Vector &y,
865  const double a = 1.0) const override;
866 
867  /** @brief Gather the degrees of freedom, i.e. goes from face E-Vector to
868  L-Vector.
869 
870  @param[in,out] x The face E-Vector degrees of freedom with the given format:
871  if L2FacesValues::DoubleValued (face_dofs x vdim x 2 x nf),
872  if L2FacesValues::SingleValued (face_dofs x vdim x nf),
873  where nf is the number of interior or boundary faces
874  requested by @a type in the constructor.
875  The face_dofs should be ordered according to the given
876  ElementDofOrdering
877  @param[in,out] y The L-vector degrees of freedom.
878 
879  @note This method is an optimization of AddMultTranspose where the @a x
880  Vector is used and modified to avoid memory allocation and memcpy. */
881  void AddMultTransposeInPlace(Vector &x, Vector &y) const override;
882 
883  /** @brief Fill the I array of SparseMatrix corresponding to the sparsity
884  pattern given by this NCL2FaceRestriction.
885 
886  @param[in,out] mat The sparse matrix for which we want to initialize the
887  row offsets.
888  @param[in] keep_nbr_block When set to true the SparseMatrix will
889  include the rows (in addition to the columns)
890  corresponding to face-neighbor dofs. The
891  default behavior is to disregard those rows.
892 
893  @warning This method is not implemented yet. */
894  void FillI(SparseMatrix &mat,
895  const bool keep_nbr_block = false) const override;
896 
897  /** @brief Fill the J and Data arrays of the SparseMatrix corresponding to
898  the sparsity pattern given by this NCL2FaceRestriction, and the values of
899  ea_data.
900 
901  @param[in] fea_data The dense matrices representing the local operators
902  on each face. The format is:
903  face_dofs x face_dofs x 2 x nf.
904  On each face the first local matrix corresponds to
905  the contribution of elem1 on elem2, and the second to
906  the contribution of elem2 on elem1.
907  @param[in,out] mat The sparse matrix that is getting filled.
908  @param[in] keep_nbr_block When set to true the SparseMatrix will
909  include the rows (in addition to the columns)
910  corresponding to face-neighbor dofs. The
911  default behavior is to disregard those rows.
912 
913  @warning This method is not implemented yet. */
914  void FillJAndData(const Vector &fea_data,
915  SparseMatrix &mat,
916  const bool keep_nbr_block = false) const override;
917 
918  /** @brief This methods adds the DG face matrices to the element matrices.
919 
920  @param[in] fea_data The dense matrices representing the local operators
921  on each face. The format is:
922  face_dofs x face_dofs x 2 x nf.
923  On each face the first and second local matrices
924  correspond to the contributions of elem1 and elem2 on
925  themselves respectively.
926  @param[in,out] ea_data The dense matrices representing the element local
927  contributions for each element to which will be
928  added the face contributions.
929  The format is: dofs x dofs x ne, where dofs is the
930  number of dofs per element and ne the number of
931  elements.
932 
933  @warning This method is not implemented yet. */
934  void AddFaceMatricesToElementMatrices(const Vector &fea_data,
935  Vector &ea_data) const override;
936 
937 private:
938  /** @brief Compute the scatter indices: L-vector to E-vector, the offsets
939  for the gathering: E-vector to L-vector, and the interpolators from
940  coarse to fine face for master non-comforming faces.
941 
942  @param[in] f_ordering Request a specific face dof ordering.
943  @param[in] type Request internal or boundary faces dofs.
944  */
945  void ComputeScatterIndicesAndOffsets(const ElementDofOrdering f_ordering,
946  const FaceType type);
947 
948  /** @brief Compute the gather indices: E-vector to L-vector.
949 
950  Note: Requires the gather offsets to be computed.
951 
952  @param[in] f_ordering Request a specific face dof ordering.
953  @param[in] type Request internal or boundary faces dofs.
954  */
955  void ComputeGatherIndices(const ElementDofOrdering f_ordering,
956  const FaceType type);
957 
958 public:
959  /** @brief Scatter the degrees of freedom, i.e. goes from L-Vector to
960  face E-Vector. Should only be used with nonconforming faces and when:
961  L2FaceValues m == L2FaceValues::DoubleValued
962 
963  @param[in] x The L-vector degrees of freedom.
964  @param[out] y The face E-Vector degrees of freedom with the given format:
965  (face_dofs x vdim x 2 x nf),
966  where nf is the number of interior or boundary faces
967  requested by @a type in the constructor.
968  The face_dofs are ordered according to the given
969  ElementDofOrdering. */
970  virtual void DoubleValuedNonconformingMult(const Vector& x, Vector& y) const;
971 
972  /** @brief Apply a change of basis from coarse element basis to fine element
973  basis for the coarse face dofs.
974 
975  @param[in,out] x The dofs vector that needs coarse dofs to be express in
976  term of the fine basis.
977  */
979 
980  /** @brief Apply a change of basis from fine element basis to coarse element
981  basis for the coarse face dofs. Should only be used when:
982  L2FaceValues m == L2FaceValues::SingleValued
983 
984  @param[in] x The dofs vector that needs coarse dofs to be express in term
985  of the coarse basis, the result is stored in x_interp.
986  */
988 
989  /** @brief Apply a change of basis from fine element basis to coarse element
990  basis for the coarse face dofs. Should only be used when:
991  L2FaceValues m == L2FaceValues::SingleValued
992 
993  @param[in,out] x The dofs vector that needs coarse dofs to be express in
994  term of the coarse basis, the result is stored in x.
995  */
997 
998  /** @brief Apply a change of basis from fine element basis to coarse element
999  basis for the coarse face dofs. Should only be used when:
1000  L2FaceValues m == L2FaceValues::DoubleValued
1001 
1002  @param[in] x The dofs vector that needs coarse dofs to be express in term
1003  of the coarse basis, the result is stored in x_interp.
1004  */
1006 
1007  /** @brief Apply a change of basis from fine element basis to coarse element
1008  basis for the coarse face dofs. Should only be used when:
1009  L2FaceValues m == L2FaceValues::DoubleValued
1010 
1011  @param[in,out] x The dofs vector that needs coarse dofs to be express in
1012  term of the coarse basis, the result is stored in
1013  x.
1014  */
1016 };
1017 
1018 
1019 /** @brief Convert a dof face index from Native ordering to lexicographic
1020  ordering for quads and hexes.
1021 
1022  @param[in] dim The dimension of the element, 2 for quad, 3 for hex
1023  @param[in] face_id The local face identifier
1024  @param[in] size1d The 1D number of degrees of freedom for each dimension
1025  @param[in] index The native index on the face
1026  @return The lexicographic index on the face
1027 */
1028 int ToLexOrdering(const int dim, const int face_id, const int size1d,
1029  const int index);
1030 
1031 /** @brief Compute the dof face index of elem2 corresponding to the given dof
1032  face index.
1033 
1034  @param[in] dim The dimension of the element, 2 for quad, 3 for hex
1035  @param[in] face_id1 The local face identifier of elem1
1036  @param[in] face_id2 The local face identifier of elem2
1037  @param[in] orientation The orientation of elem2 relative to elem1 on the
1038  face
1039  @param[in] size1d The 1D number of degrees of freedom for each dimension
1040  @param[in] index The dof index on elem1
1041  @return The dof index on elem2 facing the dof on elem1
1042 */
1043 int PermuteFaceL2(const int dim, const int face_id1,
1044  const int face_id2, const int orientation,
1045  const int size1d, const int index);
1046 
1047 }
1048 
1049 #endif // MFEM_RESTRICTION
void FillI(SparseMatrix &mat) const
Array< NCInterpConfig > nc_interp_config
const Array< NCInterpConfig > & GetNCFaceInterpConfig() const
Return an array containing the interpolation configuration for each face registered with RegisterFace...
void MultTranspose(const Vector &x, Vector &y) const override
Set the face degrees of freedom in the element degrees of freedom y to the values given in x...
void Mult(const Vector &x, Vector &y) const override=0
Extract the face degrees of freedom from x into y.
void SingleValuedNonconformingTransposeInterpolation(const Vector &x) const
Apply a change of basis from fine element basis to coarse element basis for the coarse face dofs...
void PermuteAndSetSharedFaceDofsScatterIndices2(const Mesh::FaceInformation &face, const int face_index)
Permute and set the scattering indices of elem2 for the shared face described by the face...
void FillJAndData(const Vector &ea_data, SparseMatrix &mat) const
void AddMultTranspose(const Vector &x, Vector &y, const double a=1.0) const override
Add the E-vector degrees of freedom x to the L-vector degrees of freedom y.
NCInterpConfig(int face_index, int master_side, int nc_index)
InterpConfig & operator=(const InterpConfig &rhs)=default
void Mult(const Vector &x, Vector &y) const override
Scatter the degrees of freedom, i.e. goes from L-Vector to face E-Vector.
void SetBoundaryDofsScatterIndices2(const Mesh::FaceInformation &face, const int face_index)
Set the scattering indices of elem2 for the boundary face described by the face.
std::pair< const DenseMatrix *, int > Key
virtual void AddMultTransposeUnsigned(const Vector &x, Vector &y, const double a=1.0) const
Add the face degrees of freedom x to the element degrees of freedom y ignoring the signs from DOF ori...
void AddFaceMatricesToElementMatrices(const Vector &fea_data, Vector &ea_data) const override
This methods adds the DG face matrices to the element matrices.
void SetFaceDofsScatterIndices(const Mesh::FaceInformation &face, const int face_index, const ElementDofOrdering f_ordering)
Set the scattering indices of elem1, and increment the offsets for the face described by the face...
virtual void DoubleValuedNonconformingMult(const Vector &x, Vector &y) const
Scatter the degrees of freedom, i.e. goes from L-Vector to face E-Vector. Should only be used with no...
Data type dense matrix using column-major storage.
Definition: densemat.hpp:23
Operator that extracts face degrees of freedom for L2 nonconforming spaces.
void DoubleValuedNonconformingTransposeInterpolation(const Vector &x) const
Apply a change of basis from fine element basis to coarse element basis for the coarse face dofs...
NCInterpConfig(int face_index, InterpConfig &config)
virtual void FillJAndData(const Vector &fea_data, SparseMatrix &mat, const bool keep_nbr_block=false) const
Fill the J and Data arrays of the SparseMatrix corresponding to the sparsity pattern given by this L2...
void SetFaceDofsGatherIndices1(const Mesh::FaceInformation &face, const int face_index)
Set the gathering indices of elem1 for the interior face described by the face.
FaceRestriction(int h, int w)
void AddMultTranspose(const Vector &x, Vector &y, const double a=1.0) const override=0
Add the E-vector degrees of freedom x to the L-vector degrees of freedom y.
void AddMultTransposeUnsigned(const Vector &x, Vector &y, const double a=1.0) const override
Gather the degrees of freedom, i.e. goes from face E-Vector to L-Vector not taking into account signs...
void FillI(SparseMatrix &mat, const bool keep_nbr_block=false) const override
Fill the I array of SparseMatrix corresponding to the sparsity pattern given by this NCL2FaceRestrict...
const L2FaceValues m
void AddMultTranspose(const Vector &x, Vector &y, const double a=1.0) const override
Add the E-vector degrees of freedom x to the L-vector degrees of freedom y.
Array< int > gather_indices
void PermuteAndSetFaceDofsGatherIndices2(const Mesh::FaceInformation &face, const int face_index)
Permute and set the gathering indices of elem2 for the interior face described by the face...
L2ElementRestriction(const FiniteElementSpace &)
void DoubleValuedNonconformingInterpolation(Vector &x) const
Apply a change of basis from coarse element basis to fine element basis for the coarse face dofs...
Array< InterpConfig > interp_config
InterpConfig(int master_side, int nc_index)
const Array< int > & GatherMap() const
Definition: restriction.hpp:97
void SingleValuedConformingMult(const Vector &x, Vector &y) const
Scatter the degrees of freedom, i.e. goes from L-Vector to face E-Vector. Should only be used with co...
void FillJAndData(const Vector &fea_data, SparseMatrix &mat, const bool keep_nbr_block=false) const override
Fill the J and Data arrays of the SparseMatrix corresponding to the sparsity pattern given by this NC...
void SetFaceDofsGatherIndices(const Mesh::FaceInformation &face, const int face_index, const ElementDofOrdering f_ordering)
Set the gathering indices of elem1 for the interior face described by the face.
virtual void AddFaceMatricesToElementMatrices(const Vector &fea_data, Vector &ea_data) const
This methods adds the DG face matrices to the element matrices.
FaceType
Definition: mesh.hpp:45
void RegisterFaceConformingInterpolation(const Mesh::FaceInformation &face, int face_index)
Register the face with face and index face_index as a conforming face for the interpolation of the de...
const FiniteElementSpace & fes
void MultTranspose(const Vector &x, Vector &y) const override
Action of the transpose operator: y=A^t(x). The default behavior in class Operator is to generate an ...
virtual void AddMultTransposeInPlace(Vector &x, Vector &y) const
Add the face degrees of freedom x to the element degrees of freedom y. Perform the same computation a...
void PermuteAndSetFaceDofsScatterIndices2(const Mesh::FaceInformation &face, const int face_index)
Permute and set the scattering indices of elem2, and increment the offsets for the face described by ...
Data type sparse matrix.
Definition: sparsemat.hpp:50
virtual void AddMultTranspose(const Vector &x, Vector &y, const double a=1.0) const override=0
Add the face degrees of freedom x to the element degrees of freedom y.
void Mult(const Vector &x, Vector &y) const override
Operator application: y=A(x).
Operator that extracts Face degrees of freedom for L2 spaces.
ConformingFaceRestriction(const FiniteElementSpace &fes, const ElementDofOrdering f_ordering, const FaceType type, bool build)
Construct a ConformingFaceRestriction.
void Mult(const Vector &x, Vector &y) const override
Scatter the degrees of freedom, i.e. goes from L-Vector to face E-Vector.
void LinearizeInterpolatorMapIntoVector()
Transform the interpolation matrix map into a contiguous memory structure.
const Array< int > & Offsets() const
Definition: restriction.hpp:99
void MultLeftInverse(const Vector &x, Vector &y) const
void MultTransposeUnsigned(const Vector &x, Vector &y) const
Compute MultTranspose without applying signs based on DOF orientations.
virtual void FillI(SparseMatrix &mat, const bool keep_nbr_block=false) const
Fill the I array of SparseMatrix corresponding to the sparsity pattern given by this L2FaceRestrictio...
ElementRestriction(const FiniteElementSpace &, ElementDofOrdering)
Definition: restriction.cpp:27
void Mult(const Vector &x, Vector &y) const override
Scatter the degrees of freedom, i.e. goes from L-Vector to face E-Vector.
const FiniteElementSpace & fes
This class manages the storage and computation of the interpolations from master (coarse) face to sla...
Abstract base class that defines an interface for element restrictions.
Definition: restriction.hpp:25
void AddMultTransposeInPlace(Vector &x, Vector &y) const override
Gather the degrees of freedom, i.e. goes from face E-Vector to L-Vector.
This structure is used as a human readable output format that decipheres the information contained in...
Definition: mesh.hpp:1664
const ElementDofOrdering ordering
int ToLexOrdering(const int dim, const int face_id, const int size1d, const int index)
Convert a dof face index from Native ordering to lexicographic ordering for quads and hexes...
const FiniteElementSpace & fes
Definition: restriction.hpp:45
Array< int > gather_offsets
void MultTranspose(const Vector &x, Vector &y) const override
Action of the transpose operator: y=A^t(x). The default behavior in class Operator is to generate an ...
Operator that converts FiniteElementSpace L-vectors to E-vectors.
Definition: restriction.hpp:37
Class FiniteElementSpace - responsible for providing FEM view of the mesh, mainly managing the set of...
Definition: fespace.hpp:219
virtual void DoubleValuedConformingMult(const Vector &x, Vector &y) const
Scatter the degrees of freedom, i.e. goes from L-Vector to face E-Vector. Should only be used with co...
L2FaceRestriction(const FiniteElementSpace &fes, const ElementDofOrdering f_ordering, const FaceType type, const L2FaceValues m, bool build)
Constructs an L2FaceRestriction.
const Array< InterpConfig > & GetFaceInterpConfig() const
Return an array containing the interpolation configuration for each face registered with RegisterFace...
void FillJAndData(const Vector &ea_data, SparseMatrix &mat) const
void MultUnsigned(const Vector &x, Vector &y) const
Compute Mult without applying signs based on DOF orientations.
std::map< Key, std::pair< int, const DenseMatrix * > > Map
The temporary map used to store the different interpolators.
void SingleValuedConformingAddMultTranspose(const Vector &x, Vector &y) const
Gather the degrees of freedom, i.e. goes from face E-Vector to L-Vector. Should only be used with con...
void FillSparseMatrix(const Vector &mat_ea, SparseMatrix &mat) const
Fill a Sparse Matrix with Element Matrices.
void CheckFESpace(const ElementDofOrdering f_ordering)
Verify that ConformingFaceRestriction is built from a supported finite element space.
Array< int > scatter_indices1
double a
Definition: lissajous.cpp:41
void AddMultTranspose(const Vector &x, Vector &y, const double a=1.0) const override
Gather the degrees of freedom, i.e. goes from face E-Vector to L-Vector.
Operator that extracts face degrees of freedom for H1, ND, or RT FiniteElementSpaces.
ElementDofOrdering
Constants describing the possible orderings of the DOFs in one element.
Definition: fespace.hpp:74
int dim
Definition: ex24.cpp:53
int GetNumInterpolators() const
Return the total number of interpolators.
int FillI(SparseMatrix &mat) const
int index(int i, int j, int nx, int ny)
Definition: life.cpp:235
void Mult(const Vector &x, Vector &y) const override
Operator application: y=A(x).
const FiniteElementSpace & fes
void DoubleValuedConformingAddMultTranspose(const Vector &x, Vector &y) const
Gather the degrees of freedom, i.e. goes from face E-Vector to L-Vector. Should only be used with con...
void AddMultTranspose(const Vector &x, Vector &y, const double a=1.0) const override
Gather the degrees of freedom, i.e. goes from face E-Vector to L-Vector.
void SingleValuedNonconformingTransposeInterpolationInPlace(Vector &x) const
Apply a change of basis from fine element basis to coarse element basis for the coarse face dofs...
Vector data type.
Definition: vector.hpp:58
int PermuteFaceL2(const int dim, const int face_id1, const int face_id2, const int orientation, const int size1d, const int index)
Compute the dof face index of elem2 corresponding to the given dof face index.
void SetFaceDofsScatterIndices1(const Mesh::FaceInformation &face, const int face_index)
Set the scattering indices of elem1, and increment the offsets for the face described by the face...
const Array< int > & Indices() const
Definition: restriction.hpp:98
void RegisterFaceCoarseToFineInterpolation(const Mesh::FaceInformation &face, int face_index)
Register the face with face and index face_index as a conforming face for the interpolation of the de...
Base class for operators that extracts Face degrees of freedom.
Abstract operator.
Definition: operator.hpp:24
uint32_t is_non_conforming
void CheckFESpace(const ElementDofOrdering f_ordering)
Verify that L2FaceRestriction is built from an L2 FESpace.
void AddMultTranspose(const Vector &x, Vector &y, const double a=1.0) const override
Gather the degrees of freedom, i.e. goes from face E-Vector to L-Vector.
const Vector & GetInterpolators() const
Return an mfem::Vector containing the interpolators in the following format: face_dofs x face_dofs x ...
Operator that converts L2 FiniteElementSpace L-vectors to E-vectors.
InterpolationManager interpolations
NCInterpConfig & operator=(const NCInterpConfig &rhs)=default
InterpConfig()=default
Array< int > scatter_indices2
void DoubleValuedNonconformingTransposeInterpolationInPlace(Vector &x) const
Apply a change of basis from fine element basis to coarse element basis for the coarse face dofs...
void BooleanMask(Vector &y) const
Fills the E-vector y with boolean values 0.0 and 1.0 such that each each entry of the L-vector is uni...
NCInterpConfig()=default
NCL2FaceRestriction(const FiniteElementSpace &fes, const ElementDofOrdering f_ordering, const FaceType type, const L2FaceValues m, bool build)
Constructs an NCL2FaceRestriction, this is a specialization of a L2FaceRestriction for nonconforming ...