sampledThresholdCellFaces.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 | Copyright (C) 2011-2016 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 
27 #include "dictionary.H"
28 #include "volFields.H"
29 #include "volPointInterpolation.H"
31 #include "fvMesh.H"
32 #include "thresholdCellFaces.H"
33 
34 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
35 
36 namespace Foam
37 {
38  defineTypeNameAndDebug(sampledThresholdCellFaces, 0);
40  (
41  sampledSurface,
42  sampledThresholdCellFaces,
43  word,
44  thresholdCellFaces
45  );
46 }
47 
48 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
49 
50 bool Foam::sampledThresholdCellFaces::updateGeometry() const
51 {
52  const fvMesh& fvm = static_cast<const fvMesh&>(mesh());
53 
54  // no update needed
55  if (fvm.time().timeIndex() == prevTimeIndex_)
56  {
57  return false;
58  }
59 
60  prevTimeIndex_ = fvm.time().timeIndex();
61 
62  // Optionally read volScalarField
63  autoPtr<volScalarField> readFieldPtr_;
64 
65  // 1. see if field in database
66  // 2. see if field can be read
67  const volScalarField* cellFldPtr = NULL;
68  if (fvm.foundObject<volScalarField>(fieldName_))
69  {
70  if (debug)
71  {
72  InfoInFunction<< "Lookup " << fieldName_ << endl;
73  }
74 
75  cellFldPtr = &fvm.lookupObject<volScalarField>(fieldName_);
76  }
77  else
78  {
79  // Bit of a hack. Read field and store.
80 
81  if (debug)
82  {
84  << "Reading " << fieldName_
85  << " from time " << fvm.time().timeName()
86  << endl;
87  }
88 
89  readFieldPtr_.reset
90  (
91  new volScalarField
92  (
93  IOobject
94  (
95  fieldName_,
96  fvm.time().timeName(),
97  fvm,
100  false
101  ),
102  fvm
103  )
104  );
105 
106  cellFldPtr = readFieldPtr_.operator->();
107  }
108  const volScalarField& cellFld = *cellFldPtr;
109 
110 
111  thresholdCellFaces surf
112  (
113  fvm,
114  cellFld.primitiveField(),
115  lowerThreshold_,
116  upperThreshold_,
117  triangulate_
118  );
119 
120  const_cast<sampledThresholdCellFaces&>
121  (
122  *this
124  meshCells_.transfer(surf.meshCells());
125 
126  // clear derived data
128 
129  if (debug)
130  {
131  Pout<< "sampledThresholdCellFaces::updateGeometry() : constructed"
132  << nl
133  << " field : " << fieldName_ << nl
134  << " lowerLimit : " << lowerThreshold_ << nl
135  << " upperLimit : " << upperThreshold_ << nl
136  << " point : " << points().size() << nl
137  << " faces : " << faces().size() << nl
138  << " cut cells : " << meshCells_.size() << endl;
139  }
140 
141  return true;
142 }
143 
144 
145 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
146 
148 (
149  const word& name,
150  const polyMesh& mesh,
151  const dictionary& dict
152 )
153 :
154  sampledSurface(name, mesh, dict),
155  fieldName_(dict.lookup("field")),
156  lowerThreshold_(dict.lookupOrDefault<scalar>("lowerLimit", -VGREAT)),
157  upperThreshold_(dict.lookupOrDefault<scalar>("upperLimit", VGREAT)),
158  zoneKey_(keyType::null),
159  triangulate_(dict.lookupOrDefault("triangulate", false)),
160  prevTimeIndex_(-1),
161  meshCells_(0)
162 {
163  if (!dict.found("lowerLimit") && !dict.found("upperLimit"))
164  {
166  << "require at least one of 'lowerLimit' or 'upperLimit'" << endl
167  << abort(FatalError);
168  }
169 }
170 
171 
172 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
173 
175 {}
176 
177 
178 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
179 
181 {
182  const fvMesh& fvm = static_cast<const fvMesh&>(mesh());
183 
184  return fvm.time().timeIndex() != prevTimeIndex_;
185 }
186 
187 
189 {
190  // already marked as expired
191  if (prevTimeIndex_ == -1)
192  {
193  return false;
194  }
195 
196  // force update
197  prevTimeIndex_ = -1;
198  return true;
199 }
200 
201 
203 {
204  return updateGeometry();
205 }
206 
207 
209 (
210  const volScalarField& vField
211 ) const
212 {
213  return sampleField(vField);
214 }
215 
216 
218 (
219  const volVectorField& vField
220 ) const
221 {
222  return sampleField(vField);
223 }
224 
225 
227 (
228  const volSphericalTensorField& vField
229 ) const
230 {
231  return sampleField(vField);
232 }
233 
234 
236 (
237  const volSymmTensorField& vField
238 ) const
239 {
240  return sampleField(vField);
241 }
242 
243 
245 (
246  const volTensorField& vField
247 ) const
248 {
249  return sampleField(vField);
250 }
251 
252 
254 (
255  const interpolation<scalar>& interpolator
256 ) const
257 {
258  return interpolateField(interpolator);
259 }
260 
261 
263 (
264  const interpolation<vector>& interpolator
265 ) const
266 {
267  return interpolateField(interpolator);
268 }
269 
272 (
273  const interpolation<sphericalTensor>& interpolator
274 ) const
275 {
276  return interpolateField(interpolator);
277 }
278 
279 
281 (
282  const interpolation<symmTensor>& interpolator
283 ) const
284 {
285  return interpolateField(interpolator);
286 }
287 
288 
290 (
291  const interpolation<tensor>& interpolator
292 ) const
293 {
294  return interpolateField(interpolator);
295 }
296 
297 
299 {
300  os << "sampledThresholdCellFaces: " << name() << " :"
301  << " field:" << fieldName_
302  << " lowerLimit:" << lowerThreshold_
303  << " upperLimit:" << upperThreshold_;
304  //<< " faces:" << faces().size() // possibly no geom yet
305  //<< " points:" << points().size();
306 }
307 
308 
309 // ************************************************************************* //
virtual ~sampledThresholdCellFaces()
Destructor.
label timeIndex() const
Return current time index.
Definition: TimeStateI.H:35
virtual bool update()
Update the surface as required.
virtual void clearGeom() const
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
An abstract class for surfaces with sampling.
virtual tmp< scalarField > sample(const volScalarField &) const
Sample field on surface.
virtual bool needsUpdate() const
Does the surface need an update?
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
virtual void print(Ostream &) const
Write.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
void transfer(MeshedSurface< face > &)
Transfer the contents of the argument and annul the argument.
Macros for easy insertion into run-time selection tables.
sampledThresholdCellFaces(const word &name, const polyMesh &, const dictionary &)
Construct from dictionary.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
dynamicFvMesh & mesh
const pointField & points
A class for handling words, derived from string.
Definition: word.H:59
virtual bool expire()
Mark the surface as needing an update.
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:306
errorManip< error > abort(error &err)
Definition: errorManip.H:131
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
prefixOSstream Pout(cout,"Pout")
Definition: IOstreams.H:53
static const char nl
Definition: Ostream.H:262
defineTypeNameAndDebug(combustionModel, 0)
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
addNamedToRunTimeSelectionTable(GAMGProcAgglomeration, noneGAMGProcAgglomeration, GAMGAgglomeration, none)
bool interpolate() const
Interpolation requested for surface.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
A class for managing temporary objects.
Definition: PtrList.H:54
static const keyType null
An empty keyType.
Definition: keyType.H:83
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
Namespace for OpenFOAM.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:451
#define InfoInFunction
Report an information message using Foam::Info.