← Back to team overview

yade-users team mailing list archive

Re: [Question #404070]: about adding an new InsertionSortCollider

 

Question #404070 on Yade changed:
https://answers.launchpad.net/yade/+question/404070

Haohui Xin posted a new comment:
Hello, Thanks.

The latest error is as following. It could cmakem make and make install.
But it appeared problems when I run it.

 Traceback (most recent call last):
  File "./yade-2016-10-31.git-1b22b1a", line 129, in <module>
    import yade
  File "/home/xhh/yade/install/lib/x86_64-linux-gnu/yade-2016-10-31.git-1b22b1a/py/yade/__init__.py", line 67, in <module>
    import system
RuntimeError: extension class wrapper for base class HeatCollider has not been created yet


Firstly, I define a new engine. 

#pragma once

#include<core/GlobalEngine.hpp>

class HeatEngine: public Engine {
	public :
		virtual ~HeatEngine() {};
	YADE_CLASS_BASE_DOC(HeatEngine,Engine,"Engine that used for heat simulation.");
};
REGISTER_SERIALIZABLE(HeatEngine);

Second, I define a new collider

#pragma once
#include<core/Bound.hpp>
#include<core/Interaction.hpp>
#include<pkg/common/HeatEngine.hpp>
#include<pkg/common/Dispatching.hpp>

class HeatCollider: public HeatEngine {
	public:
		static int avoidSelfInteractionMask;
		/*! Probe the Bound on a bodies presense. Returns list of body ids with which there is potential overlap. */
		virtual  vector<Body::id_t> probeBoundingVolume(const Bound&){throw;}

                static bool mayCollide(const Body*, const Body*);

                virtual void invalidatePersistentData(){}

		// ctor with functors for the integrated BoundDispatcher
		virtual void pyHandleCustomCtorArgs(boost::python::tuple& t, boost::python::dict& d);
		
		int get_avoidSelfInteractionMask(){return avoidSelfInteractionMask;}
		void set_avoidSelfInteractionMask(const int &v){avoidSelfInteractionMask = v;}
		
	YADE_CLASS_BASE_DOC_ATTRS_CTOR_PY(HeatCollider,HeatEngine,"Abstract class for finding spatial collisions between bodies. \n\n.. admonition:: Special constructor\n\n\tDerived HeatColliders (unless they override ``pyHandleCustomCtorArgs``) can be given list of :yref:`BoundFunctors <BoundFunctor>` which is used to initialize the internal :yref:`boundDispatcher <HeatCollider.boundDispatcher>` instance.",
		((shared_ptr<BoundDispatcher>,boundDispatcher,new BoundDispatcher,Attr::readonly,":yref:`BoundDispatcher` object that is used for creating :yref:`bounds <Body.bound>` on HeatCollider's request as necessary.")),
		/*ctor*/,
		.add_property("avoidSelfInteractionMask",&HeatCollider::get_avoidSelfInteractionMask,&HeatCollider::set_avoidSelfInteractionMask,"This mask is used to avoid the interactions inside a group of particles. To do so, the particles must have the same mask and this mask have to be compatible with this one.")
	);
	DECLARE_LOGGER;
};
REGISTER_SERIALIZABLE(HeatCollider);


#include<pkg/common/HeatCollider.hpp>

YADE_PLUGIN((HeatCollider));

int HeatCollider::avoidSelfInteractionMask = 0 ;

bool HeatCollider::mayCollide(const Body* b1, const Body* b2){
	return 
		// might be called with deleted bodies, i.e. NULL pointers
		(b1!=NULL && b2!=NULL) &&
		// only collide if at least one particle is standalone or they belong to different clumps
		(b1->isStandalone() || b2->isStandalone() || b1->clumpId!=b2->clumpId 
		#ifdef YADE_SPH
		// If SPH-mode enabled, we do not skip interactions between clump-members
		// to get the correct calculation of density b->rho
		|| true
		#endif
		) &&
		 // do not collide clumps, since they are just containers, never interact
		!b1->isClump() && !b2->isClump() &&
		// masks must have at least 1 bit in common
		b1->maskCompatible(b2->groupMask) &&
		// avoid contact between particles having the same mask compatible with the avoidSelfInteractionMask.
 		!( (b1->groupMask == b2->groupMask) && b1->maskCompatible(avoidSelfInteractionMask) )
	;
}

void HeatCollider::pyHandleCustomCtorArgs(boost::python::tuple& t, boost::python::dict& d){
	if(boost::python::len(t)==0) return; // nothing to do
	if(boost::python::len(t)!=1) throw invalid_argument(("HeatCollider optionally takes exactly one list of BoundFunctor's as non-keyword argument for constructor ("+boost::lexical_cast<string>(boost::python::len(t))+" non-keyword ards given instead)").c_str());
	typedef std::vector<shared_ptr<BoundFunctor> > vecBound;
	vecBound vb=boost::python::extract<vecBound>(t[0])();
	FOREACH(shared_ptr<BoundFunctor> bf, vb) this->boundDispatcher->add(bf);
	t=boost::python::tuple(); // empty the args
}


define an new sort methods

/*************************************************************************
*  Copyright (C) 2008 by Sergei Dorofeenko				 *
*  sega@xxxxxxxxxxxxxxxx                                                 *
*                                                                        *
*  This program is free software; it is licensed under the terms of the  *
*  GNU General Public License v2 or later. See file LICENSE for details. *
*************************************************************************/
#pragma once
#include <pkg/common/HeatCollider.hpp>
#include <core/InteractionContainer.hpp>
#include <vector>

class HeatSpatialQuickSortCollider : public HeatCollider {
    protected:

	struct AABBBound {
	    Vector3r min,max;
	    int id;
	};

	class xBoundComparator {
	    public:
	      bool operator() (shared_ptr<AABBBound> b1, shared_ptr<AABBBound> b2)
	      {
		 return b1->min[0] < b2->min[0];
	      }
	};

        vector<shared_ptr<AABBBound> > rank;

   public:
		virtual void action();
	YADE_CLASS_BASE_DOC(HeatSpatialQuickSortCollider,HeatCollider,"Collider using quicksort along axes at each step, using :yref:`Aabb` bounds. \n\n Its performance is lower than that of :yref:`InsertionSortCollider` (see `Colliders' performance <https://yade-dem.org/index.php/Colliders_performace>`_), but the algorithm is simple enought to make it good for checking other collider's correctness.");
};
REGISTER_SERIALIZABLE(HeatSpatialQuickSortCollider);

/*************************************************************************
*  Copyright (C) 2008 by Sergei Dorofeenko				 *
*  sega@xxxxxxxxxxxxxxxx                                                 *
*                                                                        *
*  This program is free software; it is licensed under the terms of the  *
*  GNU General Public License v2 or later. See file LICENSE for details. *
*************************************************************************/
#include "HeatSpatialQuickSortCollider.hpp"
#include <core/Scene.hpp>
#include <core/BodyContainer.hpp>
#include <cmath>
#include <algorithm>

YADE_PLUGIN((HeatSpatialQuickSortCollider));

void HeatSpatialQuickSortCollider::action()
{
	// update bounds
	boundDispatcher->scene=scene; boundDispatcher->action();

        const shared_ptr<BodyContainer>& bodies = scene->bodies;

	// This collider traverses all interactions at every step, therefore all interactions
	// that were requested for erase might be erased here and will be recreated if necessary.
	scene->interactions->eraseNonReal();

	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);
	}

	Vector3r min,max;
	int i=0;
	FOREACH(const shared_ptr<Body>& b, *bodies){
		if(!b->bound) continue;
		min = b->bound->min;
		max = b->bound->max;
		rank[i]->id = b->getId();
		rank[i]->min = min;
		rank[i]->max = max;
		i++;
	}
	
	const shared_ptr<InteractionContainer>& interactions=scene->interactions;
	scene->interactions->iterColliderLastRun=scene->iter;

        sort(rank.begin(), rank.end(), xBoundComparator()); // sorting
along X

	int id,id2; size_t j;
	shared_ptr<Interaction> interaction;
	for(int i=0,e=nbElements-1; i<e; ++i)
	{
		id  = rank[i]->id;
		min = rank[i]->min;
		max = rank[i]->max;
		j=i;
		while(++j<nbElements)
		{
			if ( rank[j]->min[0] > max[0]) break;
			if ( rank[j]->min[1] < max[1]
			&& rank[j]->max[1] > min[1]
			&& rank[j]->min[2] < max[2]
			&& rank[j]->max[2] > min[2])
			{
				id2=rank[j]->id;
				if ( (interaction = interactions->find(Body::id_t(id),Body::id_t(id2))) == 0)
				{
					interaction = shared_ptr<Interaction>(new Interaction(id,id2) );
					interactions->insert(interaction);
				}
				interaction->iterLastSeen=scene->iter; 
			}
		}
	}
}



define an new loop

#pragma once
#include<pkg/common/HeatEngine.hpp>
#include <pkg/common/Callbacks.hpp>
#include <pkg/common/Dispatching.hpp>

#ifdef USE_TIMING_DELTAS
	#define TIMING_DELTAS_CHECKPOINT(cpt) timingDeltas->checkpoint(cpt)
	#define TIMING_DELTAS_START() timingDeltas->start()
#else
	#define TIMING_DELTAS_CHECKPOINT(cpt)
	#define TIMING_DELTAS_START()
#endif

class HeatInteractionLoop: public HeatEngine {
	bool alreadyWarnedNoCollider;
	using idPair = std::pair<Body::id_t, Body::id_t>;
	// store interactions that should be deleted after loop in action, not later
	#ifdef YADE_OPENMP
		std::vector<std::list<idPair> > eraseAfterLoopIds;
		void eraseAfterLoop(Body::id_t id1,Body::id_t id2){ eraseAfterLoopIds[omp_get_thread_num()].push_back(idPair(id1,id2)); }
	#else
		list<idPair> eraseAfterLoopIds;
		void eraseAfterLoop(Body::id_t id1,Body::id_t id2){ eraseAfterLoopIds.push_back(idPair(id1,id2)); }
	#endif
	public:
		virtual void pyHandleCustomCtorArgs(boost::python::tuple& t, boost::python::dict& d);
		virtual void action();
		YADE_CLASS_BASE_DOC_ATTRS_CTOR_PY(HeatInteractionLoop,HeatEngine,"Unified dispatcher for handling interaction loop at every step, for parallel performance reasons.\n\n.. admonition:: Special constructor\n\n\tConstructs from 3 lists of :yref:`Ig2<IGeomFunctor>`, :yref:`Ip2<IPhysFunctor>`, :yref:`Law<LawFunctor>` functors respectively; they will be passed to interal dispatchers, which you might retrieve.",
			((shared_ptr<IGeomDispatcher>,geomDispatcher,new IGeomDispatcher,Attr::readonly,":yref:`IGeomDispatcher` object that is used for dispatch."))
			((shared_ptr<IPhysDispatcher>,physDispatcher,new IPhysDispatcher,Attr::readonly,":yref:`IPhysDispatcher` object used for dispatch."))
			((shared_ptr<LawDispatcher>,lawDispatcher,new LawDispatcher,Attr::readonly,":yref:`LawDispatcher` object used for dispatch."))
			((vector<shared_ptr<IntrCallback> >,callbacks,,,":yref:`Callbacks<IntrCallback>` which will be called for every :yref:`Interaction`, if activated."))
			((bool, eraseIntsInLoop, false,,"Defines if the interaction loop should erase pending interactions, else the collider takes care of that alone (depends on what collider is used)."))
			,
			/*ctor*/ alreadyWarnedNoCollider=false;
				#ifdef YADE_OPENMP
					eraseAfterLoopIds.resize(omp_get_max_threads());
				#endif
			,
			/*py*/
		);
		DECLARE_LOGGER;
};
REGISTER_SERIALIZABLE(HeatInteractionLoop);

#include"HeatInteractionLoop.hpp"

YADE_PLUGIN((HeatInteractionLoop));
CREATE_LOGGER(HeatInteractionLoop);

void HeatInteractionLoop::pyHandleCustomCtorArgs(boost::python::tuple& t, boost::python::dict& d){
	if(boost::python::len(t)==0) return; // nothing to do
	if(boost::python::len(t)!=3) throw invalid_argument("Exactly 3 lists of functors must be given");
	// parse custom arguments (3 lists) and do in-place modification of args
	using vecGeom = std::vector<shared_ptr<IGeomFunctor> >;
	using vecPhys = std::vector<shared_ptr<IPhysFunctor> >;
	using vecLaw = std::vector<shared_ptr<LawFunctor> >;
	vecGeom vg=boost::python::extract<vecGeom>(t[0])();
	vecPhys vp=boost::python::extract<vecPhys>(t[1])();
	vecLaw vl=boost::python::extract<vecLaw>(t[2])();
	for(const auto gf : vg) this->geomDispatcher->add(gf);
	for(const auto pf : vp) this->physDispatcher->add(pf);
	for(const auto cf : vl) this->lawDispatcher->add(cf);
	t=boost::python::tuple(); // empty the args; not sure if this is OK, as there is some refcounting in raw_constructor code
}

void HeatInteractionLoop::action(){
	// update Scene* of the dispatchers
	lawDispatcher->scene=scene;
	physDispatcher->scene=scene;
	geomDispatcher->scene=scene;
	
	// ask dispatchers to update Scene* of their functors
	geomDispatcher->updateScenePtr();
	physDispatcher->updateScenePtr();
	lawDispatcher->updateScenePtr();

	/*
		initialize callbacks; they return pointer (used only in this timestep) to the function to be called
		returning NULL deactivates the callback in this timestep
	*/
	// pair of callback object and pointer to the function to be called
	vector<IntrCallback::FuncPtr> callbackPtrs;
	for (const auto cb : callbacks){
		cb->scene=scene;
		callbackPtrs.push_back(cb->stepInit());
	}
	assert(callbackPtrs.size()==callbacks.size());
	const size_t callbacksSize=callbacks.size();

	// cache transformed cell size
	Matrix3r cellHsize = Matrix3r::Zero();
	if(scene->isPeriodic) {
		cellHsize=scene->cell->hSize;
	}

	// force removal of interactions that were not encountered by the collider
	// (only for some kinds of colliders; see comment for InteractionContainer::iterColliderLastRun)
	const bool removeUnseenIntrs=(scene->interactions->iterColliderLastRun>=0 && scene->interactions->iterColliderLastRun==scene->iter);

	#ifdef YADE_OPENMP
	const long size=scene->interactions->size();
	#pragma omp parallel for schedule(guided) num_threads(ompThreads>0 ? min(ompThreads,omp_get_max_threads()) : omp_get_max_threads())
	for(long i=0; i<size; i++){
		const shared_ptr<Interaction>& I=(*scene->interactions)[i];
	#else
	for (const auto & I : *scene->interactions){
	#endif
		if(removeUnseenIntrs && !I->isReal() && I->iterLastSeen<scene->iter) {
			eraseAfterLoop(I->getId1(),I->getId2());
			continue;
		}

		const shared_ptr<Body>& b1_=Body::byId(I->getId1(),scene);
		const shared_ptr<Body>& b2_=Body::byId(I->getId2(),scene);

		if(!b1_ || !b2_){
			LOG_DEBUG("Body #"<<(b1_?I->getId2():I->getId1())<<" vanished, erasing intr #"<<I->getId1()<<"+#"<<I->getId2()<<"!");
			scene->interactions->requestErase(I);
			continue;
		}
    
    // Skip interaction with clumps
    if (b1_->isClump() || b2_->isClump()) { continue; }
		
		// we know there is no geometry functor already, take the short path
		if(!I->functorCache.geomExists) { assert(!I->isReal()); continue; }
		
		// no interaction geometry for either of bodies; no interaction possible
		if(!b1_->shape || !b2_->shape) { assert(!I->isReal()); continue; }

		bool swap=false;
		// IGeomDispatcher
		if(!I->functorCache.geom){
			I->functorCache.geom=geomDispatcher->getFunctor2D(b1_->shape,b2_->shape,swap);
			// returns NULL ptr if no functor exists; remember that and shortcut
			if(!I->functorCache.geom) {
				I->functorCache.geomExists=false;
				continue;
			}
		}
		// arguments for the geom functor are in the reverse order (dispatcher would normally call goReverse).
		// we don't remember the fact that is reverse, so we swap bodies within the interaction
		// and can call go in all cases
		if(swap){I->swapOrder();}
		// body pointers must be updated, in case we swapped
		const shared_ptr<Body>& b1=swap?b2_:b1_;
		const shared_ptr<Body>& b2=swap?b1_:b2_;

		assert(I->functorCache.geom);
		
		bool wasReal=I->isReal();
		bool geomCreated;
		if(!scene->isPeriodic){
			geomCreated=I->functorCache.geom->go(b1->shape,b2->shape, *b1->state, *b2->state, Vector3r::Zero(), /*force*/false, I);
		} else {
			// handle periodicity
			Vector3r shift2=cellHsize*I->cellDist.cast<Real>();
			// in sheared cell, apply shear on the mutual position as well
			geomCreated=I->functorCache.geom->go(b1->shape,b2->shape,*b1->state,*b2->state,shift2,/*force*/false,I);
		}
		if(!geomCreated){
			if(wasReal) LOG_WARN("IGeomFunctor returned false on existing interaction!");
			if(wasReal) scene->interactions->requestErase(I); // fully created interaction without geometry is reset and perhaps erased in the next step
			continue; // in any case don't care about this one anymore
		}
		
		// IPhysDispatcher
		if(!I->functorCache.phys){
			I->functorCache.phys=physDispatcher->getFunctor2D(b1->material,b2->material,swap);
			assert(!swap); // InteractionPhysicsEngineUnits are symmetric
		}
		
		if(!I->functorCache.phys){
			throw std::runtime_error("Undefined or ambiguous IPhys dispatch for types "+b1->material->getClassName()+" and "+b2->material->getClassName()+".");
		}
		I->functorCache.phys->go(b1->material,b2->material,I);
		assert(I->phys);

                if(!wasReal) I->iterMadeReal=scene->iter; // mark the
interaction as created right now

		// LawDispatcher
		// populating constLaw cache must be done after geom and physics dispatchers have been called, since otherwise the interaction
		// would not have geom and phys yet.
		if(!I->functorCache.constLaw){
			I->functorCache.constLaw=lawDispatcher->getFunctor2D(I->geom,I->phys,swap);
			if(!I->functorCache.constLaw){
				LOG_FATAL("None of given Law2 functors can handle interaction #"<<I->getId1()<<"+"<<I->getId2()<<", types geom:"<<I->geom->getClassName()<<"="<<I->geom->getClassIndex()<<" and phys:"<<I->phys->getClassName()<<"="<<I->phys->getClassIndex()<<" (LawDispatcher::getFunctor2D returned empty functor)");
				exit(1);
			}
			assert(!swap); // reverse call would make no sense, as the arguments are of different types
		}
		assert(I->functorCache.constLaw);
		
		//If the functor return false, the interaction is reset
		if (!I->functorCache.constLaw->go(I->geom,I->phys,I.get())) scene->interactions->requestErase(I);

		// process callbacks for this interaction
		if(!I->isReal()) continue; // it is possible that Law2_ functor called requestErase, hence this check
		for(size_t i=0; i<callbacksSize; i++){
			if(callbackPtrs[i]!=NULL) (*(callbackPtrs[i]))(callbacks[i].get(),I.get());
		}
	}
}

-- 
You received this question notification because your team yade-users is
an answer contact for Yade.