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-2023 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 "RemoteData.H"
30 #include "treeBoundBox.H"
31 #include "treeDataFace.H"
33 
34 
35 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
36 
37 namespace Foam
38 {
40 
42  (
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.findIndex(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<RemoteData<scalar>> 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  forAll(probeLocations(), probei)
98  {
99  const point sample = probeLocations()[probei];
100 
101  scalar span = boundaryTree.bb().mag();
102 
103  pointIndexHit info = boundaryTree.findNearest
104  (
105  sample,
106  Foam::sqr(span)
107  );
108 
109  if (!info.hit())
110  {
111  info = boundaryTree.findNearest
112  (
113  sample,
114  Foam::sqr(great)
115  );
116  }
117 
118  label facei = boundaryTree.shapes().faceLabels()[info.index()];
119 
120  const label patchi = bm.whichPatch(facei);
121 
122  if (isA<emptyPolyPatch>(bm[patchi]))
123  {
125  << " The sample point: " << sample
126  << " belongs to " << patchi
127  << " which is an empty patch. This is not permitted. "
128  << " This sample will not be included "
129  << endl;
130  }
131  else
132  {
133  const point& fc = mesh.faceCentres()[facei];
134 
135  nearest[probei].proci = Pstream::myProcNo();
136  nearest[probei].elementi = facei;
137  nearest[probei].data = magSqr(fc-sample);
138  }
139  }
140  }
141 
142  // Find nearest.
144  (
145  nearest,
147  );
149 
150  if (debug)
151  {
152  InfoInFunction << endl;
153  forAll(nearest, sampleI)
154  {
155  Info<< " " << sampleI << " coord:" << operator[](sampleI)
156  << " found on processor:" << nearest[sampleI].proci
157  << " in local cell/face:" << nearest[sampleI].elementi
158  << endl;
159  }
160  }
161 
162  // Extract any local faces to sample
163  elementList_.setSize(nearest.size(), -1);
164  forAll(nearest, sampleI)
165  {
166  if (nearest[sampleI].proci == Pstream::myProcNo())
167  {
168  elementList_[sampleI] = nearest[sampleI].elementi;
169  }
170  }
171 }
172 
173 
174 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
175 
177 (
178  const word& name,
179  const Time& t,
180  const dictionary& dict
181 )
182 :
183  probes(name, t, dict)
184 {
185  // When constructing probes above it will have called the
186  // probes::findElements (since the virtual mechanism not yet operating).
187  // Not easy to workaround (apart from feeding through flag into constructor)
188  // so clear out any cells found for now.
190  faceList_.clear();
191 
192  read(dict);
193 }
194 
195 
196 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
197 
199 {}
200 
201 
203 {
204  if (this->size() && prepare())
205  {
206  sampleAndWrite(scalarFields_);
207  sampleAndWrite(vectorFields_);
208  sampleAndWrite(sphericalTensorFields_);
209  sampleAndWrite(symmTensorFields_);
210  sampleAndWrite(tensorFields_);
211 
212  sampleAndWriteSurfaceFields(surfaceScalarFields_);
213  sampleAndWriteSurfaceFields(surfaceVectorFields_);
214  sampleAndWriteSurfaceFields(surfaceSphericalTensorFields_);
215  sampleAndWriteSurfaceFields(surfaceSymmTensorFields_);
216  sampleAndWriteSurfaceFields(surfaceTensorFields_);
217  }
218 
219  return true;
220 }
221 
222 
224 {
225  dict.lookup("patchName") >> patchName_;
226  return probes::read(dict);
227 }
228 
229 
230 // ************************************************************************* //
Istream and Ostream manipulators taking arguments.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Macros for easy insertion into run-time selection tables.
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: List.H:91
label size() const
Return the number of elements in the UList.
Definition: ListI.H:171
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
void clear()
Clear the list, i.e. set size to zero.
Definition: ListI.H:125
void setSize(const label)
Reset size of List.
Definition: List.C:281
This class describes the interaction of (usually) a face and a point. It carries the info of a succes...
Definition: PointIndexHit.H:54
label index() const
Return index.
bool hit() const
Is there a hit.
const Field< PointType > & points() const
Return reference to global points.
static void listCombineGather(const List< commsStruct > &comms, List< T > &Value, const CombineOp &cop, const int tag, const label comm)
static void listCombineScatter(const List< commsStruct > &comms, List< T > &Value, const int tag, const label comm)
Scatter data. Reverse of combineGather.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:76
label size() const
Return the number of elements in the UList.
Definition: UListI.H:311
T & operator[](const label)
Return element of UList.
Definition: UListI.H:167
T * data()
Return a pointer to the first data element,.
Definition: UListI.H:149
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
Definition: UPstream.H:429
scalar mag() const
The magnitude of the bounding box span.
Definition: boundBoxI.H:96
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:710
Abstract base-class for Time/database functionObjects.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
Non-pointer based hierarchical recursive searching.
Definition: indexedOctree.H:72
const Type & shapes() const
Reference to shape.
const treeBoundBox & bb() const
Top bounding box.
Set of locations to sample.at patches.
Definition: patchProbes.H:59
virtual void findElements(const fvMesh &)
Find elements containing patchProbes.
Definition: patchProbes.C:51
patchProbes(const word &name, const Time &time, const dictionary &dict)
Construct from Time and dictionary.
Definition: patchProbes.C:177
virtual bool write()
Sample and write.
Definition: patchProbes.C:202
virtual ~patchProbes()
Destructor.
Definition: patchProbes.C:198
virtual bool read(const dictionary &)
Read.
Definition: patchProbes.C:223
Foam::polyBoundaryMesh.
label findIndex(const word &patchName) const
Find patch index given a name.
label whichPatch(const label faceIndex) const
Return patch index for a given face label.
const labelIOList & tetBasePtIs() const
Return the tetBasePtIs.
Definition: polyMesh.C:1023
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:404
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:70
label start() const
Return start label of this patch in the polyMesh face list.
Definition: polyPatch.H:280
const vectorField & faceCentres() const
Set of locations to sample.
Definition: probes.H:68
labelList elementList_
Definition: probes.H:126
virtual const pointField & probeLocations() const
Return locations to probe.
Definition: probes.H:210
labelList faceList_
Definition: probes.H:129
virtual bool read(const dictionary &)
Read the probes.
Definition: probes.C:294
Standard boundBox + extra functionality for use in octree.
Definition: treeBoundBox.H:90
treeBoundBox extend(const scalar s) const
Return asymmetrically extended bounding box, with guaranteed.
Encapsulation of data needed to search for faces.
Definition: treeDataFace.H:59
A class for handling words, derived from string.
Definition: word.H:62
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:334
label patchi
#define WarningInFunction
Report a warning using Foam::Warning.
#define InfoInFunction
Report an information message using Foam::Info.
Namespace for OpenFOAM.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
const doubleScalar e
Definition: doubleScalar.H:106
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
addToRunTimeSelectionTable(polyPatch, mergedCyclicPolyPatch, word)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
messageStream Info
defineTypeNameAndDebug(combustionModel, 0)
error FatalError
dimensioned< scalar > magSqr(const dimensioned< Type > &)
dictionary dict
Operator to take smallest valid value.
Definition: RemoteData.H:94