15 #include "Foundation/console.h"
16 #include "Geometry/GeometryInfo.h"
17 #include "Geometry/CircularNeighbourTable.h"
29 PackingInfo::PackingInfo(
30 const BoundingBox &bBox,
31 const BoolVector &periodicDimensions,
32 Orientation orientation,
36 m_periodicDimensions(periodicDimensions),
37 m_orientation(orientation),
38 m_minRadius(minRadius),
39 m_maxRadius(maxRadius)
41 initialiseFitPlaneVector();
44 bool PackingInfo::is3d()
const
46 return (m_bBox.getSizes().Z() > 0.0);
49 void PackingInfo::initialiseFitPlaneVector()
51 m_fitPlaneVector.clear();
52 if ((m_orientation != XZ) && (!getPeriodicDimensions()[1])) {
53 m_fitPlaneVector.push_back(
54 Plane(
Vec3(0, 1, 0), getBBox().getMinPt())
56 m_fitPlaneVector.push_back(
57 Plane(
Vec3(0, -1, 0), getBBox().getMaxPt())
60 if ((m_orientation != YZ) && (!getPeriodicDimensions()[0])) {
61 m_fitPlaneVector.push_back(
62 Plane(
Vec3( 1, 0, 0), getBBox().getMinPt())
64 m_fitPlaneVector.push_back(
65 Plane(
Vec3(-1, 0, 0), getBBox().getMaxPt())
73 (!getPeriodicDimensions()[2])
75 m_fitPlaneVector.push_back(
76 Plane(
Vec3(0, 0, 1), getBBox().getMinPt())
78 m_fitPlaneVector.push_back(
79 Plane(
Vec3(0, 0, -1), getBBox().getMaxPt())
84 const BoundingBox &PackingInfo::getBBox()
const
89 const PlaneVector &PackingInfo::getFitPlaneVector()
const
91 return m_fitPlaneVector;
94 double PackingInfo::getMinParticleRadius()
const
99 double PackingInfo::getMaxParticleRadius()
const
104 const BoolVector &PackingInfo::getPeriodicDimensions()
const
106 return m_periodicDimensions;
111 template <
typename TGrainGen>
112 GougePackingInfo<TGrainGen>::GougePackingInfo(
113 const BoundingBox &bBox,
114 const BoolVector &periodicDimensions,
115 Orientation orientation,
116 ParticleGrainGen &particleGrainGen
121 particleGrainGen.getMinParticleRadius(),
122 particleGrainGen.getMaxParticleRadius()
124 m_pParticleGrainGen(&particleGrainGen)
129 template <
typename TGrainGen>
130 typename GougePackingInfo<TGrainGen>::ParticleGrainGen &
131 GougePackingInfo<TGrainGen>::getParticleGrainGen()
const
133 return *m_pParticleGrainGen;
137 template <
typename TGrainGen>
138 const typename GougePackingInfo<TGrainGen>::ParticleGrainGen &
139 GougePackingInfo<TGrainGen>::getParticleGrainGen()
const
141 return *m_pParticleGrainGen;
145 template <
typename TGrainGen>
146 double GougePackingInfo<TGrainGen>::getMinGrainRadius()
const
148 return getParticleGrainGen().getMinGrainRadius();
151 template <
typename TGrainGen>
152 double GougePackingInfo<TGrainGen>::getMaxGrainRadius()
const
154 return getParticleGrainGen().getMaxGrainRadius();
159 ParticleRndPackPrms::ParticleRndPackPrms()
161 m_minParticleRadius(0.0),
162 m_maxParticleRadius(0.0)
166 ParticleRndPackPrms::ParticleRndPackPrms(
171 m_minParticleRadius(minRadius),
172 m_maxParticleRadius(maxRadius)
176 ParticleRndPackPrms::~ParticleRndPackPrms()
180 double ParticleRndPackPrms::getSize()
const
185 double ParticleRndPackPrms::getMinParticleRadius()
const
187 return m_minParticleRadius;
190 double ParticleRndPackPrms::getMaxParticleRadius()
const
192 return m_maxParticleRadius;
197 template <
typename TPGrainGen>
198 GrainRndPackPrms<TPGrainGen>::GrainRndPackPrms()
201 m_pParticleGrainGen(NULL)
205 template <
typename TPGrainGen>
206 GrainRndPackPrms<TPGrainGen>::GrainRndPackPrms(
208 ParticleGrainGen &particleGrainGen,
212 particleGrainGen.getMinParticleRadius(),
213 particleGrainGen.getMaxParticleRadius()
215 m_pParticleGrainGen(&particleGrainGen),
216 m_connectionTag(connectionTag)
220 template <
typename TPGrainGen>
221 typename GrainRndPackPrms<TPGrainGen>::ParticleGrainGen &
222 GrainRndPackPrms<TPGrainGen>::getParticleGrainGen()
const
224 return *m_pParticleGrainGen;
227 template <
typename TPGrainGen>
229 GrainRndPackPrms<TPGrainGen>::getConnectionTag()
const
231 return m_connectionTag;
235 template <
typename TPGrainGen>
236 const typename GrainRndPackPrms<TPGrainGen>::ParticleGrainGen &
237 GrainRndPackPrms<TPGrainGen>::getParticleGrainGen()
const
239 return *m_pParticleGrainGen;
243 template <
typename TPGrainGen>
244 double GrainRndPackPrms<TPGrainGen>::getMinGrainRadius()
246 return getParticleGrainGen().getMinGrainRadius();
249 template <
typename TPGrainGen>
250 double GrainRndPackPrms<TPGrainGen>::getMaxGrainRadius()
252 return getParticleGrainGen().getMaxGrainRadius();
257 template <
typename TPGrainGen>
258 GougeConfigPrms<TPGrainGen>::GougeConfigPrms()
264 m_periodicDimensions(3, false),
265 m_maxInsertionFailures(50),
266 m_tolerance(DBL_EPSILON*128),
267 m_connectionTolerance(DBL_EPSILON*128*10),
268 m_blockConnectionTag(0)
272 template <
typename TPGrainGen>
273 GougeConfigPrms<TPGrainGen>::GougeConfigPrms(
274 const BoundingBox &bBox,
276 Orientation orientation,
277 const ParticleRndPackPrms &faultRegionPrms,
278 const GrainRPackPrms &gougeRegionPrms,
279 const BoolVector &periodicDimensions,
280 int maxInsertionFailures,
282 double connectionTolerance,
283 int blockConnectionTag
285 m_padRadius(padRadius),
286 m_orientation(orientation),
287 m_faultPrms(faultRegionPrms),
288 m_gougePrms(gougeRegionPrms),
289 m_periodicDimensions(periodicDimensions),
290 m_maxInsertionFailures(maxInsertionFailures),
291 m_tolerance(tolerance),
292 m_connectionTolerance(connectionTolerance),
293 m_blockConnectionTag(blockConnectionTag)
295 m_bBox = GridIterator(bBox, getMaxRadius()).getSphereBBox();
298 template <
typename TPGrainGen>
299 GougeConfigPrms<TPGrainGen>::~GougeConfigPrms()
303 template <
typename TPGrainGen>
304 const BoundingBox &GougeConfigPrms<TPGrainGen>::getBBox()
const
309 template <
typename TPGrainGen>
310 int GougeConfigPrms<TPGrainGen>::getGougeConnectionTag()
const
312 return m_gougePrms.getConnectionTag();
315 template <
typename TPGrainGen>
316 int GougeConfigPrms<TPGrainGen>::getBlockConnectionTag()
const
318 return m_blockConnectionTag;
321 template <
typename TPGrainGen>
322 int GougeConfigPrms<TPGrainGen>::getMaxInsertionFailures()
const
324 return m_maxInsertionFailures;
327 template <
typename TPGrainGen>
328 const BoolVector &GougeConfigPrms<TPGrainGen>::getPeriodicDimensions()
const
330 return m_periodicDimensions;
333 template <
typename TPGrainGen>
334 Orientation GougeConfigPrms<TPGrainGen>::getOrientation()
const
336 return m_orientation;
339 template <
typename TPGrainGen>
340 double GougeConfigPrms<TPGrainGen>::getTolerance()
const
345 template <
typename TPGrainGen>
346 double GougeConfigPrms<TPGrainGen>::getConnectionTolerance()
const
348 return m_connectionTolerance;
351 template <
typename TPGrainGen>
352 int GougeConfigPrms<TPGrainGen>::getOrientationIndex()
const
355 switch (m_orientation) {
373 std::stringstream msg;
374 msg <<
"Invalid orientation: " << m_orientation;
375 throw std::runtime_error(msg.str());
381 template <
typename TPGrainGen>
382 BoundingBox GougeConfigPrms<TPGrainGen>::cutFromCentre(
double d1,
double d2)
const
384 const int idx = getOrientationIndex();
385 const BoundingBox bBox = getBBox();
387 const double cmp1 = d1 + (bBox.getMaxPt()[idx] + bBox.getMinPt()[idx])/2.0;
388 const double cmp2 = d2 + (bBox.getMaxPt()[idx] + bBox.getMinPt()[idx])/2.0;
389 Vec3 minPt = bBox.getMinPt();
390 Vec3 maxPt = bBox.getMaxPt();
392 minPt[idx] = std::min(cmp1, cmp2);
393 maxPt[idx] = std::max(cmp1, cmp2);
396 tmpPt[idx] = minPt[idx];
398 if ((tmpPt - bBox.getMinPt())[idx] > getTolerance()) {
399 const BoundingBox minBBox = GridIterator(BoundingBox(bBox.getMinPt(), tmpPt), getMaxRadius()).getSphereBBox();
401 tmpPt[idx] = minBBox.getMaxPt()[idx];
404 tmpPt = bBox.getMinPt();
406 const BoundingBox maxBBox = GridIterator(BoundingBox(bBox.getMinPt(), maxPt), getMaxRadius()).getSphereBBox();
407 BoundingBox returnBBox = BoundingBox(tmpPt, maxBBox.getMaxPt());
412 template <
typename TPGrainGen>
413 double GougeConfigPrms<TPGrainGen>::getRegularBlockRadius()
const
418 template <
typename TPGrainGen>
419 double GougeConfigPrms<TPGrainGen>::getFaultMinRadius()
const
421 return m_faultPrms.getMinParticleRadius();
424 template <
typename TPGrainGen>
425 double GougeConfigPrms<TPGrainGen>::getFaultMaxRadius()
const
427 return m_faultPrms.getMaxParticleRadius();
430 template <
typename TPGrainGen>
431 double GougeConfigPrms<TPGrainGen>::getGougeMinRadius()
const
433 return m_gougePrms.getMinParticleRadius();
436 template <
typename TPGrainGen>
437 double GougeConfigPrms<TPGrainGen>::getGougeMaxRadius()
const
439 return m_gougePrms.getMaxParticleRadius();
442 template <
typename TPGrainGen>
443 double GougeConfigPrms<TPGrainGen>::getOrientationSize()
const
447 getBBox().getMaxPt()[getOrientationIndex()]
449 getBBox().getMinPt()[getOrientationIndex()]
453 template <
typename TPGrainGen>
454 BoundingBoxVector GougeConfigPrms<TPGrainGen>::getRegularBBoxVector()
const
456 BoundingBoxVector bBoxVector;
458 (getOrientationSize() - (m_gougePrms.getSize() + 2*m_faultPrms.getSize()))
463 bBoxVector.reserve(2);
464 bBoxVector.push_back(
466 -(m_gougePrms.getSize()/2.0 + m_faultPrms.getSize()),
467 -getOrientationSize()/2.0
470 bBoxVector.push_back(
472 m_gougePrms.getSize()/2.0 + m_faultPrms.getSize(),
473 getOrientationSize()/2.0
480 template <
typename TPGrainGen>
481 typename GougeConfigPrms<TPGrainGen>::GougePackingInfoVector
482 GougeConfigPrms<TPGrainGen>::getGougePackingInfoVector()
const
484 GougePackingInfoVector infoVec;
485 if (m_gougePrms.getSize() > 0.0) {
486 Vec3 overlap = Vec3::ZERO;
487 overlap[getOrientationIndex()] = m_faultPrms.getMaxParticleRadius();
490 m_gougePrms.getSize()/2.0,
491 -m_gougePrms.getSize()/2.0
496 BoundingBox(bBox.getMinPt() - overlap, bBox.getMaxPt() + overlap),
497 getPeriodicDimensions(),
499 m_gougePrms.getParticleGrainGen()
506 template <
typename TPGrainGen>
507 PackingInfoVector GougeConfigPrms<TPGrainGen>::getFaultPackingInfoVector()
const
509 PackingInfoVector infoVec;
510 if (m_faultPrms.getSize() > 0.0)
513 (getOrientationSize() - (m_gougePrms.getSize() + 2.0*m_faultPrms.getSize()))
519 const double roughnessSize = m_faultPrms.getSize();
520 Vec3 overlap = Vec3::ZERO;
521 overlap[getOrientationIndex()] = m_padRadius;
522 const BoundingBox bBox1 =
524 -m_gougePrms.getSize()/2.0,
525 -(m_gougePrms.getSize()/2.0 + roughnessSize)
527 const BoundingBox bBox2 =
529 m_gougePrms.getSize()/2.0,
530 m_gougePrms.getSize()/2.0 + roughnessSize
535 BoundingBox(bBox1.getMinPt() - overlap, bBox1.getMaxPt()),
536 getPeriodicDimensions(),
545 BoundingBox(bBox2.getMinPt(), bBox2.getMaxPt() + overlap),
546 getPeriodicDimensions(),
555 std::stringstream msg;
557 <<
"Roughness size plus gouge size is greater than block size: "
558 <<
"2*" << m_faultPrms.getSize() <<
" + " << m_gougePrms.getSize()
559 <<
" > " << getOrientationSize();
560 throw std::runtime_error(msg.str().c_str());
566 template <
typename TPGrainGen>
567 double GougeConfigPrms<TPGrainGen>::getMaxRadius()
const
572 std::max(m_faultPrms.getMaxParticleRadius(), m_gougePrms.getMaxParticleRadius())
576 template <
typename TPGrainGen>
577 double GougeConfigPrms<TPGrainGen>::getMinRadius()
const
582 std::min(m_faultPrms.getMinParticleRadius(), m_gougePrms.getMinParticleRadius())
586 template <
typename TPGrainGen>
587 bool GougeConfigPrms<TPGrainGen>::is2d()
const
589 return (getBBox().getSizes()[2] == 0.0);
594 template <
typename TGPckr,
typename TPPckr,
typename TConn>
595 int GougeConfig<TGPckr,TPPckr,TConn>::getNumParticles()
const
597 int numParticles = 0;
599 typename GeneratorPtrVector::const_iterator it = m_genPtrVector.begin();
600 it != m_genPtrVector.end();
604 numParticles += (*it)->getNumParticles();
609 template <
typename TGPckr,
typename TPPckr,
typename TConn>
610 int GougeConfig<TGPckr,TPPckr,TConn>::getNumGrains()
const
614 typename GrainRndPackerPtrVector::const_iterator packerIt =
615 getGougeGeneratorVector().begin();
616 packerIt != getGougeGeneratorVector().end();
620 numGrains += (*packerIt)->getNumGrains();
625 template <
typename TGPckr,
typename TPPckr,
typename TConn>
626 int GougeConfig<TGPckr,TPPckr,TConn>::getNumConnections()
const
628 return m_connectionSet.size();
631 template <
typename TGPckr,
typename TPPckr,
typename TConn>
632 GougeConfig<TGPckr,TPPckr,TConn>::GougeConfig(
const GougeConfPrms &prms)
636 m_gougeGenPtrVector(),
638 m_particlePoolPtr(new ParticlePool(8*4096)),
639 m_grainPoolPtr(new GrainPool(4096)),
640 m_regularGenPtrVector(),
641 m_faultGenPtrVector()
647 const BoundingBox bBox = m_prms.getBBox();
650 m_prms.getPeriodicDimensions()[0] ? 1 : 0,
651 m_prms.getPeriodicDimensions()[1] ? 1 : 0,
652 m_prms.getPeriodicDimensions()[2] ? 1 : 0
653 )*m_prms.getMaxRadius();
654 if (m_prms.getBBox().getSizes().Z() >= 4*m_prms.getMaxRadius()) {
655 ntableAdjust +=
Vec3(m_prms.getMaxRadius(), 0, 0);
657 const BoundingBox nTableBBox(bBox.getMinPt(), bBox.getMaxPt() - ntableAdjust);
662 (4.0*m_prms.getMinRadius()),
663 m_prms.getPeriodicDimensions(),
664 2.1*m_prms.getMaxRadius()
671 template <
typename TGPckr,
typename TPPckr,
typename TConn>
672 void GougeConfig<TGPckr,TPPckr,TConn>::createRegularBlockGenerators()
674 BoundingBoxVector bBoxVector = m_prms.getRegularBBoxVector();
676 BoundingBoxVector::const_iterator it = bBoxVector.begin();
677 it != bBoxVector.end();
681 <<
"GougeConfig<TGPckr,TPPckr,TConn>::createRegularBlockGenerators:"
682 <<
"Creating RegBoxPacker in box: " << StringUtil::toString(*it)
685 GeneratorPtr genPtr =
688 RegRadiusGenPtr(
new RegRadiusGen(m_prms.getRegularBlockRadius())),
692 m_prms.getPeriodicDimensions(),
693 m_prms.getTolerance(),
694 m_prms.getRegularBlockRadius()
697 m_genPtrVector.push_back(genPtr);
698 m_regularGenPtrVector.push_back(genPtr);
702 template <
typename TGPckr,
typename TPPckr,
typename TConn>
703 void GougeConfig<TGPckr,TPPckr,TConn>::createFaultBlockGenerators()
705 PackingInfoVector infoVec = m_prms.getFaultPackingInfoVector();
707 PackingInfoVector::const_iterator it = infoVec.begin();
712 <<
"GougeConfig<TGPckr,TPPckr,TConn>::createFaultBlockGenerators:"
713 <<
"Creating RndBoxPacker in box: " << StringUtil::toString(it->getBBox())
715 GeneratorPtr genPtr =
720 it->getMinParticleRadius(),
721 it->getMaxParticleRadius()
727 it->getPeriodicDimensions(),
728 m_prms.getTolerance(),
729 it->getMaxParticleRadius(),
730 m_prms.getMaxInsertionFailures(),
731 it->getFitPlaneVector()
734 m_genPtrVector.push_back(genPtr);
735 m_faultGenPtrVector.push_back(genPtr);
739 template <
typename TGPckr,
typename TPPckr,
typename TConn>
740 void GougeConfig<TGPckr,TPPckr,TConn>::createGougeConfigGenerators()
742 GougePackingInfoVector infoVec = m_prms.getGougePackingInfoVector();
744 typename GougePackingInfoVector::const_iterator it = infoVec.begin();
749 <<
"GougeConfig<TGPckr,TPPckr,TConn>::createGougeConfigGenerators:"
750 <<
"Creating GrainRandomPacker in box: " << StringUtil::toString(it->getBBox())
752 GrainRandomPackerPtr genPtr =
753 GrainRandomPackerPtr(
754 new GrainRandomPacker(
755 typename GrainRandomPacker::ParticleGrainGenPtr(),
759 it->getPeriodicDimensions(),
760 m_prms.getTolerance(),
761 it->getParticleGrainGen().getMaxGrainRadius(),
762 m_prms.getMaxInsertionFailures(),
763 it->getFitPlaneVector(),
767 genPtr->setParticleGrainGen(it->getParticleGrainGen());
768 m_genPtrVector.push_back(genPtr);
769 m_gougeGenPtrVector.push_back(genPtr);
773 template <
typename TGPckr,
typename TPPckr,
typename TConn>
774 GougeConfig<TGPckr,TPPckr,TConn>::~GougeConfig()
778 template <
typename TGPckr,
typename TPPckr,
typename TConn>
779 void GougeConfig<TGPckr,TPPckr,TConn>::generate()
782 createRegularBlockGenerators();
783 createFaultBlockGenerators();
784 createGougeConfigGenerators();
787 console.
Info() <<
"bbox = " << m_prms.getBBox().getMinPt() <<
" " << m_prms.getBBox().getMaxPt() <<
"\n";
789 typename GeneratorPtrVector::iterator it = m_genPtrVector.begin();
790 it != m_genPtrVector.end();
796 createConnectionSet();
799 template <
typename TGPckr,
typename TPPckr,
typename TConn>
800 void GougeConfig<TGPckr,TPPckr,TConn>::writeToFile(
const std::string &fileName)
const
802 std::ofstream fStream(fileName.c_str());
806 template <
typename TGPckr,
typename TPPckr,
typename TConn>
807 const typename GougeConfig<TGPckr,TPPckr,TConn>::GrainRndPackerPtrVector &
808 GougeConfig<TGPckr,TPPckr,TConn>::getGougeGeneratorVector()
const
810 return m_gougeGenPtrVector;
813 template <
typename TGPckr,
typename TPPckr,
typename TConn>
814 typename GougeConfig<TGPckr,TPPckr,TConn>::GrainRndPackerPtrVector &
815 GougeConfig<TGPckr,TPPckr,TConn>::getGougeGeneratorVector()
817 return m_gougeGenPtrVector;
820 template <
typename TGPckr,
typename TPPckr,
typename TConn>
821 const typename GougeConfig<TGPckr,TPPckr,TConn>::GeneratorPtrVector &
822 GougeConfig<TGPckr,TPPckr,TConn>::getFaultGeneratorVector()
const
824 return m_faultGenPtrVector;
827 template <
typename TGPckr,
typename TPPckr,
typename TConn>
828 bool GougeConfig<TGPckr,TPPckr,TConn>::areInDifferentFaultBlocks(
833 const GeneratorPtrVector &generators = getFaultGeneratorVector();
834 if (generators.size() == 2) {
837 (generators[0]->contains(p1) && generators[1]->contains(p2))
839 (generators[0]->contains(p2) && generators[1]->contains(p1))
842 else if (generators.size() > 2) {
845 "GougeConfig<TGPckr,TPPckr,TConn>::areInDifferentFaultBlocks: "
846 "More than two fault blocks."
852 template <
typename TGPckr,
typename TPPckr,
typename TConn>
853 bool GougeConfig<TGPckr,TPPckr,TConn>::isGougeParticle(
const Particle &particle)
const
855 const GrainRndPackerPtrVector &generators = getGougeGeneratorVector();
857 typename GrainRndPackerPtrVector::const_iterator it = generators.begin();
858 it != generators.end();
862 if ((*it)->contains(particle)) {
869 template <
typename TGPckr,
typename TPPckr,
typename TConn>
870 void GougeConfig<TGPckr,TPPckr,TConn>::createConnectionSet()
875 typename NTable::ParticleIterator particleIt = m_nTablePtr->getParticleIterator();
876 ConnectionValidator validator = ConnectionValidator(*
this, m_prms.getConnectionTolerance());
877 while (particleIt.hasNext()) {
878 const typename NTable::Particle *pParticle = particleIt.next();
879 const typename NTable::ParticleVector neighbours =
880 m_nTablePtr->getNeighbourVector(
882 pParticle->getRad() + m_prms.getConnectionTolerance()
885 typename NTable::ParticleVector::const_iterator it = neighbours.begin();
886 it != neighbours.end();
890 if (validator.isValid(*pParticle, *(*it))) {
891 m_connectionSet.insert(
892 typename ConnectionSet::value_type(
895 m_prms.getBlockConnectionTag()
901 const int numBlockConns = m_connectionSet.size();
903 <<
"Created " << numBlockConns <<
" connections in "
904 <<
"bonded blocks.\n";
910 <<
"Prms BBox: " << StringUtil::toString(m_prms.getBBox()) <<
"\n";
912 <<
"NTbl BBox: " << StringUtil::toString(m_nTablePtr->getBBox()) <<
"\n";
914 typename GrainRndPackerPtrVector::iterator packerIt = getGougeGeneratorVector().begin();
915 packerIt != getGougeGeneratorVector().end();
919 typename GrainRandomPacker::GrainIterator grainIt = (*packerIt)->getGrainIterator();
920 while (grainIt.hasNext())
922 ConnectionFinder connFinder =
924 m_prms.getConnectionTolerance(),
925 m_prms.getGougeConnectionTag(),
926 m_nTablePtr->getBBox(),
927 m_nTablePtr->getPeriodicDimensions()
929 Grain &g = grainIt.next();
930 connFinder.create(g.getParticleIterator());
931 typename ConnectionFinder::Iterator connIt = connFinder.getIterator();
932 while (connIt.hasNext())
934 m_connectionSet.insert(connIt.next());
936 if (connFinder.getNumConnections() == 0)
939 <<
"Found no connections in grain " << g.getId()
941 typename Grain::ParticleIterator partIt = g.getParticleIterator();
942 while (partIt.hasNext())
944 console.
Info() << StringUtil::toString(partIt.next()) <<
"\n";
950 <<
"Created " << m_connectionSet.size()-numBlockConns <<
" connections in "
951 <<
"gouge region.\n";
954 template <
typename TGPckr,
typename TPPckr,
typename TConn>
955 typename GougeConfig<TGPckr,TPPckr,TConn>::ParticleCollection
956 GougeConfig<TGPckr,TPPckr,TConn>::getParticleCollection()
958 ParticleCollection pCollection(m_particlePoolPtr);
960 typename GeneratorPtrVector::iterator it = m_genPtrVector.begin();
961 it != m_genPtrVector.end();
965 ParticleIterator particleIt = (*it)->getParticleIterator();
966 while (particleIt.hasNext()) {
967 pCollection.insertRef(particleIt.next());
974 template <
typename TGPckr,
typename TPPckr,
typename TConn>
975 typename GougeConfig<TGPckr,TPPckr,TConn>::GrainCollection
976 GougeConfig<TGPckr,TPPckr,TConn>::getGrainCollection()
978 GrainCollection gCollection(m_particlePoolPtr, m_grainPoolPtr);
980 typename GrainRndPackerPtrVector::iterator packerIt =
981 getGougeGeneratorVector().begin();
982 packerIt != getGougeGeneratorVector().end();
986 GrainIterator grainIt = (*packerIt)->getGrainIterator();
987 while (grainIt.hasNext()) {
988 gCollection.insertRef(grainIt.next());
995 template <
typename TGPckr,
typename TPPckr,
typename TConn>
996 const typename GougeConfig<TGPckr,TPPckr,TConn>::ConnectionSet &
997 GougeConfig<TGPckr,TPPckr,TConn>::getConnectionSet()
const
999 return m_connectionSet;
1002 template <
typename TGPckr,
typename TPPckr,
typename TConn>
1003 void GougeConfig<TGPckr,TPPckr,TConn>::write(std::ostream &oStream)
const
1005 Vec3 minPt = m_nTablePtr->getBBox().getMinPt();
1006 Vec3 maxPt = m_nTablePtr->getBBox().getMaxPt();
1007 if (fabs(maxPt.Z() - minPt.Z()) < (2*m_prms.getMaxRadius())) {
1008 minPt.Z() = minPt.Z() - m_prms.getMaxRadius() - m_prms.getTolerance();
1009 maxPt.Z() = maxPt.Z() + m_prms.getMaxRadius() + m_prms.getTolerance();
1012 const BoundingBox geoBBox(minPt, maxPt + m_prms.getTolerance());
1019 m_prms.getPeriodicDimensions(),
1020 (m_prms.getBBox().getSizes().Z() <= 0.0)
1022 info.write(oStream);
1029 typename NTable::ParticleIterator particleIt = m_nTablePtr->getParticleIterator();
1030 typedef std::set<Particle *, IdCompare> ParticleSet;
1031 ParticleSet taggedParticleSet;
1032 while (particleIt.hasNext()) {
1033 taggedParticleSet.insert(particleIt.next());
1041 particleIt = m_nTablePtr->getParticleIterator();
1042 ParticleSet particleSet;
1043 while (particleIt.hasNext()) {
1044 Particle *pParticle = particleIt.next();
1045 if (geoBBox.contains(pParticle->getPos())) {
1046 pParticle->setTag((*(taggedParticleSet.find(pParticle)))->getTag());
1047 particleSet.insert(pParticle);
1060 << particleSet.size()
1062 const int precision = 12;
1063 GeoParticleWriter particleVisitor(oStream, precision);
1065 typename ParticleSet::const_iterator it = particleSet.begin();
1066 it != particleSet.end();
1070 particleVisitor.visitParticle(*(*it));
1073 oStream <<
"EndParticles\n" <<
"BeginConnect\n";
1074 oStream << getConnectionSet().size() <<
"\n";
1077 GeoConnectionWriter connectionVisitor(oStream);
1078 visitConnections(connectionVisitor);
1079 oStream <<
"EndConnect";
1083 template <
typename TGPckr,
typename TPPckr,
typename TConn>
1084 void GougeConfig<TGPckr,TPPckr,TConn>::tagGougeParticles(
int tag)
1087 typename GrainRndPackerPtrVector::iterator it = m_gougeGenPtrVector.begin();
1088 it != m_gougeGenPtrVector.end();
1092 typename GrainRandomPacker::ParticleIterator particleIt =
1093 (*it)->getParticleIterator();
1094 while (particleIt.hasNext()) {
1095 particleIt.next().setTag(tag);
1100 template <
typename TGPckr,
typename TPPckr,
typename TConn>
1101 void GougeConfig<TGPckr,TPPckr,TConn>::tagRndBlockParticles(
int tag)
1104 typename GeneratorPtrVector::iterator it = m_faultGenPtrVector.begin();
1105 it != m_faultGenPtrVector.end();
1109 ParticleIterator particleIt = (*it)->getParticleIterator();
1110 while (particleIt.hasNext()) {
1111 particleIt.next().setTag(tag);
1116 template <
typename TGPckr,
typename TPPckr,
typename TConn>
1117 void GougeConfig<TGPckr,TPPckr,TConn>::tagDrivingPlateParticles(
1120 double distanceFromBBoxEdge
1123 ParticleCollection particleCollection = getParticleCollection();
1124 const BoundingBox bBox = particleCollection.getParticleBBox();
1126 const int idx = this->m_prms.getOrientationIndex();
1127 const double maxLow = bBox.getMinPt()[idx] + distanceFromBBoxEdge;
1128 const double minHigh = bBox.getMaxPt()[idx] - distanceFromBBoxEdge;
1130 int lowTagCount = 0;
1131 int highTagCount = 0;
1132 typename ParticleCollection::ParticleIterator particleIt =
1133 particleCollection.getParticleIterator();
1134 while (particleIt.hasNext())
1136 Particle &particle = particleIt.next();
1137 const double dimPos = particle.getPos()[idx];
1138 const double radius = particle.getRad();
1140 if (dimPos - radius <= maxLow) {
1141 particle.setTag(lowDrivingTag);
1144 if (dimPos + radius >= minHigh) {
1145 particle.setTag(highDrivingTag);
1149 console.
Info() <<
"Tagged " << lowTagCount <<
" particles with " << lowDrivingTag <<
"\n";
1150 console.
Info() <<
"Tagged " << highTagCount <<
" particles with " << highDrivingTag <<
"\n";