MFEM v4.8.0
Finite element discretization library
Loading...
Searching...
No Matches
transfer.hpp
Go to the documentation of this file.
1// Copyright (c) 2010-2025, 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
22namespace 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{
31protected:
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 bool use_ea;
44
46
47#ifdef MFEM_USE_MPI
49#endif
50 bool Parallel() const
51 {
52#ifndef MFEM_USE_MPI
53 return false;
54#else
55 return parallel;
56#endif
57 }
58
60 FiniteElementSpace &fes_out,
61 const Operator &oper,
62 OperatorHandle &t_oper);
63
64public:
65 /** Construct a transfer algorithm between the domain, @a dom_fes_, and
66 range, @a ran_fes_, FE spaces, d_mt_ will specify memory space for
67 large data structures */
69 FiniteElementSpace &ran_fes_);
70
71 /// Virtual destructor
72 virtual ~GridTransfer() { }
73
74 /** Uses device friendly element assembly versions for L2Projection
75 transfers, L2, H1 FEM spaces currently supported */
76 void UseEA(bool use_ea_) { use_ea = use_ea_;}
77
78 /** Set memory type for large data structures */
79 void SetMemType(MemoryType d_mt_) {d_mt = d_mt_;}
80
81 /** @brief Set the desired Operator::Type for the construction of all
82 operators defined by the underlying transfer algorithm. */
83 /** The default value is Operator::ANY_TYPE which typically corresponds to a
84 matrix-free operator representation. Note that derived classes are not
85 required to support this setting and can ignore it. */
87
88 /** @brief Return an Operator that transfers GridFunction%s from the domain
89 FE space to GridFunction%s in the range FE space. */
90 virtual const Operator &ForwardOperator() = 0;
91
92 /** @brief Return an Operator that transfers GridFunction%s from the range FE
93 space back to GridFunction%s in the domain FE space. */
94 virtual const Operator &BackwardOperator() = 0;
95
96 /** @brief Return an Operator that transfers true-dof Vector%s from the
97 domain FE space to true-dof Vector%s in the range FE space. */
98 /** This method is implemented in the base class, based on ForwardOperator(),
99 however, derived classes can overload the construction, if necessary. */
104
105 /** @brief Return an Operator that transfers true-dof Vector%s from the range
106 FE space back to true-dof Vector%s in the domain FE space. */
107 /** This method is implemented in the base class, based on
108 BackwardOperator(), however, derived classes can overload the
109 construction, if necessary. */
114
115 virtual bool SupportsBackwardsOperator() const { return true; }
116};
117
118
119/** @brief Transfer data between a coarse mesh and an embedded refined mesh
120 using interpolation. */
121/** The forward, coarse-to-fine, transfer uses nodal interpolation. The
122 backward, fine-to-coarse, transfer is defined locally (on a coarse element)
123 as B = (F^t M_f F)^{-1} F^t M_f, where F is the forward transfer matrix, and
124 M_f is a mass matrix on the union of all fine elements comprising the coarse
125 element. Note that the backward transfer operator, B, is a left inverse of
126 the forward transfer operator, F, i.e. B F = I. Both F and B are defined in
127 reference space and do not depend on the actual physical shape of the mesh
128 elements.
129
130 It is assumed that both the coarse and the fine FiniteElementSpace%s use
131 compatible types of elements, e.g. finite elements with the same map-type
132 (VALUE, INTEGRAL, H_DIV, H_CURL - see class FiniteElement). Generally, the
133 FE spaces can have different orders, however, in order for the backward
134 operator to be well-defined, the (local) number of the fine dofs should not
135 be smaller than the number of coarse dofs. */
137{
138protected:
139 BilinearFormIntegrator *mass_integ; ///< Ownership depends on #own_mass_integ
140 bool own_mass_integ; ///< Ownership flag for #mass_integ
141
142 OperatorHandle F; ///< Forward, coarse-to-fine, operator
143 OperatorHandle B; ///< Backward, fine-to-coarse, operator
144
145public:
147 FiniteElementSpace &fine_fes)
148 : GridTransfer(coarse_fes, fine_fes),
149 mass_integ(NULL), own_mass_integ(false)
150 { }
151
153
154 /** @brief Assign a mass integrator to be used in the construction of the
155 backward, fine-to-coarse, transfer operator. */
157 bool own_mass_integ_ = true);
158
159 const Operator &ForwardOperator() override;
160
161 const Operator &BackwardOperator() override;
162};
163
164
165/** @brief Transfer data in L2 and H1 finite element spaces between a coarse
166 mesh and an embedded refined mesh using L2 projection. */
167/** The forward, coarse-to-fine, transfer uses L2 projection. The backward,
168 fine-to-coarse, transfer is defined as B = (F^t M_f F)^{-1} F^t M_f, where F
169 is the forward transfer matrix, and M_f is the mass matrix on the coarse
170 element. For L2 spaces, M_f is the mass matrix on the union of all fine
171 elements comprising the coarse element. For H1 spaces, M_f is a diagonal
172 (lumped) mass matrix computed through row-summation. Note that the backward
173 transfer operator, B, is a left inverse of the forward transfer operator, F,
174 i.e. B F = I. Both F and B are defined in physical space and, generally for
175 L2 spaces, vary between different mesh elements.
176
177 This class supports H1 and L2 finite element spaces. Fine meshes are a
178 uniform refinement of the coarse mesh, usually created through
179 Mesh::MakeRefined. Generally, the coarse and fine FE spaces can have
180 different orders, however, in order for the backward operator to be
181 well-defined, the number of fine dofs (in a coarse element) should not be
182 smaller than the number of coarse dofs. */
184{
185 // Must be public due to host device lambdas
186public:
187 /** Abstract class representing projection operator between a high-order
188 finite element space on a coarse mesh, and a low-order finite element
189 space on a refined mesh (LOR). We assume that the low-order space,
190 fes_lor, lives on a mesh obtained by refining the mesh of the high-order
191 space, fes_ho. */
192 class L2Projection : public Operator
193 {
194 public:
195 virtual void Prolongate(const Vector& x, Vector& y) const = 0;
196 virtual void ProlongateTranspose(const Vector& x, Vector& y) const = 0;
197 /// @brief Sets relative tolerance in preconditioned conjugate gradient
198 /// solver.
199 ///
200 /// Only used for H1 spaces.
201 virtual void SetRelTol(real_t p_rtol_) = 0;
202 /// @brief Sets absolute tolerance in preconditioned conjugate gradient
203 /// solver.
204 ///
205 /// Only used for H1 spaces.
206 virtual void SetAbsTol(real_t p_atol_) = 0;
207 protected:
210
214
215 L2Projection(const FiniteElementSpace& fes_ho_,
216 const FiniteElementSpace& fes_lor_,
218
219 void BuildHo2Lor(int nel_ho, int nel_lor,
220 const CoarseFineTransformations& cf_tr);
221
222 void ElemMixedMass(Geometry::Type geom, const FiniteElement& fe_ho,
223 const FiniteElement& fe_lor, ElementTransformation* tr_ho,
224 ElementTransformation* tr_lor,
226 DenseMatrix& M_mixed_el) const;
227
228 void ElemMixedMass(Geometry::Type geom, const FiniteElement& fe_ho,
229 const FiniteElement& fe_lor,
232 DenseMatrix& B_L, DenseMatrix& B_H) const;
233 public:
234 /* Returns the Mixed Mass M_LH via device element assembly by building the
235 basis functions and data at the quadrature points. */
236 void MixedMassEA(const FiniteElementSpace& fes_ho_,
237 const FiniteElementSpace& fes_lor_,
238 Vector &M_LH,
240 };
241
242 // Class below must be public as we now have device code
243public:
245 {
246 protected:
251 public:
253 const FiniteElementSpace* fes_lor_,
254 Table* ho2lor_, Vector* M_LH_ea_);
255 void Mult(const Vector& x, Vector& y) const;
256 void MultTranspose(const Vector& x, Vector& y) const;
257 };
258
260 {
261 protected:
264 Vector* ML_inv; // inverse of lumped M_L
265 public:
267 const FiniteElementSpace* fes_lor_,
268 Vector& ML_inv_);
269 void Mult(const Vector& x, Vector& y) const;
270 void MultTranspose(const Vector& x, Vector& y) const;
271 };
272
273 /** Class for projection operator between a L2 high-order finite element
274 space on a coarse mesh, and a L2 low-order finite element space on a
275 refined mesh (LOR). */
277 {
278 /// The restriction and prolongation operators are represented as dense
279 /// elementwise matrices (of potentially different sizes, because of mixed
280 /// meshes or p-refinement). The matrix entries are stored in the R and P
281 /// arrays. The entries of the i'th high-order element are stored at the
282 /// index given by offsets[i].
283 mutable Array<real_t> R, P;
284
285 const bool use_ea;
286
287 public:
289 const FiniteElementSpace& fes_lor_,
290 const bool use_ea_,
292
293 /*Same as above but assembles and stores R_ea, P_ea */
295
296 /// Maps <tt>x</tt>, primal field coefficients defined on a coarse mesh
297 /// with a higher order L2 finite element space, to <tt>y</tt>, primal
298 /// field coefficients defined on a refined mesh with a low order L2
299 /// finite element space. Refined mesh should be a uniform refinement of
300 /// the coarse mesh. Coefficients are computed through minimization of L2
301 /// error between the fields.
302 void Mult(const Vector& x, Vector& y) const override;
303
304 /// Perform mult on the device (same as above)
305 void EAMult(const Vector& x, Vector& y) const;
306
307 /// Maps <tt>x</tt>, dual field coefficients defined on a refined mesh
308 /// with a low order L2 finite element space, to <tt>y</tt>, dual field
309 /// coefficients defined on a coarse mesh with a higher order L2 finite
310 /// element space. Refined mesh should be a uniform refinement of the
311 /// coarse mesh. Coefficients are computed through minimization of L2
312 /// error between the primal fields. Note, if the <tt>x</tt>-coefficients
313 /// come from ProlongateTranspose, then mass is conserved.
314 void MultTranspose(const Vector& x, Vector& y) const override;
315
316 void EAMultTranspose(const Vector& x, Vector& y) const;
317
318 /// Maps <tt>x</tt>, primal field coefficients defined on a refined mesh
319 /// with a low order L2 finite element space, to <tt>y</tt>, primal field
320 /// coefficients defined on a coarse mesh with a higher order L2 finite
321 /// element space. Refined mesh should be a uniform refinement of the
322 /// coarse mesh. Coefficients are computed from the mass conservative
323 /// left-inverse prolongation operation. This functionality is also
324 /// provided as an Operator by L2Prolongation.
325 void Prolongate(const Vector& x, Vector& y) const override;
326
327 void EAProlongate(const Vector& x, Vector& y) const;
328
329 /// Maps <tt>x</tt>, dual field coefficients defined on a coarse mesh with
330 /// a higher order L2 finite element space, to <tt>y</tt>, dual field
331 /// coefficients defined on a refined mesh with a low order L2 finite
332 /// element space. Refined mesh should be a uniform refinement of the
333 /// coarse mesh. Coefficients are computed from the transpose of the mass
334 /// conservative left-inverse prolongation operation. This functionality
335 /// is also provided as an Operator by L2Prolongation.
336 void ProlongateTranspose(const Vector& x, Vector& y) const override;
337
338 void EAProlongateTranspose(const Vector& x, Vector& y) const;
339
340 void SetRelTol(real_t p_rtol_) override { } ///< No-op.
341 void SetAbsTol(real_t p_atol_) override { } ///< No-op.
342 };
343
344protected:
345
346 /// Class below must be public as we now have device code
347public:
348
349 /** Projection operator between a H1 high-order finite element space on a
350 coarse mesh, and a H1 low-order finite element space on a refined mesh
351 (LOR). */
353 {
354 const bool use_ea;
355
356 public:
358 const FiniteElementSpace &fes_lor_,
359 const bool use_ea_,
361#ifdef MFEM_USE_MPI
363 const ParFiniteElementSpace &pfes_lor_,
364 const bool use_ea_,
366#endif
367 /// Same as above but assembles action of R through 4 parts:
368 /// ( ) inv( lumped(M_L) ), which is a diagonal matrix (essentially a vector)
369 /// ( ) ElementRestrictionOperator for LOR space
370 /// ( ) mixed mass matrix M_{LH}
371 /// ( ) ElementRestrictionOperator for HO space
373
374#ifdef MFEM_USE_MPI
375 void EAL2ProjectionH1Space(const ParFiniteElementSpace &pfes_ho_,
376 const ParFiniteElementSpace &pfes_lor_);
377#endif
378 /// Maps <tt>x</tt>, primal field coefficients defined on a coarse mesh
379 /// with a higher order H1 finite element space, to <tt>y</tt>, primal
380 /// field coefficients defined on a refined mesh with a low order H1
381 /// finite element space. Refined mesh should be a uniform refinement of
382 /// the coarse mesh. Coefficients are computed through minimization of L2
383 /// error between the fields.
384 void Mult(const Vector& x, Vector& y) const override;
385
386 /// Maps <tt>x</tt>, dual field coefficients defined on a refined mesh
387 /// with a low order H1 finite element space, to <tt>y</tt>, dual field
388 /// coefficients defined on a coarse mesh with a higher order H1 finite
389 /// element space. Refined mesh should be a uniform refinement of the
390 /// coarse mesh. Coefficients are computed through minimization of L2
391 /// error between the primal fields. Note, if the <tt>x</tt>-coefficients
392 /// come from ProlongateTranspose, then mass is conserved.
393 void MultTranspose(const Vector& x, Vector& y) const override;
394
395 /// Maps <tt>x</tt>, primal field coefficients defined on a refined mesh
396 /// with a low order H1 finite element space, to <tt>y</tt>, primal field
397 /// coefficients defined on a coarse mesh with a higher order H1 finite
398 /// element space. Refined mesh should be a uniform refinement of the
399 /// coarse mesh. Coefficients are computed from the mass conservative
400 /// left-inverse prolongation operation. This functionality is also
401 /// provided as an Operator by L2Prolongation.
402 void Prolongate(const Vector& x, Vector& y) const override;
403
404 /// Maps <tt>x</tt>, dual field coefficients defined on a coarse mesh with
405 /// a higher order H1 finite element space, to <tt>y</tt>, dual field
406 /// coefficients defined on a refined mesh with a low order H1 finite
407 /// element space. Refined mesh should be a uniform refinement of the
408 /// coarse mesh. Coefficients are computed from the transpose of the mass
409 /// conservative left-inverse prolongation operation. This functionality
410 /// is also provided as an Operator by L2Prolongation.
411 void ProlongateTranspose(const Vector& x, Vector& y) const override;
412
413 /// Returns the inverse of an on-rank lumped mass matrix
414 void LumpedMassInverse(Vector& ML_inv) const;
415
416 void SetRelTol(real_t p_rtol_) override;
417 void SetAbsTol(real_t p_atol_) override;
418
419 protected:
420 /// Sets up the PCG solver (sets parameters, operator, and preconditioner)
421 void SetupPCG();
422
423 /// @brief Computes on-rank R and M_LH matrices. If true, computes mixed mass and/or
424 /// inverse lumped mass matrix error when compared to device implementation.
425 std::pair<std::unique_ptr<SparseMatrix>,
426 std::unique_ptr<SparseMatrix>> ComputeSparseRAndM_LH();
427
428 /// @brief Recovers vector of tdofs given a vector of dofs and a finite
429 /// element space
430 void GetTDofs(const FiniteElementSpace& fes, const Vector& x, Vector& X) const;
431 /// Sets dof values given a vector of tdofs and a finite element space
432 void SetFromTDofs(const FiniteElementSpace& fes,
433 const Vector& X,
434 Vector& x) const;
435 /// @brief Recovers a vector of dual field coefficients on the tdofs given
436 /// a vector of dual coefficients and a finite element space
438 const Vector& x,
439 Vector& X) const;
440 /// @brief Sets dual field coefficients given a vector of dual field
441 /// coefficients on the tdofs and a finite element space
443 const Vector& X,
444 Vector& x) const;
445 /// @brief Fills the vdofs_list array with a list of vdofs for a given
446 /// vdim and a given finite element space
447 void TDofsListByVDim(const FiniteElementSpace& fes,
448 int vdim,
449 Array<int>& vdofs_list) const;
450
451 /// @brief Computes sparsity pattern and initializes R matrix.
452 /// Based on BilinearForm::AllocMat(), except maps between coarse HO
453 /// elements and refined LOR elements.
454 std::unique_ptr<SparseMatrix> AllocR();
455
457 std::unique_ptr<Solver> precon;
458 // The restriction operator is represented as an Operator R. The
459 // prolongation operator is a dense matrix computed as the inverse of (R^T
460 // M_L R), and hence, is not stored.
461 // If element assembly is enabled
462 std::unique_ptr<Operator> R;
463 // Used to compute P = (RT*M_LH)^(-1) M_LH^T
464 std::unique_ptr<Operator> M_LH;
465 // Inverted operator in P = (RT*M_LH)^(-1) M_LH^T. Used to compute P via PCG.
466 std::unique_ptr<Operator> RTxM_LH;
467 // Lumped M_L inverse operator built via EA. Wrapped with restriction maps
468 // to multiply with scalar TDof LOR vectors.
469 std::unique_ptr<Operator> ML_inv_vea;
470 // LDof Mixed mass operator built via EA. Wrapped with restriction maps to send
471 // scalar LDof HO vectors to LDof LOR vectors.
473
474 // Scalar finite element spaces for stored Tdof-to-and-from-LDof maps.
475 std::unique_ptr<FiniteElementSpace> fes_ho_scalar;
476 std::unique_ptr<FiniteElementSpace> fes_lor_scalar;
477 // Element Assembled mixed mass
479 // Element Assembled lumped M_L inverse built via EA. Stores diagonal as a Ldof vector.
481
482#ifdef MFEM_USE_MPI
483 std::unique_ptr<ParFiniteElementSpace> pfes_ho_scalar;
484 std::unique_ptr<ParFiniteElementSpace> pfes_lor_scalar;
486#endif
487
489 };
490
491 /** Mass-conservative prolongation operator going in the opposite direction
492 as L2Projection. This operator is a left inverse to the L2Projection. */
494 {
495 const L2Projection &l2proj;
496
497 public:
499 : Operator(l2proj_.Width(), l2proj_.Height()), l2proj(l2proj_) { }
500 void Mult(const Vector &x, Vector &y) const
501 {
502 l2proj.Prolongate(x, y);
503 }
504 void MultTranspose(const Vector &x, Vector &y) const
505 {
506 l2proj.ProlongateTranspose(x, y);
507 }
508 virtual ~L2Prolongation() { }
509 };
510
511 L2Projection *F; ///< Forward, coarse-to-fine, operator
512 L2Prolongation *B; ///< Backward, fine-to-coarse, operator
514
515public:
517 FiniteElementSpace &fine_fes_,
518 bool force_l2_space_ = false,
519 MemoryType d_mt_ = Device::GetHostMemoryType()) // move to method
520 : GridTransfer(coarse_fes_, fine_fes_),
521 F(NULL), B(NULL), force_l2_space(force_l2_space_)
522 { }
524
525 const Operator &ForwardOperator() override;
526
527 const Operator &BackwardOperator() override;
528
529 bool SupportsBackwardsOperator() const override;
530private:
531 void BuildF();
532};
533
534/// Matrix-free transfer operator between finite element spaces
536{
537private:
538 Operator* opr;
539
540public:
541 /// Constructs a transfer operator from \p lFESpace to \p hFESpace.
542 /** No matrices are assembled, only the action to a vector is being computed.
543 If both spaces' FE collection pointers are pointing to the same
544 collection, we assume that the grid was refined while keeping the order
545 constant. If the FE collections are different, it is assumed that both
546 spaces are using the same mesh. If the first element of the high-order
547 space is a `TensorBasisElement`, the optimized tensor-product transfers
548 are used. If not, the general transfers used. */
549 TransferOperator(const FiniteElementSpace& lFESpace,
550 const FiniteElementSpace& hFESpace);
551
552 /// Destructor
553 virtual ~TransferOperator();
554
555 /// @brief Interpolation or prolongation of a vector \p x corresponding to
556 /// the coarse space to the vector \p y corresponding to the fine space.
557 void Mult(const Vector& x, Vector& y) const override;
558
559 /// Restriction by applying the transpose of the Mult method.
560 /** The vector \p x corresponding to the fine space is restricted to the
561 vector \p y corresponding to the coarse space. */
562 void MultTranspose(const Vector& x, Vector& y) const override;
563};
564
565/// Matrix-free transfer operator between finite element spaces on the same mesh
567{
568private:
569 const FiniteElementSpace& lFESpace;
570 const FiniteElementSpace& hFESpace;
571 bool isvar_order;
572
573public:
574 /// @brief Constructs a transfer operator from \p lFESpace to \p hFESpace
575 /// which have different FE collections.
576 /** No matrices are assembled, only the action to a vector is being computed.
577 The underlying finite elements need to implement the GetTransferMatrix
578 methods. */
580 const FiniteElementSpace& hFESpace_);
581
582 /// Destructor
584
585 /// @brief Interpolation or prolongation of a vector \p x corresponding to
586 /// the coarse space to the vector \p y corresponding to the fine space.
587 void Mult(const Vector& x, Vector& y) const override;
588
589 /// Restriction by applying the transpose of the Mult method.
590 /** The vector \p x corresponding to the fine space is restricted to the
591 vector \p y corresponding to the coarse space. */
592 void MultTranspose(const Vector& x, Vector& y) const override;
593};
594
595/// @brief Matrix-free transfer operator between finite element spaces on the
596/// same mesh exploiting the tensor product structure of the finite elements
598{
599private:
600 const FiniteElementSpace& lFESpace;
601 const FiniteElementSpace& hFESpace;
602 int dim;
603 int NE;
604 int D1D;
605 int Q1D;
607 Array<real_t> Bt;
608 const Operator* elem_restrict_lex_l;
609 const Operator* elem_restrict_lex_h;
610 Vector mask;
611 mutable Vector localL;
612 mutable Vector localH;
613
614public:
615 /// @brief Constructs a transfer operator from \p lFESpace to \p hFESpace
616 /// which have different FE collections.
617 /** No matrices are assembled, only the action to a vector is being computed.
618 The underlying finite elements need to be of type `TensorBasisElement`. It is
619 also assumed that all the elements in the spaces are of the same type. */
621 const FiniteElementSpace& lFESpace_,
622 const FiniteElementSpace& hFESpace_);
623
624 /// Destructor
626
627 /// @brief Interpolation or prolongation of a vector \p x corresponding to
628 /// the coarse space to the vector \p y corresponding to the fine space.
629 void Mult(const Vector& x, Vector& y) const override;
630
631 /// Restriction by applying the transpose of the Mult method.
632 /** The vector \p x corresponding to the fine space is restricted to the
633 vector \p y corresponding to the coarse space. */
634 void MultTranspose(const Vector& x, Vector& y) const override;
635};
636
637/// @brief Matrix-free transfer operator between finite element spaces working
638/// on true degrees of freedom
640{
641private:
642 const FiniteElementSpace& lFESpace;
643 const FiniteElementSpace& hFESpace;
644 const Operator * P = nullptr;
645 const SparseMatrix * R = nullptr;
646 TransferOperator* localTransferOperator;
647 mutable Vector tmpL;
648 mutable Vector tmpH;
649
650public:
651 /// @brief Constructs a transfer operator working on true degrees of freedom
652 /// from \p lFESpace to \p hFESpace
654 const FiniteElementSpace& hFESpace_);
655
656 /// Destructor
658
659 /// @brief Interpolation or prolongation of a true dof vector \p x to a true
660 /// dof vector \p y.
661 /** The true dof vector \p x corresponding to the coarse space is restricted
662 to the true dof vector \p y corresponding to the fine space. */
663 void Mult(const Vector& x, Vector& y) const override;
664
665 /// Restriction by applying the transpose of the Mult method.
666 /** The true dof vector \p x corresponding to the fine space is restricted to
667 the true dof vector \p y corresponding to the coarse space. */
668 void MultTranspose(const Vector& x, Vector& y) const override;
669};
670
671} // namespace mfem
672
673#endif
Abstract base class BilinearFormIntegrator.
Conjugate gradient method.
Definition solvers.hpp:538
Data type dense matrix using column-major storage.
Definition densemat.hpp:24
static MemoryType GetHostMemoryType()
Get the current Host MemoryType. This is the MemoryType used by most MFEM classes when allocating mem...
Definition device.hpp:265
Class FiniteElementSpace - responsible for providing FEM view of the mesh, mainly managing the set of...
Definition fespace.hpp:244
Abstract class for all finite elements.
Definition fe_base.hpp:244
Base class for transfer algorithms that construct transfer Operators between two finite element (FE) ...
Definition transfer.hpp:30
void UseEA(bool use_ea_)
Definition transfer.hpp:76
bool Parallel() const
Definition transfer.hpp:50
virtual ~GridTransfer()
Virtual destructor.
Definition transfer.hpp:72
FiniteElementSpace & dom_fes
Domain FE space.
Definition transfer.hpp:32
OperatorHandle fw_t_oper
Forward true-dof operator.
Definition transfer.hpp:40
FiniteElementSpace & ran_fes
Range FE space.
Definition transfer.hpp:33
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:110
OperatorHandle bw_t_oper
Backward true-dof operator.
Definition transfer.hpp:41
virtual bool SupportsBackwardsOperator() const
Definition transfer.hpp:115
void SetMemType(MemoryType d_mt_)
Definition transfer.hpp:79
MemoryType d_mt
Definition transfer.hpp:45
virtual const Operator & ForwardOperator()=0
Return an Operator that transfers GridFunctions from the domain FE space to GridFunctions in the rang...
void SetOperatorType(Operator::Type type)
Set the desired Operator::Type for the construction of all operators defined by the underlying transf...
Definition transfer.hpp:86
virtual const Operator & BackwardOperator()=0
Return an Operator that transfers GridFunctions from the range FE space back to GridFunctions in the ...
Operator::Type oper_type
Desired Operator::Type for the construction of all operators defined by the underlying transfer algor...
Definition transfer.hpp:38
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:100
const Operator & MakeTrueOperator(FiniteElementSpace &fes_in, FiniteElementSpace &fes_out, const Operator &oper, OperatorHandle &t_oper)
Definition transfer.cpp:35
GridTransfer(FiniteElementSpace &dom_fes_, FiniteElementSpace &ran_fes_)
Definition transfer.cpp:20
Transfer data between a coarse mesh and an embedded refined mesh using interpolation.
Definition transfer.hpp:137
InterpolationGridTransfer(FiniteElementSpace &coarse_fes, FiniteElementSpace &fine_fes)
Definition transfer.hpp:146
BilinearFormIntegrator * mass_integ
Ownership depends on own_mass_integ.
Definition transfer.hpp:139
OperatorHandle F
Forward, coarse-to-fine, operator.
Definition transfer.hpp:142
bool own_mass_integ
Ownership flag for mass_integ.
Definition transfer.hpp:140
const Operator & BackwardOperator() override
Return an Operator that transfers GridFunctions from the range FE space back to GridFunctions in the ...
Definition transfer.cpp:190
OperatorHandle B
Backward, fine-to-coarse, operator.
Definition transfer.hpp:143
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,...
Definition transfer.cpp:147
const Operator & ForwardOperator() override
Return an Operator that transfers GridFunctions from the domain FE space to GridFunctions in the rang...
Definition transfer.cpp:156
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 ...
H1SpaceLumpedMassOperator(const FiniteElementSpace *fes_ho_, const FiniteElementSpace *fes_lor_, Vector &ML_inv_)
void Mult(const Vector &x, Vector &y) const
Operator application: y=A(x).
void Mult(const Vector &x, Vector &y) const
Operator application: y=A(x).
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 ...
H1SpaceMixedMassOperator(const FiniteElementSpace *fes_ho_, const FiniteElementSpace *fes_lor_, Table *ho2lor_, Vector *M_LH_ea_)
Class below must be public as we now have device code.
Definition transfer.hpp:353
void Mult(const Vector &x, Vector &y) const override
void SetFromTDofsTranspose(const FiniteElementSpace &fes, const Vector &X, Vector &x) const
Sets dual field coefficients given a vector of dual field coefficients on the tdofs and a finite elem...
std::unique_ptr< FiniteElementSpace > fes_ho_scalar
Definition transfer.hpp:475
std::unique_ptr< ParFiniteElementSpace > pfes_ho_scalar
Definition transfer.hpp:483
void TDofsListByVDim(const FiniteElementSpace &fes, int vdim, Array< int > &vdofs_list) const
Fills the vdofs_list array with a list of vdofs for a given vdim and a given finite element space.
void SetupPCG()
Sets up the PCG solver (sets parameters, operator, and preconditioner)
void Prolongate(const Vector &x, Vector &y) const override
void GetTDofs(const FiniteElementSpace &fes, const Vector &x, Vector &X) const
Recovers vector of tdofs given a vector of dofs and a finite element space.
void SetFromTDofs(const FiniteElementSpace &fes, const Vector &X, Vector &x) const
Sets dof values given a vector of tdofs and a finite element space.
void GetTDofsTranspose(const FiniteElementSpace &fes, const Vector &x, Vector &X) const
Recovers a vector of dual field coefficients on the tdofs given a vector of dual coefficients and a f...
L2ProjectionH1Space(const FiniteElementSpace &fes_ho_, const FiniteElementSpace &fes_lor_, const bool use_ea_, MemoryType d_mt_=Device::GetHostMemoryType())
std::unique_ptr< SparseMatrix > AllocR()
Computes sparsity pattern and initializes R matrix. Based on BilinearForm::AllocMat(),...
void MultTranspose(const Vector &x, Vector &y) const override
void LumpedMassInverse(Vector &ML_inv) const
Returns the inverse of an on-rank lumped mass matrix.
std::unique_ptr< FiniteElementSpace > fes_lor_scalar
Definition transfer.hpp:476
void SetAbsTol(real_t p_atol_) override
Sets absolute tolerance in preconditioned conjugate gradient solver.
void ProlongateTranspose(const Vector &x, Vector &y) const override
std::pair< std::unique_ptr< SparseMatrix >, std::unique_ptr< SparseMatrix > > ComputeSparseRAndM_LH()
Computes on-rank R and M_LH matrices. If true, computes mixed mass and/or inverse lumped mass matrix ...
void SetRelTol(real_t p_rtol_) override
Sets relative tolerance in preconditioned conjugate gradient solver.
std::unique_ptr< ParFiniteElementSpace > pfes_lor_scalar
Definition transfer.hpp:484
void ProlongateTranspose(const Vector &x, Vector &y) const override
Definition transfer.cpp:960
L2ProjectionL2Space(const FiniteElementSpace &fes_ho_, const FiniteElementSpace &fes_lor_, const bool use_ea_, MemoryType d_mt_=Device::GetHostMemoryType())
Definition transfer.cpp:477
void EAMultTranspose(const Vector &x, Vector &y) const
Definition transfer.cpp:887
void EAProlongate(const Vector &x, Vector &y) const
Definition transfer.cpp:946
void MultTranspose(const Vector &x, Vector &y) const override
Definition transfer.cpp:845
void EAMult(const Vector &x, Vector &y) const
Perform mult on the device (same as above)
Definition transfer.cpp:831
void SetRelTol(real_t p_rtol_) override
No-op.
Definition transfer.hpp:340
void EAProlongateTranspose(const Vector &x, Vector &y) const
void SetAbsTol(real_t p_atol_) override
No-op.
Definition transfer.hpp:341
void Prolongate(const Vector &x, Vector &y) const override
Definition transfer.cpp:901
void Mult(const Vector &x, Vector &y) const override
Definition transfer.cpp:792
void BuildHo2Lor(int nel_ho, int nel_lor, const CoarseFineTransformations &cf_tr)
Definition transfer.cpp:239
void MixedMassEA(const FiniteElementSpace &fes_ho_, const FiniteElementSpace &fes_lor_, Vector &M_LH, MemoryType d_mt_=Device::GetHostMemoryType())
Definition transfer.cpp:326
virtual void ProlongateTranspose(const Vector &x, Vector &y) const =0
virtual void SetRelTol(real_t p_rtol_)=0
Sets relative tolerance in preconditioned conjugate gradient solver.
L2Projection(const FiniteElementSpace &fes_ho_, const FiniteElementSpace &fes_lor_, MemoryType d_mt_=Device::GetHostMemoryType())
Definition transfer.cpp:232
void ElemMixedMass(Geometry::Type geom, const FiniteElement &fe_ho, const FiniteElement &fe_lor, ElementTransformation *tr_ho, ElementTransformation *tr_lor, IntegrationPointTransformation &ip_tr, DenseMatrix &M_mixed_el) const
Definition transfer.cpp:259
virtual void SetAbsTol(real_t p_atol_)=0
Sets absolute tolerance in preconditioned conjugate gradient solver.
virtual void Prolongate(const Vector &x, Vector &y) const =0
L2Prolongation(const L2Projection &l2proj_)
Definition transfer.hpp:498
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:504
void Mult(const Vector &x, Vector &y) const
Operator application: y=A(x).
Definition transfer.hpp:500
Transfer data in L2 and H1 finite element spaces between a coarse mesh and an embedded refined mesh u...
Definition transfer.hpp:184
L2Projection * F
Forward, coarse-to-fine, operator.
Definition transfer.hpp:511
const Operator & BackwardOperator() override
Return an Operator that transfers GridFunctions from the range FE space back to GridFunctions in the ...
bool SupportsBackwardsOperator() const override
L2Prolongation * B
Backward, fine-to-coarse, operator.
Definition transfer.hpp:512
const Operator & ForwardOperator() override
Return an Operator that transfers GridFunctions from the domain FE space to GridFunctions in the rang...
L2ProjectionGridTransfer(FiniteElementSpace &coarse_fes_, FiniteElementSpace &fine_fes_, bool force_l2_space_=false, MemoryType d_mt_=Device::GetHostMemoryType())
Definition transfer.hpp:516
Pointer to an Operator of a specified type.
Definition handle.hpp:34
Abstract operator.
Definition operator.hpp:25
int Height() const
Get the height (size of output) of the Operator. Synonym with NumRows().
Definition operator.hpp:66
int Width() const
Get the width (size of input) of the Operator. Synonym with NumCols().
Definition operator.hpp:72
Type
Enumeration defining IDs for some classes derived from Operator.
Definition operator.hpp:284
Matrix-free transfer operator between finite element spaces on the same mesh.
Definition transfer.hpp:567
PRefinementTransferOperator(const FiniteElementSpace &lFESpace_, const FiniteElementSpace &hFESpace_)
Constructs a transfer operator from lFESpace to hFESpace which have different FE collections.
virtual ~PRefinementTransferOperator()
Destructor.
Definition transfer.hpp:583
void MultTranspose(const Vector &x, Vector &y) const override
Restriction by applying the transpose of the Mult method.
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...
Abstract parallel finite element space.
Definition pfespace.hpp:29
Data type sparse matrix.
Definition sparsemat.hpp:51
Matrix-free transfer operator between finite element spaces on the same mesh exploiting the tensor pr...
Definition transfer.hpp:598
virtual ~TensorProductPRefinementTransferOperator()
Destructor.
Definition transfer.hpp:625
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...
TensorProductPRefinementTransferOperator(const FiniteElementSpace &lFESpace_, const FiniteElementSpace &hFESpace_)
Constructs a transfer operator from lFESpace to hFESpace which have different FE collections.
void MultTranspose(const Vector &x, Vector &y) const override
Restriction by applying the transpose of the Mult method.
Matrix-free transfer operator between finite element spaces.
Definition transfer.hpp:536
virtual ~TransferOperator()
Destructor.
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...
void MultTranspose(const Vector &x, Vector &y) const override
Restriction by applying the transpose of the Mult method.
TransferOperator(const FiniteElementSpace &lFESpace, const FiniteElementSpace &hFESpace)
Constructs a transfer operator from lFESpace to hFESpace.
Matrix-free transfer operator between finite element spaces working on true degrees of freedom.
Definition transfer.hpp:640
~TrueTransferOperator()
Destructor.
void MultTranspose(const Vector &x, Vector &y) const override
Restriction by applying the transpose of the Mult method.
TrueTransferOperator(const FiniteElementSpace &lFESpace_, const FiniteElementSpace &hFESpace_)
Constructs a transfer operator working on true degrees of freedom from lFESpace to hFESpace.
void Mult(const Vector &x, Vector &y) const override
Interpolation or prolongation of a true dof vector x to a true dof vector y.
Vector data type.
Definition vector.hpp:82
float real_t
Definition config.hpp:43
MemoryType
Memory types supported by MFEM.
Defines the coarse-fine transformations of all fine elements.
Definition ncmesh.hpp:90