MFEM v4.7.0
Finite element discretization library
Loading...
Searching...
No Matches
mfem::Hybridization Class Reference

Auxiliary class Hybridization, used to implement BilinearForm hybridization. More...

#include <hybridization.hpp>

Collaboration diagram for mfem::Hybridization:
[legend]

Public Member Functions

 Hybridization (FiniteElementSpace *fespace, FiniteElementSpace *c_fespace)
 Constructor.
 
 ~Hybridization ()
 Destructor.
 
void SetConstraintIntegrator (BilinearFormIntegrator *c_integ)
 
void Init (const Array< int > &ess_tdof_list)
 Prepare the Hybridization object for assembly.
 
void AssembleMatrix (int el, const DenseMatrix &A)
 Assemble the element matrix A into the hybridized system matrix.
 
void AssembleBdrMatrix (int bdr_el, const DenseMatrix &A)
 Assemble the boundary element matrix A into the hybridized system matrix.
 
void Finalize ()
 Finalize the construction of the hybridized matrix.
 
SparseMatrixGetMatrix ()
 Return the serial hybridized matrix.
 
HypreParMatrixGetParallelMatrix ()
 Return the parallel hybridized matrix.
 
void GetParallelMatrix (OperatorHandle &H_h) const
 Return the parallel hybridized matrix in the format specified by SetOperatorType().
 
void SetOperatorType (Operator::Type tid)
 Set the operator type id for the parallel hybridized matrix/operator.
 
void ReduceRHS (const Vector &b, Vector &b_r) const
 
void ComputeSolution (const Vector &b, const Vector &sol_r, Vector &sol) const
 
void Reset ()
 Destroy the current hybridization matrix while preserving the computed constraint matrix and the set of essential true dofs. After Reset(), a new hybridized matrix can be assembled via AssembleMatrix() and Finalize(). The Mesh and FiniteElementSpace objects are assumed to be un-modified. If that is not the case, a new Hybridization object must be created.
 

Protected Member Functions

void ConstructC ()
 
void GetIBDofs (int el, Array< int > &i_dofs, Array< int > &b_dofs) const
 
void GetBDofs (int el, int &num_idofs, Array< int > &b_dofs) const
 
void ComputeH ()
 
void MultAfInv (const Vector &b, const Vector &lambda, Vector &bf, int mode) const
 

Protected Attributes

FiniteElementSpacefes
 
FiniteElementSpacec_fes
 
BilinearFormIntegratorc_bfi
 
SparseMatrixCt
 
SparseMatrixH
 
Array< int > hat_offsets
 
Array< int > hat_dofs_marker
 
Array< int > Af_offsets
 
Array< int > Af_f_offsets
 
real_tAf_data
 
int * Af_ipiv
 
HypreParMatrixpC
 
HypreParMatrixP_pc
 
OperatorHandle pH
 

Detailed Description

Auxiliary class Hybridization, used to implement BilinearForm hybridization.

Hybridization can be viewed as a technique for solving linear systems obtained through finite element assembly. The assembled matrix \( A \) can be written as:

\[ A = P^T \hat{A} P, \]

where \( P \) is the matrix mapping the conforming finite element space to the purely local finite element space without any inter-element constraints imposed, and \( \hat{A} \) is the block-diagonal matrix of all element matrices.

We assume that:

  • \( \hat{A} \) is invertible,
  • \( P \) has a left inverse \( R \), such that \( R P = I \),
  • a constraint matrix \( C \) can be constructed, such that \( \operatorname{Ker}(C) = \operatorname{Im}(P) \).

Under these conditions, the linear system \( A x = b \) can be solved using the following procedure:

  • solve for \( \lambda \) in the linear system:

    \[ (C \hat{A}^{-1} C^T) \lambda = C \hat{A}^{-1} R^T b \]

  • compute \( x = R \hat{A}^{-1} (R^T b - C^T \lambda) \)

Hybridization is advantageous when the matrix \( H = (C \hat{A}^{-1} C^T) \) of the hybridized system is either smaller than the original system, or is simpler to invert with a known method.

In some cases, e.g. high-order elements, the matrix \( C \) can be written as

\[ C = \begin{pmatrix} 0 & C_b \end{pmatrix}, \]

and then the hybridized matrix \( H \) can be assembled using the identity

\[ H = C_b S_b^{-1} C_b^T, \]

where \( S_b \) is the Schur complement of \( \hat{A} \) with respect to the same decomposition as the columns of \( C \):

\[ S_b = \hat{A}_b - \hat{A}_{bf} \hat{A}_{f}^{-1} \hat{A}_{fb}. \]

Hybridization can also be viewed as a discretization method for imposing (weak) continuity constraints between neighboring elements.

Definition at line 61 of file hybridization.hpp.

Constructor & Destructor Documentation

◆ Hybridization()

mfem::Hybridization::Hybridization ( FiniteElementSpace * fespace,
FiniteElementSpace * c_fespace )

Constructor.

Definition at line 30 of file hybridization.cpp.

◆ ~Hybridization()

mfem::Hybridization::~Hybridization ( )

Destructor.

Definition at line 41 of file hybridization.cpp.

Member Function Documentation

◆ AssembleBdrMatrix()

void mfem::Hybridization::AssembleBdrMatrix ( int bdr_el,
const DenseMatrix & A )

Assemble the boundary element matrix A into the hybridized system matrix.

Definition at line 482 of file hybridization.cpp.

◆ AssembleMatrix()

void mfem::Hybridization::AssembleMatrix ( int el,
const DenseMatrix & A )

Assemble the element matrix A into the hybridized system matrix.

Definition at line 442 of file hybridization.cpp.

◆ ComputeH()

void mfem::Hybridization::ComputeH ( )
protected

Definition at line 570 of file hybridization.cpp.

◆ ComputeSolution()

void mfem::Hybridization::ComputeSolution ( const Vector & b,
const Vector & sol_r,
Vector & sol ) const

Reconstruct the solution of the original system, sol, from solution of the hybridized system, sol_r, and the original r.h.s. vector, b. It is assumed that the vector sol has the right essential b.c.

Definition at line 841 of file hybridization.cpp.

◆ ConstructC()

void mfem::Hybridization::ConstructC ( )
protected

Definition at line 54 of file hybridization.cpp.

◆ Finalize()

void mfem::Hybridization::Finalize ( )

Finalize the construction of the hybridized matrix.

Definition at line 717 of file hybridization.cpp.

◆ GetBDofs()

void mfem::Hybridization::GetBDofs ( int el,
int & num_idofs,
Array< int > & b_dofs ) const
protected

Definition at line 426 of file hybridization.cpp.

◆ GetIBDofs()

void mfem::Hybridization::GetIBDofs ( int el,
Array< int > & i_dofs,
Array< int > & b_dofs ) const
protected

Definition at line 406 of file hybridization.cpp.

◆ GetMatrix()

SparseMatrix & mfem::Hybridization::GetMatrix ( )
inline

Return the serial hybridized matrix.

Definition at line 122 of file hybridization.hpp.

◆ GetParallelMatrix() [1/2]

HypreParMatrix & mfem::Hybridization::GetParallelMatrix ( )
inline

Return the parallel hybridized matrix.

Definition at line 126 of file hybridization.hpp.

◆ GetParallelMatrix() [2/2]

void mfem::Hybridization::GetParallelMatrix ( OperatorHandle & H_h) const
inline

Return the parallel hybridized matrix in the format specified by SetOperatorType().

Definition at line 130 of file hybridization.hpp.

◆ Init()

void mfem::Hybridization::Init ( const Array< int > & ess_tdof_list)

Prepare the Hybridization object for assembly.

Definition at line 215 of file hybridization.cpp.

◆ MultAfInv()

void mfem::Hybridization::MultAfInv ( const Vector & b,
const Vector & lambda,
Vector & bf,
int mode ) const
protected

Definition at line 726 of file hybridization.cpp.

◆ ReduceRHS()

void mfem::Hybridization::ReduceRHS ( const Vector & b,
Vector & b_r ) const

Perform the reduction of the given r.h.s. vector, b, to a r.h.s vector, b_r, for the hybridized system.

Definition at line 807 of file hybridization.cpp.

◆ Reset()

void mfem::Hybridization::Reset ( )

Destroy the current hybridization matrix while preserving the computed constraint matrix and the set of essential true dofs. After Reset(), a new hybridized matrix can be assembled via AssembleMatrix() and Finalize(). The Mesh and FiniteElementSpace objects are assumed to be un-modified. If that is not the case, a new Hybridization object must be created.

Definition at line 880 of file hybridization.cpp.

◆ SetConstraintIntegrator()

void mfem::Hybridization::SetConstraintIntegrator ( BilinearFormIntegrator * c_integ)
inline

Set the integrator that will be used to construct the constraint matrix C. The Hybridization object assumes ownership of the integrator, i.e. it will delete the integrator when destroyed.

Definition at line 106 of file hybridization.hpp.

◆ SetOperatorType()

void mfem::Hybridization::SetOperatorType ( Operator::Type tid)
inline

Set the operator type id for the parallel hybridized matrix/operator.

Definition at line 133 of file hybridization.hpp.

Member Data Documentation

◆ Af_data

real_t* mfem::Hybridization::Af_data
protected

Definition at line 71 of file hybridization.hpp.

◆ Af_f_offsets

Array<int> mfem::Hybridization::Af_f_offsets
protected

Definition at line 70 of file hybridization.hpp.

◆ Af_ipiv

int* mfem::Hybridization::Af_ipiv
protected

Definition at line 72 of file hybridization.hpp.

◆ Af_offsets

Array<int> mfem::Hybridization::Af_offsets
protected

Definition at line 70 of file hybridization.hpp.

◆ c_bfi

BilinearFormIntegrator* mfem::Hybridization::c_bfi
protected

Definition at line 65 of file hybridization.hpp.

◆ c_fes

FiniteElementSpace * mfem::Hybridization::c_fes
protected

Definition at line 64 of file hybridization.hpp.

◆ Ct

SparseMatrix* mfem::Hybridization::Ct
protected

Definition at line 67 of file hybridization.hpp.

◆ fes

FiniteElementSpace* mfem::Hybridization::fes
protected

Definition at line 64 of file hybridization.hpp.

◆ H

SparseMatrix * mfem::Hybridization::H
protected

Definition at line 67 of file hybridization.hpp.

◆ hat_dofs_marker

Array<int> mfem::Hybridization::hat_dofs_marker
protected

Definition at line 69 of file hybridization.hpp.

◆ hat_offsets

Array<int> mfem::Hybridization::hat_offsets
protected

Definition at line 69 of file hybridization.hpp.

◆ P_pc

HypreParMatrix * mfem::Hybridization::P_pc
protected

Definition at line 75 of file hybridization.hpp.

◆ pC

HypreParMatrix* mfem::Hybridization::pC
protected

Definition at line 75 of file hybridization.hpp.

◆ pH

OperatorHandle mfem::Hybridization::pH
protected

Definition at line 76 of file hybridization.hpp.


The documentation for this class was generated from the following files: