ESyS-Particle  2.3.2
RandomBoxPacker.hpp
Go to the documentation of this file.
1 
2 // //
3 // Copyright (c) 2003-2014 by The University of Queensland //
4 // Centre for Geoscience Computing //
5 // http://earth.uq.edu.au/centre-geoscience-computing //
6 // //
7 // Primary Business: Brisbane, Queensland, Australia //
8 // Licensed under the Open Software License version 3.0 //
9 // http://www.opensource.org/licenses/osl-3.0.php //
10 // //
12 
13 
14 #include "Foundation/console.h"
15 #include "Foundation/StringUtil.h"
18 #include "Geometry/SphereFitter.h"
19 
20 #include <algorithm>
21 #include <stdexcept>
22 #include <float.h>
23 
24 namespace esys
25 {
26  namespace lsm
27  {
28  template <typename TmplParticle>
29  class PlaneComparer
30  {
31  public:
32  typedef TmplParticle Particle;
33  PlaneComparer(const Particle &particle)
34  : m_pParticle(&particle)
35  {
36  }
37 
38  inline bool operator()(const Plane3D &plane1, const Plane3D &plane2) const
39  {
40  return
41  (
42  plane1.sep(m_pParticle->getPos())
43  <
44  plane2.sep(m_pParticle->getPos())
45  );
46  }
47  private:
49  };
50 
51  template<typename TmplTraits>
53  Packer &packer,
54  int maxInsertionFailures,
55  const PlaneVector &fitPlaneVector
56  ) :
57  m_pPacker(&packer),
58  m_fitPlaneVector(fitPlaneVector),
59  m_maxInsertionFailures(maxInsertionFailures),
60  m_lastFailCount(0),
61  m_successCount(0),
62  m_next(Fitter::getInvalidParticle()),
63  m_fitterPtrVector()
64  {
65  if (m_next.isValid())
66  {
67  throw
68  std::runtime_error(
69  std::string("isValid method returns true for INVALID particle:")
70  +
72  );
73  }
75  }
76 
77  template<typename TmplTraits>
80  {
81  return m_fitPlaneVector;
82  }
83 
84  template<typename TmplTraits>
86  {
87  return m_maxInsertionFailures;
88  }
89 
90  template<typename TmplTraits>
93  {
94  return *m_pPacker;
95  }
96 
97  template<typename TmplTraits>
100  {
101  return *m_pPacker;
102  }
103 
104  template<typename TmplTraits>
105  void
107  {
108  FitterPtrVector fitters;
109  fitters.push_back(FitterPtr(new Move2SurfaceFitter(getPacker())));
110  if (getPacker().is2d())
111  {
112  fitters.push_back(FitterPtr(new TwoDFitter(getPacker())));
113  if (getFitPlaneVector().size() > 0) {
114  fitters.push_back(FitterPtr(new TwoDPlaneFitter(getPacker())));
115  }
116  }
117  else
118  {
119  fitters.push_back(FitterPtr(new ThreeDFitter(getPacker())));
120  if (getFitPlaneVector().size() > 0) {
121  fitters.push_back(FitterPtr(new ThreeDPlaneFitter(getPacker())));
122  }
123  }
124 
125  m_fitterPtrVector = fitters;
126  }
127 
128  template<typename TmplTraits>
131  {
132  return m_fitterPtrVector;
133  }
134 
135  template<typename TmplTraits>
138  {
139  return m_fitterPtrVector;
140  }
141 
142  template<typename TmplTraits>
145  const Particle &particle
146  ) const
147  {
148  PlaneVector fitPlanes = getFitPlaneVector();
149  if (fitPlanes.size() > 0) {
150  std::sort(fitPlanes.begin(), fitPlanes.end(), PlaneComparer<Particle>(particle));
151 
152  return fitPlanes[0];
153  }
154  return Plane3D(Vec3(-1.0, 0, 0), Vec3(DBL_MAX/2, DBL_MAX/2, DBL_MAX/2));
155  }
156 
157  template<typename TmplTraits>
159  {
160  return getPacker().getRandomPoint();
161  }
162 
163  template<typename TmplTraits>
166  const Vec3 &point
167  )
168  {
169  return getPacker().getCandidateParticle(point);
170  }
171 
172  template<typename TmplTraits>
175  const Particle& particle,
176  int numClosest
177  )
178  {
179  return getPacker().getClosestNeighbours(particle, numClosest);
180  }
181 
182  template<typename TmplTraits>
185  {
186  m_next = Fitter::getInvalidParticle();
187  if (m_lastFailCount < getMaxInsertionFailures())
188  {
189  int numFails = 0;
190  //int insertCount = 0;
191  FitterPtrVector fitters = getFitterPtrVector();
192  Particle particle = Fitter::getInvalidParticle();
193  Particle fitParticle = particle;
194  Plane3D plane;
195  ParticleVector neighbours;
196  while ((!fitParticle.isValid()) && (numFails < getMaxInsertionFailures()))
197  {
198  particle = getCandidateParticle(getRandomPoint());
199  neighbours = getClosestNeighbours(particle, 4);
200  plane = getClosestFitPlane(particle);
201 
202  for (
203  typename FitterPtrVector::iterator it = getFitterPtrVector().begin();
204  (
205  (it != getFitterPtrVector().end())
206  &&
207  (!fitParticle.isValid())
208  );
209  it++
210  )
211  {
212  fitParticle = (*it)->getFitParticle(particle, neighbours, plane);
213  }
214  numFails++;
215  }
216  m_lastFailCount = numFails;
217  console.Info()
218  << "FittedParticleIterator<T>::generateNext: numFails="
219  << numFails << "\n";
220  m_successCount += ((fitParticle.isValid()) ? 1 : 0);
221  m_next = fitParticle;
222  }
223  return m_next;
224  }
225 
226  template<typename TmplTraits>
228  {
229  return (m_next.isValid() || generateNext().isValid());
230  }
231 
232  template<typename TmplTraits>
235  {
236  if (!(m_next.isValid()))
237  {
238  generateNext();
239  }
240  const Particle next = m_next;
241  m_next = Fitter::getInvalidParticle();
242  return next;
243  }
244 
245  template<typename TmplTraits>
247  {
248  console.Info()
249  << "BoundingBox: minPt = " << getPacker().getBBox().getMinPt()
250  << ", maxPt = " << getPacker().getBBox().getMaxPt() << "\n"
251  << "BoundingBox: sizes = " << getPacker().getBBox().getSizes() << "\n";
252 
253  for (
254  typename FitterPtrVector::iterator it = getFitterPtrVector().begin();
255  (it != getFitterPtrVector().end());
256  it++
257  ) {
258  console.Info() << (*(*it)).toString() << "\n";
259  }
260  console.Info() << "Total successful fits = " << m_successCount << "\n";
261  }
262 
263 
264 
265 
266 
267 
268 
269 
270 
271 
272 
273 
274 
275 
276 
277 
278 
279 
280 
281 
282 
283 
284 
285 
286 
287 
288 
289  template <typename TPartGen,template <typename TTPartGen> class TCubicPackWrap>
291  ParticleGeneratorPtr particleGeneratorPtr,
292  ParticlePoolPtr particlePoolPtr,
293  NTablePtr nTablePtr,
294  const BoundingBox &bBox,
295  const BoolVector &periodicDimensions,
296  double tolerance,
297  double cubicPackRadius,
298  int maxInsertionFailures,
299  const PlaneVector &fitPlaneVector
300  )
301  : Inherited(
302  particleGeneratorPtr,
303  particlePoolPtr,
304  nTablePtr,
305  bBox,
306  periodicDimensions,
307  tolerance,
308  cubicPackRadius
309  ),
310  m_fitPlaneVector(fitPlaneVector),
311  m_maxInsertionFailures(maxInsertionFailures)
312  {
313  }
314 
315  template <typename TPartGen,template <typename TTPartGen> class TCubicPackWrap>
317  ParticleGeneratorPtr particleGeneratorPtr,
318  ParticlePoolPtr particlePoolPtr,
319  NTablePtr nTablePtr,
320  const BoundingBox &bBox,
321  const BoolVector &periodicDimensions,
322  double tolerance,
323  double cubicPackRadius,
324  int maxInsertionFailures
325  )
326  : Inherited(
327  particleGeneratorPtr,
328  particlePoolPtr,
329  nTablePtr,
330  bBox,
331  periodicDimensions,
332  tolerance,
333  cubicPackRadius
334  ),
335  m_fitPlaneVector(),
336  m_maxInsertionFailures(maxInsertionFailures)
337  {
339  }
340 
341  template <typename TPartGen,template <typename TTPartGen> class TCubicPackWrap>
344  {
345  PlaneVector fitPlaneVector;
346  if ((!this->getPeriodicDimensions()[1])) {
347  fitPlaneVector.push_back(
348  Plane3D(Vec3(0, 1, 0), this->getBBox().getMinPt())
349  );
350  fitPlaneVector.push_back(
351  Plane3D(Vec3(0, -1, 0), this->getBBox().getMaxPt())
352  );
353  }
354  if ((!this->getPeriodicDimensions()[0])) {
355  fitPlaneVector.push_back(
356  Plane3D(Vec3( 1, 0, 0), this->getBBox().getMinPt())
357  );
358  fitPlaneVector.push_back(
359  Plane3D(Vec3(-1, 0, 0), this->getBBox().getMaxPt())
360  );
361  }
362  if (
363  !this->is2d()
364  &&
365  (!this->getPeriodicDimensions()[2])
366  ) {
367  fitPlaneVector.push_back(
368  Plane3D(Vec3(0, 0, 1), this->getBBox().getMinPt())
369  );
370  fitPlaneVector.push_back(
371  Plane3D(Vec3(0, 0, -1), this->getBBox().getMaxPt())
372  );
373  }
374  return fitPlaneVector;
375  }
376 
377  template <typename TPartGen,template <typename TTPartGen> class TCubicPackWrap>
379  {
380  }
381 
382  template <typename TPartGen,template <typename TTPartGen> class TCubicPackWrap>
383  double
385  double minVal,
386  double maxVal
387  ) const
388  {
389  return minVal + (maxVal-minVal)*rng::s_zeroOneUniform();
390  }
391 
392  template <typename TPartGen,template <typename TTPartGen> class TCubicPackWrap>
395  {
396  return m_fitPlaneVector;
397  }
398 
399 
400  template <typename TPartGen,template <typename TTPartGen> class TCubicPackWrap>
402  {
403  return
404  Vec3(
405  getRandom(this->getBBox().getMinPt().X(), this->getBBox().getMaxPt().X()),
406  getRandom(this->getBBox().getMinPt().Y(), this->getBBox().getMaxPt().Y()),
407  getRandom(this->getBBox().getMinPt().Z(), this->getBBox().getMaxPt().Z())
408  );
409  }
410 
411  template <typename TPartGen,template <typename TTPartGen> class TCubicPackWrap>
414  const Particle& particle,
415  int numClosest
416  )
417  {
418  ParticleVector neighbourVector =
419  this->getNTable().getUniqueNeighbourVector(
420  particle.getPos(),
421  this->getTolerance()
422  );
423 
424  if (static_cast<int>(neighbourVector.size()) < numClosest)
425  {
426  neighbourVector =
427  this->getNTable().getUniqueNeighbourVector(
428  particle.getPos(),
429  particle.getRad() + this->getTolerance()
430  );
431  if (static_cast<int>(neighbourVector.size()) < numClosest)
432  {
433  neighbourVector =
434  this->getNTable().getUniqueNeighbourVector(
435  particle.getPos(),
436  this->getParticleGenerator().getMaxFitRadius() + this->getTolerance()
437  );
438  }
439  }
440 
441  std::sort(
442  neighbourVector.begin(),
443  neighbourVector.end(),
445  );
446 
447  if (static_cast<int>(neighbourVector.size()) > numClosest) {
448  neighbourVector.erase(
449  neighbourVector.begin() + numClosest,
450  neighbourVector.end()
451  );
452  }
453 
454  return neighbourVector;
455  }
456 
457  template <typename TPartGen,template <typename TTPartGen> class TCubicPackWrap>
459  const Particle &particle
460  ) const
461  {
462  return
463  (
464  this->getParticleGenerator().isValidFitRadius(particle.getRad())
465  &&
466  Inherited::particleFitsInBBoxWithNeighbours(particle)
467  );
468  }
469 
470  template <typename TPartGen,template <typename TTPartGen> class TCubicPackWrap>
472  {
473  return m_maxInsertionFailures;
474  }
475 
476  template <typename TPartGen,template <typename TTPartGen> class TCubicPackWrap>
478  {
481  *this,
482  getMaxInsertionFailures(),
483  getFitPlaneVector()
484  );
485  while (it.hasNext())
486  {
487  const Particle p = it.next();
488  this->createAndInsertParticle(p);
489  }
490  it.logInfo();
491  }
492 
493  template <typename TPartGen,template <typename TTPartGen> class TCubicPackWrap>
495  {
496  this->generateCubicPacking();
497  this->generateRandomFill();
498  }
499 
500  }
501 }