faScalarMatrix.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2016-2017 Wikki Ltd
9  Copyright (C) 2019-2021 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
12  This file is part of OpenFOAM.
13 
14  OpenFOAM is free software: you can redistribute it and/or modify it
15  under the terms of the GNU General Public License as published by
16  the Free Software Foundation, either version 3 of the License, or
17  (at your option) any later version.
18 
19  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22  for more details.
23 
24  You should have received a copy of the GNU General Public License
25  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26 
27 Description
28  Finite-Area scalar matrix member functions and operators
29 
30 \*---------------------------------------------------------------------------*/
31 
32 #include "faScalarMatrix.H"
34 #include "PrecisionAdaptor.H"
35 
36 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
37 
38 template<>
40 (
41  const label patchI,
42  const label edgeI,
43  const direction,
44  const scalar value
45 )
46 {
47  const labelUList& faceLabels = psi_.mesh().boundary()[patchI].edgeFaces();
48 
49  internalCoeffs_[patchI][edgeI] += diag()[faceLabels[edgeI]];
50 
51  boundaryCoeffs_[patchI][edgeI] = value;
52 }
53 
54 
55 template<>
57 (
58  const dictionary& solverControls
59 )
60 {
62  << "solving faMatrix<scalar>"
63  << endl;
64 
65  const int logLevel =
66  solverControls.getOrDefault<int>
67  (
68  "log",
70  );
71 
72  auto& psi =
74 
75  scalarField saveDiag(diag());
76  addBoundaryDiag(diag(), 0);
77 
78  scalarField totalSource(source_);
79  addBoundarySource(totalSource, 0);
80 
81  // Solver call
83  (
84  psi_.name(),
85  *this,
86  boundaryCoeffs_,
87  internalCoeffs_,
88  psi_.boundaryField().scalarInterfaces(),
89  solverControls
90  )->solve(psi.ref(), totalSource);
91 
92  if (logLevel)
93  {
94  solverPerf.print(Info);
95  }
96 
97  diag() = saveDiag;
98 
99  psi.correctBoundaryConditions();
100 
101  psi.mesh().setSolverPerformance(psi.name(), solverPerf);
102 
103  return solverPerf;
104 }
105 
106 
107 template<>
109 {
110  scalarField boundaryDiag(psi_.size(), Zero);
111  addBoundaryDiag(boundaryDiag, 0);
112 
113  const scalarField& psif = psi_.internalField();
115  const solveScalarField& psi = tpsi();
116 
118  (
119  lduMatrix::residual
120  (
121  psi,
122  source_ - boundaryDiag*psif,
123  boundaryCoeffs_,
124  psi_.boundaryField().scalarInterfaces(),
125  0
126  )
127  );
128 
130  addBoundarySource(tres_s.ref());
131 
132  return tres_s;
133 }
134 
135 
136 template<>
138 {
140  (
141  new areaScalarField
142  (
143  IOobject
144  (
145  "H("+psi_.name()+')',
146  psi_.instance(),
147  psi_.db(),
148  IOobject::NO_READ,
149  IOobject::NO_WRITE
150  ),
151  psi_.mesh(),
152  dimensions_/dimArea,
153  zeroGradientFaPatchScalarField::typeName
154  )
155  );
156  areaScalarField& Hphi = tHphi.ref();
157 
158  Hphi.primitiveFieldRef() = (lduMatrix::H(psi_.primitiveField()) + source_);
159  addBoundarySource(Hphi.primitiveFieldRef());
160 
161  Hphi.ref() /= psi_.mesh().S();
163 
164  return tHphi;
165 }
166 
167 
168 // ************************************************************************* //
Foam::expressions::patchExpr::debug
int debug
Static debugging option.
Foam::IOobject
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:169
Foam::tmp
A class for managing temporary objects.
Definition: PtrList.H:61
Foam::Zero
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
Foam::diag
void diag(pointPatchField< vector > &, const pointPatchField< tensor > &)
Definition: pointPatchFieldFunctions.H:287
H
volScalarField H(IOobject("H", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::AUTO_WRITE), mesh, dimensionedScalar(dimLength, Zero))
PrecisionAdaptor.H
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:369
Foam::faMatrix::setComponentReference
void setComponentReference(const label patchi, const label facei, const direction cmpt, const scalar value)
Definition: faMatrixSolve.C:38
solve
CEqn solve()
Foam::SolverPerformance::print
void print(Ostream &os) const
Print summary of solver performance to the given stream.
Definition: SolverPerformance.C:93
Foam::tmp::ref
T & ref() const
Definition: tmpI.H:227
Foam::dimArea
const dimensionSet dimArea(sqr(dimLength))
Definition: dimensionSets.H:59
Foam::Field< scalar >
Foam::DimensionedField::mesh
const Mesh & mesh() const
Return mesh.
Definition: DimensionedFieldI.H:41
Foam::Info
messageStream Info
Information stream (stdout output on master, null elsewhere)
DebugInFunction
#define DebugInFunction
Report an information message using Foam::Info.
Definition: messageStream.H:388
Foam::faMatrix::solve
SolverPerformance< Type > solve()
Solve returning the solution statistics.
Definition: faMatrixSolve.C:179
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:123
Foam::faMatrix::H
tmp< GeometricField< Type, faPatchField, areaMesh > > H() const
Return the H operation source.
Definition: faMatrix.C:661
Foam::GeometricField::primitiveFieldRef
Internal::FieldType & primitiveFieldRef(const bool updateAccessTime=true)
Return a reference to the internal field.
Definition: GeometricField.C:766
zeroGradientFaPatchFields.H
Foam::GeometricField::correctBoundaryConditions
void correctBoundaryConditions()
Correct boundary field.
Definition: GeometricField.C:940
Foam::New
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
Definition: DimensionedFieldReuseFunctions.H:105
Foam::GeometricField::ref
Internal & ref(const bool updateAccessTime=true)
Return a reference to the dimensioned internal field.
Definition: GeometricField.C:749
Foam::UList< label >
Foam::ConstPrecisionAdaptor
A const Field/List wrapper with possible data conversion.
Definition: PrecisionAdaptor.H:57
Foam::direction
uint8_t direction
Definition: direction.H:52
faScalarMatrix.H
Foam::faMatrix::residual
tmp< Field< Type > > residual() const
Return the matrix residual.
Definition: faMatrixSolve.C:186
Foam::refPtr< Container< Type > >::ref
Container< Type > & ref() const
Definition: refPtrI.H:203
Foam::dictionary::getOrDefault
T getOrDefault(const word &keyword, const T &deflt, enum keyType::option matchOpt=keyType::REGEX) const
Definition: dictionaryTemplates.C:148
Foam::GeometricField< scalar, faPatchField, areaMesh >
psi
const volScalarField & psi
Definition: createFieldRefs.H:1
Foam::SolverPerformance
SolverPerformance is the class returned by the LduMatrix solver containing performance statistics.
Definition: SolverPerformance.H:52