patchProbes.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 | Copyright (C) 2011-2012 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 // For 'nearInfo' helper class only
30 #include "mappedPatchBase.H"
31 #include "treeBoundBox.H"
32 #include "treeDataFace.H"
33 
34 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
35 
36 namespace Foam
37 {
38  defineTypeNameAndDebug(patchProbes, 0);
39 }
40 
41 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
42 
44 {
45  (void)mesh.tetBasePtIs();
46 
47  const polyBoundaryMesh& bm = mesh.boundaryMesh();
48 
49  label patchI = bm.findPatchID(patchName_);
50 
51  if (patchI == -1)
52  {
54  (
55  " Foam::patchProbes::findElements(const fvMesh&)"
56  ) << " Unknown patch name "
57  << patchName_ << endl
58  << exit(FatalError);
59  }
60 
61  // All the info for nearest. Construct to miss
62  List<mappedPatchBase::nearInfo> nearest(this->size());
63 
64  const polyPatch& pp = bm[patchI];
65 
66  if (pp.size() > 0)
67  {
68  labelList bndFaces(pp.size());
69  forAll(bndFaces, i)
70  {
71  bndFaces[i] = pp.start() + i;
72  }
73 
74  treeBoundBox overallBb(pp.points());
75  Random rndGen(123456);
76  overallBb = overallBb.extend(rndGen, 1e-4);
77  overallBb.min() -= point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
78  overallBb.max() += point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
79 
80  const indexedOctree<treeDataFace> boundaryTree
81  (
82  treeDataFace // all information needed to search faces
83  (
84  false, // do not cache bb
85  mesh,
86  bndFaces // patch faces only
87  ),
88  overallBb, // overall search domain
89  8, // maxLevel
90  10, // leafsize
91  3.0 // duplicity
92  );
93 
94 
95  forAll(probeLocations(), probeI)
96  {
97  const point sample = probeLocations()[probeI];
98 
99  scalar span = boundaryTree.bb().mag();
100 
101  pointIndexHit info = boundaryTree.findNearest
102  (
103  sample,
104  Foam::sqr(span)
105  );
106 
107  if (!info.hit())
108  {
109  info = boundaryTree.findNearest
110  (
111  sample,
112  Foam::sqr(GREAT)
113  );
114  }
115 
116  label faceI = boundaryTree.shapes().faceLabels()[info.index()];
117 
118  const label patchi = bm.whichPatch(faceI);
119 
120  if (isA<emptyPolyPatch>(bm[patchi]))
121  {
122  WarningIn
123  (
124  " Foam::patchProbes::findElements(const fvMesh&)"
125  )
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  Info<< "patchProbes::findElements" << " : " << 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 
190 Foam::patchProbes::patchProbes
191 (
192  const word& name,
193  const objectRegistry& obr,
194  const dictionary& dict,
195  const bool loadFromFiles
196 )
197 :
198  probes(name, obr, dict, loadFromFiles)
199 {
200  // When constructing probes above it will have called the
201  // probes::findElements (since the virtual mechanism not yet operating).
202  // Not easy to workaround (apart from feeding through flag into constructor)
203  // so clear out any cells found for now.
204  elementList_.clear();
205  faceList_.clear();
206 
207  read(dict);
208 }
209 
210 
211 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
212 
214 {}
215 
216 
218 {
219  if (this->size() && prepare())
220  {
221  sampleAndWrite(scalarFields_);
222  sampleAndWrite(vectorFields_);
223  sampleAndWrite(sphericalTensorFields_);
224  sampleAndWrite(symmTensorFields_);
225  sampleAndWrite(tensorFields_);
226 
227  sampleAndWriteSurfaceFields(surfaceScalarFields_);
228  sampleAndWriteSurfaceFields(surfaceVectorFields_);
229  sampleAndWriteSurfaceFields(surfaceSphericalTensorFields_);
230  sampleAndWriteSurfaceFields(surfaceSymmTensorFields_);
231  sampleAndWriteSurfaceFields(surfaceTensorFields_);
232  }
233 }
234 
236 {
237  dict.lookup("patchName") >> patchName_;
238  probes::read(dict);
239 }
240 
241 
242 // ************************************************************************* //
const Type2 & second() const
Return second.
Definition: Tuple2.H:106
Simple random number generator.
Definition: Random.H:49
cachedRandom rndGen(label(0),-1)
A 2-tuple for storing two objects of different types.
Definition: Tuple2.H:47
virtual ~patchProbes()
Destructor.
Definition: patchProbes.C:213
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
This class describes the interaction of (usually) a face and a point. It carries the info of a succes...
Definition: PointIndexHit.H:53
vector point
Point is a vector.
Definition: point.H:41
scalar mag() const
The magnitude of the bounding box span.
Definition: boundBoxI.H:90
const labelList & tetBasePtIs() const
Return the tetBasePtIs.
Definition: polyMesh.C:864
dimensioned< scalar > magSqr(const dimensioned< Type > &)
virtual void read(const dictionary &)
Read.
Definition: patchProbes.C:235
const Type & shapes() const
Reference to shape.
const treeBoundBox & bb() const
Top bounding box.
A class for handling words, derived from string.
Definition: word.H:59
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
Standard boundBox + extra functionality for use in octree.
Definition: treeBoundBox.H:75
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
PointIndexHit< point > pointIndexHit
Definition: pointIndexHit.H:42
const Type1 & first() const
Return first.
Definition: Tuple2.H:94
messageStream Info
bool hit() const
Is there a hit.
dynamicFvMesh & mesh
Encapsulation of data needed to search for faces.
Definition: treeDataFace.H:58
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
virtual void read(const dictionary &)
Read the probes.
Definition: probes.C:335
T & first()
Return the first element of the list.
Definition: UListI.H:117
label index() const
Return index.
Set of locations to sample.
Definition: probes.H:66
Namespace for OpenFOAM.
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
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
void clear()
Clear the list, i.e. set size to zero.
Definition: List.C:379
const double e
Elementary charge.
Definition: doubleFloat.H:78
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
#define WarningIn(functionName)
Report a warning using Foam::Warning.
Foam::polyBoundaryMesh.
static void listCombineGather(const List< commsStruct > &comms, List< T > &Value, const CombineOp &cop, const int tag, const label comm)
const labelList & faceLabels() const
Definition: treeDataFace.H:179
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)
Definition: UList.H:421
label patchi
virtual void findElements(const fvMesh &)
Find elements containing patchProbes.
Definition: patchProbes.C:43
label whichPatch(const label faceIndex) const
Return patch index for a given face label.
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
Definition: UPstream.H:404
Istream and Ostream manipulators taking arguments.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:421
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
Non-pointer based hierarchical recursive searching.
Definition: treeDataEdge.H:47
label start() const
Return start label of this patch in the polyMesh face list.
Definition: polyPatch.H:300
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:452
error FatalError
Registry of regIOobjects.
bool read(const char *, int32_t &)
Definition: int32IO.C:87
dimensionedSymmTensor sqr(const dimensionedVector &dv)
label findPatchID(const word &patchName) const
Find patch index given a name.
virtual void write()
Public members.
Definition: patchProbes.C:217
const Field< PointType > & points() const
Return reference to global points.
defineTypeNameAndDebug(combustionModel, 0)
const vectorField & faceCentres() const