surfaceDistance.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 | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2019-2020 OpenCFD Ltd.
9 -------------------------------------------------------------------------------
10 License
11  This file is part of OpenFOAM.
12 
13  OpenFOAM is free software: you can redistribute it and/or modify it
14  under the terms of the GNU General Public License as published by
15  the Free Software Foundation, either version 3 of the License, or
16  (at your option) any later version.
17 
18  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21  for more details.
22 
23  You should have received a copy of the GNU General Public License
24  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25 
26 \*---------------------------------------------------------------------------*/
27 
28 #include "surfaceDistance.H"
30 
31 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35 namespace functionObjects
36 {
37  defineTypeNameAndDebug(surfaceDistance, 0);
38  addToRunTimeSelectionTable(functionObject, surfaceDistance, dictionary);
39 }
40 }
41 
42 
43 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
44 
46 (
47  const word& name,
48  const Time& runTime,
49  const dictionary& dict
50 )
51 :
53 {
54  read(dict);
55 
56  volScalarField* procFieldPtr
57  (
58  new volScalarField
59  (
60  IOobject
61  (
62  "surfaceDistance",
63  mesh_.time().timeName(),
64  mesh_,
65  IOobject::NO_READ,
66  IOobject::NO_WRITE
67  ),
68  mesh_,
70  )
71  );
72 
73  mesh_.objectRegistry::store(procFieldPtr);
74 }
75 
76 
77 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
78 
80 (
81  const dictionary& dict
82 )
83 {
85 
86  doCells_ = dict.getOrDefault("calculateCells", true);
87 
88  geomPtr_.reset(nullptr);
89  geomPtr_.reset
90  (
92  (
93  IOobject
94  (
95  "abc", // dummy name
96  mesh_.time().constant(), // directory
97  "triSurface", // instance
98  mesh_.time(), // registry
99  IOobject::MUST_READ,
100  IOobject::NO_WRITE
101  ),
102  dict.subDict("geometry"),
103  true // allow single-region shortcut
104  )
105  );
106 
107  return true;
108 }
109 
110 
112 {
114  (
115  "surfaceDistance"
116  );
117 
118  volScalarField::Boundary& bfld = distance.boundaryFieldRef();
119  forAll(bfld, patchi)
120  {
121  if (!polyPatch::constraintType(bfld[patchi].patch().type()))
122  {
123  const pointField& fc = mesh_.C().boundaryField()[patchi];
124 
125  labelList surfaces;
126  List<pointIndexHit> nearestInfo;
127  geomPtr_().findNearest
128  (
129  fc,
130  scalarField(fc.size(), GREAT),
131  surfaces,
132  nearestInfo
133  );
134 
135  scalarField dist(fc.size());
136  forAll(nearestInfo, i)
137  {
138  dist[i] = mag(nearestInfo[i].hitPoint()-fc[i]);
139  }
140  bfld[patchi] == dist;
141  }
142  }
143 
144  if (doCells_)
145  {
146  const pointField& cc = mesh_.C();
147 
148  labelList surfaces;
149  List<pointIndexHit> nearestInfo;
150  geomPtr_().findNearest
151  (
152  cc,
153  scalarField(cc.size(), GREAT),
154  surfaces,
155  nearestInfo
156  );
157 
158  forAll(nearestInfo, celli)
159  {
160  distance[celli] = mag(nearestInfo[celli].hitPoint()-cc[celli]);
161  }
162  }
163  distance.correctBoundaryConditions();
164 
165  return true;
166 }
167 
168 
170 {
171  Log << " functionObjects::" << type() << " " << name()
172  << " writing distance-to-surface field" << endl;
173 
174  const volScalarField& distance =
175  mesh_.lookupObject<volScalarField>("surfaceDistance");
176 
177 // volScalarField::Boundary& bfld = distance.boundaryFieldRef();
178 // forAll(bfld, patchi)
179 // {
180 // if (!polyPatch::constraintType(bfld[patchi].patch().type()))
181 // {
182 // const pointField& fc = mesh_.C().boundaryField()[patchi];
183 //
184 // labelList surfaces;
185 // List<pointIndexHit> nearestInfo;
186 // geomPtr_().findNearest
187 // (
188 // fc,
189 // scalarField(fc.size(), GREAT),
190 // surfaces,
191 // nearestInfo
192 // );
193 //
194 // scalarField dist(fc.size());
195 // forAll(nearestInfo, i)
196 // {
197 // dist[i] = mag(nearestInfo[i].hitPoint()-fc[i]);
198 // }
199 // bfld[patchi] == dist;
200 // }
201 // }
202 //
203 // if (doCells_)
204 // {
205 // const pointField& cc = mesh_.C();
206 //
207 // labelList surfaces;
208 // List<pointIndexHit> nearestInfo;
209 // geomPtr_().findNearest
210 // (
211 // cc,
212 // scalarField(cc.size(), GREAT),
213 // surfaces,
214 // nearestInfo
215 // );
216 //
217 // forAll(nearestInfo, celli)
218 // {
219 // distance[celli] = mag(nearestInfo[celli].hitPoint()-cc[celli]);
220 // }
221 // }
222 // distance.correctBoundaryConditions();
223  distance.write();
224 
225  return true;
226 }
227 
228 
229 // ************************************************************************* //
runTime
engineTime & runTime
Definition: createEngineTime.H:13
Foam::scalarField
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Definition: primitiveFieldsFwd.H:52
Foam::IOobject
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:169
Log
#define Log
Definition: PDRblock.C:35
Foam::functionObjects::surfaceDistance::write
virtual bool write()
Write the interpolated fields.
Definition: surfaceDistance.C:169
Foam::Time
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:73
Foam::word
A class for handling words, derived from Foam::string.
Definition: word.H:65
Foam::dimLength
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:52
Foam::Zero
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
surfaceDistance.H
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:369
Foam::functionObjects::fvMeshFunctionObject
Specialization of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
Definition: fvMeshFunctionObject.H:64
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
Foam::blockMeshTools::read
void read(Istream &, label &val, const dictionary &)
In-place read with dictionary lookup.
Definition: blockMeshTools.C:57
Foam::Field< vector >
Foam::functionObjects::surfaceDistance::doCells_
Switch doCells_
Switch to calculate distance-to-cells.
Definition: surfaceDistance.H:172
Foam::dimensionedScalar
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Definition: dimensionedScalarFwd.H:42
Foam::functionObjects::surfaceDistance::execute
virtual bool execute()
Calculate the interpolated fields.
Definition: surfaceDistance.C:111
Foam::polyPatch::constraintType
static bool constraintType(const word &pt)
Return true if the given type is a constraint type.
Definition: polyPatch.C:277
Foam::fvMesh::C
const volVectorField & C() const
Return cell centres as volVectorField.
Definition: fvMeshGeometry.C:341
dict
dictionary dict
Definition: searchingEngine.H:14
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:123
addToRunTimeSelectionTable.H
Macros for easy insertion into run-time selection tables.
Foam::functionObjects::surfaceDistance::surfaceDistance
surfaceDistance(const word &name, const Time &runTime, const dictionary &dict)
Construct for given objectRegistry and dictionary.
Definition: surfaceDistance.C:46
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::distance
scalar distance(const vector &p1, const vector &p2)
Definition: curveTools.C:12
Foam::objectRegistry::lookupObjectRef
Type & lookupObjectRef(const word &name, const bool recursive=false) const
Definition: objectRegistryTemplates.C:478
Foam::functionObjects::addToRunTimeSelectionTable
addToRunTimeSelectionTable(functionObject, ObukhovLength, dictionary)
Foam::foamVersion::patch
const std::string patch
OpenFOAM patch number as a std::string.
Foam::functionObject::type
virtual const word & type() const =0
Runtime type information.
Foam::List< label >
Foam::searchableSurfaces
Container for searchableSurfaces. The collection is specified as a dictionary. For example,...
Definition: searchableSurfaces.H:92
Foam::mag
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Foam::type
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
Definition: MSwindows.C:590
Foam::functionObjects::fvMeshFunctionObject::mesh_
const fvMesh & mesh_
Reference to the fvMesh.
Definition: fvMeshFunctionObject.H:73
Foam::functionObjects::surfaceDistance::geomPtr_
autoPtr< searchableSurfaces > geomPtr_
Geometry.
Definition: surfaceDistance.H:175
Foam::functionObjects::surfaceDistance::read
virtual bool read(const dictionary &)
Read the controls.
Definition: surfaceDistance.C:80
Foam::functionObjects::defineTypeNameAndDebug
defineTypeNameAndDebug(ObukhovLength, 0)
Foam::name
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
Foam::GeometricField< scalar, fvPatchField, volMesh >
Foam::GeometricField::boundaryField
const Boundary & boundaryField() const
Return const-reference to the boundary field.
Definition: GeometricFieldI.H:62