patchExprDriver.H
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) 2019-2021 OpenCFD Ltd.
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM.
12
13 OpenFOAM is free software: you can redistribute it and/or modify it
14 under the terms of the GNU General Public License as published by
15 the Free Software Foundation, either version 3 of the License, or
16 (at your option) any later version.
17
18 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 for more details.
22
23 You should have received a copy of the GNU General Public License
24 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25
26Class
27 Foam::expressions::patchExpr::parseDriver
28
29Description
30 Driver for patch expressions
31
32 In addition to the standard mathematical functions, operations and
33 logical and relational operations, the patch expressions support the
34 following driver-specific functions:
35
36 Functions
37 \table
38 Function | Description | Number of arguments |
39 pos | The face centres | 0 |
40 pts | The face points | 0 |
41 area | The face area magnitudes | 0 |
42 weightAverage| Area weighted average | 1 |
43 weightSum | Area weighted sum | 1 |
44 face | The face areaNormal vectors | 0 |
45 point | A point-field point value | 1 |
46 faceToPoint | Interpolate face values onto points | 1 |
47 pointToFace | Interpolate point values onto faces | 1 |
48 rand | Random field | 0/1 |
49 snGrad | Surface normal field | 0 |
50 internalField | Internal field next to patch | 0 |
51 neighbourField | patch field on opposite side of coupled patch | 0 |
52 \endtable
53
54Note
55 Use namespace debug switch \c patchExpr for scanner (2), parser (4)
56 or dictionary controls as per Foam::expressions::exprDriver.
57
58SourceFiles
59 patchExprDriver.C
60 patchExprDriverFields.C
61 patchExprDriverTemplates.C
62
63\*---------------------------------------------------------------------------*/
64
65#ifndef expressions_patchExprDriver_H
66#define expressions_patchExprDriver_H
67
68#include "patchExprFwd.H"
69#include "fvExprDriver.H"
71#include "Enum.H"
72#include "volFields.H"
73#include "surfaceFields.H"
74#include "pointFields.H"
76
77// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
78
79namespace Foam
80{
81namespace expressions
82{
83namespace patchExpr
84{
85
86/*---------------------------------------------------------------------------*\
87 Class parseDriver Declaration
88\*---------------------------------------------------------------------------*/
89
90class parseDriver
91:
92 public parsing::genericRagelLemonDriver,
93 public expressions::fvExprDriver
94{
95 // Private Member Functions
96
97 static const fvPatch& getFvPatch
98 (
99 const fvMesh& fvm,
100 const dictionary& dict
101 );
102
103
104protected:
105
106 // Protected Data
107
108 //- The referenced patch
109 const fvPatch& patch_;
110
111
112 // Protected Member Functions
113
114 //- Cell selections (as logical)
115 tmp<boolField> field_cellSelection
116 (
117 const word& name,
118 enum topoSetSource::sourceType setType
119 ) const;
120
121 //- Face selections (as logical)
122 tmp<boolField> field_faceSelection
123 (
124 const word& name,
125 enum topoSetSource::sourceType setType
126 ) const;
127
128
129public:
130
131 ClassName("patchExpr::driver");
132
133 // Generated Methods
134
135 // No copy copy construct
136 parseDriver(const parseDriver&) = delete;
137
138 // No copy assignment
139 void operator=(const parseDriver&) = delete;
140
141
142 // Constructors
143
144 //- Construct for specified patch, with dictionary information
145 explicit parseDriver
146 (
147 const fvPatch& p,
148 const dictionary& dict = dictionary::null
149 );
150
151 //- Construct for specified patch with copy of driver context
153 (
154 const fvPatch& p,
155 const parseDriver& driver,
156 const dictionary& dict // = dictionary::null
157 );
158
159 //- Construct with patchName for the given mesh
160 parseDriver(const word& patchName, const fvMesh& mesh);
161
162 //- Construct with "patch" (mandatory) and "region" (optional)
163 //- specified in dictionary
164 parseDriver(const dictionary& dict, const fvMesh& mesh);
165
166
167 // Not generally clonable
168
169
170 //- Destructor
171 virtual ~parseDriver() = default;
172
173
174 // Public Member Functions
175
176 //- The mesh we are attached to
177 virtual const fvMesh& mesh() const
179 return patch_.boundaryMesh().mesh();
180 }
181
182 //- The natural field size for the expression
183 virtual label size() const
184 {
185 return patch_.patch().size();
186 }
187
188 //- The point field size for the expression
189 virtual label pointSize() const
190 {
191 return patch_.patch().nPoints();
192 }
193
194 //- Field size associated with different geometric field types
195 inline label size(const FieldAssociation geoType) const;
196
197
198 // Evaluation
199
200 //- Perform parsing on (sub) string
201 using genericRagelLemonDriver::content;
202
203 //- Execute the parser.
204 // The return value currently has no meaning.
205 virtual unsigned parse
206 (
207 const std::string& expr,
208 size_t pos = 0,
209 size_t len = std::string::npos
210 );
211
212
213 // General
214
215 //- Set result
216 template<class Type>
217 void setResult(Field<Type>* ptr, bool pointVal = false)
218 {
219 result().setResult<Type>(ptr, pointVal);
220 }
221
222
223 //- Retrieve variable as field if possible.
224 // Test tmp for validity to determine success of the operation.
225 template<class Type>
226 tmp<Field<Type>> getVariableIfAvailable(const word& fldName) const;
227
228
229 // Field Retrieval
230
231 //- Return named field
232 template<class Type>
233 tmp<Field<Type>> getField(const word& fldName);
234
235 //- Retrieve field (vol field)
236 template<class Type>
237 tmp<Field<Type>> getVolField(const word& fldName);
238
239 //- Retrieve field (surface field)
240 template<class Type>
241 tmp<Field<Type>> getSurfaceField(const word& fldName);
242
243 //- Retrieve field (point field)
244 template<class Type>
245 tmp<Field<Type>> getPointField(const word& fldName);
247 //- Return internal field next to patch
248 template<class Type>
249 tmp<Field<Type>> patchInternalField(const word& fldName);
250
251 //- Return patchField on the opposite patch of a coupled patch
252 template<class Type>
254
255 //- Return surface normal field (snGrad)
256 template<class Type>
257 tmp<Field<Type>> patchNormalField(const word& fldName);
259
260 // Field "shape" conversions
261
262 //- Interpolate face to point
263 template<class Type>
265
266 //- Interpolate point to face values
267 template<class Type>
269
270
271 // Custom Field Functions
272
273 //- The area-weighted average of a field
274 template<class Type>
275 Type areaAverage(const Field<Type>& fld) const
276 {
277 return weightedAverage(patch_.magSf(), fld);
278 }
279
280 //- The area-weighted sum of a field
281 template<class Type>
282 Type areaSum(const Field<Type>& fld) const
283 {
284 return weightedSum(patch_.magSf(), fld);
285 }
287 //- The face area magnitudes [magSf] - (swak = area)
289
290 //- The face centres - (swak = pos)
292
293 //- The face areas with their vector direction [Sf] - (swak = face)
296 //- The patch point locations - (swak = pts)
298
299
300 //- Cell selection (set)
301 inline tmp<boolField> field_cellSet(const word& name) const;
303 //- Cell selection (zone)
304 inline tmp<boolField> field_cellZone(const word& name) const;
305
306 //- Face selection (set)
307 inline tmp<boolField> field_faceSet(const word& name) const;
308
309 //- Face selection (zone)
311
312
313 //- A uniform random field
314 tmp<scalarField> field_rand(label seed=0, bool gaussian=false) const;
315
316 //- A Gaussian random field
317 tmp<scalarField> field_randGaussian(label seed=0) const
319 return field_rand(seed, true);
320 }
321};
323
324// Template specializations
325
326//- Retrieve field (surface field: bool)
327template<>
328tmp<Field<bool>> parseDriver::getSurfaceField<bool>(const word& fldName);
329
330//- Retrieve field (point field: bool)
331template<>
332tmp<Field<bool>> parseDriver::getPointField<bool>(const word& fldName);
334
335// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
336
337} // End namespace patchExpr
338} // End namespace expressions
339} // End namespace Foam
340
341
342// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
343
345
346#ifdef NoRepository
348#endif
349
350// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
352#endif
353
354// ************************************************************************* //
Info<< nl<< "Wrote faMesh in vtk format: "<< writer.output().name()<< nl;}{ vtk::lineWriter writer(aMesh.points(), aMesh.edges(), fileName(aMesh.mesh().time().globalPath()/"finiteArea-edges"));writer.writeGeometry();writer.beginCellData(4);writer.writeProcIDs();{ Field< scalar > fld(faMeshTools::flattenEdgeField(aMesh.magLe(), true))
Generic templated field type.
Definition: Field.H:82
label nPoints() const
Number of points supporting patch faces.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
static const dictionary null
An empty dictionary, which is also the parent for all dictionaries.
Definition: dictionary.H:394
The field association for mesh (patch/volume) values.
static Type weightedAverage(const scalarField &weights, const Field< Type > &fld)
The (global) weighted average of a field, with stabilisation.
static Type weightedSum(const scalarField &weights, const Field< Type > &fld)
The (global) weighted sum (integral) of a field.
const exprResult & result() const noexcept
Const access to expression result.
Definition: exprDriver.H:391
const dictionary & dict() const noexcept
The dictionary with all input data/specification.
Definition: exprDriver.H:385
void setResult(Field< Type > *, bool wantPointData=false)
Set result field, taking ownership of the pointer.
Definition: exprResultI.H:359
Base driver for parsing value expressions associated with an fvMesh.
Definition: fvExprDriver.H:139
virtual label size() const
The natural field size for the expression.
tmp< Field< Type > > getVariableIfAvailable(const word &fldName) const
Retrieve variable as field if possible.
tmp< boolField > field_faceZone(const word &name) const
Face selection (zone)
tmp< Field< Type > > patchNormalField(const word &fldName)
Return surface normal field (snGrad)
tmp< Field< Type > > faceToPoint(const Field< Type > &field) const
Interpolate face to point.
tmp< Field< Type > > patchNeighbourField(const word &fldName)
Return patchField on the opposite patch of a coupled patch.
virtual label pointSize() const
The point field size for the expression.
tmp< boolField > field_cellZone(const word &name) const
Cell selection (zone)
tmp< vectorField > field_faceCentre() const
The face centres - (swak = pos)
tmp< Field< Type > > getVolField(const word &fldName)
Retrieve field (vol field)
tmp< Field< Type > > getPointField(const word &fldName)
Retrieve field (point field)
virtual unsigned parse(const std::string &expr, size_t pos=0, size_t len=std::string::npos)
Execute the parser.
void setResult(Field< Type > *ptr, bool pointVal=false)
Set result.
tmp< boolField > field_cellSet(const word &name) const
Cell selection (set)
void operator=(const parseDriver &)=delete
tmp< boolField > field_cellSelection(const word &name, enum topoSetSource::sourceType setType) const
Cell selections (as logical)
tmp< Field< Type > > pointToFace(const Field< Type > &field) const
Interpolate point to face values.
parseDriver(const parseDriver &)=delete
tmp< Field< Type > > getSurfaceField(const word &fldName)
Retrieve field (surface field)
tmp< boolField > field_faceSet(const word &name) const
Face selection (set)
tmp< scalarField > field_rand(label seed=0, bool gaussian=false) const
A uniform random field.
tmp< boolField > field_faceSelection(const word &name, enum topoSetSource::sourceType setType) const
Face selections (as logical)
tmp< Field< Type > > patchInternalField(const word &fldName)
Return internal field next to patch.
const fvPatch & patch_
The referenced patch.
Type areaAverage(const Field< Type > &fld) const
The area-weighted average of a field.
tmp< scalarField > field_faceArea() const
The face area magnitudes [magSf] - (swak = area)
tmp< vectorField > field_areaNormal() const
The face areas with their vector direction [Sf] - (swak = face)
tmp< vectorField > field_pointField() const
The patch point locations - (swak = pts)
virtual const fvMesh & mesh() const
The mesh we are attached to.
tmp< Field< Type > > getField(const word &fldName)
Return named field.
Type areaSum(const Field< Type > &fld) const
The area-weighted sum of a field.
tmp< scalarField > field_randGaussian(label seed=0) const
A Gaussian random field.
virtual ~parseDriver()=default
Destructor.
A topoSetPointSource to select all points based on usage in given faceSet(s).
Definition: faceToPoint.H:177
const fvMesh & mesh() const noexcept
Return the mesh reference.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:91
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:71
const scalarField & magSf() const
Return face area magnitudes.
Definition: fvPatch.C:156
const fvBoundaryMesh & boundaryMesh() const
Return boundaryMesh reference.
Definition: fvPatch.H:209
const polyPatch & patch() const
Return the polyPatch.
Definition: fvPatch.H:167
Generic interface code for Ragel/Lemon combination Subclasses should implement one or more process() ...
A topoSetFaceSource to select faces with any point or any edge within a given pointSet(s).
Definition: pointToFace.H:179
A class for managing temporary objects.
Definition: tmp.H:65
sourceType
Enumeration defining the types of sources.
Definition: topoSetSource.H:75
A class for handling words, derived from Foam::string.
Definition: word.H:68
#define ClassName(TypeNameString)
Add typeName information from argument TypeNameString to a class.
Definition: className.H:67
volScalarField & p
rDeltaTY field()
Namespace for OpenFOAM.
dimensionedScalar pos(const dimensionedScalar &ds)
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
Foam::surfaceFields.