MFEM  v4.0
Finite element discretization library
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
fespace.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_FESPACE
13 #define MFEM_FESPACE
14 
15 #include "../config/config.hpp"
16 #include "../linalg/sparsemat.hpp"
17 #include "../mesh/mesh.hpp"
18 #include "fe_coll.hpp"
19 #include <iostream>
20 
21 namespace mfem
22 {
23 
24 /** @brief The ordering method used when the number of unknowns per mesh node
25  (vector dimension) is bigger than 1. */
26 class Ordering
27 {
28 public:
29  /// %Ordering methods:
30  enum Type
31  {
32  byNODES, /**< loop first over the nodes (inner loop) then over the vector
33  dimension (outer loop); symbolically it can be represented
34  as: XXX...,YYY...,ZZZ... */
35  byVDIM /**< loop first over the vector dimension (inner loop) then over
36  the nodes (outer loop); symbolically it can be represented
37  as: XYZ,XYZ,XYZ,... */
38  };
39 
40  template <Type Ord>
41  static inline int Map(int ndofs, int vdim, int dof, int vd);
42 
43  template <Type Ord>
44  static void DofsToVDofs(int ndofs, int vdim, Array<int> &dofs);
45 };
46 
47 template <> inline int
48 Ordering::Map<Ordering::byNODES>(int ndofs, int vdim, int dof, int vd)
49 {
50  MFEM_ASSERT(dof < ndofs && -1-dof < ndofs && 0 <= vd && vd < vdim, "");
51  return (dof >= 0) ? dof+ndofs*vd : dof-ndofs*vd;
52 }
53 
54 template <> inline int
55 Ordering::Map<Ordering::byVDIM>(int ndofs, int vdim, int dof, int vd)
56 {
57  MFEM_ASSERT(dof < ndofs && -1-dof < ndofs && 0 <= vd && vd < vdim, "");
58  return (dof >= 0) ? vd+vdim*dof : -1-(vd+vdim*(-1-dof));
59 }
60 
61 
62 /// Constants describing the possible orderings of the DOFs in one element.
64 {
65  /// Native ordering as defined by the FiniteElement.
66  /** This ordering can be used by tensor-product elements when the
67  interpolation from the DOFs to quadrature points does not use the
68  tensor-product structure. */
69  NATIVE,
70  /// Lexicographic ordering for tensor-product FiniteElements.
71  /** This ordering can be used only with tensor-product elements. */
73 };
74 
75 
76 // Forward declarations
77 class NURBSExtension;
78 class BilinearFormIntegrator;
79 class QuadratureSpace;
80 class QuadratureInterpolator;
81 
82 
83 /** @brief Class FiniteElementSpace - responsible for providing FEM view of the
84  mesh, mainly managing the set of degrees of freedom. */
86 {
88 
89 protected:
90  /// The mesh that FE space lives on (not owned).
92 
93  /// Associated FE collection (not owned).
95 
96  /// %Vector dimension (number of unknowns per degree of freedom).
97  int vdim;
98 
99  /** Type of ordering of the vector dofs when #vdim > 1.
100  - Ordering::byNODES - first nodes, then vector dimension,
101  - Ordering::byVDIM - first vector dimension, then nodes */
103 
104  /// Number of degrees of freedom. Number of unknowns is #ndofs * #vdim.
105  int ndofs;
106 
108  int *fdofs, *bdofs;
109 
110  mutable Table *elem_dof; // if NURBS FE space, not owned; otherwise, owned.
111  Table *bdrElem_dof; // used only with NURBS FE spaces; not owned.
112 
114 
116  int own_ext;
117 
118  /** Matrix representing the prolongation from the global conforming dofs to
119  a set of intermediate partially conforming dofs, e.g. the dofs associated
120  with a "cut" space on a non-conforming mesh. */
121  mutable SparseMatrix *cP; // owned
122  /// Conforming restriction matrix such that cR.cP=I.
123  mutable SparseMatrix *cR; // owned
124  mutable bool cP_is_set;
125 
126  /// Transformation to apply to GridFunctions after space Update().
128 
129  /// The element restriction operators, see GetElementRestriction().
131 
133 
134  long sequence; // should match Mesh::GetSequence
135 
136  void UpdateNURBS();
137 
138  void Construct();
139  void Destroy();
140 
141  void BuildElementToDofTable() const;
142 
143  /// Helper to remove encoded sign from a DOF
144  static inline int DecodeDof(int dof, double& sign)
145  { return (dof >= 0) ? (sign = 1, dof) : (sign = -1, (-1 - dof)); }
146 
147  /// Helper to get vertex, edge or face DOFs (entity=0,1,2 resp.).
148  void GetEntityDofs(int entity, int index, Array<int> &dofs) const;
149 
150  /// Calculate the cP and cR matrices for a nonconforming mesh.
151  void BuildConformingInterpolation() const;
152 
153  static void AddDependencies(SparseMatrix& deps, Array<int>& master_dofs,
154  Array<int>& slave_dofs, DenseMatrix& I);
155 
156  static bool DofFinalizable(int dof, const Array<bool>& finalized,
157  const SparseMatrix& deps);
158 
159  void MakeVDimMatrix(SparseMatrix &mat) const;
160 
161  /// GridFunction interpolation operator applicable after mesh refinement.
163  {
164  const FiniteElementSpace* fespace;
166  Table* old_elem_dof; // Owned.
167 
168  public:
169  /** Construct the operator based on the elem_dof table of the original
170  (coarse) space. The class takes ownership of the table. */
171  RefinementOperator(const FiniteElementSpace* fespace,
172  Table *old_elem_dof/*takes ownership*/, int old_ndofs);
173  RefinementOperator(const FiniteElementSpace *fespace,
174  const FiniteElementSpace *coarse_fes);
175  virtual void Mult(const Vector &x, Vector &y) const;
176  virtual ~RefinementOperator();
177  };
178 
179  // Derefinement operator, used by the friend class InterpolationGridTransfer.
181  {
182  const FiniteElementSpace *fine_fes; // Not owned.
184  Table *coarse_elem_dof; // Owned.
185  Table coarse_to_fine;
186  Array<int> coarse_to_ref_type;
187  Array<Geometry::Type> ref_type_to_geom;
188  Array<int> ref_type_to_fine_elem_offset;
189 
190  public:
192  const FiniteElementSpace *c_fes,
193  BilinearFormIntegrator *mass_integ);
194  virtual void Mult(const Vector &x, Vector &y) const;
195  virtual ~DerefinementOperator();
196  };
197 
198  // This method makes the same assumptions as the method:
199  // void GetLocalRefinementMatrices(
200  // const FiniteElementSpace &coarse_fes, Geometry::Type geom,
201  // DenseTensor &localP) const
202  // which is defined below. It also assumes that the coarse fes and this have
203  // the same vector dimension, vdim.
204  SparseMatrix *RefinementMatrix_main(const int coarse_ndofs,
205  const Table &coarse_elem_dof,
206  const DenseTensor localP[]) const;
207 
209  DenseTensor &localP) const;
211  DenseTensor &localR) const;
212 
213  /** Calculate explicit GridFunction interpolation matrix (after mesh
214  refinement). NOTE: consider using the RefinementOperator class instead
215  of the fully assembled matrix, which can take a lot of memory. */
216  SparseMatrix* RefinementMatrix(int old_ndofs, const Table* old_elem_dof);
217 
218  /// Calculate GridFunction restriction matrix after mesh derefinement.
219  SparseMatrix* DerefinementMatrix(int old_ndofs, const Table* old_elem_dof);
220 
221  // This method assumes that this->mesh is a refinement of coarse_fes->mesh
222  // and that the CoarseFineTransformations of this->mesh are set accordingly.
223  // Another assumption is that the FEs of this use the same MapType as the FEs
224  // of coarse_fes. Finally, it assumes that the spaces this and coarse_fes are
225  // NOT variable-order spaces.
226  void GetLocalRefinementMatrices(const FiniteElementSpace &coarse_fes,
228  DenseTensor &localP) const;
229 
230  /// Help function for constructors + Load().
231  void Constructor(Mesh *mesh, NURBSExtension *ext,
233  int vdim = 1, int ordering = Ordering::byNODES);
234 
235 public:
236  /** @brief Default constructor: the object is invalid until initialized using
237  the method Load(). */
239 
240  /** @brief Copy constructor: deep copy all data from @a orig except the Mesh,
241  the FiniteElementCollection, ans some derived data. */
242  /** If the @a mesh or @a fec pointers are NULL (default), then the new
243  FiniteElementSpace will reuse the respective pointers from @a orig. If
244  any of these pointers is not NULL, the given pointer will be used instead
245  of the one used by @a orig.
246 
247  @note The objects pointed to by the @a mesh and @a fec parameters must be
248  either the same objects as the ones used by @a orig, or copies of them.
249  Otherwise, the behavior is undefined.
250 
251  @note Derived data objects, such as the conforming prolongation and
252  restriction matrices, and the update operator, will not be copied, even
253  if they are created in the @a orig object. */
254  FiniteElementSpace(const FiniteElementSpace &orig, Mesh *mesh = NULL,
255  const FiniteElementCollection *fec = NULL);
256 
259  int vdim = 1, int ordering = Ordering::byNODES)
260  { Constructor(mesh, NULL, fec, vdim, ordering); }
261 
262  /// Construct a NURBS FE space based on the given NURBSExtension, @a ext.
263  /** @note If the pointer @a ext is NULL, this constructor is equivalent to
264  the standard constructor with the same arguments minus the
265  NURBSExtension, @a ext. */
268  int vdim = 1, int ordering = Ordering::byNODES)
269  { Constructor(mesh, ext, fec, vdim, ordering); }
270 
271  /// Returns the mesh
272  inline Mesh *GetMesh() const { return mesh; }
273 
274  const NURBSExtension *GetNURBSext() const { return NURBSext; }
277 
278  bool Conforming() const { return mesh->Conforming(); }
279  bool Nonconforming() const { return mesh->Nonconforming(); }
280 
281  /// The returned SparseMatrix is owned by the FiniteElementSpace.
283 
284  /// The returned SparseMatrix is owned by the FiniteElementSpace.
285  const SparseMatrix *GetConformingRestriction() const;
286 
287  /// The returned Operator is owned by the FiniteElementSpace.
288  virtual const Operator *GetProlongationMatrix() const
289  { return GetConformingProlongation(); }
290 
291  /// The returned SparseMatrix is owned by the FiniteElementSpace.
292  virtual const SparseMatrix *GetRestrictionMatrix() const
293  { return GetConformingRestriction(); }
294 
295  /// Return an Operator that converts L-vectors to E-vectors.
296  /** An L-vector is a vector of size GetVSize() which is the same size as a
297  GridFunction. An E-vector represents the element-wise discontinuous
298  version of the FE space.
299 
300  The layout of the E-vector is: ND x VDIM x NE, where ND is the number of
301  degrees of freedom, VDIM is the vector dimension of the FE space, and NE
302  is the number of the mesh elements.
303 
304  The parameter @a e_ordering describes how the local DOFs in each element
305  should be ordered, see ElementDofOrdering.
306 
307  For discontinuous spaces, where the element-restriction is the identity,
308  this method will return NULL.
309 
310  The returned Operator is owned by the FiniteElementSpace. */
311  const Operator *GetElementRestriction(ElementDofOrdering e_ordering) const;
312 
313  /** @brief Return a QuadratureInterpolator that interpolates E-vectors to
314  quadrature point values and/or derivatives (Q-vectors). */
315  /** An E-vector represents the element-wise discontinuous version of the FE
316  space and can be obtained, for example, from a GridFunction using the
317  Operator returned by GetElementRestriction().
318 
319  All elements will use the same IntegrationRule, @a ir as the target
320  quadrature points. */
322  const IntegrationRule &ir) const;
323 
324  /** @brief Return a QuadratureInterpolator that interpolates E-vectors to
325  quadrature point values and/or derivatives (Q-vectors). */
326  /** An E-vector represents the element-wise discontinuous version of the FE
327  space and can be obtained, for example, from a GridFunction using the
328  Operator returned by GetElementRestriction().
329 
330  The target quadrature points in the elements are described by the given
331  QuadratureSpace, @a qs. */
333  const QuadratureSpace &qs) const;
334 
335  /// Returns vector dimension.
336  inline int GetVDim() const { return vdim; }
337 
338  /// Returns the order of the i'th finite element
339  int GetOrder(int i) const;
340  /// Returns the order of the i'th face finite element
341  int GetFaceOrder(int i) const;
342 
343  /// Returns number of degrees of freedom.
344  inline int GetNDofs() const { return ndofs; }
345 
346  /// Return the number of vector dofs, i.e. GetNDofs() x GetVDim().
347  inline int GetVSize() const { return vdim * ndofs; }
348 
349  /// Return the number of vector true (conforming) dofs.
350  virtual int GetTrueVSize() const { return GetConformingVSize(); }
351 
352  /// Returns the number of conforming ("true") degrees of freedom
353  /// (if the space is on a nonconforming mesh with hanging nodes).
354  int GetNConformingDofs() const;
355 
356  int GetConformingVSize() const { return vdim * GetNConformingDofs(); }
357 
358  /// Return the ordering method.
359  inline Ordering::Type GetOrdering() const { return ordering; }
360 
361  const FiniteElementCollection *FEColl() const { return fec; }
362 
363  /// Number of all scalar vertex dofs
364  int GetNVDofs() const { return nvdofs; }
365  /// Number of all scalar edge-interior dofs
366  int GetNEDofs() const { return nedofs; }
367  /// Number of all scalar face-interior dofs
368  int GetNFDofs() const { return nfdofs; }
369 
370  /// Returns number of vertices in the mesh.
371  inline int GetNV() const { return mesh->GetNV(); }
372 
373  /// Returns number of elements in the mesh.
374  inline int GetNE() const { return mesh->GetNE(); }
375 
376  /// Returns number of faces (i.e. co-dimension 1 entities) in the mesh.
377  /** The co-dimension 1 entities are those that have dimension 1 less than the
378  mesh dimension, e.g. for a 2D mesh, the faces are the 1D entities, i.e.
379  the edges. */
380  inline int GetNF() const { return mesh->GetNumFaces(); }
381 
382  /// Returns number of boundary elements in the mesh.
383  inline int GetNBE() const { return mesh->GetNBE(); }
384 
385  /// Returns the type of element i.
386  inline int GetElementType(int i) const
387  { return mesh->GetElementType(i); }
388 
389  /// Returns the vertices of element i.
390  inline void GetElementVertices(int i, Array<int> &vertices) const
391  { mesh->GetElementVertices(i, vertices); }
392 
393  /// Returns the type of boundary element i.
394  inline int GetBdrElementType(int i) const
395  { return mesh->GetBdrElementType(i); }
396 
397  /// Returns ElementTransformation for the @a i-th element.
399  { return mesh->GetElementTransformation(i); }
400 
401  /** @brief Returns the transformation defining the @a i-th element in the
402  user-defined variable @a ElTr. */
404  { mesh->GetElementTransformation(i, ElTr); }
405 
406  /// Returns ElementTransformation for the @a i-th boundary element.
408  { return mesh->GetBdrElementTransformation(i); }
409 
410  int GetAttribute(int i) const { return mesh->GetAttribute(i); }
411 
412  int GetBdrAttribute(int i) const { return mesh->GetBdrAttribute(i); }
413 
414  /// Returns indexes of degrees of freedom in array dofs for i'th element.
415  virtual void GetElementDofs(int i, Array<int> &dofs) const;
416 
417  /// Returns indexes of degrees of freedom for i'th boundary element.
418  virtual void GetBdrElementDofs(int i, Array<int> &dofs) const;
419 
420  /** Returns the indexes of the degrees of freedom for i'th face
421  including the dofs for the edges and the vertices of the face. */
422  virtual void GetFaceDofs(int i, Array<int> &dofs) const;
423 
424  /** Returns the indexes of the degrees of freedom for i'th edge
425  including the dofs for the vertices of the edge. */
426  void GetEdgeDofs(int i, Array<int> &dofs) const;
427 
428  void GetVertexDofs(int i, Array<int> &dofs) const;
429 
430  void GetElementInteriorDofs(int i, Array<int> &dofs) const;
431 
432  void GetFaceInteriorDofs(int i, Array<int> &dofs) const;
433 
434  int GetNumElementInteriorDofs(int i) const
436 
437  void GetEdgeInteriorDofs(int i, Array<int> &dofs) const;
438 
439  void DofsToVDofs(Array<int> &dofs, int ndofs = -1) const;
440 
441  void DofsToVDofs(int vd, Array<int> &dofs, int ndofs = -1) const;
442 
443  int DofToVDof(int dof, int vd, int ndofs = -1) const;
444 
445  int VDofToDof(int vdof) const
446  { return (ordering == Ordering::byNODES) ? (vdof%ndofs) : (vdof/vdim); }
447 
448  static void AdjustVDofs(Array<int> &vdofs);
449 
450  /// Returns indexes of degrees of freedom in array dofs for i'th element.
451  void GetElementVDofs(int i, Array<int> &vdofs) const;
452 
453  /// Returns indexes of degrees of freedom for i'th boundary element.
454  void GetBdrElementVDofs(int i, Array<int> &vdofs) const;
455 
456  /// Returns indexes of degrees of freedom for i'th face element (2D and 3D).
457  void GetFaceVDofs(int i, Array<int> &vdofs) const;
458 
459  /// Returns indexes of degrees of freedom for i'th edge.
460  void GetEdgeVDofs(int i, Array<int> &vdofs) const;
461 
462  void GetVertexVDofs(int i, Array<int> &vdofs) const;
463 
464  void GetElementInteriorVDofs(int i, Array<int> &vdofs) const;
465 
466  void GetEdgeInteriorVDofs(int i, Array<int> &vdofs) const;
467 
469 
470  /** @brief Reorder the scalar DOFs based on the element ordering.
471 
472  The new ordering is constructed as follows: 1) loop over all elements as
473  ordered in the Mesh; 2) for each element, assign new indices to all of
474  its current DOFs that are still unassigned; the new indices we assign are
475  simply the sequence `0,1,2,...`; if there are any signed DOFs their sign
476  is preserved. */
478 
479  void BuildDofToArrays();
480 
481  const Table &GetElementToDofTable() const { return *elem_dof; }
482  const Table &GetBdrElementToDofTable() const { return *bdrElem_dof; }
483 
484  int GetElementForDof(int i) const { return dof_elem_array[i]; }
485  int GetLocalDofForDof(int i) const { return dof_ldof_array[i]; }
486 
487  /// Returns pointer to the FiniteElement associated with i'th element.
488  const FiniteElement *GetFE(int i) const;
489 
490  /// Returns pointer to the FiniteElement for the i'th boundary element.
491  const FiniteElement *GetBE(int i) const;
492 
493  const FiniteElement *GetFaceElement(int i) const;
494 
495  const FiniteElement *GetEdgeElement(int i) const;
496 
497  /// Return the trace element from element 'i' to the given 'geom_type'
498  const FiniteElement *GetTraceElement(int i, Geometry::Type geom_type) const;
499 
500  /** Mark degrees of freedom associated with boundary elements with
501  the specified boundary attributes (marked in 'bdr_attr_is_ess').
502  For spaces with 'vdim' > 1, the 'component' parameter can be used
503  to restricts the marked vDOFs to the specified component. */
504  virtual void GetEssentialVDofs(const Array<int> &bdr_attr_is_ess,
505  Array<int> &ess_vdofs,
506  int component = -1) const;
507 
508  /** Get a list of essential true dofs, ess_tdof_list, corresponding to the
509  boundary attributes marked in the array bdr_attr_is_ess.
510  For spaces with 'vdim' > 1, the 'component' parameter can be used
511  to restricts the marked tDOFs to the specified component. */
512  virtual void GetEssentialTrueDofs(const Array<int> &bdr_attr_is_ess,
513  Array<int> &ess_tdof_list,
514  int component = -1);
515 
516  /// Convert a Boolean marker array to a list containing all marked indices.
517  static void MarkerToList(const Array<int> &marker, Array<int> &list);
518 
519  /** Convert an array of indices (list) to a Boolean marker array where all
520  indices in the list are marked with the given value and the rest are set
521  to zero. */
522  static void ListToMarker(const Array<int> &list, int marker_size,
523  Array<int> &marker, int mark_val = -1);
524 
525  /** For a partially conforming FE space, convert a marker array (nonzero
526  entries are true) on the partially conforming dofs to a marker array on
527  the conforming dofs. A conforming dofs is marked iff at least one of its
528  dependent dofs is marked. */
529  void ConvertToConformingVDofs(const Array<int> &dofs, Array<int> &cdofs);
530 
531  /** For a partially conforming FE space, convert a marker array (nonzero
532  entries are true) on the conforming dofs to a marker array on the
533  (partially conforming) dofs. A dof is marked iff it depends on a marked
534  conforming dofs, where dependency is defined by the ConformingRestriction
535  matrix; in other words, a dof is marked iff it corresponds to a marked
536  conforming dof. */
537  void ConvertFromConformingVDofs(const Array<int> &cdofs, Array<int> &dofs);
538 
539  /** Generate the global restriction matrix from a discontinuous
540  FE space to the continuous FE space of the same polynomial degree. */
542 
543  /** Generate the global restriction matrix from a discontinuous
544  FE space to the piecewise constant FE space. */
546 
547  /** Construct the restriction matrix from the FE space given by
548  (*this) to the lower degree FE space given by (*lfes) which
549  is defined on the same mesh. */
551 
552  /** @brief Construct and return an Operator that can be used to transfer
553  GridFunction data from @a coarse_fes, defined on a coarse mesh, to @a
554  this FE space, defined on a refined mesh. */
555  /** It is assumed that the mesh of this FE space is a refinement of the mesh
556  of @a coarse_fes and the CoarseFineTransformations returned by the method
557  Mesh::GetRefinementTransforms() of the refined mesh are set accordingly.
558  The Operator::Type of @a T can be set to request an Operator of the set
559  type. Currently, only Operator::MFEM_SPARSEMAT and Operator::ANY_TYPE
560  (matrix-free) are supported. When Operator::ANY_TYPE is requested, the
561  choice of the particular Operator sub-class is left to the method. This
562  method also works in parallel because the transfer operator is local to
563  the MPI task when the input is a synchronized ParGridFunction. */
564  void GetTransferOperator(const FiniteElementSpace &coarse_fes,
565  OperatorHandle &T) const;
566 
567  /** @brief Construct and return an Operator that can be used to transfer
568  true-dof data from @a coarse_fes, defined on a coarse mesh, to @a this FE
569  space, defined on a refined mesh.
570 
571  This method calls GetTransferOperator() and multiplies the result by the
572  prolongation operator of @a coarse_fes on the right, and by the
573  restriction operator of this FE space on the left.
574 
575  The Operator::Type of @a T can be set to request an Operator of the set
576  type. In serial, the supported types are: Operator::MFEM_SPARSEMAT and
577  Operator::ANY_TYPE (matrix-free). In parallel, the supported types are:
578  Operator::Hypre_ParCSR and Operator::ANY_TYPE. Any other type is treated
579  as Operator::ANY_TYPE: the operator representation choice is made by this
580  method. */
581  virtual void GetTrueTransferOperator(const FiniteElementSpace &coarse_fes,
582  OperatorHandle &T) const;
583 
584  /** Reflect changes in the mesh: update number of DOFs, etc. Also, calculate
585  GridFunction transformation operator (unless want_transform is false).
586  Safe to call multiple times, does nothing if space already up to date. */
587  virtual void Update(bool want_transform = true);
588 
589  /// Get the GridFunction update operator.
590  const Operator* GetUpdateOperator() { Update(); return Th.Ptr(); }
591 
592  /// Return the update operator in the given OperatorHandle, @a T.
594 
595  /** @brief Set the ownership of the update operator: if set to false, the
596  Operator returned by GetUpdateOperator() must be deleted outside the
597  FiniteElementSpace. */
598  /** The update operator ownership is automatically reset to true when a new
599  update operator is created by the Update() method. */
600  void SetUpdateOperatorOwner(bool own) { Th.SetOperatorOwner(own); }
601 
602  /// Specify the Operator::Type to be used by the update operators.
603  /** The default type is Operator::ANY_TYPE which leaves the choice to this
604  class. The other currently supported option is Operator::MFEM_SPARSEMAT
605  which is only guaranteed to be honored for a refinement update operator.
606  Any other type will be treated as Operator::ANY_TYPE.
607  @note This operation destroys the current update operator (if owned). */
609 
610  /// Free the GridFunction update operator (if any), to save memory.
611  virtual void UpdatesFinished() { Th.Clear(); }
612 
613  /// Return update counter (see Mesh::sequence)
614  long GetSequence() const { return sequence; }
615 
616  void Save(std::ostream &out) const;
617 
618  /** @brief Read a FiniteElementSpace from a stream. The returned
619  FiniteElementCollection is owned by the caller. */
620  FiniteElementCollection *Load(Mesh *m, std::istream &input);
621 
622  virtual ~FiniteElementSpace();
623 };
624 
625 
626 /// Class representing the storage layout of a QuadratureFunction.
627 /** Multiple QuadratureFunction%s can share the same QuadratureSpace. */
629 {
630 protected:
631  friend class QuadratureFunction; // Uses the element_offsets.
632 
634  int order;
635  int size;
636 
638  int *element_offsets; // scalar offsets; size = number of elements + 1
639 
640  // protected functions
641 
642  // Assuming mesh and order are set, construct the members: int_rule,
643  // element_offsets, and size.
644  void Construct();
645 
646 public:
647  /// Create a QuadratureSpace based on the global rules from #IntRules.
648  QuadratureSpace(Mesh *mesh_, int order_)
649  : mesh(mesh_), order(order_) { Construct(); }
650 
651  /// Read a QuadratureSpace from the stream @a in.
652  QuadratureSpace(Mesh *mesh_, std::istream &in);
653 
654  virtual ~QuadratureSpace() { delete [] element_offsets; }
655 
656  /// Return the total number of quadrature points.
657  int GetSize() const { return size; }
658 
659  /// Get the IntegrationRule associated with mesh element @a idx.
660  const IntegrationRule &GetElementIntRule(int idx) const
661  { return *int_rule[mesh->GetElementBaseGeometry(idx)]; }
662 
663  /// Write the QuadratureSpace to the stream @a out.
664  void Save(std::ostream &out) const;
665 };
666 
667 
668 /** @brief Base class for transfer algorithms that construct transfer Operator%s
669  between two finite element (FE) spaces. */
670 /** Generally, the two FE spaces (domain and range) can be defined on different
671  meshes. */
673 {
674 protected:
675  FiniteElementSpace &dom_fes; ///< Domain FE space
676  FiniteElementSpace &ran_fes; ///< Range FE space
677 
678  /** @brief Desired Operator::Type for the construction of all operators
679  defined by the underlying transfer algorithm. It can be ignored by
680  derived classes. */
682 
683  OperatorHandle fw_t_oper; ///< Forward true-dof operator
684  OperatorHandle bw_t_oper; ///< Backward true-dof operator
685 
686 #ifdef MFEM_USE_MPI
687  bool parallel;
688 #endif
689  bool Parallel() const
690  {
691 #ifndef MFEM_USE_MPI
692  return false;
693 #else
694  return parallel;
695 #endif
696  }
697 
699  FiniteElementSpace &fes_out,
700  const Operator &oper,
701  OperatorHandle &t_oper);
702 
703 public:
704  /** Construct a transfer algorithm between the domain, @a dom_fes_, and
705  range, @a ran_fes_, FE spaces. */
706  GridTransfer(FiniteElementSpace &dom_fes_, FiniteElementSpace &ran_fes_);
707 
708  /// Virtual destructor
709  virtual ~GridTransfer() { }
710 
711  /** @brief Set the desired Operator::Type for the construction of all
712  operators defined by the underlying transfer algorithm. */
713  /** The default value is Operator::ANY_TYPE which typically corresponds to
714  a matrix-free operator representation. Note that derived classes are not
715  required to support this setting and can ignore it. */
716  void SetOperatorType(Operator::Type type) { oper_type = type; }
717 
718  /** @brief Return an Operator that transfers GridFunction%s from the domain
719  FE space to GridFunction%s in the range FE space. */
720  virtual const Operator &ForwardOperator() = 0;
721 
722  /** @brief Return an Operator that transfers GridFunction%s from the range
723  FE space back to GridFunction%s in the domain FE space. */
724  virtual const Operator &BackwardOperator() = 0;
725 
726  /** @brief Return an Operator that transfers true-dof Vector%s from the
727  domain FE space to true-dof Vector%s in the range FE space. */
728  /** This method is implemented in the base class, based on ForwardOperator(),
729  however, derived classes can overload the construction, if necessary. */
730  virtual const Operator &TrueForwardOperator()
731  {
733  }
734 
735  /** @brief Return an Operator that transfers true-dof Vector%s from the
736  range FE space back to true-dof Vector%s in the domain FE space. */
737  /** This method is implemented in the base class, based on
738  BackwardOperator(), however, derived classes can overload the
739  construction, if necessary. */
741  {
743  }
744 };
745 
746 
747 /** @brief Transfer data between a coarse mesh and an embedded refined mesh
748  using interpolation. */
749 /** The forward, coarse-to-fine, transfer uses nodal interpolation. The
750  backward, fine-to-coarse, transfer is defined locally (on a coarse element)
751  as B = (F^t M_f F)^{-1} F^t M_f, where F is the forward transfer matrix, and
752  M_f is a mass matrix on the union of all fine elements comprising the coarse
753  element. Note that the backward transfer operator, B, is a left inverse of
754  the forward transfer operator, F, i.e. B F = I. Both F and B are defined in
755  reference space and do not depend on the actual physical shape of the mesh
756  elements.
757 
758  It is assumed that both the coarse and the fine FiniteElementSpace%s use
759  compatible types of elements, e.g. finite elements with the same map-type
760  (VALUE, INTEGRAL, H_DIV, H_CURL - see class FiniteElement). Generally, the
761  FE spaces can have different orders, however, in order for the backward
762  operator to be well-defined, the (local) number of the fine dofs should not
763  be smaller than the number of coarse dofs. */
765 {
766 protected:
767  BilinearFormIntegrator *mass_integ; ///< Ownership depends on #own_mass_integ
768  bool own_mass_integ; ///< Ownership flag for #mass_integ
769 
770  OperatorHandle F; ///< Forward, coarse-to-fine, operator
771  OperatorHandle B; ///< Backward, fine-to-coarse, operator
772 
773 public:
775  FiniteElementSpace &fine_fes)
776  : GridTransfer(coarse_fes, fine_fes),
777  mass_integ(NULL), own_mass_integ(false)
778  { }
779 
780  virtual ~InterpolationGridTransfer();
781 
782  /** @brief Assign a mass integrator to be used in the construction of the
783  backward, fine-to-coarse, transfer operator. */
784  void SetMassIntegrator(BilinearFormIntegrator *mass_integ_,
785  bool own_mass_integ_ = true);
786 
787  virtual const Operator &ForwardOperator();
788 
789  virtual const Operator &BackwardOperator();
790 };
791 
792 
793 /** @brief Transfer data between a coarse mesh and an embedded refined mesh
794  using L2 projection. */
795 /** The forward, coarse-to-fine, transfer uses L2 projection. The backward,
796  fine-to-coarse, transfer is defined locally (on a coarse element) as
797  B = (F^t M_f F)^{-1} F^t M_f, where F is the forward transfer matrix, and
798  M_f is the mass matrix on the union of all fine elements comprising the
799  coarse element. Note that the backward transfer operator, B, is a left
800  inverse of the forward transfer operator, F, i.e. B F = I. Both F and B are
801  defined in physical space and, generally, vary between different mesh
802  elements.
803 
804  This class currently only fully supports L2 finite element spaces and fine
805  meshes that are a uniform refinement of the coarse mesh. Generally, the
806  coarse and fine FE spaces can have different orders, however, in order for
807  the backward operator to be well-defined, the number of the fine dofs (in a
808  coarse element) should not be smaller than the number of coarse dofs.
809 
810  If used on H1 finite element spaces, the transfer will be performed locally,
811  and the value of shared (interface) degrees of freedom will be determined by
812  the value of the last transfer to be performed (according to the element
813  numbering in the finite element space). As a consequence, the mass
814  conservation properties for this operator from the L2 case do not carry over
815  to H1 spaces. */
817 {
818 protected:
819  /** Class representing projection operator between a high-order L2 finite
820  element space on a coarse mesh, and a low-order L2 finite element space
821  on a refined mesh (LOR). We assume that the low-order space, fes_lor,
822  lives on a mesh obtained by refining the mesh of the high-order space,
823  fes_ho. */
824  class L2Projection : public Operator
825  {
826  const FiniteElementSpace &fes_ho;
827  const FiniteElementSpace &fes_lor;
828 
829  int ndof_lor, ndof_ho, nref;
830 
831  Table ho2lor;
832 
833  DenseTensor R, P;
834 
835  public:
836  L2Projection(const FiniteElementSpace &fes_ho_,
837  const FiniteElementSpace &fes_lor_);
838  /// Perform the L2 projection onto the LOR space
839  virtual void Mult(const Vector &x, Vector &y) const;
840  /// Perform the mass conservative left-inverse prolongation operation.
841  /// This functionality is also provided as an Operator by L2Prolongation.
842  void Prolongate(const Vector &x, Vector &y) const;
843  virtual ~L2Projection() { }
844  };
845 
846  /** Mass-conservative prolongation operator going in the opposite direction
847  as L2Projection. This operator is a left inverse to the L2Projection. */
848  class L2Prolongation : public Operator
849  {
850  const L2Projection &l2proj;
851 
852  public:
853  L2Prolongation(const L2Projection &l2proj_) : l2proj(l2proj_) { }
854  void Mult(const Vector &x, Vector &y) const
855  {
856  l2proj.Prolongate(x, y);
857  }
858  virtual ~L2Prolongation() { }
859  };
860 
861  L2Projection *F; ///< Forward, coarse-to-fine, operator
862  L2Prolongation *B; ///< Backward, fine-to-coarse, operator
863 
864 public:
866  FiniteElementSpace &fine_fes)
867  : GridTransfer(coarse_fes, fine_fes),
868  F(NULL), B(NULL)
869  { }
870 
871  virtual const Operator &ForwardOperator();
872 
873  virtual const Operator &BackwardOperator();
874 };
875 
876 
877 /// Operator that converts FiniteElementSpace L-vectors to E-vectors.
878 /** Objects of this type are typically created and owned by FiniteElementSpace
879  objects, see FiniteElementSpace::GetElementRestriction(). */
881 {
882 protected:
884  const int ne;
885  const int vdim;
886  const bool byvdim;
887  const int ndofs;
888  const int dof;
889  const int nedofs;
892 
893 public:
895  void Mult(const Vector &x, Vector &y) const;
896  void MultTranspose(const Vector &x, Vector &y) const;
897 };
898 
899 
900 /** @brief A class that performs interpolation from an E-vector to quadrature
901  point values and/or derivatives (Q-vectors). */
902 /** An E-vector represents the element-wise discontinuous version of the FE
903  space and can be obtained, for example, from a GridFunction using the
904  Operator returned by FiniteElementSpace::GetElementRestriction().
905 
906  The target quadrature points in the elements can be described either by an
907  IntegrationRule (all mesh elements must be of the same type in this case) or
908  by a QuadratureSpace. */
910 {
911 protected:
912  friend class FiniteElementSpace; // Needs access to qspace and IntRule
913 
914  const FiniteElementSpace *fespace; ///< Not owned
915  const QuadratureSpace *qspace; ///< Not owned
916  const IntegrationRule *IntRule; ///< Not owned
917 
918  mutable bool use_tensor_products;
919 
920  static const int MAX_NQ2D = 100;
921  static const int MAX_ND2D = 100;
922  static const int MAX_VDIM2D = 2;
923 
924  static const int MAX_NQ3D = 1000;
925  static const int MAX_ND3D = 1000;
926  static const int MAX_VDIM3D = 3;
927 
928 public:
930  {
931  VALUES = 1 << 0, ///< Evaluate the values at quadrature points
932  DERIVATIVES = 1 << 1, ///< Evaluate the derivatives at quadrature points
933  /** @brief Assuming the derivative at quadrature points form a matrix,
934  this flag can be used to compute and store their determinants. This
935  flag can only be used in Mult(). */
936  DETERMINANTS = 1 << 2
937  };
938 
940  const IntegrationRule &ir);
941 
943  const QuadratureSpace &qs);
944 
945  /** @brief Disable the use of tensor product evaluations, for tensor-product
946  elements, e.g. quads and hexes. */
947  /** Currently, tensor product evaluations are not implemented and this method
948  has no effect. */
949  void DisableTensorProducts(bool disable = true) const
950  { use_tensor_products = !disable; }
951 
952  /// Interpolate the E-vector @a e_vec to quadrature points.
953  /** The @a eval_flags are a bitwise mask of constants from the EvalFlags
954  enumeration. When the VALUES flag is set, the values at quadrature points
955  are computed and stored in the Vector @a q_val. Similarly, when the flag
956  DERIVATIVES is set, the derivatives are computed and stored in @a q_der.
957  When the DETERMINANTS flags is set, it is assumed that the derivatives
958  form a matrix at each quadrature point (i.e. the associated
959  FiniteElementSpace is a vector space) and their determinants are computed
960  and stored in @a q_det. */
961  void Mult(const Vector &e_vec, unsigned eval_flags,
962  Vector &q_val, Vector &q_der, Vector &q_det) const;
963 
964  /// Perform the transpose operation of Mult(). (TODO)
965  void MultTranspose(unsigned eval_flags, const Vector &q_val,
966  const Vector &q_der, Vector &e_vec) const;
967 
968  // Compute kernels follow (cannot be private or protected with nvcc)
969 
970  /// Template compute kernel for 2D.
971  template<const int T_VDIM = 0, const int T_ND = 0, const int T_NQ = 0>
972  static void Eval2D(const int NE,
973  const int vdim,
974  const DofToQuad &maps,
975  const Vector &e_vec,
976  Vector &q_val,
977  Vector &q_der,
978  Vector &q_det,
979  const int eval_flags);
980 
981  /// Template compute kernel for 3D.
982  template<const int T_VDIM = 0, const int T_ND = 0, const int T_NQ = 0>
983  static void Eval3D(const int NE,
984  const int vdim,
985  const DofToQuad &maps,
986  const Vector &e_vec,
987  Vector &q_val,
988  Vector &q_der,
989  Vector &q_det,
990  const int eval_flags);
991 };
992 
993 }
994 
995 #endif
Abstract class for Finite Elements.
Definition: fe.hpp:229
const Operator & MakeTrueOperator(FiniteElementSpace &fes_in, FiniteElementSpace &fes_out, const Operator &oper, OperatorHandle &t_oper)
Definition: fespace.cpp:2259
void Mult(const Vector &x, Vector &y) const
Operator application: y=A(x).
Definition: fespace.cpp:2715
Ordering::Type GetOrdering() const
Return the ordering method.
Definition: fespace.hpp:359
int GetVSize() const
Return the number of vector dofs, i.e. GetNDofs() x GetVDim().
Definition: fespace.hpp:347
void GetVertexVDofs(int i, Array< int > &vdofs) const
Definition: fespace.cpp:196
int GetAttribute(int i) const
Definition: fespace.hpp:410
int GetNDofs() const
Returns number of degrees of freedom.
Definition: fespace.hpp:344
int ndofs
Number of degrees of freedom. Number of unknowns is ndofs * vdim.
Definition: fespace.hpp:105
Class for an integration rule - an Array of IntegrationPoint.
Definition: intrules.hpp:85
int GetBdrAttribute(int i) const
Return the attribute of boundary element i.
Definition: mesh.hpp:976
const Table & GetBdrElementToDofTable() const
Definition: fespace.hpp:482
void GetElementInteriorVDofs(int i, Array< int > &vdofs) const
Definition: fespace.cpp:202
int DofToVDof(int dof, int vd, int ndofs=-1) const
Definition: fespace.cpp:143
const SparseMatrix * GetConformingProlongation() const
The returned SparseMatrix is owned by the FiniteElementSpace.
Definition: fespace.cpp:769
virtual void Update(bool want_transform=true)
Definition: fespace.cpp:1942
static int Map(int ndofs, int vdim, int dof, int vd)
void GetElementVertices(int i, Array< int > &vertices) const
Returns the vertices of element i.
Definition: fespace.hpp:390
BilinearFormIntegrator * mass_integ
Ownership depends on own_mass_integ.
Definition: fespace.hpp:767
void BuildElementToDofTable() const
Definition: fespace.cpp:214
const IntegrationRule & GetElementIntRule(int idx) const
Get the IntegrationRule associated with mesh element idx.
Definition: fespace.hpp:660
DerefinementOperator(const FiniteElementSpace *f_fes, const FiniteElementSpace *c_fes, BilinearFormIntegrator *mass_integ)
Definition: fespace.cpp:1033
void Prolongate(const Vector &x, Vector &y) const
Definition: fespace.cpp:2596
L2Prolongation(const L2Projection &l2proj_)
Definition: fespace.hpp:853
static const int NumGeom
Definition: geom.hpp:41
static const int MAX_ND2D
Definition: fespace.hpp:921
FiniteElementSpace & dom_fes
Domain FE space.
Definition: fespace.hpp:675
const FiniteElement * GetTraceElement(int i, Geometry::Type geom_type) const
Return the trace element from element &#39;i&#39; to the given &#39;geom_type&#39;.
Definition: fespace.cpp:1835
bool Conforming() const
Definition: mesh.hpp:1135
SparseMatrix * RefinementMatrix(int old_ndofs, const Table *old_elem_dof)
Definition: fespace.cpp:926
int GetNBE() const
Returns number of boundary elements.
Definition: mesh.hpp:679
SparseMatrix * DerefinementMatrix(int old_ndofs, const Table *old_elem_dof)
Calculate GridFunction restriction matrix after mesh derefinement.
Definition: fespace.cpp:1192
void GetElementVDofs(int i, Array< int > &vdofs) const
Returns indexes of degrees of freedom in array dofs for i&#39;th element.
Definition: fespace.cpp:172
Ordering::Type ordering
Definition: fespace.hpp:102
const FiniteElementSpace * fespace
Not owned.
Definition: fespace.hpp:914
const Geometry::Type geom
Definition: ex1.cpp:40
static int DecodeDof(int dof, double &sign)
Helper to remove encoded sign from a DOF.
Definition: fespace.hpp:144
SparseMatrix * RefinementMatrix_main(const int coarse_ndofs, const Table &coarse_elem_dof, const DenseTensor localP[]) const
Definition: fespace.cpp:840
int GetSize() const
Return the total number of quadrature points.
Definition: fespace.hpp:657
Element::Type GetElementType(int i) const
Returns the type of element i.
Definition: mesh.cpp:4351
Pointer to an Operator of a specified type.
Definition: handle.hpp:33
static const int MAX_NQ3D
Definition: fespace.hpp:924
int VDofToDof(int vdof) const
Definition: fespace.hpp:445
void Mult(const Vector &x, Vector &y) const
Operator application: y=A(x).
Definition: fespace.hpp:854
Element::Type GetBdrElementType(int i) const
Returns the type of boundary element i.
Definition: mesh.cpp:4356
int GetFaceOrder(int i) const
Returns the order of the i&#39;th face finite element.
Definition: fespace.cpp:101
virtual const Operator & BackwardOperator()=0
Return an Operator that transfers GridFunctions from the range FE space back to GridFunctions in the ...
virtual void GetEssentialVDofs(const Array< int > &bdr_attr_is_ess, Array< int > &ess_vdofs, int component=-1) const
Definition: fespace.cpp:297
Data type dense matrix using column-major storage.
Definition: densemat.hpp:23
int vdim
Vector dimension (number of unknowns per degree of freedom).
Definition: fespace.hpp:97
OperatorHandle L2E_lex
Definition: fespace.hpp:130
void GetEntityDofs(int entity, int index, Array< int > &dofs) const
Helper to get vertex, edge or face DOFs (entity=0,1,2 resp.).
Definition: fespace.cpp:571
Operator::Type oper_type
Desired Operator::Type for the construction of all operators defined by the underlying transfer algor...
Definition: fespace.hpp:681
int GetNumElementInteriorDofs(int i) const
Definition: fespace.hpp:434
int GetBdrAttribute(int i) const
Definition: fespace.hpp:412
virtual void GetEssentialTrueDofs(const Array< int > &bdr_attr_is_ess, Array< int > &ess_tdof_list, int component=-1)
Definition: fespace.cpp:366
int GetNE() const
Returns number of elements.
Definition: mesh.hpp:676
int GetBdrElementType(int i) const
Returns the type of boundary element i.
Definition: fespace.hpp:394
int GetNV() const
Returns number of vertices in the mesh.
Definition: fespace.hpp:371
Evaluate the derivatives at quadrature points.
Definition: fespace.hpp:932
const FiniteElement * GetFaceElement(int i) const
Definition: fespace.cpp:1807
virtual void GetTrueTransferOperator(const FiniteElementSpace &coarse_fes, OperatorHandle &T) const
Construct and return an Operator that can be used to transfer true-dof data from coarse_fes, defined on a coarse mesh, to this FE space, defined on a refined mesh.
Definition: fespace.cpp:1901
FiniteElementCollection * Load(Mesh *m, std::istream &input)
Read a FiniteElementSpace from a stream. The returned FiniteElementCollection is owned by the caller...
Definition: fespace.cpp:2088
virtual const Operator & ForwardOperator()
Return an Operator that transfers GridFunctions from the domain FE space to GridFunctions in the rang...
Definition: fespace.cpp:2624
const IntegrationRule * int_rule[Geometry::NumGeom]
Definition: fespace.hpp:637
L2ProjectionGridTransfer(FiniteElementSpace &coarse_fes, FiniteElementSpace &fine_fes)
Definition: fespace.hpp:865
const FiniteElementCollection * fec
Associated FE collection (not owned).
Definition: fespace.hpp:94
OperatorHandle Th
Transformation to apply to GridFunctions after space Update().
Definition: fespace.hpp:127
NURBSExtension * GetNURBSext()
Definition: fespace.hpp:275
Geometry::Type GetElementBaseGeometry(int i) const
Definition: mesh.hpp:759
void GetLocalDerefinementMatrices(Geometry::Type geom, DenseTensor &localR) const
Definition: fespace.cpp:1166
int GetConformingVSize() const
Definition: fespace.hpp:356
static const int MAX_ND3D
Definition: fespace.hpp:925
Array< int > dof_elem_array
Definition: fespace.hpp:113
static void MarkerToList(const Array< int > &marker, Array< int > &list)
Convert a Boolean marker array to a list containing all marked indices.
Definition: fespace.cpp:385
Array< int > indices
Definition: fespace.hpp:891
virtual const Operator & TrueBackwardOperator()
Return an Operator that transfers true-dof Vectors from the range FE space back to true-dof Vectors i...
Definition: fespace.hpp:740
void GetFaceInteriorDofs(int i, Array< int > &dofs) const
Definition: fespace.cpp:1766
virtual int DofForGeometry(Geometry::Type GeomType) const =0
virtual void GetFaceDofs(int i, Array< int > &dofs) const
Definition: fespace.cpp:1641
virtual ~GridTransfer()
Virtual destructor.
Definition: fespace.hpp:709
ElementTransformation * GetBdrElementTransformation(int i)
Returns the transformation defining the i-th boundary element.
Definition: mesh.cpp:414
FiniteElementSpace & ran_fes
Range FE space.
Definition: fespace.hpp:676
void ConvertFromConformingVDofs(const Array< int > &cdofs, Array< int > &dofs)
Definition: fespace.cpp:422
const Operator * GetElementRestriction(ElementDofOrdering e_ordering) const
Return an Operator that converts L-vectors to E-vectors.
Definition: fespace.cpp:789
int GetNumFaces() const
Return the number of faces (3D), edges (2D) or vertices (1D).
Definition: mesh.cpp:3663
void DisableTensorProducts(bool disable=true) const
Disable the use of tensor product evaluations, for tensor-product elements, e.g. quads and hexes...
Definition: fespace.hpp:949
const SparseMatrix * GetConformingRestriction() const
The returned SparseMatrix is owned by the FiniteElementSpace.
Definition: fespace.cpp:776
virtual void Mult(const Vector &x, Vector &y) const
Operator application: y=A(x).
Definition: fespace.cpp:1135
int GetNE() const
Returns number of elements in the mesh.
Definition: fespace.hpp:374
void GetVertexDofs(int i, Array< int > &dofs) const
Definition: fespace.cpp:1729
virtual const Operator & BackwardOperator()
Return an Operator that transfers GridFunctions from the range FE space back to GridFunctions in the ...
Definition: fespace.cpp:2412
bool Nonconforming() const
Definition: mesh.hpp:1136
const FiniteElement * GetEdgeElement(int i) const
Definition: fespace.cpp:1830
virtual void Mult(const Vector &x, Vector &y) const
Operator application: y=A(x).
Definition: fespace.cpp:984
int GetNBE() const
Returns number of boundary elements in the mesh.
Definition: fespace.hpp:383
Data type sparse matrix.
Definition: sparsemat.hpp:40
Native ordering as defined by the FiniteElement.
Mesh * GetMesh() const
Returns the mesh.
Definition: fespace.hpp:272
OperatorHandle F
Forward, coarse-to-fine, operator.
Definition: fespace.hpp:770
virtual const Operator & ForwardOperator()=0
Return an Operator that transfers GridFunctions from the domain FE space to GridFunctions in the rang...
A class that performs interpolation from an E-vector to quadrature point values and/or derivatives (Q...
Definition: fespace.hpp:909
GridTransfer(FiniteElementSpace &dom_fes_, FiniteElementSpace &ran_fes_)
Definition: fespace.cpp:2244
Type
Ordering methods:
Definition: fespace.hpp:30
SparseMatrix * cR
Conforming restriction matrix such that cR.cP=I.
Definition: fespace.hpp:123
const IntegrationRule * IntRule
Not owned.
Definition: fespace.hpp:916
SparseMatrix * D2Const_GlobalRestrictionMatrix(FiniteElementSpace *cfes)
Definition: fespace.cpp:463
SparseMatrix * D2C_GlobalRestrictionMatrix(FiniteElementSpace *cfes)
Definition: fespace.cpp:431
virtual const Operator * GetProlongationMatrix() const
The returned Operator is owned by the FiniteElementSpace.
Definition: fespace.hpp:288
virtual void GetElementDofs(int i, Array< int > &dofs) const
Returns indexes of degrees of freedom in array dofs for i&#39;th element.
Definition: fespace.cpp:1439
void BuildConformingInterpolation() const
Calculate the cP and cR matrices for a nonconforming mesh.
Definition: fespace.cpp:582
NURBSExtension * StealNURBSext()
Definition: fespace.cpp:1341
int GetElementForDof(int i) const
Definition: fespace.hpp:484
ElementTransformation * GetBdrElementTransformation(int i) const
Returns ElementTransformation for the i-th boundary element.
Definition: fespace.hpp:407
L2Projection(const FiniteElementSpace &fes_ho_, const FiniteElementSpace &fes_lor_)
Definition: fespace.cpp:2454
static void AddDependencies(SparseMatrix &deps, Array< int > &master_dofs, Array< int > &slave_dofs, DenseMatrix &I)
Definition: fespace.cpp:532
void GetElementTransformation(int i, IsoparametricTransformation *ElTr)
Returns the transformation defining the i-th element in the user-defined variable ElTr...
Definition: fespace.hpp:403
Array< int > dof_ldof_array
Definition: fespace.hpp:113
InterpolationGridTransfer(FiniteElementSpace &coarse_fes, FiniteElementSpace &fine_fes)
Definition: fespace.hpp:774
virtual int GetTrueVSize() const
Return the number of vector true (conforming) dofs.
Definition: fespace.hpp:350
bool Conforming() const
Definition: fespace.hpp:278
OperatorHandle bw_t_oper
Backward true-dof operator.
Definition: fespace.hpp:684
SparseMatrix * cP
Definition: fespace.hpp:121
ElementRestriction(const FiniteElementSpace &, ElementDofOrdering)
Definition: fespace.cpp:2641
int GetVDim() const
Returns vector dimension.
Definition: fespace.hpp:336
int GetNF() const
Returns number of faces (i.e. co-dimension 1 entities) in the mesh.
Definition: fespace.hpp:380
static const int MAX_VDIM2D
Definition: fespace.hpp:922
void SetOperatorType(Operator::Type type)
Set the desired Operator::Type for the construction of all operators defined by the underlying transf...
Definition: fespace.hpp:716
Array< QuadratureInterpolator * > E2Q_array
Definition: fespace.hpp:132
void GetEdgeDofs(int i, Array< int > &dofs) const
Definition: fespace.cpp:1700
Operator * Ptr() const
Access the underlying Operator pointer.
Definition: handle.hpp:82
OperatorHandle L2E_nat
The element restriction operators, see GetElementRestriction().
Definition: fespace.hpp:130
void GetTransferOperator(const FiniteElementSpace &coarse_fes, OperatorHandle &T) const
Construct and return an Operator that can be used to transfer GridFunction data from coarse_fes...
Definition: fespace.cpp:1876
OperatorHandle B
Backward, fine-to-coarse, operator.
Definition: fespace.hpp:771
void GetElementInteriorDofs(int i, Array< int > &dofs) const
Definition: fespace.cpp:1741
void GetLocalRefinementMatrices(Geometry::Type geom, DenseTensor &localP) const
Definition: fespace.cpp:902
static void DofsToVDofs(int ndofs, int vdim, Array< int > &dofs)
virtual const Operator & BackwardOperator()
Return an Operator that transfers GridFunctions from the range FE space back to GridFunctions in the ...
Definition: fespace.cpp:2630
Base class for transfer algorithms that construct transfer Operators between two finite element (FE) ...
Definition: fespace.hpp:672
OperatorHandle fw_t_oper
Forward true-dof operator.
Definition: fespace.hpp:683
const FiniteElementSpace & fes
Definition: fespace.hpp:883
void SetMassIntegrator(BilinearFormIntegrator *mass_integ_, bool own_mass_integ_=true)
Assign a mass integrator to be used in the construction of the backward, fine-to-coarse, transfer operator.
Definition: fespace.cpp:2370
ElementTransformation * GetElementTransformation(int i) const
Returns ElementTransformation for the i-th element.
Definition: fespace.hpp:398
QuadratureInterpolator(const FiniteElementSpace &fes, const IntegrationRule &ir)
Definition: fespace.cpp:2769
Type
Enumeration defining IDs for some classes derived from Operator.
Definition: operator.hpp:135
Abstract base class BilinearFormIntegrator.
Definition: bilininteg.hpp:23
RefinementOperator(const FiniteElementSpace *fespace, Table *old_elem_dof, int old_ndofs)
Definition: fespace.cpp:943
const NURBSExtension * GetNURBSext() const
Definition: fespace.hpp:274
Operator that converts FiniteElementSpace L-vectors to E-vectors.
Definition: fespace.hpp:880
Class FiniteElementSpace - responsible for providing FEM view of the mesh, mainly managing the set of...
Definition: fespace.hpp:85
Transfer data between a coarse mesh and an embedded refined mesh using interpolation.
Definition: fespace.hpp:764
int GetLocalDofForDof(int i) const
Definition: fespace.hpp:485
virtual ~FiniteElementSpace()
Definition: fespace.cpp:1841
SparseMatrix * H2L_GlobalRestrictionMatrix(FiniteElementSpace *lfes)
Definition: fespace.cpp:494
static void Eval2D(const int NE, const int vdim, const DofToQuad &maps, const Vector &e_vec, Vector &q_val, Vector &q_der, Vector &q_det, const int eval_flags)
Template compute kernel for 2D.
Definition: fespace.cpp:2798
virtual ~QuadratureSpace()
Definition: fespace.hpp:654
void GetUpdateOperator(OperatorHandle &T)
Return the update operator in the given OperatorHandle, T.
Definition: fespace.hpp:593
L2Prolongation * B
Backward, fine-to-coarse, operator.
Definition: fespace.hpp:862
int GetOrder(int i) const
Returns the order of the i&#39;th finite element.
Definition: fespace.cpp:95
void Save(std::ostream &out) const
Write the QuadratureSpace to the stream out.
Definition: fespace.cpp:2236
void Clear()
Clear the OperatorHandle, deleting the held Operator (if owned), while leaving the type id unchanged...
Definition: handle.hpp:116
Structure representing the matrices/tensors needed to evaluate (in reference space) the values...
Definition: fe.hpp:125
void MakeVDimMatrix(SparseMatrix &mat) const
Definition: fespace.cpp:741
bool own_mass_integ
Ownership flag for mass_integ.
Definition: fespace.hpp:768
The ordering method used when the number of unknowns per mesh node (vector dimension) is bigger than ...
Definition: fespace.hpp:26
void GetEdgeInteriorVDofs(int i, Array< int > &vdofs) const
Definition: fespace.cpp:208
L2Projection * F
Forward, coarse-to-fine, operator.
Definition: fespace.hpp:861
int GetElementType(int i) const
Returns the type of element i.
Definition: fespace.hpp:386
FiniteElementSpace(Mesh *mesh, const FiniteElementCollection *fec, int vdim=1, int ordering=Ordering::byNODES)
Definition: fespace.hpp:257
bool Nonconforming() const
Definition: fespace.hpp:279
virtual void GetBdrElementDofs(int i, Array< int > &dofs) const
Returns indexes of degrees of freedom for i&#39;th boundary element.
Definition: fespace.cpp:1558
Assuming the derivative at quadrature points form a matrix, this flag can be used to compute and stor...
Definition: fespace.hpp:936
int GetNV() const
Returns number of vertices. Vertices are only at the corners of elements, where you would expect them...
Definition: mesh.hpp:673
int GetNVDofs() const
Number of all scalar vertex dofs.
Definition: fespace.hpp:364
virtual const Operator & ForwardOperator()
Return an Operator that transfers GridFunctions from the domain FE space to GridFunctions in the rang...
Definition: fespace.cpp:2379
void ConvertToConformingVDofs(const Array< int > &dofs, Array< int > &cdofs)
Definition: fespace.cpp:414
static void Eval3D(const int NE, const int vdim, const DofToQuad &maps, const Vector &e_vec, Vector &q_val, Vector &q_der, Vector &q_det, const int eval_flags)
Template compute kernel for 3D.
Definition: fespace.cpp:2886
void SetOperatorOwner(bool own=true)
Set the ownership flag for the held Operator.
Definition: handle.hpp:112
bool Parallel() const
Definition: fespace.hpp:689
virtual const Operator & TrueForwardOperator()
Return an Operator that transfers true-dof Vectors from the domain FE space to true-dof Vectors in th...
Definition: fespace.hpp:730
Mesh * mesh
The mesh that FE space lives on (not owned).
Definition: fespace.hpp:91
void GetElementTransformation(int i, IsoparametricTransformation *ElTr)
Definition: mesh.cpp:336
void ReorderElementToDofTable()
Reorder the scalar DOFs based on the element ordering.
Definition: fespace.cpp:243
GridFunction interpolation operator applicable after mesh refinement.
Definition: fespace.hpp:162
Array< int > offsets
Definition: fespace.hpp:890
int GetNConformingDofs() const
Definition: fespace.cpp:783
static const int MAX_VDIM3D
Definition: fespace.hpp:926
void GetEdgeInteriorDofs(int i, Array< int > &dofs) const
Definition: fespace.cpp:1754
void Mult(const Vector &e_vec, unsigned eval_flags, Vector &q_val, Vector &q_der, Vector &q_det) const
Interpolate the E-vector e_vec to quadrature points.
Definition: fespace.cpp:2978
const QuadratureInterpolator * GetQuadratureInterpolator(const IntegrationRule &ir) const
Return a QuadratureInterpolator that interpolates E-vectors to quadrature point values and/or derivat...
Definition: fespace.cpp:812
ElementDofOrdering
Constants describing the possible orderings of the DOFs in one element.
Definition: fespace.hpp:63
const FiniteElement * GetFE(int i) const
Returns pointer to the FiniteElement associated with i&#39;th element.
Definition: fespace.cpp:1541
const Operator * GetUpdateOperator()
Get the GridFunction update operator.
Definition: fespace.hpp:590
Transfer data between a coarse mesh and an embedded refined mesh using L2 projection.
Definition: fespace.hpp:816
Lexicographic ordering for tensor-product FiniteElements.
FiniteElementSpace()
Default constructor: the object is invalid until initialized using the method Load().
Definition: fespace.cpp:59
QuadratureSpace(Mesh *mesh_, int order_)
Create a QuadratureSpace based on the global rules from IntRules.
Definition: fespace.hpp:648
void Constructor(Mesh *mesh, NURBSExtension *ext, const FiniteElementCollection *fec, int vdim=1, int ordering=Ordering::byNODES)
Help function for constructors + Load().
Definition: fespace.cpp:1295
void SetUpdateOperatorOwner(bool own)
Set the ownership of the update operator: if set to false, the Operator returned by GetUpdateOperator...
Definition: fespace.hpp:600
void GetFaceVDofs(int i, Array< int > &vdofs) const
Returns indexes of degrees of freedom for i&#39;th face element (2D and 3D).
Definition: fespace.cpp:184
void GetEdgeVDofs(int i, Array< int > &vdofs) const
Returns indexes of degrees of freedom for i&#39;th edge.
Definition: fespace.cpp:190
static void ListToMarker(const Array< int > &list, int marker_size, Array< int > &marker, int mark_val=-1)
Definition: fespace.cpp:403
virtual void Mult(const Vector &x, Vector &y) const
Perform the L2 projection onto the LOR space.
Definition: fespace.cpp:2570
Vector data type.
Definition: vector.hpp:48
const FiniteElementCollection * FEColl() const
Definition: fespace.hpp:361
virtual void UpdatesFinished()
Free the GridFunction update operator (if any), to save memory.
Definition: fespace.hpp:611
const Table & GetElementToDofTable() const
Definition: fespace.hpp:481
static const int MAX_NQ2D
Definition: fespace.hpp:920
void SetUpdateOperatorType(Operator::Type tid)
Specify the Operator::Type to be used by the update operators.
Definition: fespace.hpp:608
static bool DofFinalizable(int dof, const Array< bool > &finalized, const SparseMatrix &deps)
Definition: fespace.cpp:556
Class representing the storage layout of a QuadratureFunction.
Definition: fespace.hpp:628
const FiniteElement * GetBE(int i) const
Returns pointer to the FiniteElement for the i&#39;th boundary element.
Definition: fespace.cpp:1781
NURBSExtension * NURBSext
Definition: fespace.hpp:115
void GetBdrElementVDofs(int i, Array< int > &vdofs) const
Returns indexes of degrees of freedom for i&#39;th boundary element.
Definition: fespace.cpp:178
OutStream out(std::cout)
Global stream used by the library for standard output. Initially it uses the same std::streambuf as s...
Definition: globals.hpp:64
Abstract operator.
Definition: operator.hpp:21
virtual const SparseMatrix * GetRestrictionMatrix() const
The returned SparseMatrix is owned by the FiniteElementSpace.
Definition: fespace.hpp:292
long GetSequence() const
Return update counter (see Mesh::sequence)
Definition: fespace.hpp:614
Rank 3 tensor (array of matrices)
Definition: densemat.hpp:656
void Save(std::ostream &out) const
Definition: fespace.cpp:2019
int GetAttribute(int i) const
Return the attribute of element i.
Definition: mesh.hpp:970
Class representing a function through its values (scalar or vector) at quadrature points...
Definition: gridfunc.hpp:494
int GetNFDofs() const
Number of all scalar face-interior dofs.
Definition: fespace.hpp:368
int GetNEDofs() const
Number of all scalar edge-interior dofs.
Definition: fespace.hpp:366
void GetElementVertices(int i, Array< int > &v) const
Returns the indices of the vertices of element i.
Definition: mesh.hpp:803
static void AdjustVDofs(Array< int > &vdofs)
Definition: fespace.cpp:159
void SetType(Operator::Type tid)
Invoke Clear() and set a new type id.
Definition: handle.hpp:124
void MultTranspose(unsigned eval_flags, const Vector &q_val, const Vector &q_der, Vector &e_vec) const
Perform the transpose operation of Mult(). (TODO)
Definition: fespace.cpp:3119
Evaluate the values at quadrature points.
Definition: fespace.hpp:931
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: fespace.cpp:2741
const QuadratureSpace * qspace
Not owned.
Definition: fespace.hpp:915
void DofsToVDofs(Array< int > &dofs, int ndofs=-1) const
Definition: fespace.cpp:107
FiniteElementSpace(Mesh *mesh, NURBSExtension *ext, const FiniteElementCollection *fec, int vdim=1, int ordering=Ordering::byNODES)
Construct a NURBS FE space based on the given NURBSExtension, ext.
Definition: fespace.hpp:266