MFEM  v4.5.1
Finite element discretization library
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
transfer.hpp
Go to the documentation of this file.
1 // Copyright (c) 2010-2022, Lawrence Livermore National Security, LLC. Produced
2 // at the Lawrence Livermore National Laboratory. All Rights reserved. See files
3 // LICENSE and NOTICE for details. LLNL-CODE-806117.
4 //
5 // This file is part of the MFEM library. For more information and source code
6 // availability visit https://mfem.org.
7 //
8 // MFEM is free software; you can redistribute it and/or modify it under the
9 // terms of the BSD-3 license. We welcome feedback and contributions, see file
10 // CONTRIBUTING.md for details.
11 
12 #ifndef MFEM_TRANSFER_HPP
13 #define MFEM_TRANSFER_HPP
14 
15 #include "../linalg/linalg.hpp"
16 #include "fespace.hpp"
17 
18 #ifdef MFEM_USE_MPI
19 #include "pfespace.hpp"
20 #endif
21 
22 namespace mfem
23 {
24 
25 /** @brief Base class for transfer algorithms that construct transfer Operator%s
26  between two finite element (FE) spaces. */
27 /** Generally, the two FE spaces (domain and range) can be defined on different
28  meshes. */
30 {
31 protected:
32  FiniteElementSpace &dom_fes; ///< Domain FE space
33  FiniteElementSpace &ran_fes; ///< Range FE space
34 
35  /** @brief Desired Operator::Type for the construction of all operators
36  defined by the underlying transfer algorithm. It can be ignored by
37  derived classes. */
39 
40  OperatorHandle fw_t_oper; ///< Forward true-dof operator
41  OperatorHandle bw_t_oper; ///< Backward true-dof operator
42 
43 #ifdef MFEM_USE_MPI
44  bool parallel;
45 #endif
46  bool Parallel() const
47  {
48 #ifndef MFEM_USE_MPI
49  return false;
50 #else
51  return parallel;
52 #endif
53  }
54 
56  FiniteElementSpace &fes_out,
57  const Operator &oper,
58  OperatorHandle &t_oper);
59 
60 public:
61  /** Construct a transfer algorithm between the domain, @a dom_fes_, and
62  range, @a ran_fes_, FE spaces. */
64 
65  /// Virtual destructor
66  virtual ~GridTransfer() { }
67 
68  /** @brief Set the desired Operator::Type for the construction of all
69  operators defined by the underlying transfer algorithm. */
70  /** The default value is Operator::ANY_TYPE which typically corresponds to a
71  matrix-free operator representation. Note that derived classes are not
72  required to support this setting and can ignore it. */
73  void SetOperatorType(Operator::Type type) { oper_type = type; }
74 
75  /** @brief Return an Operator that transfers GridFunction%s from the domain
76  FE space to GridFunction%s in the range FE space. */
77  virtual const Operator &ForwardOperator() = 0;
78 
79  /** @brief Return an Operator that transfers GridFunction%s from the range FE
80  space back to GridFunction%s in the domain FE space. */
81  virtual const Operator &BackwardOperator() = 0;
82 
83  /** @brief Return an Operator that transfers true-dof Vector%s from the
84  domain FE space to true-dof Vector%s in the range FE space. */
85  /** This method is implemented in the base class, based on ForwardOperator(),
86  however, derived classes can overload the construction, if necessary. */
87  virtual const Operator &TrueForwardOperator()
88  {
90  }
91 
92  /** @brief Return an Operator that transfers true-dof Vector%s from the range
93  FE space back to true-dof Vector%s in the domain FE space. */
94  /** This method is implemented in the base class, based on
95  BackwardOperator(), however, derived classes can overload the
96  construction, if necessary. */
97  virtual const Operator &TrueBackwardOperator()
98  {
100  }
101 };
102 
103 
104 /** @brief Transfer data between a coarse mesh and an embedded refined mesh
105  using interpolation. */
106 /** The forward, coarse-to-fine, transfer uses nodal interpolation. The
107  backward, fine-to-coarse, transfer is defined locally (on a coarse element)
108  as B = (F^t M_f F)^{-1} F^t M_f, where F is the forward transfer matrix, and
109  M_f is a mass matrix on the union of all fine elements comprising the coarse
110  element. Note that the backward transfer operator, B, is a left inverse of
111  the forward transfer operator, F, i.e. B F = I. Both F and B are defined in
112  reference space and do not depend on the actual physical shape of the mesh
113  elements.
114 
115  It is assumed that both the coarse and the fine FiniteElementSpace%s use
116  compatible types of elements, e.g. finite elements with the same map-type
117  (VALUE, INTEGRAL, H_DIV, H_CURL - see class FiniteElement). Generally, the
118  FE spaces can have different orders, however, in order for the backward
119  operator to be well-defined, the (local) number of the fine dofs should not
120  be smaller than the number of coarse dofs. */
122 {
123 protected:
124  BilinearFormIntegrator *mass_integ; ///< Ownership depends on #own_mass_integ
125  bool own_mass_integ; ///< Ownership flag for #mass_integ
126 
127  OperatorHandle F; ///< Forward, coarse-to-fine, operator
128  OperatorHandle B; ///< Backward, fine-to-coarse, operator
129 
130 public:
132  FiniteElementSpace &fine_fes)
133  : GridTransfer(coarse_fes, fine_fes),
134  mass_integ(NULL), own_mass_integ(false)
135  { }
136 
137  virtual ~InterpolationGridTransfer();
138 
139  /** @brief Assign a mass integrator to be used in the construction of the
140  backward, fine-to-coarse, transfer operator. */
141  void SetMassIntegrator(BilinearFormIntegrator *mass_integ_,
142  bool own_mass_integ_ = true);
143 
144  virtual const Operator &ForwardOperator();
145 
146  virtual const Operator &BackwardOperator();
147 };
148 
149 
150 /** @brief Transfer data in L2 and H1 finite element spaces between a coarse
151  mesh and an embedded refined mesh using L2 projection. */
152 /** The forward, coarse-to-fine, transfer uses L2 projection. The backward,
153  fine-to-coarse, transfer is defined as B = (F^t M_f F)^{-1} F^t M_f, where F
154  is the forward transfer matrix, and M_f is the mass matrix on the coarse
155  element. For L2 spaces, M_f is the mass matrix on the union of all fine
156  elements comprising the coarse element. For H1 spaces, M_f is a diagonal
157  (lumped) mass matrix computed through row-summation. Note that the backward
158  transfer operator, B, is a left inverse of the forward transfer operator, F,
159  i.e. B F = I. Both F and B are defined in physical space and, generally for
160  L2 spaces, vary between different mesh elements.
161 
162  This class supports H1 and L2 finite element spaces. Fine meshes are a
163  uniform refinement of the coarse mesh, usually created through
164  Mesh::MakeRefined. Generally, the coarse and fine FE spaces can have
165  different orders, however, in order for the backward operator to be
166  well-defined, the number of fine dofs (in a coarse element) should not be
167  smaller than the number of coarse dofs. */
169 {
170 protected:
171  /** Abstract class representing projection operator between a high-order
172  finite element space on a coarse mesh, and a low-order finite element
173  space on a refined mesh (LOR). We assume that the low-order space,
174  fes_lor, lives on a mesh obtained by refining the mesh of the high-order
175  space, fes_ho. */
176  class L2Projection : public Operator
177  {
178  public:
179  virtual void Prolongate(const Vector& x, Vector& y) const = 0;
180  virtual void ProlongateTranspose(const Vector& x, Vector& y) const = 0;
181  /// Sets relative tolerance and absolute tolerance in preconditioned
182  /// conjugate gradient solver. Only used for H1 spaces.
183  virtual void SetRelTol(double p_rtol_) = 0;
184  virtual void SetAbsTol(double p_atol_) = 0;
185  protected:
188 
190 
191  L2Projection(const FiniteElementSpace& fes_ho_,
192  const FiniteElementSpace& fes_lor_);
193 
194  void BuildHo2Lor(int nel_ho, int nel_lor,
195  const CoarseFineTransformations& cf_tr);
196 
197  void ElemMixedMass(Geometry::Type geom, const FiniteElement& fe_ho,
198  const FiniteElement& fe_lor, ElementTransformation* el_tr,
200  DenseMatrix& M_mixed_el) const;
201  };
202 
203  /** Class for projection operator between a L2 high-order finite element
204  space on a coarse mesh, and a L2 low-order finite element space on a
205  refined mesh (LOR). */
207  {
208  // The restriction and prolongation operators are represented as dense
209  // elementwise matrices (of potentially different sizes, because of mixed
210  // meshes or p-refinement). The matrix entries are stored in the R and P
211  // arrays. The entries of the i'th high-order element are stored at the
212  // index given by offsets[i].
213  mutable Array<double> R, P;
214  Array<int> offsets;
215 
216  public:
218  const FiniteElementSpace& fes_lor_);
219  /// Maps <tt>x</tt>, primal field coefficients defined on a coarse mesh
220  /// with a higher order L2 finite element space, to <tt>y</tt>, primal
221  /// field coefficients defined on a refined mesh with a low order L2
222  /// finite element space. Refined mesh should be a uniform refinement of
223  /// the coarse mesh. Coefficients are computed through minimization of L2
224  /// error between the fields.
225  virtual void Mult(const Vector& x, Vector& y) const;
226  /// Maps <tt>x</tt>, dual field coefficients defined on a refined mesh
227  /// with a low order L2 finite element space, to <tt>y</tt>, dual field
228  /// coefficients defined on a coarse mesh with a higher order L2 finite
229  /// element space. Refined mesh should be a uniform refinement of the
230  /// coarse mesh. Coefficients are computed through minimization of L2
231  /// error between the primal fields. Note, if the <tt>x</tt>-coefficients
232  /// come from ProlongateTranspose, then mass is conserved.
233  virtual void MultTranspose(const Vector& x, Vector& y) const;
234  /// Maps <tt>x</tt>, primal field coefficients defined on a refined mesh
235  /// with a low order L2 finite element space, to <tt>y</tt>, primal field
236  /// coefficients defined on a coarse mesh with a higher order L2 finite
237  /// element space. Refined mesh should be a uniform refinement of the
238  /// coarse mesh. Coefficients are computed from the mass conservative
239  /// left-inverse prolongation operation. This functionality is also
240  /// provided as an Operator by L2Prolongation.
241  virtual void Prolongate(const Vector& x, Vector& y) const;
242  /// Maps <tt>x</tt>, dual field coefficients defined on a coarse mesh with
243  /// a higher order L2 finite element space, to <tt>y</tt>, dual field
244  /// coefficients defined on a refined mesh with a low order L2 finite
245  /// element space. Refined mesh should be a uniform refinement of the
246  /// coarse mesh. Coefficients are computed from the transpose of the mass
247  /// conservative left-inverse prolongation operation. This functionality
248  /// is also provided as an Operator by L2Prolongation.
249  virtual void ProlongateTranspose(const Vector& x, Vector& y) const;
250  virtual void SetRelTol(double p_rtol_) {}
251  virtual void SetAbsTol(double p_atol_) {}
252  };
253 
254  /** Class for projection operator between a H1 high-order finite element
255  space on a coarse mesh, and a H1 low-order finite element space on a
256  refined mesh (LOR). */
258  {
259  // The restriction operator is represented as a SparseMatrix R. The
260  // prolongation operator is a dense matrix computed as the inverse of (R^T
261  // M_L R), and hence, is not stored.
262  SparseMatrix R;
263  // Used to compute P = (RTxM_LH)^(-1) M_LH^T
264  SparseMatrix M_LH;
265  SparseMatrix* RTxM_LH;
266  CGSolver pcg;
267  DSmoother Ds;
268 
269  public:
271  const FiniteElementSpace& fes_lor_);
272  virtual ~L2ProjectionH1Space();
273  /// Maps <tt>x</tt>, primal field coefficients defined on a coarse mesh
274  /// with a higher order H1 finite element space, to <tt>y</tt>, primal
275  /// field coefficients defined on a refined mesh with a low order H1
276  /// finite element space. Refined mesh should be a uniform refinement of
277  /// the coarse mesh. Coefficients are computed through minimization of L2
278  /// error between the fields.
279  virtual void Mult(const Vector& x, Vector& y) const;
280  /// Maps <tt>x</tt>, dual field coefficients defined on a refined mesh
281  /// with a low order H1 finite element space, to <tt>y</tt>, dual field
282  /// coefficients defined on a coarse mesh with a higher order H1 finite
283  /// element space. Refined mesh should be a uniform refinement of the
284  /// coarse mesh. Coefficients are computed through minimization of L2
285  /// error between the primal fields. Note, if the <tt>x</tt>-coefficients
286  /// come from ProlongateTranspose, then mass is conserved.
287  virtual void MultTranspose(const Vector& x, Vector& y) const;
288  /// Maps <tt>x</tt>, primal field coefficients defined on a refined mesh
289  /// with a low order H1 finite element space, to <tt>y</tt>, primal field
290  /// coefficients defined on a coarse mesh with a higher order H1 finite
291  /// element space. Refined mesh should be a uniform refinement of the
292  /// coarse mesh. Coefficients are computed from the mass conservative
293  /// left-inverse prolongation operation. This functionality is also
294  /// provided as an Operator by L2Prolongation.
295  virtual void Prolongate(const Vector& x, Vector& y) const;
296  /// Maps <tt>x</tt>, dual field coefficients defined on a coarse mesh with
297  /// a higher order H1 finite element space, to <tt>y</tt>, dual field
298  /// coefficients defined on a refined mesh with a low order H1 finite
299  /// element space. Refined mesh should be a uniform refinement of the
300  /// coarse mesh. Coefficients are computed from the transpose of the mass
301  /// conservative left-inverse prolongation operation. This functionality
302  /// is also provided as an Operator by L2Prolongation.
303  virtual void ProlongateTranspose(const Vector& x, Vector& y) const;
304  virtual void SetRelTol(double p_rtol_);
305  virtual void SetAbsTol(double p_atol_);
306  private:
307  /// Computes sparsity pattern and initializes R matrix. Based on
308  /// BilinearForm::AllocMat() except maps between HO elements and LOR
309  /// elements.
310  void AllocR();
311  };
312 
313  /** Mass-conservative prolongation operator going in the opposite direction
314  as L2Projection. This operator is a left inverse to the L2Projection. */
315  class L2Prolongation : public Operator
316  {
317  const L2Projection &l2proj;
318 
319  public:
320  L2Prolongation(const L2Projection &l2proj_)
321  : Operator(l2proj_.Width(), l2proj_.Height()), l2proj(l2proj_) { }
322  void Mult(const Vector &x, Vector &y) const
323  {
324  l2proj.Prolongate(x, y);
325  }
326  void MultTranspose(const Vector &x, Vector &y) const
327  {
328  l2proj.ProlongateTranspose(x, y);
329  }
330  virtual ~L2Prolongation() { }
331  };
332 
333  L2Projection *F; ///< Forward, coarse-to-fine, operator
334  L2Prolongation *B; ///< Backward, fine-to-coarse, operator
336 
337 public:
339  FiniteElementSpace &fine_fes_,
340  bool force_l2_space_ = false)
341  : GridTransfer(coarse_fes_, fine_fes_),
342  F(NULL), B(NULL), force_l2_space(force_l2_space_)
343  { }
344  virtual ~L2ProjectionGridTransfer();
345 
346  virtual const Operator &ForwardOperator();
347 
348  virtual const Operator &BackwardOperator();
349 private:
350  void BuildF();
351 };
352 
353 /// Matrix-free transfer operator between finite element spaces
355 {
356 private:
357  Operator* opr;
358 
359 public:
360  /// Constructs a transfer operator from \p lFESpace to \p hFESpace.
361  /** No matrices are assembled, only the action to a vector is being computed.
362  If both spaces' FE collection pointers are pointing to the same
363  collection we assume that the grid was refined while keeping the order
364  constant. If the FE collections are different, it is assumed that both
365  spaces have are using the same mesh. If the first element of the
366  high-order space is a `TensorBasisElement`, the optimized tensor-product
367  transfers are used. If not, the general transfers used. */
368  TransferOperator(const FiniteElementSpace& lFESpace,
369  const FiniteElementSpace& hFESpace);
370 
371  /// Destructor
372  virtual ~TransferOperator();
373 
374  /// @brief Interpolation or prolongation of a vector \p x corresponding to
375  /// the coarse space to the vector \p y corresponding to the fine space.
376  virtual void Mult(const Vector& x, Vector& y) const override;
377 
378  /// Restriction by applying the transpose of the Mult method.
379  /** The vector \p x corresponding to the fine space is restricted to the
380  vector \p y corresponding to the coarse space. */
381  virtual void MultTranspose(const Vector& x, Vector& y) const override;
382 };
383 
384 /// Matrix-free transfer operator between finite element spaces on the same mesh
386 {
387 private:
388  const FiniteElementSpace& lFESpace;
389  const FiniteElementSpace& hFESpace;
390  bool isvar_order;
391 
392 public:
393  /// @brief Constructs a transfer operator from \p lFESpace to \p hFESpace
394  /// which have different FE collections.
395  /** No matrices are assembled, only the action to a vector is being computed.
396  The underlying finite elements need to implement the GetTransferMatrix
397  methods. */
399  const FiniteElementSpace& hFESpace_);
400 
401  /// Destructor
403 
404  /// @brief Interpolation or prolongation of a vector \p x corresponding to
405  /// the coarse space to the vector \p y corresponding to the fine space.
406  virtual void Mult(const Vector& x, Vector& y) const override;
407 
408  /// Restriction by applying the transpose of the Mult method.
409  /** The vector \p x corresponding to the fine space is restricted to the
410  vector \p y corresponding to the coarse space. */
411  virtual void MultTranspose(const Vector& x, Vector& y) const override;
412 };
413 
414 /// @brief Matrix-free transfer operator between finite element spaces on the
415 /// same mesh exploiting the tensor product structure of the finite elements
417 {
418 private:
419  const FiniteElementSpace& lFESpace;
420  const FiniteElementSpace& hFESpace;
421  int dim;
422  int NE;
423  int D1D;
424  int Q1D;
425  Array<double> B;
426  Array<double> Bt;
427  const Operator* elem_restrict_lex_l;
428  const Operator* elem_restrict_lex_h;
429  Vector mask;
430  mutable Vector localL;
431  mutable Vector localH;
432 
433 public:
434  /// @brief Constructs a transfer operator from \p lFESpace to \p hFESpace
435  /// which have different FE collections.
436  /** No matrices are assembled, only the action to a vector is being computed.
437  The underlying finite elements need to be of type `TensorBasisElement`. It is
438  also assumed that all the elements in the spaces are of the same type. */
440  const FiniteElementSpace& lFESpace_,
441  const FiniteElementSpace& hFESpace_);
442 
443  /// Destructor
445 
446  /// @brief Interpolation or prolongation of a vector \p x corresponding to
447  /// the coarse space to the vector \p y corresponding to the fine space.
448  virtual void Mult(const Vector& x, Vector& y) const override;
449 
450  /// Restriction by applying the transpose of the Mult method.
451  /** The vector \p x corresponding to the fine space is restricted to the
452  vector \p y corresponding to the coarse space. */
453  virtual void MultTranspose(const Vector& x, Vector& y) const override;
454 };
455 
456 /// @brief Matrix-free transfer operator between finite element spaces working
457 /// on true degrees of freedom
459 {
460 private:
461  const FiniteElementSpace& lFESpace;
462  const FiniteElementSpace& hFESpace;
463  const Operator * P = nullptr;
464  const SparseMatrix * R = nullptr;
465  TransferOperator* localTransferOperator;
466  mutable Vector tmpL;
467  mutable Vector tmpH;
468 
469 public:
470  /// @brief Constructs a transfer operator working on true degrees of freedom
471  /// from \p lFESpace to \p hFESpace
472  TrueTransferOperator(const FiniteElementSpace& lFESpace_,
473  const FiniteElementSpace& hFESpace_);
474 
475  /// Destructor
477 
478  /// @brief Interpolation or prolongation of a true dof vector \p x to a true
479  /// dof vector \p y.
480  /** The true dof vector \p x corresponding to the coarse space is restricted
481  to the true dof vector \p y corresponding to the fine space. */
482  virtual void Mult(const Vector& x, Vector& y) const override;
483 
484  /// Restriction by applying the transpose of the Mult method.
485  /** The true dof vector \p x corresponding to the fine space is restricted to
486  the true dof vector \p y corresponding to the coarse space. */
487  virtual void MultTranspose(const Vector& x, Vector& y) const override;
488 };
489 
490 } // namespace mfem
491 
492 #endif
Abstract class for all finite elements.
Definition: fe_base.hpp:235
const Operator & MakeTrueOperator(FiniteElementSpace &fes_in, FiniteElementSpace &fes_out, const Operator &oper, OperatorHandle &t_oper)
Definition: transfer.cpp:34
Conjugate gradient method.
Definition: solvers.hpp:465
virtual void MultTranspose(const Vector &x, Vector &y) const override
Restriction by applying the transpose of the Mult method.
Definition: transfer.cpp:1406
Data type for scaled Jacobi-type smoother of sparse matrix.
virtual void Prolongate(const Vector &x, Vector &y) const
Definition: transfer.cpp:462
BilinearFormIntegrator * mass_integ
Ownership depends on own_mass_integ.
Definition: transfer.hpp:124
L2ProjectionH1Space(const FiniteElementSpace &fes_ho_, const FiniteElementSpace &fes_lor_)
Definition: transfer.cpp:530
virtual ~TensorProductPRefinementTransferOperator()
Destructor.
Definition: transfer.cpp:1376
L2Prolongation(const L2Projection &l2proj_)
Definition: transfer.hpp:320
FiniteElementSpace & dom_fes
Domain FE space.
Definition: transfer.hpp:32
Matrix-free transfer operator between finite element spaces on the same mesh.
Definition: transfer.hpp:385
int Width() const
Get the width (size of input) of the Operator. Synonym with NumCols().
Definition: operator.hpp:73
virtual void Mult(const Vector &x, Vector &y) const
Definition: transfer.cpp:677
Pointer to an Operator of a specified type.
Definition: handle.hpp:33
void Mult(const Vector &x, Vector &y) const
Operator application: y=A(x).
Definition: transfer.hpp:322
virtual void MultTranspose(const Vector &x, Vector &y) const override
Restriction by applying the transpose of the Mult method.
Definition: transfer.cpp:1486
virtual const Operator & BackwardOperator()=0
Return an Operator that transfers GridFunctions from the range FE space back to GridFunctions in the ...
Data type dense matrix using column-major storage.
Definition: densemat.hpp:23
virtual ~TransferOperator()
Destructor.
Definition: transfer.cpp:934
Operator::Type oper_type
Desired Operator::Type for the construction of all operators defined by the underlying transfer algor...
Definition: transfer.hpp:38
virtual void ProlongateTranspose(const Vector &x, Vector &y) const =0
~TrueTransferOperator()
Destructor.
Definition: transfer.cpp:1462
void ElemMixedMass(Geometry::Type geom, const FiniteElement &fe_ho, const FiniteElement &fe_lor, ElementTransformation *el_tr, IntegrationPointTransformation &ip_tr, DenseMatrix &M_mixed_el) const
Definition: transfer.cpp:258
virtual void Mult(const Vector &x, Vector &y) const override
Interpolation or prolongation of a vector x corresponding to the coarse space to the vector y corresp...
Definition: transfer.cpp:936
virtual const Operator & ForwardOperator()
Return an Operator that transfers GridFunctions from the domain FE space to GridFunctions in the rang...
Definition: transfer.cpp:872
void BuildHo2Lor(int nel_ho, int nel_lor, const CoarseFineTransformations &cf_tr)
Definition: transfer.cpp:238
virtual ~PRefinementTransferOperator()
Destructor.
Definition: transfer.cpp:955
virtual void Mult(const Vector &x, Vector &y) const
Definition: transfer.cpp:395
virtual void ProlongateTranspose(const Vector &x, Vector &y) const
Definition: transfer.cpp:497
virtual const Operator & TrueBackwardOperator()
Return an Operator that transfers true-dof Vectors from the range FE space back to true-dof Vectors i...
Definition: transfer.hpp:97
virtual ~GridTransfer()
Virtual destructor.
Definition: transfer.hpp:66
FiniteElementSpace & ran_fes
Range FE space.
Definition: transfer.hpp:33
virtual const Operator & BackwardOperator()
Return an Operator that transfers GridFunctions from the range FE space back to GridFunctions in the ...
Definition: transfer.cpp:189
virtual void MultTranspose(const Vector &x, Vector &y) const override
Restriction by applying the transpose of the Mult method.
Definition: transfer.cpp:941
Data type sparse matrix.
Definition: sparsemat.hpp:50
virtual void SetRelTol(double p_rtol_)=0
int Height() const
Get the height (size of output) of the Operator. Synonym with NumRows().
Definition: operator.hpp:67
OperatorHandle F
Forward, coarse-to-fine, operator.
Definition: transfer.hpp:127
virtual const Operator & ForwardOperator()=0
Return an Operator that transfers GridFunctions from the domain FE space to GridFunctions in the rang...
virtual void MultTranspose(const Vector &x, Vector &y) const
Definition: transfer.cpp:698
GridTransfer(FiniteElementSpace &dom_fes_, FiniteElementSpace &ran_fes_)
Definition: transfer.cpp:19
virtual void Prolongate(const Vector &x, Vector &y) const =0
TensorProductPRefinementTransferOperator(const FiniteElementSpace &lFESpace_, const FiniteElementSpace &hFESpace_)
Constructs a transfer operator from lFESpace to hFESpace which have different FE collections.
Definition: transfer.cpp:1082
L2Projection(const FiniteElementSpace &fes_ho_, const FiniteElementSpace &fes_lor_)
Definition: transfer.cpp:231
InterpolationGridTransfer(FiniteElementSpace &coarse_fes, FiniteElementSpace &fine_fes)
Definition: transfer.hpp:131
OperatorHandle bw_t_oper
Backward true-dof operator.
Definition: transfer.hpp:41
void SetOperatorType(Operator::Type type)
Set the desired Operator::Type for the construction of all operators defined by the underlying transf...
Definition: transfer.hpp:73
virtual void Prolongate(const Vector &x, Vector &y) const
Definition: transfer.cpp:719
OperatorHandle B
Backward, fine-to-coarse, operator.
Definition: transfer.hpp:128
virtual const Operator & BackwardOperator()
Return an Operator that transfers GridFunctions from the range FE space back to GridFunctions in the ...
Definition: transfer.cpp:878
Base class for transfer algorithms that construct transfer Operators between two finite element (FE) ...
Definition: transfer.hpp:29
OperatorHandle fw_t_oper
Forward true-dof operator.
Definition: transfer.hpp:40
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: transfer.cpp:146
Type
Enumeration defining IDs for some classes derived from Operator.
Definition: operator.hpp:256
Abstract base class BilinearFormIntegrator.
Definition: bilininteg.hpp:35
virtual void MultTranspose(const Vector &x, Vector &y) const override
Restriction by applying the transpose of the Mult method.
Definition: transfer.cpp:1008
Matrix-free transfer operator between finite element spaces.
Definition: transfer.hpp:354
Class FiniteElementSpace - responsible for providing FEM view of the mesh, mainly managing the set of...
Definition: fespace.hpp:96
Transfer data between a coarse mesh and an embedded refined mesh using interpolation.
Definition: transfer.hpp:121
virtual void MultTranspose(const Vector &x, Vector &y) const
Definition: transfer.cpp:427
Matrix-free transfer operator between finite element spaces on the same mesh exploiting the tensor pr...
Definition: transfer.hpp:416
virtual void Mult(const Vector &x, Vector &y) const override
Interpolation or prolongation of a vector x corresponding to the coarse space to the vector y corresp...
Definition: transfer.cpp:1380
L2Prolongation * B
Backward, fine-to-coarse, operator.
Definition: transfer.hpp:334
bool own_mass_integ
Ownership flag for mass_integ.
Definition: transfer.hpp:125
L2Projection * F
Forward, coarse-to-fine, operator.
Definition: transfer.hpp:333
virtual void SetAbsTol(double p_atol_)=0
virtual const Operator & ForwardOperator()
Return an Operator that transfers GridFunctions from the domain FE space to GridFunctions in the rang...
Definition: transfer.cpp:155
bool Parallel() const
Definition: transfer.hpp:46
virtual const Operator & TrueForwardOperator()
Return an Operator that transfers true-dof Vectors from the domain FE space to true-dof Vectors in th...
Definition: transfer.hpp:87
virtual void Mult(const Vector &x, Vector &y) const override
Interpolation or prolongation of a vector x corresponding to the coarse space to the vector y corresp...
Definition: transfer.cpp:957
virtual void ProlongateTranspose(const Vector &x, Vector &y) const
Definition: transfer.cpp:744
L2ProjectionGridTransfer(FiniteElementSpace &coarse_fes_, FiniteElementSpace &fine_fes_, bool force_l2_space_=false)
Definition: transfer.hpp:338
const FiniteElementSpace & fes_lor
Definition: transfer.hpp:187
PRefinementTransferOperator(const FiniteElementSpace &lFESpace_, const FiniteElementSpace &hFESpace_)
Constructs a transfer operator from lFESpace to hFESpace which have different FE collections.
Definition: transfer.cpp:947
TrueTransferOperator(const FiniteElementSpace &lFESpace_, const FiniteElementSpace &hFESpace_)
Constructs a transfer operator working on true degrees of freedom from lFESpace to hFESpace...
Definition: transfer.cpp:1433
virtual void Mult(const Vector &x, Vector &y) const override
Interpolation or prolongation of a true dof vector x to a true dof vector y.
Definition: transfer.cpp:1467
Transfer data in L2 and H1 finite element spaces between a coarse mesh and an embedded refined mesh u...
Definition: transfer.hpp:168
Vector data type.
Definition: vector.hpp:60
Defines the coarse-fine transformations of all fine elements.
Definition: ncmesh.hpp:70
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: transfer.hpp:326
Abstract operator.
Definition: operator.hpp:24
L2ProjectionL2Space(const FiniteElementSpace &fes_ho_, const FiniteElementSpace &fes_lor_)
Definition: transfer.cpp:285
TransferOperator(const FiniteElementSpace &lFESpace, const FiniteElementSpace &hFESpace)
Constructs a transfer operator from lFESpace to hFESpace.
Definition: transfer.cpp:902
Matrix-free transfer operator between finite element spaces working on true degrees of freedom...
Definition: transfer.hpp:458