boundaryPoints.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-2018 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 "boundaryPoints.H"
27 #include "polyMesh.H"
29 #include "treeBoundBox.H"
30 #include "treeDataFace.H"
31 #include "Time.H"
32 #include "meshTools.H"
33 #include "mappedPatchBase.H"
34 
35 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
36 
37 namespace Foam
38 {
39 namespace sampledSets
40 {
41  defineTypeNameAndDebug(boundaryPoints, 0);
42  addToRunTimeSelectionTable(sampledSet, boundaryPoints, word);
43 }
44 }
45 
46 
47 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
48 
49 void Foam::sampledSets::boundaryPoints::calcSamples
50 (
51  DynamicList<point>& samplingPts,
52  DynamicList<label>& samplingCells,
53  DynamicList<label>& samplingFaces,
54  DynamicList<label>& samplingSegments,
55  DynamicList<scalar>& samplingCurveDist
56 ) const
57 {
58  // Construct a single list of all patch faces
59  label nPatchFaces = 0;
60  forAllConstIter(labelHashSet, patches_, iter)
61  {
62  const polyPatch& pp = mesh().boundaryMesh()[iter.key()];
63  nPatchFaces += pp.size();
64  }
65  labelList patchFaces(nPatchFaces);
66  nPatchFaces = 0;
67  forAllConstIter(labelHashSet, patches_, iter)
68  {
69  const polyPatch& pp = mesh().boundaryMesh()[iter.key()];
70  forAll(pp, i)
71  {
72  patchFaces[nPatchFaces++] = pp.start()+i;
73  }
74  }
75 
76  // Construct a processor-local bound box
77  treeBoundBox patchBB(point::max, point::min);
78  forAllConstIter(labelHashSet, patches_, iter)
79  {
80  const polyPatch& pp = mesh().boundaryMesh()[iter.key()];
81  const boundBox patchBb(pp.points(), pp.meshPoints(), false);
82  patchBB.min() = min(patchBB.min(), patchBb.min());
83  patchBB.max() = max(patchBB.max(), patchBb.max());
84  }
85  patchBB = patchBB.extend(1e-4);
86 
87  // Create a search tree
88  indexedOctree<treeDataFace> patchTree
89  (
90  treeDataFace // all information needed to search faces
91  (
92  false, // do not cache the bound box
93  mesh(),
94  patchFaces
95  ),
96  patchBB, // overall search box
97  8, // maximum number of levels
98  10, // how many elements per leaf
99  3.0 // in how many leaves is a shape on average
100  );
101 
102  // Force calculation of face-diagonal decomposition
103  (void)mesh().tetBasePtIs();
104 
105  // Generate the nearest patch information for each sampling point
106  List<mappedPatchBase::nearInfo> nearest(points_.size());
107  forAll(points_, sampleI)
108  {
109  const point& sample = points_[sampleI];
110 
111  pointIndexHit& nearHit = nearest[sampleI].first();
112  scalar& nearDist = nearest[sampleI].second().first();
113  label& nearProc = nearest[sampleI].second().second();
114 
115  // Find the nearest
116  if (patchFaces.size())
117  {
118  nearHit = patchTree.findNearest(sample, sqr(maxDistance_));
119  }
120  else
121  {
122  nearHit.setMiss();
123  }
124 
125  // Fill in the information
126  if (nearHit.hit())
127  {
128  nearHit.setIndex(patchFaces[nearHit.index()]);
129  nearDist = magSqr(nearHit.hitPoint() - sample);
130  nearProc = Pstream::myProcNo();
131  }
132  else
133  {
134  nearHit.setIndex(-1);
135  nearDist = Foam::sqr(great);
136  nearProc = Pstream::myProcNo();
137  }
138  }
139 
140  // Reduce to get the nearest patch locations globally
141  Pstream::listCombineGather(nearest, mappedPatchBase::nearestEqOp());
143 
144  // Dump connecting lines from the sampling points to the hit locations
145  if (debug && Pstream::master())
146  {
147  OFstream str(mesh().time().path() / name() + "_nearest.obj");
148 
149  label verti = 0;
150 
151  forAll(nearest, i)
152  {
153  if (nearest[i].first().hit())
154  {
155  meshTools::writeOBJ(str, points_[i]);
156  verti++;
157  meshTools::writeOBJ(str, nearest[i].first().hitPoint());
158  verti++;
159  str << "l " << verti - 1 << ' ' << verti << nl;
160  }
161  }
162  }
163 
164  // Store the sampling locations on the nearest processor
165  forAll(nearest, sampleI)
166  {
167  const pointIndexHit& nearHit = nearest[sampleI].first();
168  const label& nearProc = nearest[sampleI].second().second();
169 
170  if (nearHit.hit())
171  {
172  if (nearProc == Pstream::myProcNo())
173  {
174  label facei = nearHit.index();
175 
176  samplingPts.append(nearHit.hitPoint());
177  samplingCells.append(mesh().faceOwner()[facei]);
178  samplingFaces.append(facei);
179  samplingSegments.append(0);
180  samplingCurveDist.append(sampleI);
181  }
182  }
183  else
184  {
186  << "Unable to find location on patches " << patches_
187  << " for the point " << points_[sampleI]
188  << " within a distance of " << maxDistance_ << endl;
189  }
190  }
191 }
192 
193 
194 void Foam::sampledSets::boundaryPoints::genSamples()
195 {
196  DynamicList<point> samplingPts;
197  DynamicList<label> samplingCells;
198  DynamicList<label> samplingFaces;
199  DynamicList<label> samplingSegments;
200  DynamicList<scalar> samplingCurveDist;
201 
202  calcSamples
203  (
204  samplingPts,
205  samplingCells,
206  samplingFaces,
207  samplingSegments,
208  samplingCurveDist
209  );
210 
211  samplingPts.shrink();
212  samplingCells.shrink();
213  samplingFaces.shrink();
214  samplingSegments.shrink();
215  samplingCurveDist.shrink();
216 
217  setSamples
218  (
219  samplingPts,
220  samplingCells,
221  samplingFaces,
222  samplingSegments,
223  samplingCurveDist
224  );
225 }
226 
227 
228 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
229 
231 (
232  const word& name,
233  const polyMesh& mesh,
234  const meshSearch& searchEngine,
235  const dictionary& dict
236 )
237 :
238  sampledSet(name, mesh, searchEngine, dict),
239  points_(dict.lookup("points")),
240  patches_
241  (
242  mesh.boundaryMesh().patchSet
243  (
244  wordReList(dict.lookup("patches"))
245  )
246  ),
247  maxDistance_(readScalar(dict.lookup("maxDistance")))
248 {
249  genSamples();
250 
251  if (debug)
252  {
253  write(Info);
254  }
255 }
256 
257 
258 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
259 
261 {}
262 
263 
264 // ************************************************************************* //
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:424
Various (local, not parallel) searches on polyMesh; uses (demand driven) octree to search...
Definition: meshSearch.H:57
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:428
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
boundaryPoints(const word &name, const polyMesh &mesh, const meshSearch &searchEngine, const dictionary &dict)
Construct from dictionary.
const word & name() const
Return the name of this functionObject.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
Definition: UPstream.H:427
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:421
PointIndexHit< point > pointIndexHit
Definition: pointIndexHit.H:42
void writeOBJ(Ostream &os, const point &pt)
Write obj representation of point.
Definition: meshTools.C:203
labelHashSet patchSet(const UList< wordRe > &patchNames, const bool warnNotFound=true, const bool usePatchGroups=true) const
Return the set of patch IDs corresponding to the given names.
static void listCombineGather(const List< commsStruct > &comms, List< T > &Value, const CombineOp &cop, const int tag, const label comm)
Macros for easy insertion into run-time selection tables.
HashSet< label, Hash< label > > labelHashSet
A HashSet with label keys.
Definition: HashSet.H:210
dynamicFvMesh & mesh
Holds list of sampling points which is filled at construction time. Various implementations of this b...
Definition: sampledSet.H:61
A class for handling words, derived from string.
Definition: word.H:59
List< label > labelList
A List of labels.
Definition: labelList.H:56
bool readScalar(const char *buf, doubleScalar &s)
Read whole of buf as a scalar. Return true if successful.
Definition: doubleScalar.H:68
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
addToRunTimeSelectionTable(sampledSet, arcUniform, word)
dimensioned< scalar > magSqr(const dimensioned< Type > &)
virtual ~boundaryPoints()
Destructor.
static const char nl
Definition: Ostream.H:265
fileName path(UMean.rootPath()/UMean.caseName()/"graphs"/UMean.instance())
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
vector point
Point is a vector.
Definition: point.H:41
#define WarningInFunction
Report a warning using Foam::Warning.
virtual bool write()
Sample and write.
Definition: sampledSets.C:241
messageStream Info
const doubleScalar e
Elementary charge.
Definition: doubleScalar.H:98
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
List< wordRe > wordReList
A List of wordRe (word or regular expression)
Definition: wordReList.H:50
defineTypeNameAndDebug(arcUniform, 0)
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:576
T & first()
Return reference to the first element of the list.
Definition: UPtrListI.H:43