45 if (psi_.needReference())
47 if (Pstream::master())
49 internalCoeffs_[patchi][facei].component(cmpt) +=
50 diag()[psi_.mesh().boundary()[patchi].faceCells()[facei]];
52 boundaryCoeffs_[patchi][facei].component(cmpt) +=
53 diag()[psi_.mesh().boundary()[patchi].faceCells()[facei]]
63 const dictionary& solverControls
67 if (psi_.mesh().name() != polyMesh::defaultRegion)
76 <<
"fvMatrix<Type>::solveSegregatedOrCoupled" 77 "(const dictionary& solverControls) : " 78 "solving fvMatrix<Type>" 83 if (solverControls.getOrDefault<label>(
"maxIter", -1) == 0)
85 return SolverPerformance<Type>();
88 word
type(solverControls.getOrDefault<word>(
"type",
"segregated"));
90 if (
type ==
"segregated")
92 return solveSegregated(solverControls);
94 else if (
type ==
"coupled")
96 return solveCoupled(solverControls);
101 <<
"Unknown type " <<
type 102 <<
"; currently supported solver types are segregated and coupled" 113 const dictionary& solverControls
119 <<
"Implicit option is not allowed for type: " << Type::typeName
126 <<
"fvMatrix<Type>::solveSegregated" 127 "(const dictionary& solverControls) : " 128 "solving fvMatrix<Type>" 133 solverControls.getOrDefault<
int>
140 const_cast<GeometricField<Type, fvPatchField, volMesh>&
>(psi_);
142 SolverPerformance<Type> solverPerfVec
144 "fvMatrix<Type>::solveSegregated",
150 Field<Type> source(source_);
155 addBoundarySource(source);
157 typename Type::labelType validComponents
159 psi.mesh().template validComponents<Type>()
164 if (validComponents[cmpt] == -1)
continue;
169 addBoundaryDiag(
diag(), cmpt);
173 FieldField<Field, scalar> bouCoeffsCmpt
175 boundaryCoeffs_.component(cmpt)
178 FieldField<Field, scalar> intCoeffsCmpt
180 internalCoeffs_.component(cmpt)
184 psi.boundaryField().scalarInterfaces();
190 PrecisionAdaptor<solveScalar, scalar> sourceCmpt_ss(sourceCmpt);
191 ConstPrecisionAdaptor<solveScalar, scalar> psiCmpt_ss(psiCmpt);
193 const label startRequest = UPstream::nRequests();
205 updateMatrixInterfaces
222 psi.name() + pTraits<Type>::componentNames[cmpt],
228 )->
solve(psiCmpt, sourceCmpt, cmpt);
235 solverPerfVec.replace(cmpt, solverPerf);
236 solverPerfVec.solverName() = solverPerf.solverName();
238 psi.primitiveFieldRef().replace(cmpt, psiCmpt);
242 psi.correctBoundaryConditions();
244 psi.mesh().data().setSolverPerformance(
psi.name(), solverPerfVec);
246 return solverPerfVec;
259 <<
"fvMatrix<Type>::solveCoupled" 260 "(const dictionary& solverControls) : " 261 "solving fvMatrix<Type>" 273 const_cast<GeometricField<Type, fvPatchField, volMesh>&
>(psi_);
275 LduMatrix<Type, scalar, scalar> coupledMatrix(
psi.mesh());
276 coupledMatrix.diag() =
diag();
277 coupledMatrix.upper() =
upper();
278 coupledMatrix.lower() =
lower();
279 coupledMatrix.source() = source();
281 addBoundaryDiag(coupledMatrix.diag(), 0);
282 addBoundarySource(coupledMatrix.source(),
false);
284 coupledMatrix.interfaces() =
psi.boundaryFieldRef().interfaces();
285 coupledMatrix.interfacesUpper() = boundaryCoeffs().component(0);
286 coupledMatrix.interfacesLower() = internalCoeffs().component(0);
288 autoPtr<typename LduMatrix<Type, scalar, scalar>::solver>
299 SolverPerformance<Type> solverPerf
301 coupledMatrixSolver->solve(
psi)
309 psi.correctBoundaryConditions();
311 psi.mesh().data().setSolverPerformance(
psi.name(), solverPerf);
320 return this->solveSegregatedOrCoupled(solverDict());
330 return psi_.mesh().solve(*
this, solverControls);
338 return solver(solverDict());
345 return solve(fvMat_.solverDict());
367 auto& res = tres.ref();
381 boundaryCoeffs_.component(cmpt)
390 res.component(cmpt) - boundaryDiagCmpt*psiCmpt,
UPtrList< const lduInterfaceField > lduInterfaceFieldPtrsList
List of coupled interface fields to be used in coupling.
void size(const label n)
Older name for setAddressableSize.
const Internal::FieldType & primitiveField() const noexcept
Return a const-reference to the internal field values.
void addBoundaryDiag(scalarField &diag, const direction cmpt) const
errorManipArg< error, int > exit(error &err, const int errNo=1)
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
A list of keyword definitions, which are a keyword followed by a number of values (eg...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
void residual(solveScalarField &rA, const solveScalarField &psi, const scalarField &source, const FieldField< Field, scalar > &interfaceBouCoeffs, const lduInterfaceFieldPtrsList &interfaces, const direction cmpt) const
lduInterfaceFieldPtrsList scalarInterfaces() const
Return a list of pointers for each patch field with only those pointing to interfaces being set...
Ostream & endl(Ostream &os)
Add newline and flush stream.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tf1, const word &name, const dimensionSet &dimensions, const bool initCopy=false)
Global function forwards to reuseTmpDimensionedField::New.
#define addProfiling(Name,...)
Define profiling trigger with specified name and description string. The description is generated by ...
string upper(const std::string &s)
Return string copy transformed with std::toupper on each character.
::Foam::direction nComponents(const expressions::valueTypeCode) noexcept
The number of components associated with given valueTypeCode.
A field of fields is a PtrList of fields with reference counting.
SolverPerformance< Type > solveCoupled(const dictionary &)
Solve coupled returning the solution statistics.
void diag(pointPatchField< vector > &, const pointPatchField< tensor > &)
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for expressions::valueTypeCode::INVALID.
SolverPerformance< Type > solve()
Solve returning the solution statistics.
A class for handling words, derived from Foam::string.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Abstract base-class for lduMatrix solvers.
tmp< Field< cmptType > > component(const direction) const
Return a component field of the field.
SolverPerformance< Type > solveSegregatedOrCoupled()
Solve segregated or coupled returning the solution statistics.
int debug
Static debugging option.
void replace(const direction, const FieldField< Field, cmptType > &)
Replace a component field of the field.
OSstream & masterStream(const label communicator)
Return OSstream for output operations on the master process only, Snull on other processes.
string lower(const std::string &s)
Return string copy transformed with std::tolower on each character.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Foam::word regionName(args.getOrDefault< word >("region", Foam::polyMesh::defaultRegion))
autoPtr< fvSolver > solver()
Construct and return the solver.
messageStream Info
Information stream (stdout output on master, null elsewhere)
void addBoundarySource(Field< Type > &source, const bool couples=true) const
Pointer management similar to std::unique_ptr, with some additional methods and type checking...
const volScalarField & psi
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.
SolverPerformance< Type > solveSegregated(const dictionary &)
Solve segregated returning the solution statistics.
T getOrDefault(const word &keyword, const T &deflt, enum keyType::option matchOpt=keyType::REGEX) const
Find and return a T, or return the given default value. FatalIOError if it is found and the number of...
A class for managing temporary objects.
tmp< Field< Type > > residual() const
Return the matrix residual.
SolverPerformance< Type > solve()
Solve returning the solution statistics.
const Boundary & boundaryField() const noexcept
Return const-reference to the boundary field.
SolverPerformance< scalar > solverPerformance
SolverPerformance instantiated for a scalar.
const dictionary & solverDict() const
Return the solver dictionary for psi, taking into account finalIteration.
IOerror FatalIOError
Error stream (stdout output on all processes), with additional 'FOAM FATAL IO ERROR' header text and ...
static constexpr const zero Zero
Global zero (0)