FreeFOAM The Cross-Platform CFD Toolkit
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) 1991-2010 OpenCFD Ltd.
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 
28 #include <OpenFOAM/dictionary.H>
29 #include <finiteVolume/volFields.H>
32 #include <finiteVolume/fvMesh.H>
33 #include "thresholdCellFaces.H"
34 
35 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
36 
37 namespace Foam
38 {
39  defineTypeNameAndDebug(sampledThresholdCellFaces, 0);
41  (
42  sampledSurface,
43  sampledThresholdCellFaces,
44  word,
45  thresholdCellFaces
46  );
47 }
48 
49 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
50 
51 bool Foam::sampledThresholdCellFaces::updateGeometry() const
52 {
53  const fvMesh& fvm = static_cast<const fvMesh&>(mesh());
54 
55  // no update needed
56  if (fvm.time().timeIndex() == prevTimeIndex_)
57  {
58  return false;
59  }
60 
61  prevTimeIndex_ = fvm.time().timeIndex();
62 
63  // Optionally read volScalarField
64  autoPtr<volScalarField> readFieldPtr_;
65 
66  // 1. see if field in database
67  // 2. see if field can be read
68  const volScalarField* cellFldPtr = NULL;
69  if (fvm.foundObject<volScalarField>(fieldName_))
70  {
71  if (debug)
72  {
73  Info<< "sampledThresholdCellFaces::updateGeometry() : lookup "
74  << fieldName_ << endl;
75  }
76 
77  cellFldPtr = &fvm.lookupObject<volScalarField>(fieldName_);
78  }
79  else
80  {
81  // Bit of a hack. Read field and store.
82 
83  if (debug)
84  {
85  Info<< "sampledThresholdCellFaces::updateGeometry() : reading "
86  << fieldName_ << " from time " << fvm.time().timeName()
87  << endl;
88  }
89 
90  readFieldPtr_.reset
91  (
92  new volScalarField
93  (
94  IOobject
95  (
96  fieldName_,
97  fvm.time().timeName(),
98  fvm,
101  false
102  ),
103  fvm
104  )
105  );
106 
107  cellFldPtr = readFieldPtr_.operator->();
108  }
109  const volScalarField& cellFld = *cellFldPtr;
110 
111 
112  thresholdCellFaces surf
113  (
114  fvm,
115  cellFld.internalField(),
116  lowerThreshold_,
117  upperThreshold_,
118  triangulate_
119  );
120 
121  const_cast<sampledThresholdCellFaces&>
122  (
123  *this
125  meshCells_.transfer(surf.meshCells());
126 
127  // clear derived data
129 
130  if (debug)
131  {
132  Pout<< "sampledThresholdCellFaces::updateGeometry() : constructed"
133  << nl
134  << " field : " << fieldName_ << nl
135  << " lowerLimit : " << lowerThreshold_ << nl
136  << " upperLimit : " << upperThreshold_ << nl
137  << " point : " << points().size() << nl
138  << " faces : " << faces().size() << nl
139  << " cut cells : " << meshCells_.size() << endl;
140  }
141 
142  return true;
143 }
144 
145 
146 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
147 
149 (
150  const word& name,
151  const polyMesh& mesh,
152  const dictionary& dict
153 )
154 :
155  sampledSurface(name, mesh, dict),
156  fieldName_(dict.lookup("field")),
157  lowerThreshold_(dict.lookupOrDefault<scalar>("lowerLimit", -VGREAT)),
158  upperThreshold_(dict.lookupOrDefault<scalar>("upperLimit", VGREAT)),
159  zoneName_(word::null),
160  triangulate_(dict.lookupOrDefault("triangulate", false)),
161  prevTimeIndex_(-1),
162  meshCells_(0)
163 {
164  if (!dict.found("lowerLimit") && !dict.found("upperLimit"))
165  {
167  (
168  "sampledThresholdCellFaces::sampledThresholdCellFaces(..)"
169  )
170  << "require at least one of 'lowerLimit' or 'upperLimit'" << endl
171  << abort(FatalError);
172  }
173 
174 
175 // dict.readIfPresent("zone", zoneName_);
176 //
177 // if (debug && zoneName_.size())
178 // {
179 // if (mesh.cellZones().findZoneID(zoneName_) < 0)
180 // {
181 // Info<< "cellZone \"" << zoneName_
182 // << "\" not found - using entire mesh" << endl;
183 // }
184 // }
185 }
186 
187 
188 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
189 
191 {}
192 
193 
194 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
195 
197 {
198  const fvMesh& fvm = static_cast<const fvMesh&>(mesh());
199 
200  return fvm.time().timeIndex() != prevTimeIndex_;
201 }
202 
203 
205 {
206  // already marked as expired
207  if (prevTimeIndex_ == -1)
208  {
209  return false;
210  }
211 
212  // force update
213  prevTimeIndex_ = -1;
214  return true;
215 }
216 
217 
219 {
220  return updateGeometry();
221 }
222 
223 
226 (
227  const volScalarField& vField
228 ) const
229 {
230  return sampleField(vField);
231 }
232 
233 
236 (
237  const volVectorField& vField
238 ) const
239 {
240  return sampleField(vField);
241 }
242 
243 
246 (
247  const volSphericalTensorField& vField
248 ) const
249 {
250  return sampleField(vField);
251 }
252 
253 
256 (
257  const volSymmTensorField& vField
258 ) const
259 {
260  return sampleField(vField);
261 }
262 
263 
266 (
267  const volTensorField& vField
268 ) const
269 {
270  return sampleField(vField);
271 }
272 
273 
276 (
277  const interpolation<scalar>& interpolator
278 ) const
279 {
280  return interpolateField(interpolator);
281 }
282 
283 
286 (
287  const interpolation<vector>& interpolator
288 ) const
289 {
290  return interpolateField(interpolator);
291 }
292 
295 (
296  const interpolation<sphericalTensor>& interpolator
297 ) const
298 {
299  return interpolateField(interpolator);
300 }
301 
302 
305 (
306  const interpolation<symmTensor>& interpolator
307 ) const
308 {
309  return interpolateField(interpolator);
310 }
311 
312 
315 (
316  const interpolation<tensor>& interpolator
317 ) const
318 {
319  return interpolateField(interpolator);
320 }
321 
322 
324 {
325  os << "sampledThresholdCellFaces: " << name() << " :"
326  << " field:" << fieldName_
327  << " lowerLimit:" << lowerThreshold_
328  << " upperLimit:" << upperThreshold_;
329  //<< " faces:" << faces().size() // possibly no geom yet
330  //<< " points:" << points().size();
331 }
332 
333 
334 // ************************ vim: set sw=4 sts=4 et: ************************ //