snGradScheme.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) 2011-2016 OpenFOAM Foundation
9 Copyright (C) 2019-2021 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
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\*---------------------------------------------------------------------------*/
28
29#include "fv.H"
30#include "snGradScheme.H"
31#include "volFields.H"
32#include "surfaceFields.H"
33#include "HashTable.H"
34
35// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
36
37namespace Foam
38{
39
40// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
41
42namespace fv
43{
44
45// * * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * //
46
47template<class Type>
49(
50 const fvMesh& mesh,
51 Istream& schemeData
52)
53{
54 if (fv::debug)
55 {
56 InfoInFunction << "Constructing snGradScheme<Type>" << endl;
57 }
58
59 if (schemeData.eof())
60 {
61 FatalIOErrorInFunction(schemeData)
62 << "Discretisation scheme not specified"
63 << nl << nl
64 << "Valid schemes are :" << nl
65 << MeshConstructorTablePtr_->sortedToc()
67 }
68
69 const word schemeName(schemeData);
70
71 auto* ctorPtr = MeshConstructorTable(schemeName);
72
73 if (!ctorPtr)
74 {
76 (
77 schemeData,
78 "discretisation",
79 schemeName,
80 *MeshConstructorTablePtr_
81 ) << exit(FatalIOError);
82 }
83
84 return ctorPtr(mesh, schemeData);
85}
86
87
88// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
89
90template<class Type>
93(
95 const tmp<surfaceScalarField>& tdeltaCoeffs,
96 const word& snGradName
97)
98{
99 const fvMesh& mesh = vf.mesh();
100
101 // construct GeometricField<Type, fvsPatchField, surfaceMesh>
103 (
105 (
107 (
108 snGradName + "("+vf.name()+')',
109 vf.instance(),
110 vf.mesh(),
113 ),
114 mesh,
115 vf.dimensions()*tdeltaCoeffs().dimensions()
116 )
117 );
119 ssf.setOriented();
120
121 // set reference to difference factors array
122 const scalarField& deltaCoeffs = tdeltaCoeffs();
123
124 // owner/neighbour addressing
125 const labelUList& owner = mesh.owner();
126 const labelUList& neighbour = mesh.neighbour();
127
128 forAll(owner, facei)
129 {
130 ssf[facei] =
131 deltaCoeffs[facei]*(vf[neighbour[facei]] - vf[owner[facei]]);
132 }
133
135 Boundary& ssfbf = ssf.boundaryFieldRef();
136
137 forAll(vf.boundaryField(), patchi)
138 {
139 const fvPatchField<Type>& pvf = vf.boundaryField()[patchi];
140
141 if (pvf.coupled())
142 {
143 ssfbf[patchi] = pvf.snGrad(tdeltaCoeffs().boundaryField()[patchi]);
144 }
145 else
146 {
147 ssfbf[patchi] = pvf.snGrad();
148 }
149 }
150
151 return tsf;
152}
153
154
155template<class Type>
158(
160 const word& sndGradName
161)
162{
163 return snGrad(vf, vf.mesh().nonOrthDeltaCoeffs(), sndGradName);
164}
165
166
167template<class Type>
170(
172) const
173{
175 (
176 snGrad(vf, deltaCoeffs(vf))
177 );
178
179 if (corrected())
180 {
181 tsf.ref() += correction(vf);
182 }
183
184 return tsf;
185}
186
187
188template<class Type>
191(
193) const
194{
196 (
197 snGrad(tvf())
198 );
199
200 tsf.clear();
201 return tsf;
202}
203
204
205// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
206
207} // End namespace fv
208
209// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
210
211} // End namespace Foam
212
213// ************************************************************************* //
const dimensionSet & dimensions() const
Return dimensions.
const Mesh & mesh() const
Return mesh.
void setOriented(const bool oriented=true) noexcept
Set the oriented flag.
Generic GeometricField class.
Boundary & boundaryFieldRef(const bool updateAccessTime=true)
Return a reference to the boundary field.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:170
const word & name() const noexcept
Return the object name.
Definition: IOobjectI.H:65
const fileName & instance() const noexcept
Read access to instance path component.
Definition: IOobjectI.H:196
bool eof() const noexcept
True if end of input seen.
Definition: IOstream.H:239
An Istream is an abstract base class for all input systems (streams, files, token lists etc)....
Definition: Istream.H:64
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
Definition: Time.C:717
virtual tmp< Field< Type > > snGrad() const
Return gradient at boundary.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:91
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:82
virtual bool coupled() const
Return true if this patch field is coupled.
Definition: fvPatchField.H:350
virtual tmp< Field< Type > > snGrad() const
Return patch-normal gradient.
Definition: fvPatchField.C:229
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > sndGrad(const GeometricField< Type, fvPatchField, volMesh > &, const word &snGradName="sndGrad")
Return the sndGrad of the given cell field.
Definition: snGradScheme.C:158
A class for managing temporary objects.
Definition: tmp.H:65
void clear() const noexcept
Definition: tmpI.H:287
T & ref() const
Definition: tmpI.H:227
A class for handling words, derived from Foam::string.
Definition: word.H:68
dynamicFvMesh & mesh
#define FatalIOErrorInLookup(ios, lookupTag, lookupName, lookupTable)
Report an error message using Foam::FatalIOError.
Definition: error.H:478
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:473
#define InfoInFunction
Report an information message using Foam::Info.
Namespace for OpenFOAM.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
IOerror FatalIOError
tmp< fvMatrix< Type > > correction(const fvMatrix< Type > &)
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
labelList fv(nPoints)
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333
Foam::surfaceFields.