triSurfaceRegionSearch.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 "triSurfaceRegionSearch.H"
27 #include "indexedOctree.H"
28 #include "triSurface.H"
29 #include "PatchTools.H"
30 
31 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
32 
34 :
35  triSurfaceSearch(surface),
36  indirectRegionPatches_(),
37  treeByRegion_()
38 {}
39 
40 
42 (
43  const triSurface& surface,
44  const dictionary& dict
45 )
46 :
47  triSurfaceSearch(surface, dict),
48  indirectRegionPatches_(),
49  treeByRegion_()
50 {}
51 
52 
53 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
54 
56 {
57  clearOut();
58 }
59 
60 
62 {
64  treeByRegion_.clear();
65 }
66 
67 
68 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
69 
72 {
73  if (treeByRegion_.empty())
74  {
75  Map<label> regionSizes;
76  forAll(surface(), fI)
77  {
78  const label regionI = surface()[fI].region();
79 
80  regionSizes(regionI)++;
81  }
82 
83  label nRegions = regionSizes.size();
84 
85  indirectRegionPatches_.setSize(nRegions);
86  treeByRegion_.setSize(nRegions);
87 
88  labelListList regionsAddressing(nRegions);
89 
90  forAll(regionsAddressing, regionI)
91  {
92  regionsAddressing[regionI] = labelList(regionSizes[regionI], -1);
93  }
94 
95  labelList nFacesInRegions(nRegions, 0);
96 
97  forAll(surface(), fI)
98  {
99  const label regionI = surface()[fI].region();
100 
101  regionsAddressing[regionI][nFacesInRegions[regionI]++] = fI;
102  }
103 
104  forAll(regionsAddressing, regionI)
105  {
106  scalar oldTol = treeType::perturbTol();
108 
109  indirectRegionPatches_.set
110  (
111  regionI,
113  (
115  (
116  surface(),
117  regionsAddressing[regionI]
118  ),
119  surface().points()
120  )
121  );
122 
123  // Calculate bb without constructing local point numbering.
124  treeBoundBox bb(Zero, Zero);
125 
126  if (indirectRegionPatches_[regionI].size())
127  {
128  label nPoints;
130  (
131  indirectRegionPatches_[regionI],
132  bb,
133  nPoints
134  );
135 
136  // if (nPoints != surface().points().size())
137  // {
138  // WarningInFunction
139  // << "Surface does not have compact point numbering. "
140  // << "Of " << surface().points().size()
141  // << " only " << nPoints
142  // << " are used."
143  // << " This might give problems in some routines."
144  // << endl;
145  // }
146 
147  // Slightly extended bb. Slightly off-centred just so
148  // on symmetric geometry there are fewer face/edge
149  // aligned items.
150  bb = bb.extend(1e-4);
151  }
152 
153  treeByRegion_.set
154  (
155  regionI,
156  new treeType
157  (
159  (
160  true,
161  indirectRegionPatches_[regionI],
162  tolerance()
163  ),
164  bb,
165  maxTreeDepth(), // maxLevel
166  10, // leafsize
167  3.0 // duplicity
168  )
169  );
170 
171  treeType::perturbTol() = oldTol;
172  }
173  }
174 
175  return treeByRegion_;
176 }
177 
178 
180 (
181  const pointField& samples,
182  const scalarField& nearestDistSqr,
183  const labelList& regionIndices,
184  List<pointIndexHit>& info
185 ) const
186 {
187  if (regionIndices.empty())
188  {
189  triSurfaceSearch::findNearest(samples, nearestDistSqr, info);
190  }
191  else
192  {
193  scalar oldTol = treeType::perturbTol();
195 
196  const PtrList<treeType>& octrees = treeByRegion();
197 
198  info.setSize(samples.size());
199 
200  forAll(octrees, treeI)
201  {
202  if (findIndex(regionIndices, treeI) == -1)
203  {
204  continue;
205  }
206 
207  const treeType& octree = octrees[treeI];
208 
209  forAll(samples, i)
210  {
211 // if (!octree.bb().contains(samples[i]))
212 // {
213 // continue;
214 // }
215 
216  pointIndexHit currentRegionHit = octree.findNearest
217  (
218  samples[i],
219  nearestDistSqr[i],
221  );
222 
223  if
224  (
225  currentRegionHit.hit()
226  &&
227  (
228  !info[i].hit()
229  ||
230  (
231  magSqr(currentRegionHit.hitPoint() - samples[i])
232  < magSqr(info[i].hitPoint() - samples[i])
233  )
234  )
235  )
236  {
237  info[i] = currentRegionHit;
238  }
239  }
240  }
241 
242  treeType::perturbTol() = oldTol;
243  }
244 }
245 
246 
247 // ************************************************************************* //
label maxTreeDepth() const
Return max tree depth of octree.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
bool empty() const
Return true if the UList is empty (ie, size() is zero)
Definition: UListI.H:325
void findNearest(const pointField &samples, const scalarField &nearestDistSqr, const labelList &regionIndices, List< pointIndexHit > &info) const
Find the nearest point on the surface out of the regions.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
triSurfaceRegionSearch(const triSurface &)
Construct from surface. Holds reference to surface!
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
static scalar & perturbTol()
Get the perturbation tolerance.
This class describes the interaction of (usually) a face and a point. It carries the info of a succes...
Definition: PointIndexHit.H:53
label size() const
Return number of elements in table.
Definition: HashTableI.H:65
void findNearest(const pointField &samples, const scalarField &nearestDistSqr, List< pointIndexHit > &info) const
Helper class to search on triSurface.
A list of faces which address into the list of points.
const Point & hitPoint() const
Return hit point.
static void calcBounds(const PrimitivePatch< FaceList, PointField > &p, boundBox &bb, label &nPoints)
const pointField & points
label nPoints
triSurfaceSearch(const triSurface &)
Construct from surface. Holds reference to surface!
void clearOut()
Clear storage.
Encapsulation of data needed to search on PrimitivePatches.
bool hit() const
Is there a hit.
List< label > labelList
A List of labels.
Definition: labelList.H:56
static const zero Zero
Definition: zero.H:97
dimensioned< scalar > magSqr(const dimensioned< Type > &)
const PtrList< treeType > & treeByRegion() const
Demand driven construction of octree for each region.
label findIndex(const ListType &, typename ListType::const_reference, const label start=0)
Find first occurrence of given element and return index,.
const triSurface & surface() const
Return reference to the surface.
void setSize(const label)
Reset size of List.
Definition: List.C:281
Non-pointer based hierarchical recursive searching.
Definition: treeDataEdge.H:47
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:70
scalar tolerance() const
Return tolerance to use in searches.
Standard boundBox + extra functionality for use in octree.
Definition: treeBoundBox.H:87
const doubleScalar e
Elementary charge.
Definition: doubleScalar.H:105
Triangulated surface description with patch information.
Definition: triSurface.H:66
A List with indirect addressing.
Definition: IndirectList.H:101
treeBoundBox extend(const scalar s) const
Return asymetrically extended bounding box, with guaranteed.