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...
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 dimensionSet & | dimensions () 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... | |
surfaceTypeFieldPtr & | faceFluxCorrectionPtr () |
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< fvSolver > | solver (const dictionary &) |
Construct and return the solver. More... | |
autoPtr< fvSolver > | solver () |
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< scalarField > | D () const |
Return the matrix scalar diagonal. More... | |
tmp< Field< Type > > | DD () const |
Return the matrix Type diagonal. More... | |
tmp< volScalarField > | A () const |
Return the central coefficient. More... | |
tmp< GeometricField< Type, fvPatchField, volMesh > > | H () const |
Return the H operation source. More... | |
tmp< volScalarField > | H1 () 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 >::fvSolver > | solver (const dictionary &solverControls) |
Foam::solverPerformance | solveSegregated (const dictionary &solverControls) |
Foam::tmp< Foam::scalarField > | residual () const |
Foam::tmp< Foam::volScalarField > | H () const |
Foam::tmp< Foam::volScalarField > | H1 () const |
void | setComponentReference (const label patchi, const label facei, const direction, const scalar value) |
autoPtr< fvMatrix< scalar >::fvSolver > | solver (const dictionary &) |
solverPerformance | solveSegregated (const dictionary &) |
tmp< scalarField > | residual () const |
tmp< volScalarField > | H () const |
tmp< volScalarField > | H1 () 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 > > &) |
Ostream & | operator (Ostream &, const fvMatrix< Type > &) |
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.
Definition at line 68 of file fvPatchField.H.
typedef GeometricField<Type, fvsPatchField, surfaceMesh>* surfaceTypeFieldPtr |
Declare return type of the faceFluxCorrectionPtr() function.
Definition at line 318 of file fvMatrix.H.
fvMatrix | ( | const GeometricField< Type, fvPatchField, volMesh > & | , |
const dimensionSet & | |||
) |
Construct given a field to solve for.
Construct as copy of tmp<fvMatrix<Type> > deleting argument.
fvMatrix | ( | const GeometricField< Type, fvPatchField, volMesh > & | , |
Istream & | |||
) |
Construct from Istream given field to solve for.
|
virtual |
Destructor.
Definition at line 462 of file fvMatrix.C.
|
protected |
Add patch contribution to internal field.
Definition at line 38 of file fvMatrix.C.
|
protected |
Definition at line 61 of file fvMatrix.C.
|
protected |
Subtract patch contribution from internal field.
Definition at line 75 of file fvMatrix.C.
|
protected |
Definition at line 98 of file fvMatrix.C.
|
protected |
Definition at line 111 of file fvMatrix.C.
Referenced by fvMatrix< Type >::residual().
|
protected |
Definition at line 129 of file fvMatrix.C.
|
protected |
Definition at line 145 of file fvMatrix.C.
Referenced by fvMatrix< Type >::residual().
|
protected |
Set solution in given cells to the specified values.
Definition at line 178 of file fvMatrix.C.
ClassName | ( | "fvMatrix" | ) |
|
inline |
Definition at line 281 of file fvMatrix.H.
Referenced by velocityDampingConstraint::addDamping(), tabulatedAccelerationSource::addSup(), explicitPorositySource::addSup(), interRegionExplicitPorositySource::addSup(), canopySource::addSup(), radialActuationDiskSource::addSup(), interRegionHeatTransferModel::addSup(), SemiImplicitSource::addSup(), svenssonHaggkvistCanopySource::addSup(), actuationDiskSource::addSup(), dalpeMassonCanopySource::addSup(), solidificationMeltingSource::addSup(), rotorDiskSource::addSup(), solidificationMeltingSource::apply(), Foam::checkMethod(), optionList::constrain(), and immersedBoundaryFvPatchField< Type >::manipulateMatrix().
|
inline |
Definition at line 286 of file fvMatrix.H.
Referenced by explicitPorositySource::addSup(), interRegionExplicitPorositySource::addSup(), meanVelocityForce::addSup(), SemiImplicitSource::addSup(), directionalPressureGradientExplicitSource::addSup(), rotorDiskSource::addSup(), and Foam::checkMethod().
|
inline |
Definition at line 291 of file fvMatrix.H.
Referenced by interRegionExplicitPorositySource::addSup(), radialActuationDiskSource::addSup(), actuationDiskSource::addSup(), solidificationMeltingSource::addSup(), effectivenessHeatExchangerSource::addSup(), immersedBoundaryFvPatchField< Type >::correctOffDiag(), EulerD2dt2Scheme< Type >::fvmD2dt2(), EulerDdtScheme< Type >::fvmDdt(), backwardDdtScheme< Type >::fvmDdt(), CoEulerDdtScheme< Type >::fvmDdt(), SLTSDdtScheme< Type >::fvmDdt(), localEulerDdtScheme< Type >::fvmDdt(), CrankNicolsonDdtScheme< Type >::fvmDdt(), gaussLaplacianScheme< Type, GType >::fvmLaplacian(), ThermoCloud< CloudType >::Sh(), ReactingCloud< CloudType >::Srho(), KinematicCloud< CloudType >::SU(), and ReactingCloud< CloudType >::SYi().
|
inline |
Definition at line 296 of file fvMatrix.H.
|
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().
|
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().
|
inline |
Return pointer to face-flux non-orthogonal correction field.
Definition at line 321 of file fvMatrix.H.
Referenced by gaussLaplacianScheme< Type, GType >::fvmLaplacian().
void setValues | ( | const labelUList & | cells, |
const UList< Type > & | values | ||
) |
Set solution in given cells to the specified values.
and eliminate the corresponding equations from the matrix.
Referenced by ExplicitSetValue< Type >::constrain(), fixedTemperatureConstraint::constrain(), fixedInternalValueFvPatchField< Type >::manipulateMatrix(), epsilonWallFunctionFvPatchScalarField::manipulateMatrix(), omegaWallFunctionFvPatchScalarField::manipulateMatrix(), and immersedBoundaryFvPatchField< Type >::manipulateMatrix().
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.
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().
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.
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().
void relax | ( | ) |
Relax matrix (for steady-state solution).
alpha is read from controlDict
Definition at line 671 of file fvMatrix.C.
void boundaryManipulate | ( | typename GeometricField< Type, fvPatchField, volMesh >::GeometricBoundaryField & | values | ) |
Manipulate based on a boundary field.
Definition at line 688 of file fvMatrix.C.
autoPtr<fvSolver> solver | ( | const dictionary & | ) |
Construct and return the solver.
Use the given solver controls
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.
SolverPerformance<Type> solve | ( | const dictionary & | ) |
Solve segregated or coupled returning the solution statistics.
Use the given solver controls
Referenced by dynamicLagrangian< BasicTurbulenceModel >::correct(), thixotropicViscosity::correct(), advectionDiffusion::correct(), Foam::CorrectPhi(), scalarTransport::execute(), main(), surfaceAlignedSBRStressFvMotionSolver::solve(), reactingOneDim::solveContinuity(), thermalBaffle::solveEnergy(), reactingOneDim::solveEnergy(), reactingOneDim::solveSpeciesMass(), and kinematicSingleLayer::solveThickness().
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.
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.
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().
Foam::tmp< Foam::Field< Type > > residual | ( | ) | const |
Return the matrix residual.
Definition at line 330 of file fvMatrixSolve.C.
Foam::tmp< Foam::scalarField > D | ( | ) | const |
Return the matrix scalar diagonal.
Definition at line 701 of file fvMatrix.C.
Foam::tmp< Foam::Field< Type > > DD | ( | ) | const |
Return the matrix Type diagonal.
Definition at line 710 of file fvMatrix.C.
Foam::tmp< Foam::volScalarField > A | ( | ) | const |
Return the central coefficient.
Definition at line 734 of file fvMatrix.C.
Referenced by meanVelocityForce::constrain(), and directionalPressureGradientExplicitSource::constrain().
Foam::tmp< Foam::GeometricField< Type, Foam::fvPatchField, Foam::volMesh > > H | ( | ) | const |
Return the H operation source.
Definition at line 763 of file fvMatrix.C.
Foam::tmp< Foam::volScalarField > H1 | ( | ) | const |
Return H(1)
Definition at line 825 of file fvMatrix.C.
Foam::tmp< Foam::GeometricField< Type, Foam::fvsPatchField, Foam::surfaceMesh > > flux | ( | ) | const |
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().
void negate | ( | ) |
Definition at line 1001 of file fvMatrix.C.
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 > & | ) |
Definition at line 1169 of file fvMatrix.C.
Definition at line 1177 of file fvMatrix.C.
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 | ||
) |
Definition at line 33 of file fvScalarMatrix.C.
Foam::autoPtr< Foam::fvMatrix< Foam::scalar >::fvSolver > solver | ( | const dictionary & | solverControls | ) |
Definition at line 58 of file fvScalarMatrix.C.
Foam::solverPerformance solveSegregated | ( | const dictionary & | solverControls | ) |
Definition at line 138 of file fvScalarMatrix.C.
Foam::tmp< Foam::scalarField > residual | ( | ) | const |
Definition at line 187 of file fvScalarMatrix.C.
Foam::tmp< Foam::volScalarField > H | ( | ) | const |
Definition at line 211 of file fvScalarMatrix.C.
Foam::tmp< Foam::volScalarField > H1 | ( | ) | const |
Definition at line 243 of file fvScalarMatrix.C.
void setComponentReference | ( | const label | patchi, |
const label | facei, | ||
const | direction, | ||
const scalar | value | ||
) |
autoPtr< fvMatrix< scalar >::fvSolver > solver | ( | const dictionary & | ) |
solverPerformance solveSegregated | ( | const dictionary & | ) |
tmp< scalarField > residual | ( | ) | const |
tmp< volScalarField > H | ( | ) | const |
tmp< volScalarField > H1 | ( | ) | const |
|
friend |
Declare friendship with the fvSolver class.
Definition at line 148 of file fvMatrix.H.
|
friend |
|
friend |
|
friend |
|
friend |
|
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().
|
private |
Dimension set.
Definition at line 126 of file fvMatrix.H.
Referenced by fvMatrix< Type >::dimensions().
|
private |
Source term.
Definition at line 129 of file fvMatrix.H.
Referenced by fvMatrix< Type >::residual(), and fvMatrix< Type >::source().
|
private |
Boundary scalar field containing pseudo-matrix coeffs.
for internal cells
Definition at line 133 of file fvMatrix.H.
|
private |
Boundary scalar field containing pseudo-matrix coeffs.
for boundary cells
Definition at line 137 of file fvMatrix.H.
Referenced by fvMatrix< Type >::residual().
|
mutableprivate |
Face flux field for non-orthogonal correction.
Definition at line 142 of file fvMatrix.H.
Referenced by fvMatrix< Type >::faceFluxCorrectionPtr().
Copyright © 2011-2018 OpenFOAM | OPENFOAM® is a registered trademark of OpenCFD Ltd.