All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
sampledCuttingPlane.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2019 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "sampledCuttingPlane.H"
28 
29 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
30 
31 namespace Foam
32 {
33 namespace sampledSurfaces
34 {
35  defineTypeNameAndDebug(cuttingPlane, 0);
36  addToRunTimeSelectionTable(sampledSurface, cuttingPlane, word);
37 }
38 }
39 
40 
41 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
42 
43 void Foam::sampledSurfaces::cuttingPlane::createGeometry()
44 {
45  if (debug)
46  {
47  Pout<< "cuttingPlane::createGeometry :updating geometry."
48  << endl;
49  }
50 
51  // Clear any stored topologies
52  facesPtr_.clear();
53  isoSurfPtr_.ptr();
54  pointDistance_.clear();
55  cellDistancePtr_.clear();
56 
57  // Clear derived data
58  clearGeom();
59 
60  // Get any subMesh
61  if (zoneID_.index() != -1 && !subMeshPtr_.valid())
62  {
63  const polyBoundaryMesh& patches = mesh().boundaryMesh();
64 
65  // Patch to put exposed internal faces into
66  const label exposedPatchi = patches.findPatchID(exposedPatchName_);
67 
68  DebugInfo
69  << "Allocating subset of size "
70  << mesh().cellZones()[zoneID_.index()].size()
71  << " with exposed faces into patch "
72  << patches[exposedPatchi].name() << endl;
73 
74  subMeshPtr_.reset
75  (
76  new fvMeshSubset(static_cast<const fvMesh&>(mesh()))
77  );
78  subMeshPtr_().setLargeCellSubset
79  (
80  labelHashSet(mesh().cellZones()[zoneID_.index()]),
81  exposedPatchi
82  );
83  }
84 
85 
86  // Select either the submesh or the underlying mesh
87  const fvMesh& mesh =
88  (
89  subMeshPtr_.valid()
90  ? subMeshPtr_().subMesh()
91  : static_cast<const fvMesh&>(this->mesh())
92  );
93 
94 
95  // Distance to cell centres
96  // ~~~~~~~~~~~~~~~~~~~~~~~~
97 
98  cellDistancePtr_.reset
99  (
100  new volScalarField
101  (
102  IOobject
103  (
104  "cellDistance",
105  mesh.time().timeName(),
106  mesh.time(),
109  false
110  ),
111  mesh,
113  )
114  );
115  volScalarField& cellDistance = cellDistancePtr_();
116 
117  // Internal field
118  {
119  const pointField& cc = mesh.cellCentres();
120  scalarField& fld = cellDistance.primitiveFieldRef();
121 
122  forAll(cc, i)
123  {
124  // Signed distance
125  fld[i] = (cc[i] - plane_.refPoint()) & plane_.normal();
126  }
127  }
128 
129  volScalarField::Boundary& cellDistanceBf =
130  cellDistance.boundaryFieldRef();
131 
132  // Patch fields
133  {
134  forAll(cellDistanceBf, patchi)
135  {
136  if
137  (
138  isA<emptyFvPatchScalarField>
139  (
140  cellDistanceBf[patchi]
141  )
142  )
143  {
144  cellDistanceBf.set
145  (
146  patchi,
147  new calculatedFvPatchScalarField
148  (
149  mesh.boundary()[patchi],
150  cellDistance
151  )
152  );
153 
154  const polyPatch& pp = mesh.boundary()[patchi].patch();
155  pointField::subField cc = pp.patchSlice(mesh.faceCentres());
156 
157  fvPatchScalarField& fld = cellDistanceBf[patchi];
158  fld.setSize(pp.size());
159  forAll(fld, i)
160  {
161  fld[i] = (cc[i] - plane_.refPoint()) & plane_.normal();
162  }
163  }
164  else
165  {
166  const pointField& cc = mesh.C().boundaryField()[patchi];
167  fvPatchScalarField& fld = cellDistanceBf[patchi];
168 
169  forAll(fld, i)
170  {
171  fld[i] = (cc[i] - plane_.refPoint()) & plane_.normal();
172  }
173  }
174  }
175  }
176 
177 
178  // On processor patches the mesh.C() will already be the cell centre
179  // on the opposite side so no need to swap cellDistance.
180 
181 
182  // Distance to points
183  pointDistance_.setSize(mesh.nPoints());
184  {
185  const pointField& pts = mesh.points();
186 
187  forAll(pointDistance_, i)
188  {
189  pointDistance_[i] = (pts[i] - plane_.refPoint()) & plane_.normal();
190  }
191  }
192 
193 
194  if (debug)
195  {
196  Pout<< "Writing cell distance:" << cellDistance.objectPath() << endl;
197  cellDistance.write();
198  pointScalarField pDist
199  (
200  IOobject
201  (
202  "pointDistance",
203  mesh.time().timeName(),
204  mesh.time(),
207  false
208  ),
209  pointMesh::New(mesh),
211  );
212  pDist.primitiveFieldRef() = pointDistance_;
213 
214  Pout<< "Writing point distance:" << pDist.objectPath() << endl;
215  pDist.write();
216  }
217 
218 
219  //- Direct from cell field and point field.
220  isoSurfPtr_.reset
221  (
222  new isoSurface
223  (
224  mesh,
225  cellDistance,
226  pointDistance_,
227  0,
228  filter_
229  )
230  );
231 
232  if (debug)
233  {
234  print(Pout);
235  Pout<< endl;
236  }
237 }
238 
239 
240 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
241 
243 (
244  const word& name,
245  const polyMesh& mesh,
246  const dictionary& dict
247 )
248 :
249  sampledSurface(name, mesh, dict),
250  plane_(dict),
251  filter_
252  (
253  dict.found("filtering")
254  ? isoSurface::filterTypeNames_.read(dict.lookup("filtering"))
256  ),
257  average_(dict.lookupOrDefault("average", false)),
258  zoneID_(dict.lookupOrDefault("zone", word::null), mesh.cellZones()),
259  exposedPatchName_(word::null),
260  needsUpdate_(true),
261  subMeshPtr_(nullptr),
262  cellDistancePtr_(nullptr),
263  isoSurfPtr_(nullptr)
264 {
265  if (zoneID_.index() != -1)
266  {
267  dict.lookup("exposedPatchName") >> exposedPatchName_;
268 
269  if (mesh.boundaryMesh().findPatchID(exposedPatchName_) == -1)
270  {
272  << "Cannot find patch " << exposedPatchName_
273  << " in which to put exposed faces." << endl
274  << "Valid patches are " << mesh.boundaryMesh().names()
275  << exit(FatalError);
276  }
277 
278  if (debug && zoneID_.index() != -1)
279  {
280  Info<< "Restricting to cellZone " << zoneID_.name()
281  << " with exposed internal faces into patch "
282  << exposedPatchName_ << endl;
283  }
284  }
285 }
286 
287 
288 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
289 
291 {}
292 
293 
294 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
295 
297 {
298  return needsUpdate_;
299 }
300 
301 
303 {
304  if (debug)
305  {
306  Pout<< "cuttingPlane::expire :"
307  << " needsUpdate_:" << needsUpdate_ << endl;
308  }
309 
310  // Clear derived data
311  clearGeom();
312 
313  // Already marked as expired
314  if (needsUpdate_)
315  {
316  return false;
317  }
318 
319  needsUpdate_ = true;
320  return true;
321 }
322 
323 
325 {
326  if (debug)
327  {
328  Pout<< "cuttingPlane::update :"
329  << " needsUpdate_:" << needsUpdate_ << endl;
330  }
331 
332  if (!needsUpdate_)
333  {
334  return false;
335  }
336 
337  createGeometry();
338 
339  needsUpdate_ = false;
340  return true;
341 }
342 
343 
346 (
347  const volScalarField& vField
348 ) const
349 {
350  return sampleField(vField);
351 }
352 
353 
356 (
357  const volVectorField& vField
358 ) const
359 {
360  return sampleField(vField);
361 }
362 
363 
366 (
367  const volSphericalTensorField& vField
368 ) const
369 {
370  return sampleField(vField);
371 }
372 
373 
376 (
377  const volSymmTensorField& vField
378 ) const
379 {
380  return sampleField(vField);
381 }
382 
383 
386 (
387  const volTensorField& vField
388 ) const
389 {
390  return sampleField(vField);
391 }
392 
393 
396 (
397  const interpolation<scalar>& interpolator
398 ) const
399 {
400  return interpolateField(interpolator);
401 }
402 
403 
406 (
407  const interpolation<vector>& interpolator
408 ) const
409 {
410  return interpolateField(interpolator);
411 }
412 
415 (
416  const interpolation<sphericalTensor>& interpolator
417 ) const
418 {
419  return interpolateField(interpolator);
420 }
421 
422 
425 (
426  const interpolation<symmTensor>& interpolator
427 ) const
428 {
429  return interpolateField(interpolator);
430 }
431 
432 
435 (
436  const interpolation<tensor>& interpolator
437 ) const
438 {
439  return interpolateField(interpolator);
440 }
441 
442 
444 {
445  os << "cuttingPlane: " << name() << " :"
446  << " plane:" << plane_
447  << " faces:" << faces().size()
448  << " points:" << points().size();
449 }
450 
451 
452 // ************************************************************************* //
GeometricField< scalar, pointPatchField, pointMesh > pointScalarField
Definition: pointFields.H:49
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:434
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:667
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
virtual Ostream & write(const char)=0
Write character.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
cuttingPlane(const word &name, const polyMesh &mesh, const dictionary &dict)
Construct from dictionary.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
const word & name() const
Return the name of this functionObject.
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
An abstract class for surfaces with sampling.
virtual void print(Ostream &) const
Write.
defineTypeNameAndDebug(distanceSurface, 0)
virtual bool update()
Update the surface as required.
bool interpolate() const
Interpolation requested for surface.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
label findPatchID(const word &patchName) const
Find patch index given a name.
Macros for easy insertion into run-time selection tables.
static const pointMesh & New(const polyMesh &mesh)
Definition: MeshObject.C:44
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:57
HashSet< label, Hash< label > > labelHashSet
A HashSet with label keys.
Definition: HashSet.H:211
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
dynamicFvMesh & mesh
const pointField & points
fvPatchField< scalar > fvPatchScalarField
A class for handling words, derived from string.
Definition: word.H:59
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
const cellZoneMesh & cellZones() const
Return cell zone mesh.
Definition: polyMesh.H:482
wordList names() const
Return a list of patch names.
static const NamedEnum< filterType, 3 > filterTypeNames_
Definition: isoSurface.H:73
static const word null
An empty word.
Definition: word.H:77
virtual tmp< scalarField > sample(const volScalarField &) const
Sample field on surface.
#define DebugInfo
Report an information message using Foam::Info.
void reset(const label nPoints, const label nInternalFaces, const label nFaces, const label nCells)
Reset this primitiveMesh given the primitive array sizes.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
Internal::FieldType & primitiveFieldRef()
Return a reference to the internal field.
addToRunTimeSelectionTable(sampledSurface, distanceSurface, word)
virtual bool needsUpdate() const
Does the surface need an update?
label patchi
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:50
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
prefixOSstream Pout(cout, "Pout")
Definition: IOstreams.H:53
messageStream Info
SubField< Type > subField
Declare type of subField.
Definition: Field.H:100
Enum read(Istream &) const
Read a word from Istream and return the corresponding.
Definition: NamedEnum.C:61
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
A class for managing temporary objects.
Definition: PtrList.H:53
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:812
virtual bool expire()
Mark the surface as needing an update.