patchProbes.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-2021 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 "patchProbes.H"
27 #include "volFields.H"
28 #include "IOmanip.H"
29 #include "mappedPatchBase.H"
30 #include "treeBoundBox.H"
31 #include "treeDataFace.H"
33 
34 
35 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
36 
37 namespace Foam
38 {
39  defineTypeNameAndDebug(patchProbes, 0);
40 
42  (
43  functionObject,
44  patchProbes,
45  dictionary
46  );
47 }
48 
49 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
50 
52 {
53  (void)mesh.tetBasePtIs();
54 
55  const polyBoundaryMesh& bm = mesh.boundaryMesh();
56 
57  label patchi = bm.findPatchID(patchName_);
58 
59  if (patchi == -1)
60  {
62  << " Unknown patch name "
63  << patchName_ << endl
64  << exit(FatalError);
65  }
66 
67  // All the info for nearest. Construct to miss
68  List<mappedPatchBase::nearInfo> nearest(this->size());
69 
70  const polyPatch& pp = bm[patchi];
71 
72  if (pp.size() > 0)
73  {
74  labelList bndFaces(pp.size());
75  forAll(bndFaces, i)
76  {
77  bndFaces[i] = pp.start() + i;
78  }
79 
80  treeBoundBox overallBb(pp.points());
81  overallBb = overallBb.extend(1e-4);
82 
83  const indexedOctree<treeDataFace> boundaryTree
84  (
85  treeDataFace // all information needed to search faces
86  (
87  false, // do not cache bb
88  mesh,
89  bndFaces // patch faces only
90  ),
91  overallBb, // overall search domain
92  8, // maxLevel
93  10, // leafsize
94  3.0 // duplicity
95  );
96 
97 
98  forAll(probeLocations(), probei)
99  {
100  const point sample = probeLocations()[probei];
101 
102  scalar span = boundaryTree.bb().mag();
103 
104  pointIndexHit info = boundaryTree.findNearest
105  (
106  sample,
107  Foam::sqr(span)
108  );
109 
110  if (!info.hit())
111  {
112  info = boundaryTree.findNearest
113  (
114  sample,
115  Foam::sqr(great)
116  );
117  }
118 
119  label facei = boundaryTree.shapes().faceLabels()[info.index()];
120 
121  const label patchi = bm.whichPatch(facei);
122 
123  if (isA<emptyPolyPatch>(bm[patchi]))
124  {
126  << " The sample point: " << sample
127  << " belongs to " << patchi
128  << " which is an empty patch. This is not permitted. "
129  << " This sample will not be included "
130  << endl;
131  }
132  else
133  {
134  const point& fc = mesh.faceCentres()[facei];
135 
136  mappedPatchBase::nearInfo sampleInfo;
137 
138  sampleInfo.first() = pointIndexHit
139  (
140  true,
141  fc,
142  facei
143  );
144 
145  sampleInfo.second().first() = magSqr(fc-sample);
146  sampleInfo.second().second() = Pstream::myProcNo();
147 
148  nearest[probei]= sampleInfo;
149  }
150  }
151  }
152 
153 
154  // Find nearest.
157 
158  if (debug)
159  {
160  InfoInFunction << endl;
161  forAll(nearest, sampleI)
162  {
163  label proci = nearest[sampleI].second().second();
164  label localI = nearest[sampleI].first().index();
165 
166  Info<< " " << sampleI << " coord:"<< operator[](sampleI)
167  << " found on processor:" << proci
168  << " in local cell/face:" << localI
169  << " with fc:" << nearest[sampleI].first().rawPoint() << endl;
170  }
171  }
172 
173 
174  // Extract any local faces to sample
175  elementList_.setSize(nearest.size(), -1);
176 
177  forAll(nearest, sampleI)
178  {
179  if (nearest[sampleI].second().second() == Pstream::myProcNo())
180  {
181  // Store the face to sample
182  elementList_[sampleI] = nearest[sampleI].first().index();
183  }
184  }
185 }
186 
187 
188 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
189 
191 (
192  const word& name,
193  const Time& t,
194  const dictionary& dict
195 )
196 :
197  probes(name, t, dict)
198 {
199  // When constructing probes above it will have called the
200  // probes::findElements (since the virtual mechanism not yet operating).
201  // Not easy to workaround (apart from feeding through flag into constructor)
202  // so clear out any cells found for now.
203  elementList_.clear();
204  faceList_.clear();
205 
206  read(dict);
207 }
208 
209 
210 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
211 
213 {}
214 
215 
217 {
218  if (this->size() && prepare())
219  {
220  sampleAndWrite(scalarFields_);
221  sampleAndWrite(vectorFields_);
222  sampleAndWrite(sphericalTensorFields_);
223  sampleAndWrite(symmTensorFields_);
224  sampleAndWrite(tensorFields_);
225 
226  sampleAndWriteSurfaceFields(surfaceScalarFields_);
227  sampleAndWriteSurfaceFields(surfaceVectorFields_);
228  sampleAndWriteSurfaceFields(surfaceSphericalTensorFields_);
229  sampleAndWriteSurfaceFields(surfaceSymmTensorFields_);
230  sampleAndWriteSurfaceFields(surfaceTensorFields_);
231  }
232 
233  return true;
234 }
235 
236 
238 {
239  dict.lookup("patchName") >> patchName_;
240  return probes::read(dict);
241 }
242 
243 
244 // ************************************************************************* //
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:453
static void listCombineScatter(const List< commsStruct > &comms, List< T > &Value, const int tag, const label comm)
Scatter data. Reverse of combineGather.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:306
A 2-tuple for storing two objects of different types.
Definition: HashTable.H:65
const labelIOList & tetBasePtIs() const
Return the tetBasePtIs.
Definition: polyMesh.C:928
Encapsulation of data needed to search for faces.
Definition: treeDataFace.H:58
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
const Type1 & first() const
Return first.
Definition: Tuple2.H:116
dimensionedSymmTensor sqr(const dimensionedVector &dv)
virtual ~patchProbes()
Destructor.
Definition: patchProbes.C:212
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
Definition: UPstream.H:429
Set of locations to sample.
Definition: probes.H:64
virtual bool read(const dictionary &)
Read.
Definition: patchProbes.C:237
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
label findPatchID(const word &patchName) const
Find patch index given a name.
patchProbes(const word &name, const Time &time, const dictionary &dict)
Construct from Time and dictionary.
Definition: patchProbes.C:191
This class describes the interaction of (usually) a face and a point. It carries the info of a succes...
Definition: PointIndexHit.H:53
PointIndexHit< point > pointIndexHit
Definition: pointIndexHit.H:42
fvMesh & mesh
T & first()
Return the first element of the list.
Definition: UListI.H:114
static void listCombineGather(const List< commsStruct > &comms, List< T > &Value, const CombineOp &cop, const int tag, const label comm)
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:69
Macros for easy insertion into run-time selection tables.
bool read(const char *, int32_t &)
Definition: int32IO.C:85
void clear()
Clear the list, i.e. set size to zero.
Definition: ListI.H:125
A class for handling words, derived from string.
Definition: word.H:59
scalar mag() const
The magnitude of the bounding box span.
Definition: boundBoxI.H:90
bool hit() const
Is there a hit.
const Field< PointType > & points() const
Return reference to global points.
const Type & shapes() const
Reference to shape.
const treeBoundBox & bb() const
Top bounding box.
Foam::polyBoundaryMesh.
dimensioned< scalar > magSqr(const dimensioned< Type > &)
Istream and Ostream manipulators taking arguments.
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
defineTypeNameAndDebug(combustionModel, 0)
labelList second(const UList< labelPair > &p)
Definition: patchToPatch.C:48
const vectorField & faceCentres() const
virtual void findElements(const fvMesh &)
Find elements containing patchProbes.
Definition: patchProbes.C:51
label patchi
Non-pointer based hierarchical recursive searching.
Definition: treeDataEdge.H:47
#define WarningInFunction
Report a warning using Foam::Warning.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:95
label start() const
Return start label of this patch in the polyMesh face list.
Definition: polyPatch.H:306
virtual bool write()
Sample and write.
Definition: patchProbes.C:216
virtual bool read(const dictionary &)
Read the probes.
Definition: probes.C:298
Standard boundBox + extra functionality for use in octree.
Definition: treeBoundBox.H:87
messageStream Info
const Type2 & second() const
Return second.
Definition: Tuple2.H:128
const doubleScalar e
Elementary charge.
Definition: doubleScalar.H:105
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
label whichPatch(const label faceIndex) const
Return patch index for a given face label.
label index() const
Return index.
const labelList & faceLabels() const
Definition: treeDataFace.H:179
treeBoundBox extend(const scalar s) const
Return asymetrically extended bounding box, with guaranteed.
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:864
#define InfoInFunction
Report an information message using Foam::Info.