distanceSurface.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-2019 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 "distanceSurface.H"
28 
29 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
30 
31 namespace Foam
32 {
33 namespace sampledSurfaces
34 {
37 }
38 }
39 
40 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
41 
42 void Foam::sampledSurfaces::distanceSurface::createGeometry()
43 {
44  if (debug)
45  {
46  Pout<< "distanceSurface::createGeometry :updating geometry." << endl;
47  }
48 
49  // Clear any stored topologies
50  isoSurfPtr_.clear();
51 
52  // Clear derived data
53  clearGeom();
54 
55  const fvMesh& mesh = static_cast<const fvMesh&>(this->mesh());
56 
57  // Distance to cell centres
58  // ~~~~~~~~~~~~~~~~~~~~~~~~
59 
60  cellDistancePtr_.reset
61  (
62  new volScalarField
63  (
64  IOobject
65  (
66  "cellDistance",
67  mesh.time().timeName(),
68  mesh.time(),
71  false
72  ),
73  mesh,
75  )
76  );
77  volScalarField& cellDistance = cellDistancePtr_();
78 
79  // Internal field
80  {
81  const pointField& cc = mesh.C();
82  scalarField& fld = cellDistance.primitiveFieldRef();
83 
84  List<pointIndexHit> nearest;
85  surfPtr_().findNearest
86  (
87  cc,
88  scalarField(cc.size(), great),
89  nearest
90  );
91 
92  if (signed_)
93  {
94  List<volumeType> volType;
95 
96  surfPtr_().getVolumeType(cc, volType);
97 
98  forAll(volType, i)
99  {
100  volumeType vT = volType[i];
101 
102  if (vT == volumeType::outside)
103  {
104  fld[i] = Foam::mag(cc[i] - nearest[i].hitPoint());
105  }
106  else if (vT == volumeType::inside)
107  {
108  fld[i] = -Foam::mag(cc[i] - nearest[i].hitPoint());
109  }
110  else
111  {
113  << "getVolumeType failure, neither INSIDE or OUTSIDE"
114  << exit(FatalError);
115  }
116  }
117  }
118  else
119  {
120  forAll(nearest, i)
121  {
122  fld[i] = Foam::mag(cc[i] - nearest[i].hitPoint());
123  }
124  }
125  }
126 
127  volScalarField::Boundary& cellDistanceBf =
128  cellDistance.boundaryFieldRef();
129 
130  // Patch fields
131  {
132  forAll(mesh.C().boundaryField(), patchi)
133  {
134  const pointField& cc = mesh.C().boundaryField()[patchi];
135  fvPatchScalarField& fld = cellDistanceBf[patchi];
136 
137  List<pointIndexHit> nearest;
138  surfPtr_().findNearest
139  (
140  cc,
141  scalarField(cc.size(), great),
142  nearest
143  );
144 
145  if (signed_)
146  {
147  List<volumeType> volType;
148 
149  surfPtr_().getVolumeType(cc, volType);
150 
151  forAll(volType, i)
152  {
153  volumeType vT = volType[i];
154 
155  if (vT == volumeType::outside)
156  {
157  fld[i] = Foam::mag(cc[i] - nearest[i].hitPoint());
158  }
159  else if (vT == volumeType::inside)
160  {
161  fld[i] = -Foam::mag(cc[i] - nearest[i].hitPoint());
162  }
163  else
164  {
166  << "getVolumeType failure, "
167  << "neither INSIDE or OUTSIDE"
168  << exit(FatalError);
169  }
170  }
171  }
172  else
173  {
174  forAll(nearest, i)
175  {
176  fld[i] = Foam::mag(cc[i] - nearest[i].hitPoint());
177  }
178  }
179  }
180  }
181 
182 
183  // On processor patches the mesh.C() will already be the cell centre
184  // on the opposite side so no need to swap cellDistance.
185 
186 
187  // Distance to points
188  pointDistance_.setSize(mesh.nPoints());
189  {
190  const pointField& pts = mesh.points();
191 
192  List<pointIndexHit> nearest;
193  surfPtr_().findNearest
194  (
195  pts,
196  scalarField(pts.size(), great),
197  nearest
198  );
199 
200  if (signed_)
201  {
202  List<volumeType> volType;
203 
204  surfPtr_().getVolumeType(pts, volType);
205 
206  forAll(volType, i)
207  {
208  volumeType vT = volType[i];
209 
210  if (vT == volumeType::outside)
211  {
212  pointDistance_[i] =
213  Foam::mag(pts[i] - nearest[i].hitPoint());
214  }
215  else if (vT == volumeType::inside)
216  {
217  pointDistance_[i] =
218  -Foam::mag(pts[i] - nearest[i].hitPoint());
219  }
220  else
221  {
223  << "getVolumeType failure, neither INSIDE or OUTSIDE"
224  << exit(FatalError);
225  }
226  }
227  }
228  else
229  {
230  forAll(nearest, i)
231  {
232  pointDistance_[i] = Foam::mag(pts[i]-nearest[i].hitPoint());
233  }
234  }
235  }
236 
237 
238  if (debug)
239  {
240  Pout<< "Writing cell distance:" << cellDistance.objectPath() << endl;
241  cellDistance.write();
242  pointScalarField pDist
243  (
244  IOobject
245  (
246  "pointDistance",
247  mesh.time().timeName(),
248  mesh.time(),
251  false
252  ),
253  pointMesh::New(mesh),
255  );
256  pDist.primitiveFieldRef() = pointDistance_;
257 
258  Pout<< "Writing point distance:" << pDist.objectPath() << endl;
259  pDist.write();
260  }
261 
262 
263  //- Direct from cell field and point field.
264  isoSurfPtr_.reset
265  (
266  new isoSurface
267  (
268  mesh,
269  cellDistance,
270  pointDistance_,
271  distance_,
272  filter_
273  )
274  );
275 
276  if (debug)
277  {
278  print(Pout);
279  Pout<< endl;
280  }
281 }
282 
283 
284 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
285 
287 (
288  const word& name,
289  const polyMesh& mesh,
290  const dictionary& dict
291 )
292 :
293  sampledSurface(name, mesh, dict),
294  surfPtr_
295  (
297  (
298  dict.lookup("surfaceType"),
299  IOobject
300  (
301  dict.lookupOrDefault("surfaceName", name), // name
302  mesh.time().constant(), // directory
303  "triSurface", // instance
304  mesh.time(), // registry
307  ),
308  dict
309  )
310  ),
311  distance_(dict.lookup<scalar>("distance")),
312  signed_(readBool(dict.lookup("signed"))),
313  filter_
314  (
315  dict.found("filtering")
316  ? isoSurface::filterTypeNames_.read(dict.lookup("filtering"))
318  ),
319  average_(dict.lookupOrDefault("average", false)),
320  zoneKey_(keyType::null),
321  needsUpdate_(true),
322  isoSurfPtr_(nullptr)
323 {}
324 
325 
327 (
328  const word& name,
329  const polyMesh& mesh,
330  const bool interpolate,
331  const word& surfaceType,
332  const word& surfaceName,
333  const scalar distance,
334  const bool signedDistance,
335  const isoSurface::filterType filter,
336  const Switch average
337 )
338 :
339  sampledSurface(name, mesh, interpolate),
340  surfPtr_
341  (
343  (
344  surfaceType,
345  IOobject
346  (
347  surfaceName, // name
348  mesh.time().constant(), // directory
349  "triSurface", // instance
350  mesh.time(), // registry
353  ),
354  dictionary()
355  )
356  ),
357  distance_(distance),
358  signed_(signedDistance),
359  filter_(filter),
360  average_(average),
361  zoneKey_(keyType::null),
362  needsUpdate_(true),
363  isoSurfPtr_(nullptr)
364 {}
365 
366 
367 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
368 
370 {}
371 
372 
373 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
374 
376 {
377  return needsUpdate_;
378 }
379 
380 
382 {
383  if (debug)
384  {
385  Pout<< "distanceSurface::expire :"
386  << " needsUpdate_:" << needsUpdate_ << endl;
387  }
388 
389  // Clear derived data
390  clearGeom();
391 
392  // already marked as expired
393  if (needsUpdate_)
394  {
395  return false;
396  }
397 
398  needsUpdate_ = true;
399  return true;
400 }
401 
402 
404 {
405  if (debug)
406  {
407  Pout<< "distanceSurface::update :"
408  << " needsUpdate_:" << needsUpdate_ << endl;
409  }
410 
411  if (!needsUpdate_)
412  {
413  return false;
414  }
415 
416  createGeometry();
417 
418  needsUpdate_ = false;
419  return true;
420 }
421 
422 
425 (
426  const volScalarField& vField
427 ) const
428 {
429  return sampleField(vField);
430 }
431 
432 
435 (
436  const volVectorField& vField
437 ) const
438 {
439  return sampleField(vField);
440 }
441 
442 
445 (
446  const volSphericalTensorField& vField
447 ) const
448 {
449  return sampleField(vField);
450 }
451 
452 
455 (
456  const volSymmTensorField& vField
457 ) const
458 {
459  return sampleField(vField);
460 }
461 
462 
465 (
466  const volTensorField& vField
467 ) const
468 {
469  return sampleField(vField);
470 }
471 
472 
475 (
476  const interpolation<scalar>& interpolator
477 ) const
478 {
479  return interpolateField(interpolator);
480 }
481 
482 
485 (
486  const interpolation<vector>& interpolator
487 ) const
488 {
489  return interpolateField(interpolator);
490 }
491 
494 (
495  const interpolation<sphericalTensor>& interpolator
496 ) const
497 {
498  return interpolateField(interpolator);
499 }
500 
501 
504 (
505  const interpolation<symmTensor>& interpolator
506 ) const
507 {
508  return interpolateField(interpolator);
509 }
510 
511 
514 (
515  const interpolation<tensor>& interpolator
516 ) const
517 {
518  return interpolateField(interpolator);
519 }
520 
521 
523 {
524  os << "distanceSurface: " << name() << " :"
525  << " surface:" << surfPtr_().name()
526  << " distance:" << distance_
527  << " faces:" << faces().size()
528  << " points:" << points().size();
529 }
530 
531 
532 // ************************************************************************* //
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:667
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
virtual Ostream & write(const char)=0
Write character.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
const word & name() const
Return the name of this functionObject.
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
An abstract class for surfaces with sampling.
defineTypeNameAndDebug(distanceSurface, 0)
const Boundary & boundaryField() const
Return const-reference to the boundary field.
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: HashTable.H:59
virtual bool update()
Update the surface as required.
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
bool interpolate() const
Interpolation requested for surface.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
A simple wrapper around bool so that it can be read as a word: true/false, on/off, yes/no, y/n, t/f, or none/any.
Definition: Switch.H:60
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:622
bool readBool(Istream &)
Definition: boolIO.C:60
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:239
A sampledSurface defined by a distance to a surface.
Macros for easy insertion into run-time selection tables.
static const pointMesh & New(const polyMesh &mesh)
Definition: MeshObject.C:44
virtual const pointField & points() const
Return raw points.
Definition: polyMesh.C:1131
virtual void print(Ostream &) const
Write.
dynamicFvMesh & mesh
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().x()<< ' ';}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().y()<< ' ';}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().z()<< ' ';}gmvFile<< nl;forAll(lagrangianScalarNames, i){ const word &name=lagrangianScalarNames[i];IOField< scalar > fld(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
const pointField & points
A class for handling words, derived from string.
Definition: word.H:59
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
const word & constant() const
Return constant name.
Definition: TimePaths.H:124
static const NamedEnum< filterType, 3 > filterTypeNames_
Definition: isoSurface.H:73
virtual bool needsUpdate() const
Does the surface need an update?
virtual bool expire()
Mark the surface as needing an update.
A sampledSurface defined by a surface of iso value. To be used in sampleSurfaces / functionObjects...
void reset(const label nPoints, const label nInternalFaces, const label nFaces, const label nCells)
Reset this primitiveMesh given the primitive array sizes.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
Internal::FieldType & primitiveFieldRef()
Return a reference to the internal field.
const Time & time() const
Return time.
addToRunTimeSelectionTable(sampledSurface, distanceSurface, word)
distanceSurface(const word &name, const polyMesh &mesh, const dictionary &dict)
Construct from dictionary.
label patchi
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:50
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
prefixOSstream Pout(cout, "Pout")
Definition: IOstreams.H:53
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
virtual tmp< scalarField > sample(const volScalarField &) const
Sample field on surface.
static autoPtr< searchableSurface > New(const word &surfaceType, const IOobject &io, const dictionary &dict)
Return a reference to the selected searchableSurface.
Set of surfaces to sample.
dimensioned< scalar > mag(const dimensioned< Type > &)
label nPoints() const
Enum read(Istream &) const
Read a word from Istream and return the corresponding.
Definition: NamedEnum.C:61
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
const volVectorField & C() const
Return cell centres as volVectorField.
virtual bool write(const bool write=true) const
Write using setting from DB.
A class for managing temporary objects.
Definition: PtrList.H:53
static const keyType null
An empty keyType.
Definition: keyType.H:92
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
fileName objectPath() const
Return complete path + object name.
Definition: IOobject.H:419
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:812