Data Structures | Public Types | Public Member Functions | Protected Member Functions | Private Attributes | Friends
fvMatrix Class Reference

A special matrix type and solver, designed for finite volume solutions of scalar equations. Face addressing is used to make all matrix assembly and solution loops vectorise. More...

Inheritance diagram for fvMatrix:
Inheritance graph
[legend]
Collaboration diagram for fvMatrix:
Collaboration graph
[legend]

Data Structures

class  fvSolver
 Solver class returned by the solver function. More...
 

Public Types

typedef GeometricField< Type, fvsPatchField, surfaceMesh > * surfaceTypeFieldPtr
 Declare return type of the faceFluxCorrectionPtr() function. More...
 

Public Member Functions

 ClassName ("fvMatrix")
 
 fvMatrix (const GeometricField< Type, fvPatchField, volMesh > &, const dimensionSet &)
 Construct given a field to solve for. More...
 
 fvMatrix (const fvMatrix< Type > &)
 Construct as copy. More...
 
 fvMatrix (const tmp< fvMatrix< Type > > &)
 Construct as copy of tmp<fvMatrix<Type> > deleting argument. More...
 
 fvMatrix (const GeometricField< Type, fvPatchField, volMesh > &, Istream &)
 Construct from Istream given field to solve for. More...
 
virtual ~fvMatrix ()
 Destructor. More...
 
const GeometricField< Type, fvPatchField, volMesh > & psi () const
 
const dimensionSetdimensions () const
 
Field< Type > & source ()
 
const Field< Type > & source () const
 
FieldField< Field, Type > & internalCoeffs ()
 fvBoundary scalar field containing pseudo-matrix coeffs More...
 
FieldField< Field, Type > & boundaryCoeffs ()
 fvBoundary scalar field containing pseudo-matrix coeffs More...
 
surfaceTypeFieldPtrfaceFluxCorrectionPtr ()
 Return pointer to face-flux non-orthogonal correction field. More...
 
void setValues (const labelUList &cells, const UList< Type > &values)
 Set solution in given cells to the specified values. More...
 
void setValues (const labelUList &cells, const UIndirectList< Type > &values)
 Set solution in given cells to the specified values. More...
 
void setReference (const label celli, const Type &value, const bool forceReference=false)
 Set reference level for solution. More...
 
void setComponentReference (const label patchi, const label facei, const direction cmpt, const scalar value)
 Set reference level for a component of the solution. More...
 
void relax (const scalar alpha)
 Relax matrix (for steady-state solution). More...
 
void relax ()
 Relax matrix (for steady-state solution). More...
 
void boundaryManipulate (typename GeometricField< Type, fvPatchField, volMesh >::GeometricBoundaryField &values)
 Manipulate based on a boundary field. More...
 
autoPtr< fvSolversolver (const dictionary &)
 Construct and return the solver. More...
 
autoPtr< fvSolversolver ()
 Construct and return the solver. More...
 
SolverPerformance< Type > solve (const dictionary &)
 Solve segregated or coupled returning the solution statistics. More...
 
SolverPerformance< Type > solveSegregated (const dictionary &)
 Solve segregated returning the solution statistics. More...
 
SolverPerformance< Type > solveCoupled (const dictionary &)
 Solve coupled returning the solution statistics. More...
 
SolverPerformance< Type > solve ()
 Solve returning the solution statistics. More...
 
tmp< Field< Type > > residual () const
 Return the matrix residual. More...
 
tmp< scalarFieldD () const
 Return the matrix scalar diagonal. More...
 
tmp< Field< Type > > DD () const
 Return the matrix Type diagonal. More...
 
tmp< volScalarFieldA () const
 Return the central coefficient. More...
 
tmp< GeometricField< Type, fvPatchField, volMesh > > H () const
 Return the H operation source. More...
 
tmp< volScalarFieldH1 () const
 Return H(1) More...
 
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > flux () const
 Return the face-flux field from the matrix. More...
 
void operator= (const fvMatrix< Type > &)
 
void operator= (const tmp< fvMatrix< Type > > &)
 
void negate ()
 
void operator+= (const fvMatrix< Type > &)
 
void operator+= (const tmp< fvMatrix< Type > > &)
 
void operator-= (const fvMatrix< Type > &)
 
void operator-= (const tmp< fvMatrix< Type > > &)
 
void operator+= (const DimensionedField< Type, volMesh > &)
 
void operator+= (const tmp< DimensionedField< Type, volMesh > > &)
 
void operator+= (const tmp< GeometricField< Type, fvPatchField, volMesh > > &)
 
void operator-= (const DimensionedField< Type, volMesh > &)
 
void operator-= (const tmp< DimensionedField< Type, volMesh > > &)
 
void operator-= (const tmp< GeometricField< Type, fvPatchField, volMesh > > &)
 
void operator+= (const dimensioned< Type > &)
 
void operator-= (const dimensioned< Type > &)
 
void operator+= (const zero &)
 
void operator-= (const zero &)
 
void operator*= (const DimensionedField< scalar, volMesh > &)
 
void operator*= (const tmp< DimensionedField< scalar, volMesh > > &)
 
void operator*= (const tmp< volScalarField > &)
 
void operator*= (const dimensioned< scalar > &)
 
void setComponentReference (const label patchi, const label facei, const direction, const scalar value)
 
Foam::autoPtr< Foam::fvMatrix< Foam::scalar >::fvSolversolver (const dictionary &solverControls)
 
Foam::solverPerformance solveSegregated (const dictionary &solverControls)
 
Foam::tmp< Foam::scalarFieldresidual () const
 
Foam::tmp< Foam::volScalarFieldH () const
 
Foam::tmp< Foam::volScalarFieldH1 () const
 
void setComponentReference (const label patchi, const label facei, const direction, const scalar value)
 
autoPtr< fvMatrix< scalar >::fvSolversolver (const dictionary &)
 
solverPerformance solveSegregated (const dictionary &)
 
tmp< scalarFieldresidual () const
 
tmp< volScalarFieldH () const
 
tmp< volScalarFieldH1 () const
 

Protected Member Functions

template<class Type2 >
void addToInternalField (const labelUList &addr, const Field< Type2 > &pf, Field< Type2 > &intf) const
 Add patch contribution to internal field. More...
 
template<class Type2 >
void addToInternalField (const labelUList &addr, const tmp< Field< Type2 > > &tpf, Field< Type2 > &intf) const
 
template<class Type2 >
void subtractFromInternalField (const labelUList &addr, const Field< Type2 > &pf, Field< Type2 > &intf) const
 Subtract patch contribution from internal field. More...
 
template<class Type2 >
void subtractFromInternalField (const labelUList &addr, const tmp< Field< Type2 > > &tpf, Field< Type2 > &intf) const
 
void addBoundaryDiag (scalarField &diag, const direction cmpt) const
 
void addCmptAvBoundaryDiag (scalarField &diag) const
 
void addBoundarySource (Field< Type > &source, const bool couples=true) const
 
template<template< class > class ListType>
void setValuesFromList (const labelUList &cells, const ListType< Type > &values)
 Set solution in given cells to the specified values. More...
 

Private Attributes

const GeometricField< Type, fvPatchField, volMesh > & psi_
 Const reference to GeometricField<Type, fvPatchField, volMesh> More...
 
dimensionSet dimensions_
 Dimension set. More...
 
Field< Type > source_
 Source term. More...
 
FieldField< Field, Type > internalCoeffs_
 Boundary scalar field containing pseudo-matrix coeffs. More...
 
FieldField< Field, Type > boundaryCoeffs_
 Boundary scalar field containing pseudo-matrix coeffs. More...
 
GeometricField< Type, fvsPatchField, surfaceMesh > * faceFluxCorrectionPtr_
 Face flux field for non-orthogonal correction. More...
 

Friends

class fvSolver
 Declare friendship with the fvSolver class. More...
 
tmp< GeometricField< Type, fvPatchField, volMesh > > operator& (const fvMatrix< Type > &, const DimensionedField< Type, volMesh > &)
 
tmp< GeometricField< Type, fvPatchField, volMesh > > operator& (const fvMatrix< Type > &, const tmp< GeometricField< Type, fvPatchField, volMesh > > &)
 
tmp< GeometricField< Type, fvPatchField, volMesh > > operator& (const tmp< fvMatrix< Type > > &, const DimensionedField< Type, volMesh > &)
 
tmp< GeometricField< Type, fvPatchField, volMesh > > operator& (const tmp< fvMatrix< Type > > &, const tmp< GeometricField< Type, fvPatchField, volMesh > > &)
 
Ostreamoperator (Ostream &, const fvMatrix< Type > &)
 

Detailed Description

A special matrix type and solver, designed for finite volume solutions of scalar equations. Face addressing is used to make all matrix assembly and solution loops vectorise.

Source files

Definition at line 68 of file fvPatchField.H.

Member Typedef Documentation

◆ surfaceTypeFieldPtr

Declare return type of the faceFluxCorrectionPtr() function.

Definition at line 318 of file fvMatrix.H.

Constructor & Destructor Documentation

◆ fvMatrix() [1/4]

fvMatrix ( const GeometricField< Type, fvPatchField, volMesh > &  ,
const dimensionSet  
)

Construct given a field to solve for.

◆ fvMatrix() [2/4]

fvMatrix ( const fvMatrix< Type > &  )

Construct as copy.

◆ fvMatrix() [3/4]

fvMatrix ( const tmp< fvMatrix< Type > > &  )

Construct as copy of tmp<fvMatrix<Type> > deleting argument.

◆ fvMatrix() [4/4]

fvMatrix ( const GeometricField< Type, fvPatchField, volMesh > &  ,
Istream  
)

Construct from Istream given field to solve for.

◆ ~fvMatrix()

~fvMatrix ( )
virtual

Destructor.

Definition at line 462 of file fvMatrix.C.

Member Function Documentation

◆ addToInternalField() [1/2]

void addToInternalField ( const labelUList addr,
const Field< Type2 > &  pf,
Field< Type2 > &  intf 
) const
protected

Add patch contribution to internal field.

Definition at line 38 of file fvMatrix.C.

◆ addToInternalField() [2/2]

void addToInternalField ( const labelUList addr,
const tmp< Field< Type2 > > &  tpf,
Field< Type2 > &  intf 
) const
protected

Definition at line 61 of file fvMatrix.C.

◆ subtractFromInternalField() [1/2]

void subtractFromInternalField ( const labelUList addr,
const Field< Type2 > &  pf,
Field< Type2 > &  intf 
) const
protected

Subtract patch contribution from internal field.

Definition at line 75 of file fvMatrix.C.

◆ subtractFromInternalField() [2/2]

void subtractFromInternalField ( const labelUList addr,
const tmp< Field< Type2 > > &  tpf,
Field< Type2 > &  intf 
) const
protected

Definition at line 98 of file fvMatrix.C.

◆ addBoundaryDiag()

void addBoundaryDiag ( scalarField diag,
const direction  cmpt 
) const
protected

Definition at line 111 of file fvMatrix.C.

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

Here is the caller graph for this function:

◆ addCmptAvBoundaryDiag()

void addCmptAvBoundaryDiag ( scalarField diag) const
protected

Definition at line 129 of file fvMatrix.C.

◆ addBoundarySource()

void addBoundarySource ( Field< Type > &  source,
const bool  couples = true 
) const
protected

Definition at line 145 of file fvMatrix.C.

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

Here is the caller graph for this function:

◆ setValuesFromList()

void setValuesFromList ( const labelUList cells,
const ListType< Type > &  values 
)
protected

Set solution in given cells to the specified values.

Definition at line 178 of file fvMatrix.C.

◆ ClassName()

ClassName ( "fvMatrix"  )

◆ psi()

const GeometricField<Type, fvPatchField, volMesh>& psi ( ) const
inline

◆ dimensions()

const dimensionSet& dimensions ( ) const
inline

◆ source() [1/2]

Field<Type>& source ( )
inline

◆ source() [2/2]

const Field<Type>& source ( ) const
inline

Definition at line 296 of file fvMatrix.H.

◆ internalCoeffs()

FieldField<Field, Type>& internalCoeffs ( )
inline

fvBoundary scalar field containing pseudo-matrix coeffs

for internal cells

Definition at line 303 of file fvMatrix.H.

Referenced by immersedBoundaryFvPatchField< Type >::correctOffDiag(), gaussConvectionScheme< Type >::fvmDiv(), and gaussLaplacianScheme< Type, GType >::fvmLaplacianUncorrected().

Here is the caller graph for this function:

◆ boundaryCoeffs()

FieldField<Field, Type>& boundaryCoeffs ( )
inline

fvBoundary scalar field containing pseudo-matrix coeffs

for boundary cells

Definition at line 310 of file fvMatrix.H.

Referenced by immersedBoundaryFvPatchField< Type >::correctOffDiag(), gaussConvectionScheme< Type >::fvmDiv(), and gaussLaplacianScheme< Type, GType >::fvmLaplacianUncorrected().

Here is the caller graph for this function:

◆ faceFluxCorrectionPtr()

surfaceTypeFieldPtr& faceFluxCorrectionPtr ( )
inline

Return pointer to face-flux non-orthogonal correction field.

Definition at line 321 of file fvMatrix.H.

Referenced by gaussLaplacianScheme< Type, GType >::fvmLaplacian().

Here is the caller graph for this function:

◆ setValues() [1/2]

void setValues ( const labelUList cells,
const UList< Type > &  values 
)

◆ setValues() [2/2]

void setValues ( const labelUList cells,
const UIndirectList< Type > &  values 
)

Set solution in given cells to the specified values.

and eliminate the corresponding equations from the matrix.

◆ setReference()

void setReference ( const label  celli,
const Type &  value,
const bool  forceReference = false 
)

Set reference level for solution.

Definition at line 504 of file fvMatrix.C.

Referenced by Foam::CorrectPhi(), and main().

Here is the caller graph for this function:

◆ setComponentReference() [1/3]

void setComponentReference ( const label  patchi,
const label  facei,
const direction  cmpt,
const scalar  value 
)

Set reference level for a component of the solution.

on a given patch face

Definition at line 33 of file fvMatrixSolve.C.

◆ relax() [1/2]

void relax ( const scalar  alpha)

Relax matrix (for steady-state solution).

alpha = 1 : diagonally equal alpha < 1 : diagonally dominant alpha = 0 : do nothing Note: Requires positive diagonal.

Definition at line 519 of file fvMatrix.C.

Referenced by dynamicLagrangian< BasicTurbulenceModel >::correct(), thixotropicViscosity::correct(), advectionDiffusion::correct(), scalarTransport::execute(), main(), thermalBaffle::solveEnergy(), and reactingOneDim::solveEnergy().

Here is the caller graph for this function:

◆ relax() [2/2]

void relax ( )

Relax matrix (for steady-state solution).

alpha is read from controlDict

Definition at line 671 of file fvMatrix.C.

◆ boundaryManipulate()

void boundaryManipulate ( typename GeometricField< Type, fvPatchField, volMesh >::GeometricBoundaryField &  values)

Manipulate based on a boundary field.

Definition at line 688 of file fvMatrix.C.

◆ solver() [1/4]

autoPtr<fvSolver> solver ( const dictionary )

Construct and return the solver.

Use the given solver controls

◆ solver() [2/4]

Foam::autoPtr< typename Foam::fvMatrix< Type >::fvSolver > solver ( )

Construct and return the solver.

Solver controls read from fvSolution

Definition at line 279 of file fvMatrixSolve.C.

◆ solve() [1/2]

SolverPerformance<Type> solve ( const dictionary )

◆ solveSegregated() [1/3]

Foam::SolverPerformance< Type > solveSegregated ( const dictionary solverControls)

Solve segregated returning the solution statistics.

Use the given solver controls

Definition at line 104 of file fvMatrixSolve.C.

◆ solveCoupled()

Foam::SolverPerformance< Type > solveCoupled ( const dictionary solverControls)

Solve coupled returning the solution statistics.

Use the given solver controls

Definition at line 219 of file fvMatrixSolve.C.

◆ solve() [2/2]

Foam::SolverPerformance< Type > solve ( )

Solve returning the solution statistics.

Solver controls read from fvSolution

Definition at line 313 of file fvMatrixSolve.C.

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

Here is the caller graph for this function:

◆ residual() [1/3]

Foam::tmp< Foam::Field< Type > > residual ( ) const

Return the matrix residual.

Definition at line 330 of file fvMatrixSolve.C.

◆ D()

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

Return the matrix scalar diagonal.

Definition at line 701 of file fvMatrix.C.

◆ DD()

Foam::tmp< Foam::Field< Type > > DD ( ) const

Return the matrix Type diagonal.

Definition at line 710 of file fvMatrix.C.

◆ A()

Return the central coefficient.

Definition at line 734 of file fvMatrix.C.

Referenced by meanVelocityForce::constrain(), and directionalPressureGradientExplicitSource::constrain().

Here is the caller graph for this function:

◆ H() [1/3]

Return the H operation source.

Definition at line 763 of file fvMatrix.C.

◆ H1() [1/3]

Return H(1)

Definition at line 825 of file fvMatrix.C.

◆ flux()

Return the face-flux field from the matrix.

Definition at line 873 of file fvMatrix.C.

Referenced by Foam::CorrectPhi(), Foam::MULES::implicitSolve(), main(), and kinematicSingleLayer::solveThickness().

Here is the caller graph for this function:

◆ operator=() [1/2]

void operator= ( const fvMatrix< Type > &  )

◆ operator=() [2/2]

void operator= ( const tmp< fvMatrix< Type > > &  )

◆ negate()

void negate ( )

Definition at line 1001 of file fvMatrix.C.

◆ operator+=() [1/7]

void operator+= ( const fvMatrix< Type > &  )

◆ operator+=() [2/7]

void operator+= ( const tmp< fvMatrix< Type > > &  )

◆ operator-=() [1/7]

void operator-= ( const fvMatrix< Type > &  )

◆ operator-=() [2/7]

void operator-= ( const tmp< fvMatrix< Type > > &  )

◆ operator+=() [3/7]

void operator+= ( const DimensionedField< Type, volMesh > &  )

◆ operator+=() [4/7]

void operator+= ( const tmp< DimensionedField< Type, volMesh > > &  )

◆ operator+=() [5/7]

void operator+= ( const tmp< GeometricField< Type, fvPatchField, volMesh > > &  )

◆ operator-=() [3/7]

void operator-= ( const DimensionedField< Type, volMesh > &  )

◆ operator-=() [4/7]

void operator-= ( const tmp< DimensionedField< Type, volMesh > > &  )

◆ operator-=() [5/7]

void operator-= ( const tmp< GeometricField< Type, fvPatchField, volMesh > > &  )

◆ operator+=() [6/7]

void operator+= ( const dimensioned< Type > &  )

◆ operator-=() [6/7]

void operator-= ( const dimensioned< Type > &  )

◆ operator+=() [7/7]

void operator+= ( const zero )

Definition at line 1169 of file fvMatrix.C.

◆ operator-=() [7/7]

void operator-= ( const zero )

Definition at line 1177 of file fvMatrix.C.

◆ operator*=() [1/4]

void operator*= ( const DimensionedField< scalar, volMesh > &  )

◆ operator*=() [2/4]

void operator*= ( const tmp< DimensionedField< scalar, volMesh > > &  )

◆ operator*=() [3/4]

void operator*= ( const tmp< volScalarField > &  )

◆ operator*=() [4/4]

void operator*= ( const dimensioned< scalar > &  )

◆ setComponentReference() [2/3]

void setComponentReference ( const label  patchi,
const label  facei,
const  direction,
const scalar  value 
)

Definition at line 33 of file fvScalarMatrix.C.

◆ solver() [3/4]

Foam::autoPtr< Foam::fvMatrix< Foam::scalar >::fvSolver > solver ( const dictionary solverControls)

Definition at line 58 of file fvScalarMatrix.C.

◆ solveSegregated() [2/3]

Foam::solverPerformance solveSegregated ( const dictionary solverControls)

Definition at line 138 of file fvScalarMatrix.C.

◆ residual() [2/3]

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

Definition at line 187 of file fvScalarMatrix.C.

◆ H() [2/3]

Definition at line 211 of file fvScalarMatrix.C.

◆ H1() [2/3]

Definition at line 243 of file fvScalarMatrix.C.

◆ setComponentReference() [3/3]

void setComponentReference ( const label  patchi,
const label  facei,
const  direction,
const scalar  value 
)

◆ solver() [4/4]

autoPtr< fvMatrix< scalar >::fvSolver > solver ( const dictionary )

◆ solveSegregated() [3/3]

solverPerformance solveSegregated ( const dictionary )

◆ residual() [3/3]

tmp< scalarField > residual ( ) const

◆ H() [3/3]

tmp< volScalarField > H ( ) const

◆ H1() [3/3]

tmp< volScalarField > H1 ( ) const

Friends And Related Function Documentation

◆ fvSolver

friend class fvSolver
friend

Declare friendship with the fvSolver class.

Definition at line 148 of file fvMatrix.H.

◆ operator& [1/4]

tmp<GeometricField<Type, fvPatchField, volMesh> > operator& ( const fvMatrix< Type > &  ,
const DimensionedField< Type, volMesh > &   
)
friend

◆ operator& [2/4]

tmp<GeometricField<Type, fvPatchField, volMesh> > operator& ( const fvMatrix< Type > &  ,
const tmp< GeometricField< Type, fvPatchField, volMesh > > &   
)
friend

◆ operator& [3/4]

tmp<GeometricField<Type, fvPatchField, volMesh> > operator& ( const tmp< fvMatrix< Type > > &  ,
const DimensionedField< Type, volMesh > &   
)
friend

◆ operator& [4/4]

tmp<GeometricField<Type, fvPatchField, volMesh> > operator& ( const tmp< fvMatrix< Type > > &  ,
const tmp< GeometricField< Type, fvPatchField, volMesh > > &   
)
friend

◆ operator

Ostream& operator ( Ostream ,
const fvMatrix< Type > &   
)
friend

Field Documentation

◆ psi_

const GeometricField<Type, fvPatchField, volMesh>& psi_
private

Const reference to GeometricField<Type, fvPatchField, volMesh>

Converted into a non-const reference at the point of solution.

Definition at line 123 of file fvMatrix.H.

Referenced by fvMatrix< Type >::psi(), fvMatrix< Type >::residual(), fvMatrix::fvSolver::solve(), and fvMatrix< Type >::solve().

◆ dimensions_

dimensionSet dimensions_
private

Dimension set.

Definition at line 126 of file fvMatrix.H.

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

◆ source_

Field<Type> source_
private

Source term.

Definition at line 129 of file fvMatrix.H.

Referenced by fvMatrix< Type >::residual(), and fvMatrix< Type >::source().

◆ internalCoeffs_

FieldField<Field, Type> internalCoeffs_
private

Boundary scalar field containing pseudo-matrix coeffs.

for internal cells

Definition at line 133 of file fvMatrix.H.

◆ boundaryCoeffs_

FieldField<Field, Type> boundaryCoeffs_
private

Boundary scalar field containing pseudo-matrix coeffs.

for boundary cells

Definition at line 137 of file fvMatrix.H.

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

◆ faceFluxCorrectionPtr_

GeometricField<Type, fvsPatchField, surfaceMesh>* faceFluxCorrectionPtr_
mutableprivate

Face flux field for non-orthogonal correction.

Definition at line 142 of file fvMatrix.H.

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


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