FreeFOAM The Cross-Platform CFD Toolkit
treeBoundBoxI.H
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) 1991-2010 OpenCFD Ltd.
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 "treeBoundBox.H"
27 #include <OpenFOAM/Random.H>
28 
29 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
30 
32 :
33  boundBox()
34 {}
35 
36 
38 :
39  boundBox(min, max)
40 {}
41 
42 
44 :
45  boundBox(bb)
46 {}
47 
48 
49 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
50 
51 inline Foam::scalar Foam::treeBoundBox::typDim() const
52 {
53  return avgDim();
54 }
55 
56 
58 {
59  return point
60  (
61  (octant & RIGHTHALF) ? max().x() : min().x(),
62  (octant & TOPHALF) ? max().y() : min().y(),
63  (octant & FRONTHALF) ? max().z() : min().z()
64  );
65 }
66 
67 
68 // Returns octant in which point resides. Reverse of subBbox.
70 {
71  return subOctant(midpoint(), pt);
72 }
73 
74 
75 // Returns octant in which point resides. Reverse of subBbox.
76 // Precalculated midpoint
78 (
79  const point& mid,
80  const point& pt
81 )
82 {
83  direction octant = 0;
84 
85  if (pt.x() > mid.x())
86  {
87  octant |= treeBoundBox::RIGHTHALF;
88  }
89 
90  if (pt.y() > mid.y())
91  {
92  octant |= treeBoundBox::TOPHALF;
93  }
94 
95  if (pt.z() > mid.z())
96  {
97  octant |= treeBoundBox::FRONTHALF;
98  }
99 
100  return octant;
101 }
102 
103 
104 // Returns octant in which point resides. Reverse of subBbox.
105 // Flags point exactly on edge.
107 (
108  const point& pt,
109  bool& onEdge
110 ) const
111 {
112  return subOctant(midpoint(), pt, onEdge);
113 }
114 
115 
116 // Returns octant in which point resides. Reverse of subBbox.
117 // Precalculated midpoint
119 (
120  const point& mid,
121  const point& pt,
122  bool& onEdge
123 )
124 {
125  direction octant = 0;
126  onEdge = false;
127 
128  if (pt.x() > mid.x())
129  {
130  octant |= treeBoundBox::RIGHTHALF;
131  }
132  else if (pt.x() == mid.x())
133  {
134  onEdge = true;
135  }
136 
137  if (pt.y() > mid.y())
138  {
139  octant |= treeBoundBox::TOPHALF;
140  }
141  else if (pt.y() == mid.y())
142  {
143  onEdge = true;
144  }
145 
146  if (pt.z() > mid.z())
147  {
148  octant |= treeBoundBox::FRONTHALF;
149  }
150  else if (pt.z() == mid.z())
151  {
152  onEdge = true;
153  }
154 
155  return octant;
156 }
157 
158 
159 // Returns octant in which intersection resides.
160 // Precalculated midpoint. If the point is on the dividing line between
161 // the octants the direction vector determines which octant to use
162 // (i.e. in which octant the point would be if it were moved along dir)
164 (
165  const point& mid,
166  const vector& dir,
167  const point& pt,
168  bool& onEdge
169 )
170 {
171  direction octant = 0;
172  onEdge = false;
173 
174  if (pt.x() > mid.x())
175  {
176  octant |= treeBoundBox::RIGHTHALF;
177  }
178  else if (pt.x() == mid.x())
179  {
180  onEdge = true;
181  if (dir.x() > 0)
182  {
183  octant |= treeBoundBox::RIGHTHALF;
184  }
185  }
186 
187  if (pt.y() > mid.y())
188  {
189  octant |= treeBoundBox::TOPHALF;
190  }
191  else if (pt.y() == mid.y())
192  {
193  onEdge = true;
194  if (dir.y() > 0)
195  {
196  octant |= treeBoundBox::TOPHALF;
197  }
198  }
199 
200  if (pt.z() > mid.z())
201  {
202  octant |= treeBoundBox::FRONTHALF;
203  }
204  else if (pt.z() == mid.z())
205  {
206  onEdge = true;
207  if (dir.z() > 0)
208  {
209  octant |= treeBoundBox::FRONTHALF;
210  }
211  }
212 
213  return octant;
214 }
215 
216 
217 // Returns reference to octantOrder which defines the
218 // order to do the search.
220 (
221  const point& pt,
222  FixedList<direction,8>& octantOrder
223 ) const
224 {
225  vector dist = midpoint() - pt;
226 
227  direction octant = 0;
228 
229  if (dist.x() < 0)
230  {
231  octant |= treeBoundBox::RIGHTHALF;
232  dist.x() *= -1;
233  }
234 
235  if (dist.y() < 0)
236  {
237  octant |= treeBoundBox::TOPHALF;
238  dist.y() *= -1;
239  }
240 
241  if (dist.z() < 0)
242  {
243  octant |= treeBoundBox::FRONTHALF;
244  dist.z() *= -1;
245  }
246 
247  direction min = 0;
248  direction mid = 0;
249  direction max = 0;
250 
251  if (dist.x() < dist.y())
252  {
253  if (dist.y() < dist.z())
254  {
256  mid = treeBoundBox::TOPHALF;
258  }
259  else if (dist.z() < dist.x())
260  {
263  max = treeBoundBox::TOPHALF;
264  }
265  else
266  {
269  max = treeBoundBox::TOPHALF;
270  }
271  }
272  else
273  {
274  if (dist.z() < dist.y())
275  {
277  mid = treeBoundBox::TOPHALF;
279  }
280  else if (dist.x() < dist.z())
281  {
282  min = treeBoundBox::TOPHALF;
285  }
286  else
287  {
288  min = treeBoundBox::TOPHALF;
291  }
292  }
293 
294  // Primary subOctant
295  octantOrder[0] = octant;
296  // subOctants joined to the primary by faces.
297  octantOrder[1] = octant ^ min;
298  octantOrder[2] = octant ^ mid;
299  octantOrder[3] = octant ^ max;
300  // subOctants joined to the primary by edges.
301  octantOrder[4] = octantOrder[1] ^ mid;
302  octantOrder[5] = octantOrder[1] ^ max;
303  octantOrder[6] = octantOrder[2] ^ max;
304  // subOctants joined to the primary by corner.
305  octantOrder[7] = octantOrder[4] ^ max;
306 }
307 
308 
309 // true if bb's intersect or overlap.
310 // Note: <= to make sure we catch all.
311 inline bool Foam::treeBoundBox::overlaps(const treeBoundBox& bb) const
312 {
313  return boundBox::overlaps(bb);
314 }
315 
316 
317 inline bool Foam::treeBoundBox::contains(const point& pt) const
318 {
319  return boundBox::contains(pt);
320 }
321 
322 
323 //- Return slightly wider bounding box
325 (
326  Random& rndGen,
327  const scalar s
328 ) const
329 {
330  treeBoundBox bb(*this);
331 
332  vector newSpan = bb.span();
333 
334  // Make 3D
335  scalar minSpan = s * Foam::mag(newSpan);
336 
337  for (direction dir = 0; dir < vector::nComponents; dir++)
338  {
339  newSpan[dir] = Foam::max(newSpan[dir], minSpan);
340  }
341 
342  bb.min() -= cmptMultiply(s * rndGen.vector01(), newSpan);
343  bb.max() += cmptMultiply(s * rndGen.vector01(), newSpan);
344 
345  return bb;
346 }
347 
348 
349 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
350 
351 // ************************ vim: set sw=4 sts=4 et: ************************ //