43#ifndef Foam_faMatrix_H
44#define Foam_faMatrix_H
61template<
class Type>
class faMatrix;
62template<
class T>
class UIndirectList;
130 template<
class Type2>
138 template<
class Type2>
147 template<
class Type2>
155 template<
class Type2>
177 const bool couples =
true
184 template<
template<
class>
class ListType>
188 const ListType<Type>& values
209 solver_(std::move(sol))
295 return internalCoeffs_;
302 return internalCoeffs_;
309 return boundaryCoeffs_;
316 return boundaryCoeffs_;
326 return faceFluxCorrectionPtr_;
332 return bool(faceFluxCorrectionPtr_);
367 const bool forceReference =
false
375 const bool forceReference =
false
383 const bool forceReference =
false
474 friend Ostream& operator<< <Type>
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
A field of fields is a PtrList of fields with reference counting.
Generic templated field type.
Generic GeometricField class.
An Istream is an abstract base class for all input systems (streams, files, token lists etc)....
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
A List with indirect addressing. Like IndirectList but does not store addressing.
Mesh data needed to do the Finite Area discretisation.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Dimension set for the base types, which can be used to implement rigorous dimension checking for alge...
Generic dimensioned Type class.
SolverPerformance< Type > solve()
Solve returning the solution statistics.
faSolver(faMatrix< Type > &faMat, autoPtr< lduMatrix::solver > &&sol)
A special matrix type and solver, designed for finite area solutions of scalar equations....
void addToInternalField(const labelUList &addr, const Field< Type2 > &pf, Field< Type2 > &intf) const
Add patch contribution to internal field.
const Field< Type > & source() const noexcept
void operator=(const faMatrix< Type > &)
GeometricField< Type, faePatchField, edgeMesh > faceFluxFieldType
Field type for face flux (for non-orthogonal correction)
const FieldField< Field, Type > & internalCoeffs() const noexcept
void relax()
Relax matrix (for steady-state solution).
SolverPerformance< Type > solve()
Solve returning the solution statistics.
tmp< faMatrix< Type > > clone() const
Construct and return a clone.
faceFluxFieldPtrType & faceFluxCorrectionPtr()
Return pointer to face-flux non-orthogonal correction field.
FieldField< Field, Type > & internalCoeffs() noexcept
tmp< GeometricField< Type, faePatchField, edgeMesh > > flux() const
Return the face-flux field from the matrix.
void setValues(const labelUList &faceLabels, const Type &value)
const FieldField< Field, Type > & boundaryCoeffs() const noexcept
void setReference(const label facei, const Type &value, const bool forceReference=false)
Set reference level for solution.
void setReferences(const labelUList &faceLabels, const Type &value, const bool forceReference=false)
Set reference level for solution.
void operator+=(const faMatrix< Type > &)
const GeometricField< Type, faPatchField, areaMesh > & psi() const
GeometricField< Type, faePatchField, edgeMesh > * faceFluxFieldPtrType
Declare return type of the faceFluxCorrectionPtr() function.
const dimensionSet & dimensions() const noexcept
tmp< GeometricField< Type, faPatchField, areaMesh > > H() const
Return the H operation source.
void addCmptAvBoundaryDiag(scalarField &diag) const
void addBoundarySource(Field< Type > &source, const bool couples=true) const
void setComponentReference(const label patchi, const label facei, const direction cmpt, const scalar value)
GeometricField< Type, faPatchField, areaMesh > psiFieldType
The geometric field type for psi.
virtual ~faMatrix()
Destructor.
tmp< areaScalarField > A() const
Return the central coefficient.
void operator-=(const faMatrix< Type > &)
tmp< Field< Type > > residual() const
Return the matrix residual.
void addBoundaryDiag(scalarField &diag, const direction cmpt) const
void negate()
Inplace negate.
void setValuesFromList(const labelUList &faceLabels, const ListType< Type > &values)
Set solution in given faces to the specified values.
tmp< scalarField > D() const
Return the matrix diagonal.
void subtractFromInternalField(const labelUList &addr, const Field< Type2 > &pf, Field< Type2 > &intf) const
Subtract patch contribution from internal field.
bool hasFaceFluxCorrection() const noexcept
True if face-flux non-orthogonal correction field exists.
FieldField< Field, Type > & boundaryCoeffs() noexcept
Field< Type > & source() noexcept
void operator*=(const areaScalarField::Internal &)
faPatchField<Type> abstract base class. This class gives a fat-interface to all derived classes cover...
lduMatrix is a general matrix class in which the coefficients are stored as three arrays,...
Reference counter for various OpenFOAM components.
A class for managing temporary objects.
A class representing the concept of 0 (zero) that can be used to avoid manipulating objects known to ...
Macro definitions for declaring ClassName(), NamespaceName(), etc.
#define ClassName(TypeNameString)
Add typeName information from argument TypeNameString to a class.
void checkMethod(const faMatrix< Type > &, const faMatrix< Type > &, const char *)
Ostream & operator<<(Ostream &, const boundaryPatch &p)
Write boundaryPatch as dictionary entries (without surrounding braces)
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh > > &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
#define FOAM_DEPRECATED_FOR(since, replacement)