Bilinear form for the stiffness matrix compression with aglomerated wavelets. More...

#include <bform.hh>

Inheritance diagram for aglowav2::AglowavBF00< F >:
concepts::BilinearForm< concepts::Real > concepts::Cloneable concepts::OutputOperator

Classes

struct  Cell
 Entry in the cache to compute wavelet entries. More...
 
struct  Nfld
 Structure to hold the one scale near field one each level. More...
 

Public Member Functions

 AglowavBF00 (concepts::BilinearForm< F > &bf, F00< F > &XFY)
 Constructor. More...
 
virtual BilinearForm * clone () const=0
 Virtual constructor. More...
 
uint memory () const
 
void operator() (const concepts::Element< F > &elmX, const concepts::Element< F > &elmY, concepts::ElementMatrix< F > &em)
 
virtual void operator() (const Element< typename Realtype< concepts::Real >::type > &elmX, const Element< typename Realtype< concepts::Real >::type > &elmY, ElementMatrix< concepts::Real > &em) const=0
 Evaluates the bilinear form for all shape functions on elmX and elmY and stores the result in the matrix em. More...
 
virtual void operator() (const Element< typename Realtype< concepts::Real >::type > &elmX, const Element< typename Realtype< concepts::Real >::type > &elmY, ElementMatrix< concepts::Real > &em, const ElementPair< typename Realtype< concepts::Real >::type > &ep) const
 Evaluates the bilinear form for all shape functions on elmX and elmY and stores the result in the matrix em. More...
 
void operator() (const Haar3d000< F > &elmX, const Haar3d000< F > &elmY, concepts::ElementMatrix< F > &em)
 
 ~AglowavBF00 ()
 

Protected Member Functions

virtual std::ostream & info (std::ostream &os) const
 Returns information in an output stream. More...
 

Private Types

typedef struct aglowav2::AglowavBF00::Nfld Nfld
 Structure to hold the one scale near field one each level. More...
 

Private Member Functions

virtual AglowavBF00clone () const
 Virtual constructor. More...
 
void nearfield_ (const Haar3d000< F > &elmX, const Haar3d000< F > &elmY, F *valXY, uint ldXY, F *valYX, uint ldYX)
 Computation of the one scale near field on each level (used by the constructor). More...
 
void updateChld_ (const Haar3d000< F > &tst, const Haar3d000< F > &trl, F *valtsttrl, uint ldtsttrl, F *valtrltst, uint ldtrltst)
 Computation of one scale entries. More...
 
void updateElm_ (const Haar3d000< F > &tst, const Haar3d000< F > &trl, F *valtsttrl, uint ldtsttrl, F *valtrltst, uint ldtrltst)
 Computation of one scale entries. More...
 

Private Attributes

concepts::BilinearForm< F > & bf_
 Bilinear form to compute the one scale matrix entries. More...
 
Cellbuf_
 Buffer to store wavelet entries "for a short time". More...
 
Cellcache_
 Cache to store wavelet entries for later usage. More...
 
F * foo_
 
uint mem_
 memory More...
 
Nfld ** nfld_
 Array to hold the near field for each element. More...
 
Nfld *** nfldLast_
 Array to hold pointers to the last list entry of the near field (used only by nearfield_) More...
 
bool reuse_
 Entry reused and thus saved in cache_ otherwise stored in buf_. More...
 
F * valXY_
 Auxiliary vectors to hold intermediate values. More...
 
F * valYX_
 
F * vXY_
 
F * vYX_
 
F00< F > & XFY_
 Cluster approximation of the far field. More...
 

Detailed Description

template<class F = concepts::Real>
class aglowav2::AglowavBF00< F >

Bilinear form for the stiffness matrix compression with aglomerated wavelets.

After computation of the compressed stiffness matrix delete this bilinear form to free the memory used by the one scale near field. Unsymmetric kernels are allowed by the sparsity pattern (Delta) has to be symmetric. Symmetry of the kernel is not exploited.

Parameters
FField (Real or Cmplx)

Definition at line 29 of file bform.hh.

Member Typedef Documentation

◆ Nfld

template<class F = concepts::Real>
typedef struct aglowav2::AglowavBF00::Nfld aglowav2::AglowavBF00< F >::Nfld
private

Structure to hold the one scale near field one each level.

Constructor & Destructor Documentation

◆ AglowavBF00()

template<class F = concepts::Real>
aglowav2::AglowavBF00< F >::AglowavBF00 ( concepts::BilinearForm< F > &  bf,
F00< F > &  XFY 
)

Constructor.

Computation of the one scale near field on each level.

Parameters
bfBilinear form to compute the one scale matrix entries.
XFYCluster approximation of the far field.

◆ ~AglowavBF00()

template<class F = concepts::Real>
aglowav2::AglowavBF00< F >::~AglowavBF00 ( )

Member Function Documentation

◆ clone() [1/2]

template<class F = concepts::Real>
virtual AglowavBF00* aglowav2::AglowavBF00< F >::clone ( ) const
inlineprivatevirtual

Virtual constructor.

Returns a pointer to a copy of itself. The caller is responsible to destroy this copy.

Implements concepts::Cloneable.

Definition at line 53 of file bform.hh.

◆ clone() [2/2]

virtual BilinearForm* concepts::BilinearForm< concepts::Real , typename Realtype<concepts::Real >::type >::clone
pure virtualinherited

Virtual constructor.

Returns a pointer to a copy of itself. The caller is responsible to destroy this copy.

◆ info()

virtual std::ostream& concepts::BilinearForm< concepts::Real , typename Realtype<concepts::Real >::type >::info ( std::ostream &  os) const
protectedvirtualinherited

Returns information in an output stream.

Reimplemented from concepts::OutputOperator.

◆ memory()

template<class F = concepts::Real>
uint aglowav2::AglowavBF00< F >::memory ( ) const
inline

Definition at line 69 of file bform.hh.

◆ nearfield_()

template<class F = concepts::Real>
void aglowav2::AglowavBF00< F >::nearfield_ ( const Haar3d000< F > &  elmX,
const Haar3d000< F > &  elmY,
F *  valXY,
uint  ldXY,
F *  valYX,
uint  ldYX 
)
private

Computation of the one scale near field on each level (used by the constructor).

Call only for elmX.key() >= elmY.key() !!!

Parameters
elmXSpace element
elmYSpace element
valXYNear field of the $\phi$ of elmX with elmY (row wise)
ldXYLeading dimension of valXY
valYXNear field of the $\phi$ of elmY with elmX (row wise)
ldYXLeading dimension of valYX

◆ operator()() [1/4]

template<class F = concepts::Real>
void aglowav2::AglowavBF00< F >::operator() ( const concepts::Element< F > &  elmX,
const concepts::Element< F > &  elmY,
concepts::ElementMatrix< F > &  em 
)

◆ operator()() [2/4]

virtual void concepts::BilinearForm< concepts::Real , typename Realtype<concepts::Real >::type >::operator() ( const Element< typename Realtype<concepts::Real >::type > &  elmX,
const Element< typename Realtype<concepts::Real >::type > &  elmY,
ElementMatrix< concepts::Real > &  em 
) const
pure virtualinherited

Evaluates the bilinear form for all shape functions on elmX and elmY and stores the result in the matrix em.

Postcondition
The returned matrix em has the correct size.
Parameters
elmXLeft element (test functions)
elmYRight element (trial functions)
emReturn element matrix

◆ operator()() [3/4]

virtual void concepts::BilinearForm< concepts::Real , typename Realtype<concepts::Real >::type >::operator() ( const Element< typename Realtype<concepts::Real >::type > &  elmX,
const Element< typename Realtype<concepts::Real >::type > &  elmY,
ElementMatrix< concepts::Real > &  em,
const ElementPair< typename Realtype<concepts::Real >::type > &  ep 
) const
inlinevirtualinherited

Evaluates the bilinear form for all shape functions on elmX and elmY and stores the result in the matrix em.

If this method is not reimplemented in a derived class, the default behaviour is to call the application operator without ep.

Postcondition
The returned matrix em has the correct size.
Parameters
elmXLeft element
elmYRight element
emReturn element matrix
epElement pair holding more information on the pair elmX and elmY

Definition at line 57 of file bilinearForm.hh.

◆ operator()() [4/4]

template<class F = concepts::Real>
void aglowav2::AglowavBF00< F >::operator() ( const Haar3d000< F > &  elmX,
const Haar3d000< F > &  elmY,
concepts::ElementMatrix< F > &  em 
)

◆ updateChld_()

template<class F = concepts::Real>
void aglowav2::AglowavBF00< F >::updateChld_ ( const Haar3d000< F > &  tst,
const Haar3d000< F > &  trl,
F *  valtsttrl,
uint  ldtsttrl,
F *  valtrltst,
uint  ldtrltst 
)
private

Computation of one scale entries.

Call only for tst.key() >= trl.key() !!!

Parameters
tstTest element
trlTrial element
valtsttrlOne scale entries of tst with trl (row wise)
ldtsttrlLeading dimension of valtsttrl
valtrltstOne scale entries of trl with tst (row wise)
ldtrltstLeading dimension of valtrltst

◆ updateElm_()

template<class F = concepts::Real>
void aglowav2::AglowavBF00< F >::updateElm_ ( const Haar3d000< F > &  tst,
const Haar3d000< F > &  trl,
F *  valtsttrl,
uint  ldtsttrl,
F *  valtrltst,
uint  ldtrltst 
)
private

Computation of one scale entries.

On trl only the elements are considered.

Parameters
tstTest element
trlTrial element
valtsttrlOne scale entries of tst with trl elements (row wise)
ldtsttrlLeading dimension of valtsttrl
valtrltstOne scale entries of trl elements with tst (row wise)
ldtrltstLeading dimension of valtrltst

Member Data Documentation

◆ bf_

template<class F = concepts::Real>
concepts::BilinearForm<F>& aglowav2::AglowavBF00< F >::bf_
private

Bilinear form to compute the one scale matrix entries.

Definition at line 109 of file bform.hh.

◆ buf_

template<class F = concepts::Real>
Cell* aglowav2::AglowavBF00< F >::buf_
private

Buffer to store wavelet entries "for a short time".

Definition at line 132 of file bform.hh.

◆ cache_

template<class F = concepts::Real>
Cell* aglowav2::AglowavBF00< F >::cache_
private

Cache to store wavelet entries for later usage.

Definition at line 130 of file bform.hh.

◆ foo_

template<class F = concepts::Real>
F* aglowav2::AglowavBF00< F >::foo_
private

Definition at line 127 of file bform.hh.

◆ mem_

template<class F = concepts::Real>
uint aglowav2::AglowavBF00< F >::mem_
private

memory

Definition at line 115 of file bform.hh.

◆ nfld_

template<class F = concepts::Real>
Nfld** aglowav2::AglowavBF00< F >::nfld_
private

Array to hold the near field for each element.

Definition at line 113 of file bform.hh.

◆ nfldLast_

template<class F = concepts::Real>
Nfld*** aglowav2::AglowavBF00< F >::nfldLast_
private

Array to hold pointers to the last list entry of the near field (used only by nearfield_)

Definition at line 120 of file bform.hh.

◆ reuse_

template<class F = concepts::Real>
bool aglowav2::AglowavBF00< F >::reuse_
private

Entry reused and thus saved in cache_ otherwise stored in buf_.

Definition at line 134 of file bform.hh.

◆ valXY_

template<class F = concepts::Real>
F* aglowav2::AglowavBF00< F >::valXY_
private

Auxiliary vectors to hold intermediate values.

Definition at line 123 of file bform.hh.

◆ valYX_

template<class F = concepts::Real>
F* aglowav2::AglowavBF00< F >::valYX_
private

Definition at line 125 of file bform.hh.

◆ vXY_

template<class F = concepts::Real>
F* aglowav2::AglowavBF00< F >::vXY_
private

Definition at line 124 of file bform.hh.

◆ vYX_

template<class F = concepts::Real>
F* aglowav2::AglowavBF00< F >::vYX_
private

Definition at line 126 of file bform.hh.

◆ XFY_

template<class F = concepts::Real>
F00<F>& aglowav2::AglowavBF00< F >::XFY_
private

Cluster approximation of the far field.

Definition at line 111 of file bform.hh.


The documentation for this class was generated from the following file:
Page URL: http://wiki.math.ethz.ch/bin/view/Concepts/WebHome
21 August 2020
© 2020 Eidgenössische Technische Hochschule Zürich