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...
Data Structures | |
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... | |
Public Member Functions | |
ClassName ("lduMatrix") | |
lduMatrix (const lduMesh &) | |
Construct given an LDU addressed mesh. More... | |
lduMatrix (const lduMatrix &) | |
Construct as copy. More... | |
lduMatrix (lduMatrix &, bool reUse) | |
Construct as copy or re-use as specified. More... | |
lduMatrix (const lduMesh &, Istream &) | |
Construct given an LDU addressed mesh and an Istream. More... | |
~lduMatrix () | |
Destructor. More... | |
const lduMesh & | mesh () const |
Return the LDU mesh from which the addressing is obtained. More... | |
const lduAddressing & | lduAddr () const |
Return the LDU addressing. More... | |
const lduSchedule & | patchSchedule () const |
Return the patch evaluation schedule. More... | |
scalarField & | lower () |
scalarField & | diag () |
scalarField & | upper () |
scalarField & | lower (const label size) |
scalarField & | diag (const label nCoeffs) |
scalarField & | upper (const label nCoeffs) |
const scalarField & | lower () const |
const scalarField & | diag () const |
const scalarField & | upper () 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. More... | |
void | Tmul (scalarField &, const tmp< scalarField > &, const FieldField< Field, scalar > &, const lduInterfaceFieldPtrsList &, const direction cmpt) const |
Matrix transpose multiplication with updated interfaces. More... | |
void | sumA (scalarField &, const FieldField< Field, scalar > &, const lduInterfaceFieldPtrsList &) const |
Sum the coefficients on each row of the matrix. More... | |
void | residual (scalarField &rA, const scalarField &psi, const scalarField &source, const FieldField< Field, scalar > &interfaceBouCoeffs, const lduInterfaceFieldPtrsList &interfaces, const direction cmpt) const |
tmp< scalarField > | residual (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. More... | |
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. More... | |
template<class Type > | |
tmp< Field< Type > > | H (const Field< Type > &) const |
template<class Type > | |
tmp< Field< Type > > | H (const tmp< Field< Type > > &) const |
tmp< scalarField > | H1 () const |
template<class Type > | |
tmp< Field< Type > > | faceH (const Field< Type > &) const |
template<class Type > | |
tmp< Field< Type > > | faceH (const tmp< Field< Type > > &) const |
InfoProxy< lduMatrix > | info () const |
Return info proxy. More... | |
void | operator= (const lduMatrix &) |
void | negate () |
void | operator+= (const lduMatrix &) |
void | operator-= (const lduMatrix &) |
void | operator*= (const scalarField &) |
void | operator*= (scalar) |
template<class Type > | |
Foam::tmp< Foam::Field< Type > > | H (const Field< Type > &psi) const |
template<class Type > | |
Foam::tmp< Foam::Field< Type > > | H (const tmp< Field< Type > > &tpsi) const |
template<class Type > | |
Foam::tmp< Foam::Field< Type > > | faceH (const Field< Type > &psi) const |
template<class Type > | |
Foam::tmp< Foam::Field< Type > > | faceH (const tmp< Field< Type > > &tpsi) const |
Private Attributes | |
const lduMesh & | lduMesh_ |
LDU mesh reference. More... | |
scalarField * | lowerPtr_ |
Coefficients (not including interfaces) More... | |
scalarField * | diagPtr_ |
scalarField * | upperPtr_ |
Friends | |
Ostream & | operator<< (Ostream &, const lduMatrix &) |
Ostream & | operator<< (Ostream &, const InfoProxy< lduMatrix > &) |
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.
Definition at line 77 of file lduMatrix.H.
Construct given an LDU addressed mesh.
The coefficients are initially empty for subsequent setting.
Definition at line 40 of file lduMatrix.C.
Construct as copy.
Definition at line 49 of file lduMatrix.C.
References A(), lduMatrix::diagPtr_, lduMatrix::lowerPtr_, and lduMatrix::upperPtr_.
Construct as copy or re-use as specified.
Definition at line 73 of file lduMatrix.C.
References A(), lduMatrix::diagPtr_, lduMatrix::lowerPtr_, and lduMatrix::upperPtr_.
Construct given an LDU addressed mesh and an Istream.
from which the coefficients are read
Definition at line 120 of file lduMatrix.C.
References lduMatrix::diagPtr_, lduMatrix::hasDiag(), lduMatrix::lowerPtr_, and lduMatrix::upperPtr_.
~lduMatrix | ( | ) |
Destructor.
Definition at line 146 of file lduMatrix.C.
ClassName | ( | "lduMatrix" | ) |
|
inline |
Return the LDU mesh from which the addressing is obtained.
Definition at line 539 of file lduMatrix.H.
References lduMatrix::lduMesh_.
Referenced by GAMGSolver::agglomerateMatrix(), algebraicPairGAMGAgglomeration::algebraicPairGAMGAgglomeration(), GAMGAgglomeration::New(), and GAMGSolver::procAgglomerateMatrix().
|
inline |
Return the LDU addressing.
Definition at line 545 of file lduMatrix.H.
References lduMesh::lduAddr(), and lduMatrix::lduMesh_.
Referenced by DICPreconditioner::calcReciprocalD(), DILUPreconditioner::calcReciprocalD(), LUscalarMatrix::convert(), lduMatrix::H(), lduMatrix::H1(), GAMGSolver::interpolate(), LUscalarMatrix::LUscalarMatrix(), lduMatrix::patchSchedule(), symGaussSeidelSmoother::smooth(), GaussSeidelSmoother::smooth(), nonBlockingGaussSeidelSmoother::smooth(), and lduMatrix::sumDiag().
|
inline |
Return the patch evaluation schedule.
Definition at line 551 of file lduMatrix.H.
References lduMatrix::lduAddr(), and lduAddressing::patchSchedule().
const Foam::scalarField & lower | ( | ) |
Definition at line 165 of file lduMatrix.C.
Referenced by GAMGSolver::agglomerateMatrix(), DILUPreconditioner::calcReciprocalD(), LUscalarMatrix::convert(), lduMatrix::faceH(), GAMGSolver::GAMGSolver(), lduMatrix::H(), lduMatrix::H1(), GAMGSolver::interpolate(), lduMatrix::negSumDiag(), Foam::operator<<(), GAMGSolver::procAgglomerateMatrix(), GaussSeidelSmoother::smooth(), symGaussSeidelSmoother::smooth(), nonBlockingGaussSeidelSmoother::smooth(), lduMatrix::sumDiag(), and lduMatrix::sumMagOffDiag().
const Foam::scalarField & diag | ( | ) |
Definition at line 183 of file lduMatrix.C.
Referenced by GAMGSolver::agglomerateMatrix(), LUscalarMatrix::convert(), GAMGSolver::GAMGSolver(), GAMGSolver::initVcycle(), GAMGSolver::interpolate(), nonBlockingGaussSeidelSmoother::nonBlockingGaussSeidelSmoother(), Foam::operator<<(), GAMGSolver::procAgglomerateMatrix(), GaussSeidelSmoother::smooth(), symGaussSeidelSmoother::smooth(), nonBlockingGaussSeidelSmoother::smooth(), and lduMatrix::sumDiag().
const Foam::scalarField & upper | ( | ) |
Definition at line 194 of file lduMatrix.C.
Referenced by GAMGSolver::agglomerateMatrix(), algebraicPairGAMGAgglomeration::algebraicPairGAMGAgglomeration(), DILUPreconditioner::calcReciprocalD(), DICPreconditioner::calcReciprocalD(), LUscalarMatrix::convert(), lduMatrix::faceH(), lduMatrix::H(), lduMatrix::H1(), GAMGSolver::interpolate(), lduMatrix::negSumDiag(), Foam::operator<<(), GAMGSolver::procAgglomerateMatrix(), symGaussSeidelSmoother::smooth(), GaussSeidelSmoother::smooth(), nonBlockingGaussSeidelSmoother::smooth(), lduMatrix::sumDiag(), and lduMatrix::sumMagOffDiag().
Foam::scalarField & lower | ( | const label | size | ) |
Definition at line 212 of file lduMatrix.C.
Foam::scalarField & diag | ( | const label | nCoeffs | ) |
Definition at line 230 of file lduMatrix.C.
Foam::scalarField & upper | ( | const label | nCoeffs | ) |
Definition at line 241 of file lduMatrix.C.
const scalarField& lower | ( | ) | const |
const scalarField& diag | ( | ) | const |
const scalarField& upper | ( | ) | const |
|
inline |
Definition at line 575 of file lduMatrix.H.
References lduMatrix::diagPtr_.
Referenced by lduMatrix::lduMatrix(), Foam::operator<<(), and GAMGSolver::procAgglomerateMatrix().
|
inline |
Definition at line 580 of file lduMatrix.H.
References lduMatrix::upperPtr_.
Referenced by Foam::operator<<(), and GAMGSolver::procAgglomerateMatrix().
|
inline |
Definition at line 585 of file lduMatrix.H.
References lduMatrix::lowerPtr_.
Referenced by GAMGSolver::agglomerateMatrix(), Foam::operator<<(), and GAMGSolver::procAgglomerateMatrix().
|
inline |
Definition at line 590 of file lduMatrix.H.
References lduMatrix::diagPtr_, lduMatrix::lowerPtr_, and lduMatrix::upperPtr_.
Referenced by lduMatrix::solver::New().
|
inline |
Definition at line 595 of file lduMatrix.H.
References lduMatrix::diagPtr_, lduMatrix::lowerPtr_, and lduMatrix::upperPtr_.
Referenced by lduMatrix::solver::New(), lduMatrix::smoother::New(), and lduMatrix::preconditioner::New().
|
inline |
Definition at line 600 of file lduMatrix.H.
References lduMatrix::diagPtr_, lduMatrix::lowerPtr_, and lduMatrix::upperPtr_.
Referenced by lduMatrix::solver::New(), lduMatrix::smoother::New(), and lduMatrix::preconditioner::New().
void sumDiag | ( | ) |
Definition at line 33 of file lduMatrixOperations.C.
References lduMatrix::diag(), lduMatrix::lduAddr(), lduMatrix::lower(), lduAddressing::lowerAddr(), UList::size(), lduMatrix::upper(), and lduAddressing::upperAddr().
void negSumDiag | ( | ) |
Definition at line 50 of file lduMatrixOperations.C.
References Foam::diag(), lduMatrix::lower(), UList::size(), and lduMatrix::upper().
void sumMagOffDiag | ( | scalarField & | sumOff | ) | const |
Definition at line 68 of file lduMatrixOperations.C.
References lduMatrix::lower(), Foam::mag(), UList::size(), and lduMatrix::upper().
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 tmp::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 tmp::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 Foam::diag(), forAll, and UPtrList::set().
void residual | ( | scalarField & | rA, |
const scalarField & | psi, | ||
const scalarField & | source, | ||
const FieldField< Field, scalar > & | interfaceBouCoeffs, | ||
const lduInterfaceFieldPtrsList & | interfaces, | ||
const direction | cmpt | ||
) | const |
Definition at line 204 of file lduMatrixATmul.C.
References Foam::diag(), forAll, patchi, psi, and UPtrList::set().
Referenced by fvMatrix< Type >::residual().
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 psi, and scalarField().
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, FatalErrorInFunction, forAll, UPtrList::set(), List::size(), and UPtrList::size().
Referenced by GAMGSolver::interpolate(), GaussSeidelSmoother::smooth(), symGaussSeidelSmoother::smooth(), and nonBlockingGaussSeidelSmoother::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 97 of file lduMatrixUpdateMatrixInterfaces.C.
References Foam::exit(), Foam::FatalError, FatalErrorInFunction, forAll, UPtrList::set(), List::size(), and UPtrList::size().
Referenced by GAMGSolver::interpolate(), GaussSeidelSmoother::smooth(), symGaussSeidelSmoother::smooth(), and nonBlockingGaussSeidelSmoother::smooth().
Foam::tmp< Foam::scalarField > H1 | ( | ) | const |
Definition at line 298 of file lduMatrixATmul.C.
References UList::begin(), lduMatrix::lduAddr(), lduMatrix::lower(), lduAddressing::lowerAddr(), lduMatrix::lowerPtr_, lduMatrix::upper(), lduAddressing::upperAddr(), and lduMatrix::upperPtr_.
Return info proxy.
Used to print matrix information to a stream
Definition at line 704 of file lduMatrix.H.
void operator= | ( | const lduMatrix & | A | ) |
Definition at line 88 of file lduMatrixOperations.C.
References A(), Foam::abort(), Foam::diag(), and Foam::FatalError.
void negate | ( | ) |
Definition at line 125 of file lduMatrixOperations.C.
void operator+= | ( | const lduMatrix & | A | ) |
Definition at line 144 of file lduMatrixOperations.C.
References A(), Foam::diag(), Foam::endl(), Foam::nl, and WarningInFunction.
void operator-= | ( | const lduMatrix & | A | ) |
Definition at line 223 of file lduMatrixOperations.C.
References A(), Foam::diag(), Foam::endl(), Foam::nl, and WarningInFunction.
void operator*= | ( | const scalarField & | sf | ) |
Definition at line 302 of file lduMatrixOperations.C.
References sf().
void operator*= | ( | scalar | s | ) |
Definition at line 332 of file lduMatrixOperations.C.
References s().
Foam::tmp<Foam::Field<Type> > H | ( | const Field< Type > & | psi | ) | const |
Definition at line 34 of file lduMatrixTemplates.C.
References UList::begin(), lduMatrix::lduAddr(), lduMatrix::lower(), lduAddressing::lowerAddr(), lduMatrix::lowerPtr_, psi, lduMatrix::upper(), lduAddressing::upperAddr(), and lduMatrix::upperPtr_.
Foam::tmp<Foam::Field<Type> > H | ( | const tmp< Field< Type > > & | tpsi | ) | const |
Definition at line 69 of file lduMatrixTemplates.C.
Foam::tmp<Foam::Field<Type> > faceH | ( | const Field< Type > & | psi | ) | const |
Definition at line 79 of file lduMatrixTemplates.C.
References Foam::exit(), Foam::FatalError, FatalErrorInFunction, lduMatrix::lower(), psi, UList::size(), and lduMatrix::upper().
Foam::tmp<Foam::Field<Type> > faceH | ( | const tmp< Field< Type > > & | tpsi | ) | const |
Definition at line 115 of file lduMatrixTemplates.C.
|
private |
LDU mesh reference.
Definition at line 82 of file lduMatrix.H.
Referenced by lduMatrix::lduAddr(), and lduMatrix::mesh().
|
private |
Coefficients (not including interfaces)
Definition at line 85 of file lduMatrix.H.
Referenced by lduMatrix::asymmetric(), lduMatrix::diagonal(), lduMatrix::H(), lduMatrix::H1(), lduMatrix::hasLower(), lduMatrix::lduMatrix(), and lduMatrix::symmetric().
|
private |
Definition at line 85 of file lduMatrix.H.
Referenced by lduMatrix::asymmetric(), lduMatrix::diagonal(), lduMatrix::hasDiag(), lduMatrix::lduMatrix(), and lduMatrix::symmetric().
|
private |
Definition at line 85 of file lduMatrix.H.
Referenced by lduMatrix::asymmetric(), lduMatrix::diagonal(), lduMatrix::H(), lduMatrix::H1(), lduMatrix::hasUpper(), lduMatrix::lduMatrix(), and lduMatrix::symmetric().
Copyright © 2011-2018 OpenFOAM | OPENFOAM® is a registered trademark of OpenCFD Ltd.