volumeExprDriverTemplates.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) 2019-2022 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
26\*---------------------------------------------------------------------------*/
27
28#include "exprOps.H"
29#include "FieldOps.H"
30#include "surfaceInterpolate.H"
33
34// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
35
36template<class Type>
38(
39 const Field<Type>& fld
40)
41{
42 if (isLogical_)
43 {
44 // Eg, volScalarField -> volLogicalField
45 resultType_.replace("Scalar", "Logical");
46
47 Field<bool> bools(fld.size());
49
50 this->result().setResult(std::move(bools), this->isPointData());
51 }
52 else
53 {
54 // Deep copy
55 this->result().setResult(fld, this->isPointData());
56 }
57}
58
59
60// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
61
62template<class Type>
64(
66 bool logical
67)
68{
69 resultField_.reset(nullptr);
70
71 // Characteristics
72 resultType_ = VolumeField<Type>::typeName;
73 isLogical_ = logical;
74 fieldGeoType_ = VOLUME_DATA;
75
76 // Assign dimensions
77 if (hasDimensions_ && !logical)
78 {
79 ptr->dimensions().reset(resultDimensions_);
80 }
81
82 setInternalFieldResult(ptr->primitiveField());
83
84 // Take ownership
85 resultField_.reset(ptr);
86}
87
88
89template<class Type>
91(
93 bool logical
94)
95{
96 resultField_.reset(nullptr);
97
98 // Characteristics
99 resultType_ = SurfaceField<Type>::typeName;
100 isLogical_ = logical;
101 fieldGeoType_ = FACE_DATA;
102
103 // Assign dimensions
104 if (hasDimensions_ && !logical)
105 {
106 ptr->dimensions().reset(resultDimensions_);
107 }
108
109 setInternalFieldResult(ptr->primitiveField());
110
111 // Take ownership
112 resultField_.reset(ptr);
113}
114
115
116template<class Type>
118(
119 PointField<Type>* ptr,
120 bool logical
121)
122{
123 resultField_.reset(nullptr);
124
125 // Characteristics
126 resultType_ = PointField<Type>::typeName;
127 isLogical_ = logical;
128 fieldGeoType_ = POINT_DATA;
129
130 // Assign dimensions
131 if (hasDimensions_ && !logical)
132 {
133 ptr->dimensions().reset(resultDimensions_);
134 }
135
136 setInternalFieldResult(ptr->primitiveField());
137
138 // Take ownership
139 resultField_.reset(ptr);
140}
141
142
143template<class GeomField>
144const GeomField*
146{
147 return dynamic_cast<const GeomField*>(resultField_.get());
148}
149
150
151template<class GeomField>
152const GeomField*
154(
155 bool logical,
156 bool dieOnNull
157) const
158{
159 const regIOobject* ptr = resultField_.get();
160
161 if (dieOnNull && !ptr)
162 {
164 << "No result available. Requested "
166 << exit(FatalError);
167 }
168
169 if (isLogical_ == logical)
170 {
171 return dynamic_cast<const GeomField*>(ptr);
172 }
173
174 return nullptr;
175}
176
177
178template<class Type>
181(
182 const word& fldName,
183 bool getOldTime
184)
185{
186 return this->getOrReadField<VolumeField<Type>>
187 (
188 fldName,
189 true, // mandatory
190 getOldTime
191 );
192}
193
194
195template<class Type>
198(
199 const word& fldName,
200 bool getOldTime
201)
202{
203 return this->getOrReadField<SurfaceField<Type>>
204 (
205 fldName,
206 true, // mandatory
207 getOldTime
208 );
209}
210
211
212template<class Type>
215(
216 const word& fldName,
217 bool getOldTime
218)
219{
220 return this->getOrReadPointField<PointField<Type>>
221 (
222 fldName,
223 true, // mandatory
224 getOldTime
225 );
226}
227
228
229template<class Type>
232(
233 const Type& val
234) const
235{
237 (
238 word("constant.") + word(pTraits<Type>::typeName),
239 mesh(),
240 dimensioned<Type>("", dimless, val)
241 );
242}
243
244
245template<class Type>
248(
249 const Type& val
250) const
251{
253 (
254 word("constant.") + word(pTraits<Type>::typeName),
255 mesh(),
256 dimensioned<Type>("", dimless, val)
257 );
258}
259
260
261template<class Type>
264(
265 const Type& val
266) const
267{
269 (
270 word("constant.") + word(pTraits<Type>::typeName),
272 dimensioned<Type>("", dimless, val)
273 );
274}
275
276
277template<class Type>
280(
282) const
283{
284 return fvc::interpolate(field);
285}
286
287
288template<class Type>
291(
293) const
294{
295 volPointInterpolation interp(this->mesh());
296 return interp.interpolate(field);
297}
298
299
300template<class Type>
303(
305) const
306{
307 auto tresult = newVolField<Type>();
308 auto& result = tresult.ref();
309
310 forAll(result,celli)
311 {
312 result[celli] = interpolatePointToCell(field, celli);
313 }
314
315 return tresult;
316}
317
318
319// ************************************************************************* //
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))
const dimensionSet & dimensions() const
Return dimensions.
Generic templated field type.
Definition: Field.H:82
Generic GeometricField class.
const Internal::FieldType & primitiveField() const
Return a const-reference to the internal field.
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
Definition: Time.C:717
A topoSetFaceSource to select all the faces from given cellSet(s).
Definition: cellToFace.H:177
A topoSetPointSource to select all the points from given cellSet(s).
Definition: cellToPoint.H:176
void reset(const dimensionSet &ds)
Copy assign the exponents from the dimensionSet.
Definition: dimensionSet.C:149
Generic dimensioned Type class.
const exprResult & result() const noexcept
Const access to expression result.
Definition: exprDriver.H:391
void setResult(Field< Type > *, bool wantPointData=false)
Set result field, taking ownership of the pointer.
Definition: exprResultI.H:359
tmp< VolumeField< Type > > newVolField(const Type &val=pTraits< Type >::zero) const
Return a new volume field with the mesh size.
bool isPointData() const noexcept
A point field.
tmp< PointField< Type > > getPointField(const word &fldName, bool getOldTime=false)
Retrieve field (surface field)
tmp< PointField< Type > > newPointField(const Type &val=pTraits< Type >::zero) const
Return a new point field with the mesh nPoints size.
const GeoField * isResultType() const
Test if stored result pointer is the specified type.
tmp< VolumeField< Type > > getVolField(const word &fldName, bool getOldTime=false)
Retrieve field (vol field)
tmp< SurfaceField< Type > > newSurfaceField(const Type &val=pTraits< Type >::zero) const
Return a new surface field with the mesh nInternalFaces size.
tmp< SurfaceField< Type > > getSurfaceField(const word &fldName, bool getOldTime=false)
Retrieve field (surface field)
void setInternalFieldResult(const Field< Type > &fld)
Deep-copy the internalField as a result.
bool isLogical_
A logical (bool-like) field (but actually a scalar)
void setResult(VolumeField< Type > *ptr, bool logical=false)
Set result (vol field)
A traits class, which is primarily used for primitives.
Definition: pTraits.H:59
A topoSetCellSource to select cells with any point or any edge within a given pointSet(s).
Definition: pointToCell.H:178
regIOobject is an abstract class derived from IOobject to handle automatic object registration with t...
Definition: regIOobject.H:76
string & replace(const std::string &s1, const std::string &s2, size_type pos=0)
Definition: string.C:108
A class for managing temporary objects.
Definition: tmp.H:65
Interpolate from cell centres to points (vertices) using inverse distance weighting.
tmp< GeometricField< Type, pointPatchField, pointMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &) const
Interpolate volField using inverse distance weighting.
A class for handling words, derived from Foam::string.
Definition: word.H:68
rDeltaTY field()
dynamicFvMesh & mesh
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
Operations involving expressions.
Interpolates (averages) the vertex values to the cell center.
void assign(Field< Tout > &result, const Field< T1 > &a, const UnaryOp &op)
Populate a field as the result of a unary operation on an input.
Definition: FieldOps.C:36
@ VOLUME_DATA
Volume data.
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
const dimensionSet dimless
Dimensionless.
Type interpolatePointToCell(const GeometricField< Type, pointPatchField, pointMesh > &ptf, const label celli)
error FatalError
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
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333
Convert [0-1] values (usually scalars) as false/true.
Definition: exprOps.H:57