FreeFOAM The Cross-Platform CFD Toolkit
sampledIsoSurfaceCell.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 
26 #include "sampledIsoSurfaceCell.H"
27 #include <OpenFOAM/dictionary.H>
28 #include <finiteVolume/volFields.H>
31 #include <finiteVolume/fvMesh.H>
32 #include "isoSurfaceCell.H"
33 
34 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
35 
36 namespace Foam
37 {
38  defineTypeNameAndDebug(sampledIsoSurfaceCell, 0);
39  addNamedToRunTimeSelectionTable(sampledSurface, sampledIsoSurfaceCell, word, isoSurfaceCell);
40 }
41 
42 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
43 
44 bool Foam::sampledIsoSurfaceCell::updateGeometry() const
45 {
46  const fvMesh& fvm = static_cast<const fvMesh&>(mesh());
47 
48  // no update needed
49  if (fvm.time().timeIndex() == prevTimeIndex_)
50  {
51  return false;
52  }
53 
54  prevTimeIndex_ = fvm.time().timeIndex();
55 
56  // Clear any stored topo
57  facesPtr_.clear();
58 
59  // Clear derived data
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>(isoField_))
69  {
70  if (debug)
71  {
72  Info<< "sampledIsoSurfaceCell::updateGeometry() : lookup "
73  << isoField_ << endl;
74  }
75 
76  cellFldPtr = &fvm.lookupObject<volScalarField>(isoField_);
77  }
78  else
79  {
80  // Bit of a hack. Read field and store.
81 
82  if (debug)
83  {
84  Info<< "sampledIsoSurfaceCell::updateGeometry() : reading "
85  << isoField_ << " from time " <<fvm.time().timeName()
86  << endl;
87  }
88 
89  readFieldPtr_.reset
90  (
91  new volScalarField
92  (
93  IOobject
94  (
95  isoField_,
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  tmp<pointScalarField> pointFld
111  (
113  );
114 
115  if (average_)
116  {
117  //- From point field and interpolated cell.
118  scalarField cellAvg(fvm.nCells(), scalar(0.0));
119  labelField nPointCells(fvm.nCells(), 0);
120  {
121  for (label pointI = 0; pointI < fvm.nPoints(); pointI++)
122  {
123  const labelList& pCells = fvm.pointCells(pointI);
124 
125  forAll(pCells, i)
126  {
127  label cellI = pCells[i];
128 
129  cellAvg[cellI] += pointFld().internalField()[pointI];
130  nPointCells[cellI]++;
131  }
132  }
133  }
134  forAll(cellAvg, cellI)
135  {
136  cellAvg[cellI] /= nPointCells[cellI];
137  }
138 
139  const isoSurfaceCell iso
140  (
141  fvm,
142  cellAvg,
143  pointFld().internalField(),
144  isoVal_,
145  regularise_
146  );
147 
148  const_cast<sampledIsoSurfaceCell&>
149  (
150  *this
151  ).triSurface::operator=(iso);
152  meshCells_ = iso.meshCells();
153  }
154  else
155  {
156  //- Direct from cell field and point field. Gives bad continuity.
157  const isoSurfaceCell iso
158  (
159  fvm,
160  cellFld.internalField(),
161  pointFld().internalField(),
162  isoVal_,
163  regularise_
164  );
165 
166  const_cast<sampledIsoSurfaceCell&>
167  (
168  *this
169  ).triSurface::operator=(iso);
170  meshCells_ = iso.meshCells();
171  }
172 
173 
174  if (debug)
175  {
176  Pout<< "sampledIsoSurfaceCell::updateGeometry() : constructed iso:"
177  << nl
178  << " regularise : " << regularise_ << nl
179  << " average : " << average_ << nl
180  << " isoField : " << isoField_ << nl
181  << " isoValue : " << isoVal_ << nl
182  << " points : " << points().size() << nl
183  << " tris : " << triSurface::size() << nl
184  << " cut cells : " << meshCells_.size() << endl;
185  }
186 
187  return true;
188 }
189 
190 
191 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
192 
194 (
195  const word& name,
196  const polyMesh& mesh,
197  const dictionary& dict
198 )
199 :
200  sampledSurface(name, mesh, dict),
201  isoField_(dict.lookup("isoField")),
202  isoVal_(readScalar(dict.lookup("isoValue"))),
203  regularise_(dict.lookupOrDefault("regularise", true)),
204  average_(dict.lookupOrDefault("average", true)),
205  zoneName_(word::null),
206  facesPtr_(NULL),
207  prevTimeIndex_(-1),
208  meshCells_(0)
209 {
210 // dict.readIfPresent("zone", zoneName_);
211 //
212 // if (debug && zoneName_.size())
213 // {
214 // if (mesh.cellZones().findZoneID(zoneName_) < 0)
215 // {
216 // Info<< "cellZone \"" << zoneName_
217 // << "\" not found - using entire mesh" << endl;
218 // }
219 // }
220 }
221 
222 
223 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
224 
226 {}
227 
228 
229 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
230 
232 {
233  const fvMesh& fvm = static_cast<const fvMesh&>(mesh());
234 
235  return fvm.time().timeIndex() != prevTimeIndex_;
236 }
237 
238 
240 {
241  facesPtr_.clear();
242 
243  // Clear derived data
245 
246  // already marked as expired
247  if (prevTimeIndex_ == -1)
248  {
249  return false;
250  }
251 
252  // force update
253  prevTimeIndex_ = -1;
254  return true;
255 }
256 
257 
259 {
260  return updateGeometry();
261 }
262 
263 
266 (
267  const volScalarField& vField
268 ) const
269 {
270  return sampleField(vField);
271 }
272 
273 
276 (
277  const volVectorField& vField
278 ) const
279 {
280  return sampleField(vField);
281 }
282 
283 
286 (
287  const volSphericalTensorField& vField
288 ) const
289 {
290  return sampleField(vField);
291 }
292 
293 
296 (
297  const volSymmTensorField& vField
298 ) const
299 {
300  return sampleField(vField);
301 }
302 
303 
306 (
307  const volTensorField& vField
308 ) const
309 {
310  return sampleField(vField);
311 }
312 
313 
316 (
317  const interpolation<scalar>& interpolator
318 ) const
319 {
320  return interpolateField(interpolator);
321 }
322 
323 
326 (
327  const interpolation<vector>& interpolator
328 ) const
329 {
330  return interpolateField(interpolator);
331 }
332 
335 (
336  const interpolation<sphericalTensor>& interpolator
337 ) const
338 {
339  return interpolateField(interpolator);
340 }
341 
342 
345 (
346  const interpolation<symmTensor>& interpolator
347 ) const
348 {
349  return interpolateField(interpolator);
350 }
351 
352 
355 (
356  const interpolation<tensor>& interpolator
357 ) const
358 {
359  return interpolateField(interpolator);
360 }
361 
362 
364 {
365  os << "sampledIsoSurfaceCell: " << name() << " :"
366  << " field:" << isoField_
367  << " value:" << isoVal_;
368  //<< " faces:" << faces().size() // possibly no geom yet
369  //<< " points:" << points().size();
370 }
371 
372 
373 // ************************ vim: set sw=4 sts=4 et: ************************ //