All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
sampledPlane.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-2020 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 "sampledPlane.H"
28 
29 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
30 
31 namespace Foam
32 {
33 namespace sampledSurfaces
34 {
35  defineTypeNameAndDebug(plane, 0);
36  addToRunTimeSelectionTable(sampledSurface, plane, word);
37 }
38 }
39 
40 
41 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
42 
44 (
45  const word& name,
46  const polyMesh& mesh,
47  const dictionary& dict
48 )
49 :
50  sampledSurface(name, mesh, dict),
52  zoneKey_(keyType::null),
53  triangulate_(dict.lookupOrDefault("triangulate", true)),
54  needsUpdate_(true)
55 {
56  // Make plane relative to the coordinateSystem (Cartesian)
57  // allow lookup from global coordinate systems
58  if (dict.found("coordinateSystem"))
59  {
60  coordinateSystem cs(mesh, dict.subDict("coordinateSystem"));
61 
62  point base = cs.globalPosition(planeDesc().refPoint());
63  vector norm = cs.globalVector(planeDesc().normal());
64 
65  // Assign the plane description
66  static_cast<Foam::plane&>(*this) = Foam::plane(base, norm);
67  }
68 
69  dict.readIfPresent("zone", zoneKey_);
70 
71  if (debug && zoneKey_.size() && mesh.cellZones().findIndex(zoneKey_) < 0)
72  {
73  Info<< "cellZone " << zoneKey_
74  << " not found - using entire mesh" << endl;
75  }
76 }
77 
78 
79 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
80 
82 {}
83 
84 
85 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
86 
88 {
89  return needsUpdate_;
90 }
91 
92 
94 {
95  // Already marked as expired
96  if (needsUpdate_)
97  {
98  return false;
99  }
100 
102 
103  needsUpdate_ = true;
104  return true;
105 }
106 
107 
109 {
110  if (!needsUpdate_)
111  {
112  return false;
113  }
114 
116 
117  labelList selectedCells = mesh().cellZones().findMatching(zoneKey_).used();
118 
119  if (selectedCells.empty())
120  {
121  reCut(mesh(), triangulate_);
122  }
123  else
124  {
125  reCut(mesh(), triangulate_, selectedCells);
126  }
127 
128  if (debug)
129  {
130  print(Pout);
131  Pout<< endl;
132  }
133 
134  needsUpdate_ = false;
135  return true;
136 }
137 
138 
140 (
141  const volScalarField& vField
142 ) const
143 {
144  return sampleField(vField);
145 }
146 
147 
149 (
150  const volVectorField& vField
151 ) const
152 {
153  return sampleField(vField);
154 }
155 
156 
158 (
159  const volSphericalTensorField& vField
160 ) const
161 {
162  return sampleField(vField);
163 }
164 
165 
167 (
168  const volSymmTensorField& vField
169 ) const
170 {
171  return sampleField(vField);
172 }
173 
174 
176 (
177  const volTensorField& vField
178 ) const
179 {
180  return sampleField(vField);
181 }
182 
183 
185 (
186  const interpolation<scalar>& interpolator
187 ) const
188 {
189  return interpolateField(interpolator);
190 }
191 
192 
194 (
195  const interpolation<vector>& interpolator
196 ) const
197 {
198  return interpolateField(interpolator);
199 }
200 
202 (
203  const interpolation<sphericalTensor>& interpolator
204 ) const
205 {
206  return interpolateField(interpolator);
207 }
208 
209 
211 (
212  const interpolation<symmTensor>& interpolator
213 ) const
214 {
215  return interpolateField(interpolator);
216 }
217 
218 
220 (
221  const interpolation<tensor>& interpolator
222 ) const
223 {
224  return interpolateField(interpolator);
225 }
226 
227 
229 {
230  os << "plane: " << name() << " :"
231  << " base:" << refPoint()
232  << " normal:" << normal()
233  << " triangulate:" << triangulate_
234  << " faces:" << faces().size()
235  << " points:" << points().size();
236 }
237 
238 
239 // ************************************************************************* //
label findIndex(const keyType &) const
Return zone index for the first match, return -1 if not found.
Definition: ZoneMesh.C:306
Base class for other coordinate system specifications.
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:667
bool empty() const
Return true if the UList is empty (ie, size() is zero)
Definition: UListI.H:313
const word & name() const
Return the name of this functionObject.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
An abstract class for surfaces with sampling.
defineTypeNameAndDebug(distanceSurface, 0)
virtual bool needsUpdate() const
Does the surface need an update?
Definition: sampledPlane.C:87
bool interpolate() const
Interpolation requested for surface.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
Geometric class that creates a 2D plane and can return the intersection point between a line and the ...
Definition: plane.H:60
virtual bool expire()
Mark the surface as needing an update.
Definition: sampledPlane.C:93
Macros for easy insertion into run-time selection tables.
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:934
A sampledSurface defined by a plane using the iso-surface algorithm to &#39;cut&#39; the mesh.
virtual ~plane()
Destructor.
Definition: sampledPlane.C:81
dynamicFvMesh & mesh
plane(const word &name, const polyMesh &mesh, const dictionary &dict)
Construct from dictionary.
Definition: sampledPlane.C:44
virtual void clearGeom() const
const pointField & points
A class for handling words, derived from string.
Definition: word.H:59
const cellZoneMesh & cellZones() const
Return cell zone mesh.
Definition: polyMesh.H:482
bool readIfPresent(const word &, T &, bool recursive=false, bool patternMatch=true) const
Find an entry if present, and assign to T.
virtual bool update()
Update the surface as required.
Definition: sampledPlane.C:108
virtual void print(Ostream &) const
Write.
Definition: sampledPlane.C:228
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
virtual tmp< scalarField > sample(const volScalarField &) const
Sample field on surface.
Definition: sampledPlane.C:140
addToRunTimeSelectionTable(sampledSurface, distanceSurface, word)
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
prefixOSstream Pout(cout, "Pout")
Definition: IOstreams.H:53
messageStream Info
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
A class for managing temporary objects.
Definition: PtrList.H:53
static const keyType null
An empty keyType.
Definition: keyType.H:92
Namespace for OpenFOAM.