Author: eudoxos
Date: 2009-05-24 11:19:26 +0200 (Sun, 24 May 2009)
New Revision: 1775

Modified:
   trunk/SConstruct
   trunk/core/Interaction.cpp
   trunk/core/Interaction.hpp
   trunk/core/InteractionContainer.cpp
   trunk/core/InteractionContainer.hpp
   trunk/examples/collider-perf/mkGraph.py
   trunk/extra/Brefcom.cpp
   trunk/extra/Brefcom.hpp
   trunk/gui/py/utils.py
   trunk/pkg/common/Engine/StandAloneEngine/InsertionSortCollider.cpp
   trunk/pkg/common/Engine/StandAloneEngine/InsertionSortCollider.hpp
   trunk/pkg/common/Engine/StandAloneEngine/SpatialQuickSortCollider.cpp
   trunk/pkg/dem/Engine/StandAloneEngine/ElasticContactLaw.cpp
   trunk/scripts/test/insertion-sort-collider.py
Log:
1. Add InteractionContainer::requestErase to hint colliders that wouldn't 
otherwise see !isReal interactions.
2. Use this logic in InsertionSortCollider, PersistentSAPCollider, 
SpatialQuickSortCollider (noop in the last one)
3. Add InsertionSortCollider::sortThenCollide to make it behave as 
non-persistent collider (for debugging only?)
4. Add Interaction::reset() that has common initialization code.
5. Assign zero inertia to utils.facet (better than uninitialized binary garbage)
6. Fix contact logic in Brefcom, finally I get the same results with 
InsertionSortCollider as with SpatialQuickSortCollider on large simulation 
(more fixes to come)
7. Do not install pkg-config; fixes compilation error reported by the build bot.



Modified: trunk/SConstruct
===================================================================
--- trunk/SConstruct    2009-05-23 10:29:03 UTC (rev 1774)
+++ trunk/SConstruct    2009-05-24 09:19:26 UTC (rev 1775)
@@ -536,6 +536,7 @@
        if 0: # do not install headers, nor make pkg-config (was never used, I 
think)
                installHeaders(env.subst('$PREFIX')) # install to $PREFIX if 
specifically requested: like "scons /usr/local/include"
                makePkgConfig('$buildDir/yade${SUFFIX}.pc')
+               env.Install(pcDir,'$buildDir/yade${SUFFIX}.pc')
        if not env['haveForeach']:
                boostDir=buildDir+'/include/yade-'+env['version']+'/boost'
                foreachLink=boostDir+'/foreach.hpp'
@@ -545,7 +546,6 @@
                        if lexists(foreachLink): os.remove(foreachLink) # 
broken symlink: remove it
                        
os.symlink(relpath(foreachLink,foreachTarget),foreachLink)
                
env.InstallAs(env['PREFIX']+'/include/yade-'+env['version']+'/boost/foreach.hpp',foreachTarget)
-       env.Install(pcDir,'$buildDir/yade${SUFFIX}.pc')
        installAlias=env.Alias('install',instDirs) # build and install 
everything that should go to instDirs, which are $PREFIX/{bin,lib} (uses scons' 
Install); include pkgconfig stuff
        env.Default([installAlias,'$PREFIX'])
 

Modified: trunk/core/Interaction.cpp
===================================================================
--- trunk/core/Interaction.cpp  2009-05-23 10:29:03 UTC (rev 1774)
+++ trunk/core/Interaction.cpp  2009-05-24 09:19:26 UTC (rev 1775)
@@ -10,26 +10,15 @@
 
 #include "Interaction.hpp"
 
-Interaction::Interaction ()
-{
-       // FIXME : -1
-       id1 = 0; 
-       id2 = 0;
-       isNew = true;
-       isReal = false; // maybe we can remove this, and check if 
InteractingGeometry, and InteractionPhysics are empty?
-       
-       isNeighbor = true;//NOTE : TriangulationCollider needs that
+Interaction::Interaction(): id1(0), id2(0){ reset(); }
+Interaction::Interaction(body_id_t newId1,body_id_t newId2): id1(newId1), 
id2(newId2){ reset(); }
 
-}
-
-
-Interaction::Interaction(body_id_t newId1,body_id_t newId2) : id1(newId1) , 
id2(newId2)
-{      
-       isNew = true;
-       isReal = false;
+void Interaction::reset(){
+       isNew=true;
+       isReal=false;
        isNeighbor = true;//NOTE : TriangulationCollider needs that
-
        functorCache.geomExists=true;
+       //functorCache.geom=shared_ptr<InteractionGeometryEngineUnit>(); 
functorCache.phys=shared_ptr<InteractionPhysicsEngineUnit>(); 
functorCache.constLaw=shared_ptr<ConstitutiveLaw>();
 }
 
 

Modified: trunk/core/Interaction.hpp
===================================================================
--- trunk/core/Interaction.hpp  2009-05-23 10:29:03 UTC (rev 1774)
+++ trunk/core/Interaction.hpp  2009-05-24 09:19:26 UTC (rev 1775)
@@ -21,13 +21,16 @@
 class Interaction : public Serializable
 {
        private :
-               body_id_t id1,id2;      // FIXME  this should be 
vector<body_id_t> ids;
-
+               body_id_t id1,id2;
        public :
-               bool isNew;             // FIXME : better to test if 
InteractionPhysics==0 and remove this flag; we can also remove this flag, if we 
make another container for PotetntialInteraction with only ids
-               bool isReal;            // maybe we can remove this, and check 
if InteractingGeometry, and InteractionPhysics are empty?
-               bool cycle; // phase flag to mark (for example, 
SpatialQuickSortCollider mark by it the stale interactions) 
-               bool isNeighbor;        // Has a meaning only with 
triangulationCollider atm NOTE : TriangulationCollider needs that
+               // FIXME : test if InteractionPhysics==0 and remove this flag; 
we can also remove this flag, if we make another container for 
PotetntialInteraction with only ids
+               bool isNew;             
+               // maybe we can remove this, and check if InteractingGeometry, 
and InteractionPhysics are empty?
+               bool isReal;            
+               //! phase flag to mark (for example, SpatialQuickSortCollider 
mark by it the stale interactions) 
+               bool cycle;      
+               //! NOTE : TriangulationCollider needs this (nothing else)
+               bool isNeighbor;        
 
                shared_ptr<InteractionGeometry> interactionGeometry;
                shared_ptr<InteractionPhysics> interactionPhysics;
@@ -44,22 +47,17 @@
                //! cache functors that are called for this interaction. 
Currently used by InteractionDispatchers.
                struct {
                        // Whether geometry dispatcher exists at all; this is 
different from !geom, since that can mean we haven't populated the cache yet.
-                       // Therefore, geomExists must be initialized to true 
first (done in Interaction ctor).
+                       // Therefore, geomExists must be initialized to true 
first (done in Interaction::reset() called from ctor).
                        bool geomExists;
                        // shared_ptr's are initialized to NULLs automagically
                        shared_ptr<InteractionGeometryEngineUnit> geom;
                        shared_ptr<InteractionPhysicsEngineUnit> phys;
                        shared_ptr<ConstitutiveLaw> constLaw;
                } functorCache;
-                       
 
-               #if 0
-                       //! Whether both bodies involved in interaction 
satisfies given mask; provide rootBody for faster lookup
-                       bool maskBothOK(int mask, MetaBody* 
rootBody=NULL){return (mask==0) || (Body::byId(id1,rootBody)->maskOK(mask) && 
Body::byId(id2,rootBody)->maskOK(mask));}
-                       //! Whether at least one body in interaction satisfies 
given mask; provide rootBody for faster lookup
-                       bool maskAnyOK(int mask, MetaBody* 
rootBody=NULL){return (mask==0) || Body::byId(id1,rootBody)->maskOK(mask) || 
Body::byId(id2,rootBody)->maskOK(mask);}
-               #endif
-
+               //! Reset interaction to the intial state (keep only body ids)
+               void reset();
+                       
        REGISTER_ATTRIBUTES(/*no base*/,
                (id1)
                (id2)

Modified: trunk/core/InteractionContainer.cpp
===================================================================
--- trunk/core/InteractionContainer.cpp 2009-05-23 10:29:03 UTC (rev 1774)
+++ trunk/core/InteractionContainer.cpp 2009-05-24 09:19:26 UTC (rev 1775)
@@ -13,6 +13,8 @@
 
 
 
+void InteractionContainer::requestErase(body_id_t id1, body_id_t id2){ 
find(id1,id2)->reset(); bodyIdPair v(0,2); v.push_back(id1); v.push_back(id2); 
pendingErase.push_back(v); }
+
 void InteractionContainer::preProcessAttributes(bool deserializing)
 {
        if(deserializing)

Modified: trunk/core/InteractionContainer.hpp
===================================================================
--- trunk/core/InteractionContainer.hpp 2009-05-23 10:29:03 UTC (rev 1774)
+++ trunk/core/InteractionContainer.hpp 2009-05-24 09:19:26 UTC (rev 1775)
@@ -83,7 +83,7 @@
        public :
                boost::mutex    drawloopmutex; // FIXME a hack, containers have 
to be rewritten lock-free.
 
-               InteractionContainer() { interaction.clear(); };
+               InteractionContainer() { };
                virtual ~InteractionContainer() {};
 
                virtual bool insert(body_id_t /*id1*/,body_id_t /*id2*/)        
                        {throw;};
@@ -101,13 +101,19 @@
                virtual shared_ptr<Interaction>& operator[] (unsigned int) 
{throw;};
                virtual const shared_ptr<Interaction>& operator[] (unsigned 
int) const {throw;};
 
+               // std::pair is not handle by yade::serialization, use 
vector<body_id_t> instead
+               typedef vector<body_id_t> bodyIdPair;
+               // Ask for erasing the interaction given (from the constitutive 
law); this resets the interaction (to the initial=potential state)
+               // and collider should traverse pendingErase to decide whether 
to delete the interaction completely or keep it potential
+               void requestErase(body_id_t id1, body_id_t id2);
+               list<bodyIdPair> pendingErase;
        private :
+               // used only during serialization/deserialization
                vector<shared_ptr<Interaction> > interaction;
-
        protected :
                virtual void preProcessAttributes(bool deserializing);
                virtual void postProcessAttributes(bool deserializing);
-       REGISTER_ATTRIBUTES(/*no base*/,(interaction));
+       REGISTER_ATTRIBUTES(/*no base*/,(interaction)(pendingErase));
        REGISTER_CLASS_AND_BASE(InteractionContainer,Serializable);
 };
 

Modified: trunk/examples/collider-perf/mkGraph.py
===================================================================
--- trunk/examples/collider-perf/mkGraph.py     2009-05-23 10:29:03 UTC (rev 
1774)
+++ trunk/examples/collider-perf/mkGraph.py     2009-05-24 09:19:26 UTC (rev 
1775)
@@ -31,7 +31,7 @@
 legend(('SAP init','IS init'),'upper left')
 
 ax2=twinx()
-plot(SAP_N,SAPstep,'r-',IS_N,ISstep,'g-',QS_N,QSstep,'g-',QS_N,QSinit,'b-')
+plot(SAP_N,SAPstep,'r-',IS_N,ISstep,'k-',QS_N,QSstep,'g-',QS_N,QSinit,'b-')
 ylabel(u"Linear time per 1 step [s]")
 legend(('SAP step','InsertionSort step','QuickSort step','QuickSort 
init'),'right')
 grid()

Modified: trunk/extra/Brefcom.cpp
===================================================================
--- trunk/extra/Brefcom.cpp     2009-05-23 10:29:03 UTC (rev 1774)
+++ trunk/extra/Brefcom.cpp     2009-05-24 09:19:26 UTC (rev 1775)
@@ -73,7 +73,7 @@
 
        if(!interaction->isNew && interaction->interactionPhysics){ /* relax */ 
} 
        else {
-               interaction->isNew; // just in case
+               interaction->isNew=false; // just in case
 
                const shared_ptr<BodyMacroParameters>& 
elast1=static_pointer_cast<BodyMacroParameters>(pp1);
                const shared_ptr<BodyMacroParameters>& 
elast2=static_pointer_cast<BodyMacroParameters>(pp2);
@@ -233,15 +233,17 @@
        NNAN(sigmaN); NNANV(sigmaT); NNAN(crossSection);
        //timingDeltas->checkpoint("material");
 
-       const int watch1=6300, watch2=6299;
-       #define SHOW(a) if((I->getId1()==watch1 && I->getId2()==watch2) || 
(I->getId2()==watch1 && I->getId1()==watch2)) cerr<<__FILE__<<":"<<__LINE__<<" 
"<<a<<endl;
-       SHOW("epsN"<<epsN);
+       //const int watch1=6300, watch2=6299;
+       //#define SHOW(a) if((I->getId1()==watch1 && I->getId2()==watch2) || 
(I->getId2()==watch1 && I->getId1()==watch2)) cerr<<__FILE__<<":"<<__LINE__<<" 
"<<a<<endl;
+       //SHOW("epsN"<<epsN);
        if(epsN>0. && ((isCohesive && omega>omegaThreshold) || !isCohesive)){
-               I->isReal=false;
-               const shared_ptr<Body>& body1=Body::byId(I->getId1(),rootBody), 
body2=Body::byId(I->getId2(),rootBody); assert(body1); assert(body2);
-               const shared_ptr<BrefcomPhysParams>& 
rbp1=YADE_PTR_CAST<BrefcomPhysParams>(body1->physicalParameters), 
rbp2=YADE_PTR_CAST<BrefcomPhysParams>(body2->physicalParameters);
-               if(BC->isCohesive){rbp1->numBrokenCohesive+=1; 
rbp2->numBrokenCohesive+=1; rbp1->epsPlBroken+=epsPlSum; 
rbp2->epsPlBroken+=epsPlSum;}
-               LOG_DEBUG("Contact #"<<I->getId1()<<"=#"<<I->getId2()<<" is 
damaged over thershold ("<<omega<<">"<<omegaThreshold<<") and has been deleted 
(isReal="<<I->isReal<<")");
+               rootBody->interactions->requestErase(I->getId1(),I->getId2());
+               if(isCohesive){
+                       const shared_ptr<Body>& 
body1=Body::byId(I->getId1(),rootBody), body2=Body::byId(I->getId2(),rootBody); 
assert(body1); assert(body2);
+                       const shared_ptr<BrefcomPhysParams>& 
rbp1=YADE_PTR_CAST<BrefcomPhysParams>(body1->physicalParameters), 
rbp2=YADE_PTR_CAST<BrefcomPhysParams>(body2->physicalParameters);
+                       if(BC->isCohesive){rbp1->numBrokenCohesive+=1; 
rbp2->numBrokenCohesive+=1; rbp1->epsPlBroken+=epsPlSum; 
rbp2->epsPlBroken+=epsPlSum;}
+                       LOG_DEBUG("Contact 
#"<<I->getId1()<<"=#"<<I->getId2()<<" is damaged over thershold 
("<<omega<<">"<<omegaThreshold<<") and will be deleted.");
+               }
                return;
        }
 

Modified: trunk/extra/Brefcom.hpp
===================================================================
--- trunk/extra/Brefcom.hpp     2009-05-23 10:29:03 UTC (rev 1774)
+++ trunk/extra/Brefcom.hpp     2009-05-24 09:19:26 UTC (rev 1775)
@@ -106,7 +106,7 @@
 
 
 
-               BrefcomContact(): NormalShearInteraction(),E(0), G(0), 
tanFrictionAngle(0), undamagedCohesion(0), crossSection(0), xiShear(0), 
dmgTau(-1), dmgRateExp(0), dmgStrain(0), plTau(-1), plRateExp(0), 
isoPrestress(0.), epsPlSum(0.) { createIndex(); epsT=Vector3r::ZERO; kappaD=0; 
isCohesive=false; neverDamage=false; omega=0; Fn=0; Fs=Vector3r::ZERO; 
epsPlSum=0; dmgOverstress=0; dmgStrain=0; }
+               BrefcomContact(): NormalShearInteraction(),E(0), G(0), 
tanFrictionAngle(0), undamagedCohesion(0), crossSection(0), xiShear(0), 
dmgTau(-1), dmgRateExp(0), dmgStrain(0), plTau(-1), plRateExp(0), 
isoPrestress(0.), kappaD(0.), epsTrans(0.), epsPlSum(0.) { createIndex(); 
epsT=Vector3r::ZERO; isCohesive=false; neverDamage=false; omega=0; Fn=0; 
Fs=Vector3r::ZERO; epsPlSum=0; dmgOverstress=0; }
                virtual ~BrefcomContact();
 
                REGISTER_ATTRIBUTES(NormalShearInteraction,

Modified: trunk/gui/py/utils.py
===================================================================
--- trunk/gui/py/utils.py       2009-05-23 10:29:03 UTC (rev 1774)
+++ trunk/gui/py/utils.py       2009-05-24 09:19:26 UTC (rev 1775)
@@ -124,7 +124,7 @@
        vStr='['+' '.join(['{%g %g %g}'%(v[0],v[1],v[2]) for v in vertices])+']'
        b.shape.setRaw('vertices',vStr)
        b.mold.setRaw('vertices',vStr)
-       
pp={'se3':[center[0],center[1],center[2],1,0,0,0],'refSe3':[center[0],center[1],center[2],1,0,0,0],'young':young,'poisson':poisson,'frictionAngle':frictionAngle}
+       
pp={'se3':[center[0],center[1],center[2],1,0,0,0],'refSe3':[center[0],center[1],center[2],1,0,0,0],'young':young,'poisson':poisson,'frictionAngle':frictionAngle,'inertia':[0,0,0]}
        pp.update(physParamsAttr)
        b.phys=PhysicalParameters(physParamsClass)
        for k in [attr for attr in pp.keys() if attr in b.phys.keys()]:

Modified: trunk/pkg/common/Engine/StandAloneEngine/InsertionSortCollider.cpp
===================================================================
--- trunk/pkg/common/Engine/StandAloneEngine/InsertionSortCollider.cpp  
2009-05-23 10:29:03 UTC (rev 1774)
+++ trunk/pkg/common/Engine/StandAloneEngine/InsertionSortCollider.cpp  
2009-05-24 09:19:26 UTC (rev 1775)
@@ -31,7 +31,8 @@
        const shared_ptr<Interaction>& I=interactions->find(id1,id2);
        bool hasInter=(bool)I;
        // interaction doesn't exist and shouldn't, or it exists and should
-       if((!overlap && !hasInter) || (overlap && hasInter)) return;
+       if(!overlap && !hasInter) return;
+       if(overlap && hasInter){ /* FIXME: should check I->isNew and I->isReal; 
etc */ return; }
        // create interaction if not yet existing
        if(overlap && !hasInter){ // second condition only for readability
                
if(!Collider::mayCollide(Body::byId(id1,rb).get(),Body::byId(id2,rb).get())) 
return;
@@ -40,21 +41,18 @@
                interactions->insert(newI);
                return;
        }
-       /* Note: this doesn't cover all disappearing interactions, only those 
that broke in the sortAxis direction;
-               it is only a minor optimization (to be verified) to have it 
here.
-               The rest of interaction will be deleted at the end of action. */
-       if(!overlap && hasInter){ if(!I->isReal) interactions->erase(id1,id2); 
return; }
+       if(!overlap && hasInter){ if(!I->isReal && I->isNew) 
interactions->erase(id1,id2); return; }
        assert(false); // unreachable
 }
 
-void InsertionSortCollider::insertionSort(vector<Bound>& v, 
InteractionContainer* interactions, MetaBody* rb){
+void InsertionSortCollider::insertionSort(vector<Bound>& v, 
InteractionContainer* interactions, MetaBody* rb, bool doCollide){
        long size=v.size();
        for(long i=0; i<size; i++){
                Bound viInit=v[i]; long j=i-1; /* cache hasBB(); otherwise 1% 
overall performance hit */ bool viInitBB=viInit.hasBB();
                while(j>=0 && v[j]>viInit){
                        v[j+1]=v[j];
                        // no collisions without bounding boxes
-                       if(viInitBB && v[j].hasBB()) 
handleBoundInversion(viInit.id,v[j].id,interactions,rb);
+                       if(doCollide && viInitBB && v[j].hasBB()) 
handleBoundInversion(viInit.id,v[j].id,interactions,rb);
                        j--;
                }
                v[j+1]=viInit;
@@ -112,19 +110,32 @@
                }
 
        //timingDeltas->checkpoint("copy");
+
+       // process interactions that the constitutive law asked to be erased
+       FOREACH(const InteractionContainer::bodyIdPair& p, 
interactions->pendingErase){
+               // remove those that do not overlap spatially anymore
+               if(!spatialOverlap(p[0],p[1])){ interactions->erase(p[0],p[1]); 
LOG_TRACE("Deleted interaction #"<<p[0]<<"+#"<<p[1]); }
+               else
+               {
+                       const shared_ptr<Interaction>& 
I=interactions->find(p[0],p[1]);
+                       if(!I){ LOG_FATAL("Requested deletion of a non-existent 
interaction #"<<p[0]<<"+#"<<p[1]<<"?!"); throw; }
+                       I->reset();
+               }
+       }
+       interactions->pendingErase.clear();
        
 
        // sort
-               if(!doInitSort){
+               if(!doInitSort && !sortThenCollide){
                        /* each inversion in insertionSort calls 
handleBoundInversion, which in turns may add/remove interaction */
-                       insertionSort(XX,interactions,rb);
-                       insertionSort(YY,interactions,rb);
-                       insertionSort(ZZ,interactions,rb);
+                       insertionSort(XX,interactions,rb); 
insertionSort(YY,interactions,rb); insertionSort(ZZ,interactions,rb);
                }
                else {
-                       std::sort(XX.begin(),XX.end());
-                       std::sort(YY.begin(),YY.end());
-                       std::sort(ZZ.begin(),ZZ.end());
+                       if(doInitSort){
+                               std::sort(XX.begin(),XX.end()); 
std::sort(YY.begin(),YY.end()); std::sort(ZZ.begin(),ZZ.end());
+                       } else {
+                               insertionSort(XX,interactions,rb,false); 
insertionSort(YY,interactions,rb,false); 
insertionSort(ZZ,interactions,rb,false);
+                       }
                        // traverse the container along requested axis
                        assert(sortAxis==0 || sortAxis==1 || sortAxis==2);
                        vector<Bound>& V=(sortAxis==0?XX:(sortAxis==1?YY:ZZ));
@@ -137,17 +148,18 @@
                                // TRVAR3(i,iid,V[i].coord);
                                // go up until we meet the upper bound
                                for(size_t j=i+1; V[j].id!=iid; j++){
-                                       // skip bodies with smaller (arbitrary, 
could be greater as well) id,
-                                       // since they will detect us when their 
turn comes
                                        const body_id_t& jid=V[j].id;
-                                       if(jid<iid) { /* LOG_TRACE("Skip 
#"<<V[j].id<<(V[j].isMin()?"(min)":"(max)")<<" with "<<iid<<" (smaller id)"); 
*/ continue; }
+                                       /// FIXME: not sure why this doesn't 
work. If this condition is commented out, we have exact same interactions as 
from SpatialQuickSort. Otherwise some interactions are missing!
+                                       // skip bodies with smaller (arbitrary, 
could be greater as well) id, since they will detect us when their turn comes
+                                       // if(jid<iid) { /* LOG_TRACE("Skip 
#"<<V[j].id<<(V[j].isMin()?"(min)":"(max)")<<" with "<<iid<<" (smaller id)"); 
*/ continue; }
                                        /* abuse the same function here; since 
it does spatial overlap check first, it is OK to use it */
                                        
handleBoundInversion(iid,jid,interactions,rb);
                                }
                        }
                }
        //timingDeltas->checkpoint("sort&collide");
-
+       
+#if 0
        // garbage collection once in a while: for interactions that were still 
real when the bounding boxes separated
        // the collider would never get to see them again otherwise
        if(iter%1000==0){
@@ -159,4 +171,5 @@
                FOREACH(const bodyIdPair& p, toBeDeleted){ 
interactions->erase(p.first,p.second); LOG_TRACE("Deleted interaction 
#"<<p.first<<"+#"<<p.second); }
        }
        //timingDeltas->checkpoint("stale");
+#endif
 }

Modified: trunk/pkg/common/Engine/StandAloneEngine/InsertionSortCollider.hpp
===================================================================
--- trunk/pkg/common/Engine/StandAloneEngine/InsertionSortCollider.hpp  
2009-05-23 10:29:03 UTC (rev 1774)
+++ trunk/pkg/common/Engine/StandAloneEngine/InsertionSortCollider.hpp  
2009-05-24 09:19:26 UTC (rev 1775)
@@ -44,18 +44,20 @@
        /*! sorting routine; insertion sort is very fast for strongly 
pre-sorted lists, which is our case
            http://en.wikipedia.org/wiki/Insertion_sort has the algorithm and 
other details
        */
-       void insertionSort(std::vector<Bound>& 
v,InteractionContainer*,MetaBody*);
+       void insertionSort(std::vector<Bound>& 
v,InteractionContainer*,MetaBody*,bool doCollide=true);
        void 
handleBoundInversion(body_id_t,body_id_t,InteractionContainer*,MetaBody*);
        bool spatialOverlap(body_id_t,body_id_t);
 
        public:
        //! axis for the initial sort
        int sortAxis;
+       //! if true, separate sorting and colliding phase; MUCH slower, but 
processes all interactions at every step
+       bool sortThenCollide;
 
-       InsertionSortCollider(): sortAxis(0){ /* 
timingDeltas=shared_ptr<TimingDeltas>(new TimingDeltas);*/ }
+       InsertionSortCollider(): sortAxis(0), sortThenCollide(false){ /* 
timingDeltas=shared_ptr<TimingDeltas>(new TimingDeltas);*/ }
        virtual void action(MetaBody*);
        REGISTER_CLASS_AND_BASE(InsertionSortCollider,Collider);
-       REGISTER_ATTRIBUTES(Collider,(sortAxis));
+       REGISTER_ATTRIBUTES(Collider,(sortAxis)(sortThenCollide));
        DECLARE_LOGGER;
 };
 REGISTER_SERIALIZABLE(InsertionSortCollider);

Modified: trunk/pkg/common/Engine/StandAloneEngine/SpatialQuickSortCollider.cpp
===================================================================
--- trunk/pkg/common/Engine/StandAloneEngine/SpatialQuickSortCollider.cpp       
2009-05-23 10:29:03 UTC (rev 1774)
+++ trunk/pkg/common/Engine/StandAloneEngine/SpatialQuickSortCollider.cpp       
2009-05-24 09:19:26 UTC (rev 1775)
@@ -11,7 +11,8 @@
 #include <math.h>
 #include <algorithm>
 
-
+YADE_PLUGIN("SpatialQuickSortCollider");
+ 
 SpatialQuickSortCollider::SpatialQuickSortCollider() : Collider()
 {
        haveDistantTransient=false;
@@ -24,16 +25,19 @@
 
 void SpatialQuickSortCollider::action(MetaBody* ncb)
 {
-     shared_ptr<BodyContainer> bodies = ncb->bodies;
+       const shared_ptr<BodyContainer>& bodies = ncb->bodies;
 
+       // This collider traverses all interactions at every step, therefore 
interactions that were reset() will be deleted automatically as needed
+       ncb->interactions->pendingErase.clear();
+
        size_t nbElements=bodies->size();
-       if (nbElements!=rank.size())
-       {
-          size_t n = rank.size();
-          rank.resize(nbElements);
-          for (; n<nbElements; ++n)
-              rank[n] = shared_ptr<AABBBound>(new AABBBound);
-       }
+       if (nbElements!=rank.size())
+       {
+               size_t n = rank.size();
+               rank.resize(nbElements);
+               for (; n<nbElements; ++n)
+                       rank[n] = shared_ptr<AABBBound>(new AABBBound);
+       }
 
        Vector3r min,max;
        shared_ptr<Body> b;

Modified: trunk/pkg/dem/Engine/StandAloneEngine/ElasticContactLaw.cpp
===================================================================
--- trunk/pkg/dem/Engine/StandAloneEngine/ElasticContactLaw.cpp 2009-05-23 
10:29:03 UTC (rev 1774)
+++ trunk/pkg/dem/Engine/StandAloneEngine/ElasticContactLaw.cpp 2009-05-24 
09:19:26 UTC (rev 1775)
@@ -70,7 +70,7 @@
                        SpheresContactGeometry*    currentContactGeometry= 
static_cast<SpheresContactGeometry*>(ig.get());
                        ElasticContactInteraction* currentContactPhysics = 
static_cast<ElasticContactInteraction*>(ip.get());
                        // delete interaction where spheres don't touch
-                       if(currentContactGeometry->penetrationDepth<0){ 
contact->isReal=false; return; }
+                       if(currentContactGeometry->penetrationDepth<0){ 
ncb->interactions->requestErase(id1,id2); return; }
        
                        BodyMacroParameters* de1                                
= 
YADE_CAST<BodyMacroParameters*>(Body::byId(id1,ncb)->physicalParameters.get());
                        BodyMacroParameters* de2                                
= 
YADE_CAST<BodyMacroParameters*>(Body::byId(id2,ncb)->physicalParameters.get());
@@ -116,7 +116,7 @@
        Dem3DofGeom* geom=static_cast<Dem3DofGeom*>(ig.get());
        ElasticContactInteraction* 
phys=static_cast<ElasticContactInteraction*>(ip.get());
        Real displN=geom->displacementN();
-       if(displN>0){contact->isReal=false; return; }
+       
if(displN>0){rootBody->interactions->requestErase(contact->getId1(),contact->getId2());
 return; }
        phys->normalForce=phys->kn*displN*geom->normal;
        Real 
maxFsSq=phys->normalForce.SquaredLength()*pow(phys->tangensOfFrictionAngle,2);
        Vector3r trialFs=phys->ks*geom->displacementT();

Modified: trunk/scripts/test/insertion-sort-collider.py
===================================================================
--- trunk/scripts/test/insertion-sort-collider.py       2009-05-23 10:29:03 UTC 
(rev 1774)
+++ trunk/scripts/test/insertion-sort-collider.py       2009-05-24 09:19:26 UTC 
(rev 1775)
@@ -4,21 +4,22 @@
        BexResetter(),
        
BoundingVolumeMetaEngine([InteractingSphere2AABB(),InteractingBox2AABB(),InteractingFacet2AABB(),MetaInteractingGeometry2AABB()]),
        InsertionSortCollider(),
-       
InteractionDispatchers([ef2_Facet_Sphere_Dem3DofGeom()],[SimpleElasticRelationships()],[ef2_Dem3Dof_Elastic_ElasticLaw()],),
+       
#InteractionDispatchers([ef2_Facet_Sphere_Dem3DofGeom()],[SimpleElasticRelationships()],[ef2_Dem3Dof_Elastic_ElasticLaw()],),
+       
InteractionDispatchers([ef2_Facet_Sphere_Dem3DofGeom()],[BrefcomMakeContact(cohesiveThresholdIter=0)],[ef2_Spheres_Brefcom_BrefcomLaw()],),
        GravityEngine(gravity=[0,0,-10]),
        NewtonsDampedLaw(damping=0.01),
 ]
 
 O.bodies.append([
-       
utils.facet([[-1,-1,0],[1,-1,0],[0,1,0]],dynamic=False,color=[1,0,0],young=1e3),
-       utils.facet([[1,-1,0],[0,1,0,],[1,.5,.5]],dynamic=False,young=1e3)
+       
utils.facet([[-1,-1,0],[1,-1,0],[0,1,0]],dynamic=False,color=[1,0,0],young=1e3,physParamsClass='BrefcomPhysParams'),
+       
utils.facet([[1,-1,0],[0,1,0,],[1,.5,.5]],dynamic=False,young=1e3,physParamsClass='BrefcomPhysParams')
 ])
 import random
 if 1:
        for i in range(0,100):
-               
O.bodies.append(utils.sphere([random.gauss(0,1),random.gauss(0,1),random.uniform(1,2)],random.uniform(.02,.05),velocity=[random.gauss(0,.1),random.gauss(0,.1),random.gauss(0,.1)]))
+               
O.bodies.append(utils.sphere([random.gauss(0,1),random.gauss(0,1),random.uniform(1,2)],random.uniform(.02,.05),velocity=[random.gauss(0,.1),random.gauss(0,.1),random.gauss(0,.1)],,physParamsClass='BrefcomPhysParams'))
 else:
-       O.bodies.append(utils.sphere([0,0,.6],.5))
+       
O.bodies.append(utils.sphere([0,0,.6],.5,,physParamsClass='BrefcomPhysParams'))
 O.dt=1e-4
 O.saveTmp('init')
 import yade.log


_______________________________________________
Mailing list: https://launchpad.net/~yade-dev
Post to     : [email protected]
Unsubscribe : https://launchpad.net/~yade-dev
More help   : https://help.launchpad.net/ListHelp
_______________________________________________
yade-dev mailing list
[email protected]
https://lists.berlios.de/mailman/listinfo/yade-dev

Reply via email to