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

Class representing the DPG weak formulation for complex valued systems (see the class DPGWeakForm). More...

#include <complexweakform.hpp>

Inheritance diagram for mfem::ComplexDPGWeakForm:
[legend]
Collaboration diagram for mfem::ComplexDPGWeakForm:
[legend]

Public Member Functions

 ComplexDPGWeakForm ()
 
 ComplexDPGWeakForm (Array< FiniteElementSpace * > &fes_, Array< FiniteElementCollection * > &fecol_)
 Creates bilinear form associated with FE spaces fes_.
 
void SetTestFECollVdim (int test_fec, int vdim)
 
void SetSpaces (Array< FiniteElementSpace * > &fes_, Array< FiniteElementCollection * > &fecol_)
 
int Size () const
 
void AllocateMatrix ()
 
void Finalize (int skip_zeros=1)
 Finalizes the matrix initialization.
 
BlockMatrixBlockMat_r ()
 Returns a reference to the BlockMatrix: \( M_r \).
 
BlockMatrixBlockMat_i ()
 Returns a reference to the BlockMatrix: \( M_i \).
 
BlockMatrixBlockMatElim_r ()
 Returns a reference to the BlockMatrix of eliminated b.c.: \( M_e_r \).
 
BlockMatrixBlockMatElim_i ()
 Returns a reference to the BlockMatrix of eliminated b.c.: \( M_e_i \).
 
void AddTrialIntegrator (BilinearFormIntegrator *bfi_r, BilinearFormIntegrator *bfi_i, int n, int m)
 Adds new Domain BF Integrator. Assumes ownership of bfi.
 
void AddTestIntegrator (BilinearFormIntegrator *bfi_r, BilinearFormIntegrator *bfi_i, int n, int m)
 Adds new Test Integrator. Assumes ownership of bfi_r and bfi_i.
 
void AddDomainLFIntegrator (LinearFormIntegrator *lfi_r, LinearFormIntegrator *lfi_i, int n)
 Adds new Domain LF Integrator. Assumes ownership of lfi_r and lfi_i.
 
void Assemble (int skip_zeros=1)
 Assembles the form i.e. sums over all integrators.
 
virtual void FormLinearSystem (const Array< int > &ess_tdof_list, Vector &x, OperatorHandle &A, Vector &X, Vector &B, int copy_interior=0)
 
template<typename OpType >
void FormLinearSystem (const Array< int > &ess_tdof_list, Vector &x, OpType &A, Vector &X, Vector &B, int copy_interior=0)
 
virtual void FormSystemMatrix (const Array< int > &ess_tdof_list, OperatorHandle &A)
 
template<typename OpType >
void FormSystemMatrix (const Array< int > &ess_tdof_list, OpType &A)
 
void EliminateVDofs (const Array< int > &vdofs, Operator::DiagonalPolicy dpolicy=Operator::DIAG_ONE)
 
void EliminateVDofsInRHS (const Array< int > &vdofs, const Vector &x_r, const Vector &x_i, Vector &b_r, Vector &b_i)
 
virtual void RecoverFEMSolution (const Vector &X, Vector &x)
 
void SetDiagonalPolicy (Operator::DiagonalPolicy policy)
 Sets diagonal policy used upon construction of the linear system.
 
virtual void Update ()
 
void StoreMatrices (bool store_matrices_=true)
 
void EnableStaticCondensation ()
 
VectorComputeResidual (const Vector &x)
 
virtual ~ComplexDPGWeakForm ()
 Destroys bilinear form.
 

Protected Member Functions

void Init ()
 
void ReleaseInitMemory ()
 
void AllocMat ()
 
void ConformingAssemble ()
 
void ComputeOffsets ()
 
virtual void BuildProlongation ()
 

Protected Attributes

ComplexBlockStaticCondensationstatic_cond
 Owned.
 
bool initialized = false
 
Meshmesh = nullptr
 
int height
 
int width
 
int nblocks
 
Array< int > dof_offsets
 
Array< int > tdof_offsets
 
BlockMatrixmat_r = nullptr
 Block matrix \( M \) to be associated with the real/imag Block bilinear form. Owned.
 
BlockMatrixmat_i = nullptr
 
ComplexOperatormat = nullptr
 
BlockVectory_r = nullptr
 BlockVectors to be associated with the real/imag Block linear form.
 
BlockVectory_i = nullptr
 
Vectory = nullptr
 
BlockMatrixmat_e_r = nullptr
 Block Matrix \( M_e \) used to store the eliminations from the b.c. Owned. \( M + M_e = M_{original} \).
 
BlockMatrixmat_e_i = nullptr
 
Array< FiniteElementSpace * > trial_fes
 Trial FE spaces.
 
Array< int > IsTraceFes
 Flags to determine if a FiniteElementSpace is Trace.
 
Array< FiniteElementCollection * > test_fecols
 Test FE Collections (Broken)
 
Array< int > test_fecols_vdims
 
Array2D< Array< BilinearFormIntegrator * > * > trial_integs_r
 Set of Trial Integrators to be applied for matrix B.
 
Array2D< Array< BilinearFormIntegrator * > * > trial_integs_i
 
Array2D< Array< BilinearFormIntegrator * > * > test_integs_r
 Set of Test Space (broken) Integrators to be applied for matrix G.
 
Array2D< Array< BilinearFormIntegrator * > * > test_integs_i
 
Array< Array< LinearFormIntegrator * > * > lfis_r
 Set of LinearForm Integrators to be applied.
 
Array< Array< LinearFormIntegrator * > * > lfis_i
 
BlockMatrixP = nullptr
 Block Prolongation.
 
BlockMatrixR = nullptr
 Block Restriction.
 
mfem::Operator::DiagonalPolicy diag_policy
 
bool store_matrices = false
 
Array< ComplexDenseMatrix * > Bmat
 
Array< Vector * > fvec
 
Vector residuals
 

Detailed Description

Class representing the DPG weak formulation for complex valued systems (see the class DPGWeakForm).

Definition at line 23 of file complexweakform.hpp.

Constructor & Destructor Documentation

◆ ComplexDPGWeakForm() [1/2]

mfem::ComplexDPGWeakForm::ComplexDPGWeakForm ( )
inline

Definition at line 107 of file complexweakform.hpp.

◆ ComplexDPGWeakForm() [2/2]

mfem::ComplexDPGWeakForm::ComplexDPGWeakForm ( Array< FiniteElementSpace * > & fes_,
Array< FiniteElementCollection * > & fecol_ )
inline

Creates bilinear form associated with FE spaces fes_.

Definition at line 114 of file complexweakform.hpp.

◆ ~ComplexDPGWeakForm()

mfem::ComplexDPGWeakForm::~ComplexDPGWeakForm ( )
virtual

Destroys bilinear form.

Definition at line 1024 of file complexweakform.cpp.

Member Function Documentation

◆ AddDomainLFIntegrator()

void mfem::ComplexDPGWeakForm::AddDomainLFIntegrator ( LinearFormIntegrator * lfi_r,
LinearFormIntegrator * lfi_i,
int n )

Adds new Domain LF Integrator. Assumes ownership of lfi_r and lfi_i.

Adds new Domain LF Integrator. Assumes ownership of bfi.

Definition at line 162 of file complexweakform.cpp.

◆ AddTestIntegrator()

void mfem::ComplexDPGWeakForm::AddTestIntegrator ( BilinearFormIntegrator * bfi_r,
BilinearFormIntegrator * bfi_i,
int n,
int m )

Adds new Test Integrator. Assumes ownership of bfi_r and bfi_i.

Adds new Domain BF Integrator. Assumes ownership of bfi.

Definition at line 144 of file complexweakform.cpp.

◆ AddTrialIntegrator()

void mfem::ComplexDPGWeakForm::AddTrialIntegrator ( BilinearFormIntegrator * bfi_r,
BilinearFormIntegrator * bfi_i,
int n,
int m )

Adds new Domain BF Integrator. Assumes ownership of bfi.

Adds new Trial Integrator. Assumes ownership of bfi_r and bfi_i. n and m correspond to the trial FESpace and test FEColl respectively

Definition at line 124 of file complexweakform.cpp.

◆ AllocateMatrix()

void mfem::ComplexDPGWeakForm::AllocateMatrix ( )
inline

Definition at line 152 of file complexweakform.hpp.

◆ AllocMat()

void mfem::ComplexDPGWeakForm::AllocMat ( )
protected

Definition at line 83 of file complexweakform.cpp.

◆ Assemble()

void mfem::ComplexDPGWeakForm::Assemble ( int skip_zeros = 1)

Assembles the form i.e. sums over all integrators.

Assembles the form i.e. sums over all domain integrators.

Definition at line 318 of file complexweakform.cpp.

◆ BlockMat_i()

BlockMatrix & mfem::ComplexDPGWeakForm::BlockMat_i ( )
inline

Returns a reference to the BlockMatrix: \( M_i \).

Definition at line 164 of file complexweakform.hpp.

◆ BlockMat_r()

BlockMatrix & mfem::ComplexDPGWeakForm::BlockMat_r ( )
inline

Returns a reference to the BlockMatrix: \( M_r \).

Definition at line 158 of file complexweakform.hpp.

◆ BlockMatElim_i()

BlockMatrix & mfem::ComplexDPGWeakForm::BlockMatElim_i ( )
inline

Returns a reference to the BlockMatrix of eliminated b.c.: \( M_e_i \).

Definition at line 178 of file complexweakform.hpp.

◆ BlockMatElim_r()

BlockMatrix & mfem::ComplexDPGWeakForm::BlockMatElim_r ( )
inline

Returns a reference to the BlockMatrix of eliminated b.c.: \( M_e_r \).

Definition at line 171 of file complexweakform.hpp.

◆ BuildProlongation()

void mfem::ComplexDPGWeakForm::BuildProlongation ( )
protectedvirtual

Reimplemented in mfem::ParComplexDPGWeakForm.

Definition at line 178 of file complexweakform.cpp.

◆ ComputeOffsets()

void mfem::ComplexDPGWeakForm::ComputeOffsets ( )
protected

Definition at line 67 of file complexweakform.cpp.

◆ ComputeResidual()

Vector & mfem::ComplexDPGWeakForm::ComputeResidual ( const Vector & x)

Definition at line 924 of file complexweakform.cpp.

◆ ConformingAssemble()

void mfem::ComplexDPGWeakForm::ConformingAssemble ( )
protected

Definition at line 196 of file complexweakform.cpp.

◆ EliminateVDofs()

void mfem::ComplexDPGWeakForm::EliminateVDofs ( const Array< int > & vdofs,
Operator::DiagonalPolicy dpolicy = Operator::DIAG_ONE )

Definition at line 740 of file complexweakform.cpp.

◆ EliminateVDofsInRHS()

void mfem::ComplexDPGWeakForm::EliminateVDofsInRHS ( const Array< int > & vdofs,
const Vector & x_r,
const Vector & x_i,
Vector & b_r,
Vector & b_i )

Definition at line 727 of file complexweakform.cpp.

◆ EnableStaticCondensation()

void mfem::ComplexDPGWeakForm::EnableStaticCondensation ( )

Definition at line 918 of file complexweakform.cpp.

◆ Finalize()

void mfem::ComplexDPGWeakForm::Finalize ( int skip_zeros = 1)

Finalizes the matrix initialization.

Definition at line 108 of file complexweakform.cpp.

◆ FormLinearSystem() [1/2]

void mfem::ComplexDPGWeakForm::FormLinearSystem ( const Array< int > & ess_tdof_list,
Vector & x,
OperatorHandle & A,
Vector & X,
Vector & B,
int copy_interior = 0 )
virtual

Reimplemented in mfem::ParComplexDPGWeakForm.

Definition at line 613 of file complexweakform.cpp.

◆ FormLinearSystem() [2/2]

template<typename OpType >
void mfem::ComplexDPGWeakForm::FormLinearSystem ( const Array< int > & ess_tdof_list,
Vector & x,
OpType & A,
Vector & X,
Vector & B,
int copy_interior = 0 )
inline

Definition at line 210 of file complexweakform.hpp.

◆ FormSystemMatrix() [1/2]

void mfem::ComplexDPGWeakForm::FormSystemMatrix ( const Array< int > & ess_tdof_list,
OperatorHandle & A )
virtual

Reimplemented in mfem::ParComplexDPGWeakForm.

Definition at line 690 of file complexweakform.cpp.

◆ FormSystemMatrix() [2/2]

template<typename OpType >
void mfem::ComplexDPGWeakForm::FormSystemMatrix ( const Array< int > & ess_tdof_list,
OpType & A )
inline

Definition at line 226 of file complexweakform.hpp.

◆ Init()

void mfem::ComplexDPGWeakForm::Init ( )
protected

Definition at line 17 of file complexweakform.cpp.

◆ RecoverFEMSolution()

void mfem::ComplexDPGWeakForm::RecoverFEMSolution ( const Vector & X,
Vector & x )
virtual

Reimplemented in mfem::ParComplexDPGWeakForm.

Definition at line 768 of file complexweakform.cpp.

◆ ReleaseInitMemory()

void mfem::ComplexDPGWeakForm::ReleaseInitMemory ( )
protected

Definition at line 807 of file complexweakform.cpp.

◆ SetDiagonalPolicy()

void mfem::ComplexDPGWeakForm::SetDiagonalPolicy ( Operator::DiagonalPolicy policy)
inline

Sets diagonal policy used upon construction of the linear system.

Policies include:

  • DIAG_ZERO (Set the diagonal values to zero)
  • DIAG_ONE (Set the diagonal values to one)
  • DIAG_KEEP (Keep the diagonal values)

Definition at line 250 of file complexweakform.hpp.

◆ SetSpaces()

void mfem::ComplexDPGWeakForm::SetSpaces ( Array< FiniteElementSpace * > & fes_,
Array< FiniteElementCollection * > & fecol_ )
inline

Definition at line 125 of file complexweakform.hpp.

◆ SetTestFECollVdim()

void mfem::ComplexDPGWeakForm::SetTestFECollVdim ( int test_fec,
int vdim )
inline

Definition at line 120 of file complexweakform.hpp.

◆ Size()

int mfem::ComplexDPGWeakForm::Size ( ) const
inline

Definition at line 149 of file complexweakform.hpp.

◆ StoreMatrices()

void mfem::ComplexDPGWeakForm::StoreMatrices ( bool store_matrices_ = true)
inline

Definition at line 257 of file complexweakform.hpp.

◆ Update()

void mfem::ComplexDPGWeakForm::Update ( )
virtual

Reimplemented in mfem::ParComplexDPGWeakForm.

Definition at line 867 of file complexweakform.cpp.

Member Data Documentation

◆ Bmat

Array<ComplexDenseMatrix * > mfem::ComplexDPGWeakForm::Bmat
protected

Store the matrix L^-1 B and Vector L^-1 l where G = L L^t

Definition at line 99 of file complexweakform.hpp.

◆ diag_policy

mfem::Operator::DiagonalPolicy mfem::ComplexDPGWeakForm::diag_policy
protected

Definition at line 81 of file complexweakform.hpp.

◆ dof_offsets

Array<int> mfem::ComplexDPGWeakForm::dof_offsets
protected

Definition at line 35 of file complexweakform.hpp.

◆ fvec

Array<Vector * > mfem::ComplexDPGWeakForm::fvec
protected

Definition at line 100 of file complexweakform.hpp.

◆ height

int mfem::ComplexDPGWeakForm::height
protected

Definition at line 33 of file complexweakform.hpp.

◆ initialized

bool mfem::ComplexDPGWeakForm::initialized = false
protected

Definition at line 30 of file complexweakform.hpp.

◆ IsTraceFes

Array<int> mfem::ComplexDPGWeakForm::IsTraceFes
protected

Flags to determine if a FiniteElementSpace is Trace.

Definition at line 58 of file complexweakform.hpp.

◆ lfis_i

Array<Array<LinearFormIntegrator * > * > mfem::ComplexDPGWeakForm::lfis_i
protected

Definition at line 74 of file complexweakform.hpp.

◆ lfis_r

Array<Array<LinearFormIntegrator * > * > mfem::ComplexDPGWeakForm::lfis_r
protected

Set of LinearForm Integrators to be applied.

Definition at line 73 of file complexweakform.hpp.

◆ mat

ComplexOperator* mfem::ComplexDPGWeakForm::mat = nullptr
protected

Definition at line 41 of file complexweakform.hpp.

◆ mat_e_i

BlockMatrix* mfem::ComplexDPGWeakForm::mat_e_i = nullptr
protected

Definition at line 52 of file complexweakform.hpp.

◆ mat_e_r

BlockMatrix* mfem::ComplexDPGWeakForm::mat_e_r = nullptr
protected

Block Matrix \( M_e \) used to store the eliminations from the b.c. Owned. \( M + M_e = M_{original} \).

Definition at line 51 of file complexweakform.hpp.

◆ mat_i

BlockMatrix* mfem::ComplexDPGWeakForm::mat_i = nullptr
protected

Definition at line 40 of file complexweakform.hpp.

◆ mat_r

BlockMatrix* mfem::ComplexDPGWeakForm::mat_r = nullptr
protected

Block matrix \( M \) to be associated with the real/imag Block bilinear form. Owned.

Definition at line 39 of file complexweakform.hpp.

◆ mesh

Mesh* mfem::ComplexDPGWeakForm::mesh = nullptr
protected

Definition at line 32 of file complexweakform.hpp.

◆ nblocks

int mfem::ComplexDPGWeakForm::nblocks
protected

Definition at line 34 of file complexweakform.hpp.

◆ P

BlockMatrix* mfem::ComplexDPGWeakForm::P = nullptr
protected

Block Prolongation.

Definition at line 77 of file complexweakform.hpp.

◆ R

BlockMatrix* mfem::ComplexDPGWeakForm::R = nullptr
protected

Block Restriction.

Definition at line 79 of file complexweakform.hpp.

◆ residuals

Vector mfem::ComplexDPGWeakForm::residuals
protected

Definition at line 101 of file complexweakform.hpp.

◆ static_cond

ComplexBlockStaticCondensation* mfem::ComplexDPGWeakForm::static_cond
protected

Owned.

Definition at line 28 of file complexweakform.hpp.

◆ store_matrices

bool mfem::ComplexDPGWeakForm::store_matrices = false
protected

Definition at line 95 of file complexweakform.hpp.

◆ tdof_offsets

Array<int> mfem::ComplexDPGWeakForm::tdof_offsets
protected

Definition at line 36 of file complexweakform.hpp.

◆ test_fecols

Array<FiniteElementCollection *> mfem::ComplexDPGWeakForm::test_fecols
protected

Test FE Collections (Broken)

Definition at line 61 of file complexweakform.hpp.

◆ test_fecols_vdims

Array<int> mfem::ComplexDPGWeakForm::test_fecols_vdims
protected

Definition at line 62 of file complexweakform.hpp.

◆ test_integs_i

Array2D<Array<BilinearFormIntegrator * > * > mfem::ComplexDPGWeakForm::test_integs_i
protected

Definition at line 70 of file complexweakform.hpp.

◆ test_integs_r

Array2D<Array<BilinearFormIntegrator * > * > mfem::ComplexDPGWeakForm::test_integs_r
protected

Set of Test Space (broken) Integrators to be applied for matrix G.

Definition at line 69 of file complexweakform.hpp.

◆ trial_fes

Array<FiniteElementSpace * > mfem::ComplexDPGWeakForm::trial_fes
protected

Trial FE spaces.

Definition at line 55 of file complexweakform.hpp.

◆ trial_integs_i

Array2D<Array<BilinearFormIntegrator * > * > mfem::ComplexDPGWeakForm::trial_integs_i
protected

Definition at line 66 of file complexweakform.hpp.

◆ trial_integs_r

Array2D<Array<BilinearFormIntegrator * > * > mfem::ComplexDPGWeakForm::trial_integs_r
protected

Set of Trial Integrators to be applied for matrix B.

Definition at line 65 of file complexweakform.hpp.

◆ width

int mfem::ComplexDPGWeakForm::width
protected

Definition at line 33 of file complexweakform.hpp.

◆ y

Vector* mfem::ComplexDPGWeakForm::y = nullptr
protected

Definition at line 46 of file complexweakform.hpp.

◆ y_i

BlockVector* mfem::ComplexDPGWeakForm::y_i = nullptr
protected

Definition at line 45 of file complexweakform.hpp.

◆ y_r

BlockVector* mfem::ComplexDPGWeakForm::y_r = nullptr
protected

BlockVectors to be associated with the real/imag Block linear form.

Definition at line 44 of file complexweakform.hpp.


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