FreeFOAM The Cross-Platform CFD Toolkit
lduMatrix Class Reference

lduMatrix is a general matrix class in which the coefficients are stored as three arrays, one for the upper triangle, one for the lower triangle and a third for the diagonal. More...

#include <OpenFOAM/lduMatrix.H>


Detailed Description

lduMatrix is a general matrix class in which the coefficients are stored as three arrays, one for the upper triangle, one for the lower triangle and a third for the diagonal.

Addressing arrays must be supplied for the upper and lower triangles.

It might be better if this class were organised as a hierachy starting from an empty matrix, then deriving diagonal, symmetric and asymmetric matrices.

Source files

Definition at line 76 of file lduMatrix.H.

+ Inheritance diagram for lduMatrix:

List of all members.

Classes

class  preconditioner
 Abstract base-class for lduMatrix preconditioners. More...
class  smoother
 Abstract base-class for lduMatrix smoothers. More...
class  solver
 Abstract base-class for lduMatrix solvers. More...
class  solverPerformance
 Class returned by the solver, containing performance statistics. More...

Public Member Functions

 ClassName ("lduMatrix")
 lduMatrix (const lduMesh &)
 Construct given an LDU addressed mesh.
 lduMatrix (const lduMatrix &)
 Construct as copy.
 lduMatrix (lduMatrix &, bool reUse)
 Construct as copy or re-use as specified.
 lduMatrix (const lduMesh &, Istream &)
 Construct given an LDU addressed mesh and an Istream.
 ~lduMatrix ()
const lduMeshmesh () const
 Return the LDU mesh from which the addressing is obtained.
const lduAddressinglduAddr () const
 Return the LDU addressing.
const lduSchedulepatchSchedule () const
 Return the patch evaluation schedule.
scalarFieldlower ()
scalarFielddiag ()
scalarFieldupper ()
const scalarFieldlower () const
const scalarFielddiag () const
const scalarFieldupper () const
bool hasDiag () const
bool hasUpper () const
bool hasLower () const
bool diagonal () const
bool symmetric () const
bool asymmetric () const
void sumDiag ()
void negSumDiag ()
void sumMagOffDiag (scalarField &sumOff) const
void Amul (scalarField &, const tmp< scalarField > &, const FieldField< Field, scalar > &, const lduInterfaceFieldPtrsList &, const direction cmpt) const
 Matrix multiplication with updated interfaces.
void Tmul (scalarField &, const tmp< scalarField > &, const FieldField< Field, scalar > &, const lduInterfaceFieldPtrsList &, const direction cmpt) const
 Matrix transpose multiplication with updated interfaces.
void sumA (scalarField &, const FieldField< Field, scalar > &, const lduInterfaceFieldPtrsList &) const
 Sum the coefficients on each row of the matrix.
void residual (scalarField &rA, const scalarField &psi, const scalarField &source, const FieldField< Field, scalar > &interfaceBouCoeffs, const lduInterfaceFieldPtrsList &interfaces, const direction cmpt) const
tmp< scalarFieldresidual (const scalarField &psi, const scalarField &source, const FieldField< Field, scalar > &interfaceBouCoeffs, const lduInterfaceFieldPtrsList &interfaces, const direction cmpt) const
void initMatrixInterfaces (const FieldField< Field, scalar > &interfaceCoeffs, const lduInterfaceFieldPtrsList &interfaces, const scalarField &psiif, scalarField &result, const direction cmpt) const
 Initialise the update of interfaced interfaces.
void updateMatrixInterfaces (const FieldField< Field, scalar > &interfaceCoeffs, const lduInterfaceFieldPtrsList &interfaces, const scalarField &psiif, scalarField &result, const direction cmpt) const
 Update interfaced interfaces for matrix operations.
template<class Type >
tmp< Field< Type > > H (const Field< Type > &) const
template<class Type >
tmp< Field< Type > > H (const tmp< Field< Type > > &) const
tmp< scalarFieldH1 () const
template<class Type >
tmp< Field< Type > > faceH (const Field< Type > &) const
template<class Type >
tmp< Field< Type > > faceH (const tmp< Field< Type > > &) const
void operator= (const lduMatrix &)
void negate ()
void operator+= (const lduMatrix &)
void operator-= (const lduMatrix &)
void operator*= (const scalarField &)
void operator*= (scalar)

Static Public Attributes

static const scalar great_ = 1.0e+20
 Large scalar for the use in solvers.
static const scalar small_ = 1.0e-20
 Small scalar for the use in solvers.

Friends

Ostreamoperator<< (Ostream &, const lduMatrix &)

Constructor & Destructor Documentation

lduMatrix ( const lduMesh mesh)

Construct given an LDU addressed mesh.

The coefficients are initially empty for subsequent setting.

Definition at line 42 of file lduMatrix.C.

lduMatrix ( const lduMatrix A)

Construct as copy.

Definition at line 51 of file lduMatrix.C.

lduMatrix ( lduMatrix A,
bool  reUse 
)

Construct as copy or re-use as specified.

Definition at line 75 of file lduMatrix.C.

lduMatrix ( const lduMesh mesh,
Istream is 
)

Construct given an LDU addressed mesh and an Istream.

from which the coefficients are read

Definition at line 123 of file lduMatrix.C.

~lduMatrix ( )

Definition at line 135 of file lduMatrix.C.


Member Function Documentation

ClassName ( "lduMatrix"  )
const lduMesh& mesh ( ) const
inline

Return the LDU mesh from which the addressing is obtained.

Definition at line 671 of file lduMatrix.H.

Referenced by algebraicPairGAMGAgglomeration::algebraicPairGAMGAgglomeration(), and GAMGAgglomeration::New().

const lduSchedule& patchSchedule ( ) const
inline

Return the patch evaluation schedule.

Definition at line 683 of file lduMatrix.H.

References lduMatrix::lduAddr(), and lduAddressing::patchSchedule().

const Foam::scalarField & lower ( ) const

Definition at line 201 of file lduMatrix.C.

References Foam::abort(), Foam::FatalError, and FatalErrorIn.

const Foam::scalarField & diag ( ) const

Definition at line 221 of file lduMatrix.C.

References Foam::abort(), Foam::FatalError, and FatalErrorIn.

const Foam::scalarField & upper ( ) const

Definition at line 234 of file lduMatrix.C.

References Foam::abort(), Foam::FatalError, and FatalErrorIn.

bool hasDiag ( ) const
inline

Definition at line 699 of file lduMatrix.H.

bool hasUpper ( ) const
inline

Definition at line 704 of file lduMatrix.H.

Referenced by Foam::correction().

bool hasLower ( ) const
inline

Definition at line 709 of file lduMatrix.H.

Referenced by Foam::correction().

bool diagonal ( ) const
inline
void negSumDiag ( )
void sumMagOffDiag ( scalarField sumOff) const
void Amul ( scalarField Apsi,
const tmp< scalarField > &  tpsi,
const FieldField< Field, scalar > &  interfaceBouCoeffs,
const lduInterfaceFieldPtrsList interfaces,
const direction  cmpt 
) const

Matrix multiplication with updated interfaces.

Definition at line 35 of file lduMatrixATmul.C.

References UList< T >::begin(), tmp< T >::clear(), Foam::diag(), and psi.

void Tmul ( scalarField Tpsi,
const tmp< scalarField > &  tpsi,
const FieldField< Field, scalar > &  interfaceIntCoeffs,
const lduInterfaceFieldPtrsList interfaces,
const direction  cmpt 
) const

Matrix transpose multiplication with updated interfaces.

Definition at line 96 of file lduMatrixATmul.C.

References UList< T >::begin(), tmp< T >::clear(), Foam::diag(), and psi.

void sumA ( scalarField sumA,
const FieldField< Field, scalar > &  interfaceBouCoeffs,
const lduInterfaceFieldPtrsList interfaces 
) const

Sum the coefficients on each row of the matrix.

Definition at line 155 of file lduMatrixATmul.C.

References UList< T >::begin(), Foam::diag(), forAll, and UPtrList< T >::set().

void residual ( scalarField rA,
const scalarField psi,
const scalarField source,
const FieldField< Field, scalar > &  interfaceBouCoeffs,
const lduInterfaceFieldPtrsList interfaces,
const direction  cmpt 
) const
Foam::tmp< Foam::scalarField > residual ( const scalarField psi,
const scalarField source,
const FieldField< Field, scalar > &  interfaceBouCoeffs,
const lduInterfaceFieldPtrsList interfaces,
const direction  cmpt 
) const

Definition at line 284 of file lduMatrixATmul.C.

References scalarField(), and List< T >::size().

void initMatrixInterfaces ( const FieldField< Field, scalar > &  interfaceCoeffs,
const lduInterfaceFieldPtrsList interfaces,
const scalarField psiif,
scalarField result,
const direction  cmpt 
) const

Initialise the update of interfaced interfaces.

for matrix operations

Definition at line 31 of file lduMatrixUpdateMatrixInterfaces.C.

References Foam::exit(), Foam::FatalError, FatalErrorIn, forAll, UPtrList< T >::set(), List< T >::size(), and UPtrList< T >::size().

Referenced by GaussSeidelSmoother::smooth().

void updateMatrixInterfaces ( const FieldField< Field, scalar > &  interfaceCoeffs,
const lduInterfaceFieldPtrsList interfaces,
const scalarField psiif,
scalarField result,
const direction  cmpt 
) const

Update interfaced interfaces for matrix operations.

Definition at line 99 of file lduMatrixUpdateMatrixInterfaces.C.

References Foam::exit(), Foam::FatalError, FatalErrorIn, forAll, UPtrList< T >::set(), List< T >::size(), and UPtrList< T >::size().

Referenced by GaussSeidelSmoother::smooth().

Foam::tmp< Foam::Field< Type > > H ( const tmp< Field< Type > > &  tpsi) const

Definition at line 69 of file lduMatrixTemplates.C.

Foam::tmp< Foam::scalarField > H1 ( ) const

Reimplemented in fvMatrix< Type >, fvMatrix< Type >, and fvMatrix< Type >.

Definition at line 332 of file lduMatrixOperations.C.

References UList< T >::begin().

Referenced by fvMatrix< Type >::H1().

Foam::tmp< Foam::Field< Type > > faceH ( const Field< Type > &  psi) const
Foam::tmp< Foam::Field< Type > > faceH ( const tmp< Field< Type > > &  tpsi) const

Definition at line 113 of file lduMatrixTemplates.C.

void operator= ( const lduMatrix A)
void negate ( )

Reimplemented in fvMatrix< Type >.

Definition at line 125 of file lduMatrixOperations.C.

Referenced by fvMatrix< Type >::negate().

void operator*= ( const scalarField sf)

Definition at line 280 of file lduMatrixOperations.C.

References sf(), and List< T >::size().

Referenced by fvMatrix< Type >::operator*=().

void operator*= ( scalar  s)

Definition at line 313 of file lduMatrixOperations.C.


Friends And Related Function Documentation

Ostream& operator<< ( Ostream ,
const lduMatrix  
)
friend

Member Data Documentation

const Foam::scalar great_ = 1.0e+20
static

Large scalar for the use in solvers.

Definition at line 638 of file lduMatrix.H.

const Foam::scalar small_ = 1.0e-20
static

Small scalar for the use in solvers.

Definition at line 641 of file lduMatrix.H.


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