sampledIsoSurfaceCell.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-2018 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 "sampledIsoSurfaceCell.H"
27 #include "dictionary.H"
28 #include "volFields.H"
29 #include "volPointInterpolation.H"
31 #include "fvMesh.H"
32 #include "isoSurfaceCell.H"
33 
34 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
35 
36 namespace Foam
37 {
38  defineTypeNameAndDebug(sampledIsoSurfaceCell, 0);
40  (
41  sampledSurface,
42  sampledIsoSurfaceCell,
43  word,
44  isoSurfaceCell
45  );
46 }
47 
48 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
49 
50 bool Foam::sampledIsoSurfaceCell::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  // Clear any stored topo
63  facesPtr_.clear();
64 
65  // Clear derived data
67 
68  // Optionally read volScalarField
69  autoPtr<volScalarField> readFieldPtr_;
70 
71  // 1. see if field in database
72  // 2. see if field can be read
73  const volScalarField* cellFldPtr = nullptr;
74  if (fvm.foundObject<volScalarField>(isoField_))
75  {
76  if (debug)
77  {
78  InfoInFunction << "Lookup " << isoField_ << endl;
79  }
80 
81  cellFldPtr = &fvm.lookupObject<volScalarField>(isoField_);
82  }
83  else
84  {
85  // Bit of a hack. Read field and store.
86 
87  if (debug)
88  {
90  << "Reading " << isoField_
91  << " from time " <<fvm.time().timeName()
92  << endl;
93  }
94 
95  readFieldPtr_.reset
96  (
97  new volScalarField
98  (
99  IOobject
100  (
101  isoField_,
102  fvm.time().timeName(),
103  fvm,
106  false
107  ),
108  fvm
109  )
110  );
111 
112  cellFldPtr = readFieldPtr_.operator->();
113  }
114  const volScalarField& cellFld = *cellFldPtr;
115 
116  tmp<pointScalarField> pointFld
117  (
119  );
120 
121  if (average_)
122  {
123  //- From point field and interpolated cell.
124  scalarField cellAvg(fvm.nCells(), scalar(0));
125  labelField nPointCells(fvm.nCells(), 0);
126  {
127  for (label pointi = 0; pointi < fvm.nPoints(); pointi++)
128  {
129  const labelList& pCells = fvm.pointCells(pointi);
130 
131  forAll(pCells, i)
132  {
133  label celli = pCells[i];
134 
135  cellAvg[celli] += pointFld().primitiveField()[pointi];
136  nPointCells[celli]++;
137  }
138  }
139  }
140  forAll(cellAvg, celli)
141  {
142  cellAvg[celli] /= nPointCells[celli];
143  }
144 
145  const isoSurfaceCell iso
146  (
147  fvm,
148  cellAvg,
149  pointFld().primitiveField(),
150  isoVal_,
151  regularise_
152  );
153 
154  const_cast<sampledIsoSurfaceCell&>
155  (
156  *this
157  ).triSurface::operator=(iso);
158  meshCells_ = iso.meshCells();
159  }
160  else
161  {
162  //- Direct from cell field and point field. Gives bad continuity.
163  const isoSurfaceCell iso
164  (
165  fvm,
166  cellFld.primitiveField(),
167  pointFld().primitiveField(),
168  isoVal_,
169  regularise_
170  );
171 
172  const_cast<sampledIsoSurfaceCell&>
173  (
174  *this
175  ).triSurface::operator=(iso);
176  meshCells_ = iso.meshCells();
177  }
178 
179 
180  if (debug)
181  {
182  Pout<< "sampledIsoSurfaceCell::updateGeometry() : constructed iso:"
183  << nl
184  << " regularise : " << regularise_ << nl
185  << " average : " << average_ << nl
186  << " isoField : " << isoField_ << nl
187  << " isoValue : " << isoVal_ << nl
188  << " points : " << points().size() << nl
189  << " tris : " << triSurface::size() << nl
190  << " cut cells : " << meshCells_.size() << endl;
191  }
192 
193  return true;
194 }
195 
196 
197 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
198 
200 (
201  const word& name,
202  const polyMesh& mesh,
203  const dictionary& dict
204 )
205 :
206  sampledSurface(name, mesh, dict),
207  isoField_(dict.lookup("isoField")),
208  isoVal_(readScalar(dict.lookup("isoValue"))),
209  regularise_(dict.lookupOrDefault("regularise", true)),
210  average_(dict.lookupOrDefault("average", true)),
211  zoneKey_(keyType::null),
212  facesPtr_(nullptr),
213  prevTimeIndex_(-1),
214  meshCells_(0)
215 {}
216 
217 
218 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
219 
221 {}
222 
223 
224 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
225 
227 {
228  const fvMesh& fvm = static_cast<const fvMesh&>(mesh());
229 
230  return fvm.time().timeIndex() != prevTimeIndex_;
231 }
232 
233 
235 {
236  facesPtr_.clear();
237 
238  // Clear derived data
240 
241  // already marked as expired
242  if (prevTimeIndex_ == -1)
243  {
244  return false;
245  }
246 
247  // force update
248  prevTimeIndex_ = -1;
249  return true;
250 }
251 
252 
254 {
255  return updateGeometry();
256 }
257 
258 
261 (
262  const volScalarField& vField
263 ) const
264 {
265  return sampleField(vField);
266 }
267 
268 
271 (
272  const volVectorField& vField
273 ) const
274 {
275  return sampleField(vField);
276 }
277 
278 
281 (
282  const volSphericalTensorField& vField
283 ) const
284 {
285  return sampleField(vField);
286 }
287 
288 
291 (
292  const volSymmTensorField& vField
293 ) const
294 {
295  return sampleField(vField);
296 }
297 
298 
301 (
302  const volTensorField& vField
303 ) const
304 {
305  return sampleField(vField);
306 }
307 
308 
311 (
312  const interpolation<scalar>& interpolator
313 ) const
314 {
315  return interpolateField(interpolator);
316 }
317 
318 
321 (
322  const interpolation<vector>& interpolator
323 ) const
324 {
325  return interpolateField(interpolator);
326 }
327 
330 (
331  const interpolation<sphericalTensor>& interpolator
332 ) const
333 {
334  return interpolateField(interpolator);
335 }
336 
337 
340 (
341  const interpolation<symmTensor>& interpolator
342 ) const
343 {
344  return interpolateField(interpolator);
345 }
346 
347 
350 (
351  const interpolation<tensor>& interpolator
352 ) const
353 {
354  return interpolateField(interpolator);
355 }
356 
357 
359 {
360  os << "sampledIsoSurfaceCell: " << name() << " :"
361  << " field:" << isoField_
362  << " value:" << isoVal_;
363  //<< " faces:" << faces().size() // possibly no geom yet
364  //<< " points:" << points().size();
365 }
366 
367 
368 // ************************************************************************* //
Field< label > labelField
Specialisation of Field<T> for label.
Definition: labelField.H:49
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
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
friend Ostream & operator(Ostream &, const UList< T > &)
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
An abstract class for surfaces with sampling.
bool interpolate(const vector &p1, const vector &p2, const vector &o, vector &n, scalar l)
Definition: curveTools.C:75
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
virtual void print(Ostream &) const
Write.
bool interpolate() const
Interpolation requested for surface.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243
virtual bool update()
Update the surface as required.
Macros for easy insertion into run-time selection tables.
static const volPointInterpolation & New(const fvMesh &mesh)
Definition: MeshObject.C:44
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
virtual tmp< scalarField > sample(const volScalarField &) const
Sample field on surface.
dynamicFvMesh & mesh
virtual void clearGeom() const
const pointField & points
A class for handling words, derived from string.
Definition: word.H:59
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
virtual ~sampledIsoSurfaceCell()
Destructor.
List< label > labelList
A List of labels.
Definition: labelList.H:56
bool readScalar(const char *buf, doubleScalar &s)
Read whole of buf as a scalar. Return true if successful.
Definition: doubleScalar.H:68
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
static const char nl
Definition: Ostream.H:265
defineTypeNameAndDebug(combustionModel, 0)
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
addNamedToRunTimeSelectionTable(GAMGProcAgglomeration, noneGAMGProcAgglomeration, GAMGAgglomeration, none)
label timeIndex() const
Return current time index.
Definition: TimeStateI.H:35
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
sampledIsoSurfaceCell(const word &name, const polyMesh &mesh, const dictionary &dict)
Construct from dictionary.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
virtual bool expire()
Mark the surface as needing an update.
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:83
label size() const
Return the number of elements in the UList.
Definition: ListI.H:170
virtual bool needsUpdate() const
Does the surface need an update?
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:576
#define InfoInFunction
Report an information message using Foam::Info.