yade-dev team mailing list archive
-
yade-dev team
-
Mailing list archive
-
Message #02304
[Branch ~yade-dev/yade/trunk] Rev 1814: 1. Fix uninitialized memebrs in Material (thans Anton)
------------------------------------------------------------
revno: 1814
committer: Václav Šmilauer <eudoxos@xxxxxxxx>
branch nick: trunk
timestamp: Mon 2009-11-23 15:56:10 +0100
message:
1. Fix uninitialized memebrs in Material (thans Anton)
2. Add equivalent of dispSe3 to OpenGLRenderingEngine directly
3. Fix pack.py
modified:
core/Material.hpp
core/MetaBody.cpp
pkg/common/RenderingEngine/OpenGLRenderingEngine.cpp
pkg/common/RenderingEngine/OpenGLRenderingEngine.hpp
pkg/dem/Engine/StandAloneEngine/VTKRecorder.cpp
pkg/dem/meta/ConcretePM.cpp
pkg/dem/meta/ConcretePM.hpp
py/pack.py
--
lp:yade
https://code.launchpad.net/~yade-dev/yade/trunk
Your team Yade developers is subscribed to branch lp:yade.
To unsubscribe from this branch go to https://code.launchpad.net/~yade-dev/yade/trunk/+edit-subscription.
=== modified file 'core/Material.hpp'
--- core/Material.hpp 2009-11-23 12:32:14 +0000
+++ core/Material.hpp 2009-11-23 14:56:10 +0000
@@ -12,7 +12,7 @@
*/
class Material: public Serializable, public Indexable{
public:
- Material(){ createIndex(); }
+ Material(): id(-1), density(-1){ createIndex(); }
~Material();
//! global id of the material; if >= 0, the material is shared and can be found under this index in MetaBody::materials
//! (necessary since yade::serialization doesn't track shared pointers)
=== modified file 'core/MetaBody.cpp'
--- core/MetaBody.cpp 2009-11-23 12:32:14 +0000
+++ core/MetaBody.cpp 2009-11-23 14:56:10 +0000
@@ -70,9 +70,10 @@
void MetaBody::postProcessAttributes(bool deserializing){
/* since yade::serialization doesn't properly handle shared pointers, iterate over all bodies and make materials shared again, if id>=0 */
+ int numMaterials=materials.size();
FOREACH(const shared_ptr<Body>& b, *bodies){
if(b->material->id<0) continue; // not a shared material
- assert((size_t)b->material->id < materials.size());
+ assert(b->material->id < numMaterials);
b->material=materials[b->material->id];
}
}
=== modified file 'pkg/common/RenderingEngine/OpenGLRenderingEngine.cpp'
--- pkg/common/RenderingEngine/OpenGLRenderingEngine.cpp 2009-11-23 12:32:14 +0000
+++ pkg/common/RenderingEngine/OpenGLRenderingEngine.cpp 2009-11-23 14:56:10 +0000
@@ -46,7 +46,6 @@
scaleDisplacements=false; scaleRotations=false;
displacementScale=Vector3r(1,1,1); rotationScale=1;
- numBodiesWhenRefSe3LastSet=0;
for(int i=0; i<clipPlaneNum; i++){clipPlaneSe3.push_back(Se3r(Vector3r::ZERO,Quaternionr::IDENTITY)); clipPlaneActive.push_back(false); clipPlaneNormals.push_back(Vector3r(1,0,0));}
@@ -75,6 +74,12 @@
glutInitDone=true;
}
+void OpenGLRenderingEngine::setBodiesRefSe3(const shared_ptr<MetaBody>& rootBody){
+ LOG_DEBUG("(re)initializing reference positions and orientations.");
+ FOREACH(const shared_ptr<Body>& b, *rootBody->bodies) if(b && b->state) { b->state->refPos=b->state->pos; b->state->refOri=b->state->ori; }
+}
+
+
void OpenGLRenderingEngine::initgl(){
LOG_INFO("(re)initializing GL for gldraw methods.\n");
BOOST_FOREACH(vector<string>& s,stateFunctorNames)
@@ -116,12 +121,6 @@
return false;
}
-void OpenGLRenderingEngine::setBodiesRefSe3(const shared_ptr<MetaBody>& rootBody){
- LOG_DEBUG("(re)initializing reference positions and orientations.");
- FOREACH(const shared_ptr<Body>& b, *rootBody->bodies) if(b && b->state) { b->state->refPos=b->state->pos; b->state->refOri=b->state->ori; }
- numBodiesWhenRefSe3LastSet=rootBody->bodies->size();
- numIterWhenRefSe3LastSet=Omega::instance().getCurrentIteration();
-}
/* mostly copied from PeriodicInsertionSortCollider
FIXME: common implementation somewhere */
@@ -134,27 +133,29 @@
return Vector3r(wrapCell(pt[0],rb->cellMin[0],rb->cellMax[0]),wrapCell(pt[1],rb->cellMin[1],rb->cellMax[1]),wrapCell(pt[2],rb->cellMin[2],rb->cellMax[2]));
}
-void OpenGLRenderingEngine::setBodiesDispSe3(const shared_ptr<MetaBody>& rootBody){
- #ifdef YADE_PHYSPAR
- FOREACH(const shared_ptr<Body>& b, *rootBody->bodies){
- // FIXME: dispSe3 must be put somewhere; probably not in Body::state directly, however
- if(!b || !b->physicalParameters) continue;
- const Se3r& se3=b->physicalParameters->se3; const Se3r& refSe3=b->physicalParameters->refSe3; Se3r& dispSe3=b->physicalParameters->dispSe3;
- Vector3r posCell=wrapCellPt(se3.position,rootBody.get());
- b->physicalParameters->isDisplayed=!pointClipped(posCell);
- // if no scaling, return quickly
- if(!(scaleDisplacements||scaleRotations||rootBody->isPeriodic)){ b->physicalParameters->dispSe3=b->physicalParameters->se3; continue; }
- // apply scaling
- dispSe3.position=(scaleDisplacements ? diagMult(displacementScale,se3.position-refSe3.position)+wrapCellPt(refSe3.position,rootBody.get()) : posCell );
- if(scaleRotations){
- Quaternionr relRot=refSe3.orientation.Conjugate()*se3.orientation;
- Vector3r axis; Real angle; relRot.ToAxisAngle(axis,angle);
- angle*=rotationScale;
- dispSe3.orientation=refSe3.orientation*Quaternionr(axis,angle);
- } else {dispSe3.orientation=se3.orientation;}
+void OpenGLRenderingEngine::setBodiesDispInfo(const shared_ptr<MetaBody>& rootBody){
+ if(rootBody->bodies->size()!=bodyDisp.size()) bodyDisp.resize(rootBody->bodies->size());
+ FOREACH(const shared_ptr<Body>& b, *rootBody->bodies){
+ if(!b || !b->state) continue;
+ size_t id=b->getId();
+ const Vector3r& pos=b->state->pos; const Vector3r& refPos=b->state->refPos;
+ const Quaternionr& ori=b->state->ori; const Quaternionr& refOri=b->state->refOri;
+ Vector3r posCell=wrapCellPt(pos,rootBody.get());
+ bodyDisp[id].isDisplayed=!pointClipped(posCell);
+ // if no scaling and no periodic, return quickly
+ if(!(scaleDisplacements||scaleRotations||rootBody->isPeriodic)){ bodyDisp[id].pos=pos; bodyDisp[id].ori=ori; continue; }
+ // apply scaling
+ bodyDisp[id].pos=(scaleDisplacements ? diagMult(displacementScale,pos-refPos)+wrapCellPt(refPos,rootBody.get()) : posCell );
+ if(!scaleRotations) bodyDisp[id].ori=ori;
+ else{
+ Quaternionr relRot=refOri.Conjugate()*ori;
+ Vector3r axis; Real angle; relRot.ToAxisAngle(axis,angle);
+ angle*=rotationScale;
+ bodyDisp[id].ori=refOri*Quaternionr(axis,angle);
}
- #endif
+ }
}
+
// draw periodic cell, if active
void OpenGLRenderingEngine::drawPeriodicCell(MetaBody* rootBody){
if(!rootBody->isPeriodic) return;
@@ -208,49 +209,12 @@
/* glBegin(GL_LINES);glVertex3v(clipPlaneSe3[i].position);glVertex3v(clipPlaneSe3[i].position+clipPlaneNormals[i]);glEnd(); */
}
- // if scaling positions or orientations, save reference values if not already done; if # of bodies changes, we have to reset those; remember iteration when this was done to detect (at least in most cases) that the simulation was reloaded
- if((scaleDisplacements || scaleRotations) && (rootBody->bodies->size()!=numBodiesWhenRefSe3LastSet||Omega::instance().getCurrentIteration()<=numIterWhenRefSe3LastSet)){setBodiesRefSe3(rootBody);}
-
// set displayed Se3 of body (scaling) and isDisplayed (clipping)
- setBodiesDispSe3(rootBody);
+ setBodiesDispInfo(rootBody);
drawPeriodicCell(rootBody.get());
if (Show_DOF || Show_ID) renderDOF_ID(rootBody);
- #ifdef YADE_SHAPE
- if (Body_geometrical_model){
- if (Cast_shadows){
- if (Fast_shadow_volume) renderSceneUsingFastShadowVolumes(rootBody,Light_position);
- else renderSceneUsingShadowVolumes(rootBody,Light_position);
- // draw transparent shadow volume
- if (Shadow_volumes) {
- glAlphaFunc(GL_GREATER, 1.0f/255.0f);
- glEnable(GL_ALPHA_TEST);
- glBlendFunc(GL_SRC_ALPHA, GL_ONE_MINUS_SRC_ALPHA);
- glEnable(GL_BLEND);
-
- glColor4f(0.86,0.058,0.9,0.3);
- glEnable(GL_LIGHTING);
-
- glEnable(GL_CULL_FACE);
-
- glCullFace(GL_FRONT);
- renderShadowVolumes(rootBody,Light_position);
-
- glCullFace(GL_BACK);
- renderShadowVolumes(rootBody,Light_position);
-
- glEnable(GL_DEPTH_TEST);
- glDisable(GL_BLEND);
- glDisable(GL_ALPHA_TEST);
- }
- } else{
- glEnable(GL_CULL_FACE);
- glEnable(GL_NORMALIZE);
- renderGeometricalModel(rootBody);
- }
- }
- #endif
#ifdef YADE_PHYSPAR
if (Body_state) renderState(rootBody);
#endif
@@ -264,155 +228,6 @@
if (Interaction_physics) renderInteractionPhysics(rootBody);
}
-
-#ifdef YADE_SHAPE
-void OpenGLRenderingEngine::renderSceneUsingShadowVolumes(const shared_ptr<MetaBody>& rootBody,Vector3r Light_position)
-{
- glColorMask(GL_FALSE, GL_FALSE, GL_FALSE, GL_FALSE);
- glEnable(GL_CULL_FACE);
- glCullFace(GL_BACK);
- glEnable(GL_NORMALIZE);
- renderGeometricalModel(rootBody);
-
- glClear(GL_STENCIL_BUFFER_BIT);
- glEnable(GL_STENCIL_TEST);
- glDepthMask(GL_FALSE);
- glStencilFunc(GL_ALWAYS, 0, 0);
-
- glStencilOp(GL_KEEP, GL_KEEP, GL_INCR);
- glCullFace(GL_BACK); /* increment using front face of shadow volume */
- renderShadowVolumes(rootBody,Light_position);
-
- glStencilOp(GL_KEEP, GL_KEEP, GL_DECR);
- glCullFace(GL_FRONT); /* increment using front face of shadow volume */
- renderShadowVolumes(rootBody,Light_position);
-
- glDepthMask(GL_TRUE);
- glColorMask(GL_TRUE, GL_TRUE, GL_TRUE, GL_TRUE);
- glCullFace(GL_BACK);
- glDepthFunc(GL_LEQUAL);
- glStencilOp(GL_KEEP, GL_KEEP, GL_KEEP);
-
- glStencilFunc(GL_EQUAL, 1, 1); /* draw shadowed part */
- glStencilFunc(GL_NOTEQUAL, 0, (GLuint)(-1));
- glDisable(GL_LIGHT0);
- glEnable(GL_NORMALIZE);
- renderGeometricalModel(rootBody);
-
- glStencilFunc(GL_EQUAL, 0, 1); /* draw lit part */
- glStencilFunc(GL_EQUAL, 0, (GLuint)(-1));
- glEnable(GL_LIGHT0);
- glEnable(GL_NORMALIZE);
- renderGeometricalModel(rootBody);
-
- glDepthFunc(GL_LESS);
- glDisable(GL_STENCIL_TEST);
-
-}
-void OpenGLRenderingEngine::renderSceneUsingFastShadowVolumes(const shared_ptr<MetaBody>& rootBody,Vector3r Light_position)
-{
- glEnable(GL_CULL_FACE);
- glEnable(GL_NORMALIZE);
- renderGeometricalModel(rootBody);
-
- glClear(GL_STENCIL_BUFFER_BIT);
- glEnable(GL_STENCIL_TEST);
- glDepthMask(GL_FALSE);
- glStencilFunc(GL_ALWAYS, 0, 0);
-
- glColorMask(GL_FALSE, GL_FALSE, GL_FALSE, GL_FALSE);
-
- glStencilOp(GL_KEEP, GL_KEEP, GL_INCR);
- glCullFace(GL_BACK); /* increment using front face of shadow volume */
- renderShadowVolumes(rootBody,Light_position);
-
- glStencilOp(GL_KEEP, GL_KEEP, GL_DECR);
- glCullFace(GL_FRONT); /* increment using front face of shadow volume */
- renderShadowVolumes(rootBody,Light_position);
-
- // Need to do that to remove shadow that are not on object but if glClear is 0
-/* glStencilOp(GL_KEEP, GL_KEEP, GL_INCR);
- glCullFace(GL_BACK);
- glDepthMask(GL_TRUE);
- Real clearDepthValue=0;
- glGetDoublev(GL_DEPTH_CLEAR_VALUE,&clearDepthValue);
- glDepthFunc(GL_EQUAL);
- glMatrixMode(GL_PROJECTION);
- glPushMatrix();
- glLoadIdentity();
- glOrtho(0, 1, 1, 0, 0.0, -clearDepthValue);
- glMatrixMode(GL_MODELVIEW);
- glPushMatrix();
- glLoadIdentity();
-
- glColor3f(1,0,0);
- glBegin(GL_QUADS);
- glVertex3f(0,0,clearDepthValue);
- glVertex3f(0,1,clearDepthValue);
- glVertex3f(1,1,clearDepthValue);
- glVertex3f(1,0,clearDepthValue);
- glEnd();
-
- glMatrixMode(GL_PROJECTION);
- glPopMatrix();
- glMatrixMode(GL_MODELVIEW);
- glPopMatrix(); */
-
- //glDepthMask(GL_TRUE);
- glColorMask(GL_TRUE, GL_TRUE, GL_TRUE, GL_TRUE);
- glCullFace(GL_BACK);
- glDepthFunc(GL_LEQUAL);
- glStencilOp(GL_KEEP, GL_KEEP, GL_KEEP);
-
-
- glStencilFunc(GL_NOTEQUAL, 0, (GLuint)(-1));
- glMatrixMode(GL_PROJECTION);
- glPushMatrix();
- glLoadIdentity();
- glOrtho(0, 1, 1, 0, 0.0, -1.0);
- glMatrixMode(GL_MODELVIEW);
- glPushMatrix();
- glLoadIdentity();
- glDisable(GL_CULL_FACE);
- glAlphaFunc(GL_GREATER, 1.0f/255.0f);
- glEnable(GL_ALPHA_TEST);
- glBlendFunc(GL_SRC_ALPHA, GL_ONE_MINUS_SRC_ALPHA);
- glEnable(GL_BLEND);
- glDisable(GL_LIGHTING);
- glColor4f(0,0,0,0.5);
- glBegin(GL_QUADS);
- glVertex2f(0,0);
- glVertex2f(0,1);
- glVertex2f(1,1);
- glVertex2f(1,0);
- glEnd();
- glEnable(GL_DEPTH_TEST);
- glDisable(GL_BLEND);
- glDisable(GL_ALPHA_TEST);
- glEnable(GL_CULL_FACE);
- glMatrixMode(GL_PROJECTION);
- glPopMatrix();
- glMatrixMode(GL_MODELVIEW);
- glPopMatrix();
-
- glDepthMask(GL_TRUE);
- glDepthFunc(GL_LESS);
- glDisable(GL_STENCIL_TEST);
-}
-
-
-void OpenGLRenderingEngine::renderShadowVolumes(const shared_ptr<MetaBody>& rootBody,Vector3r Light_position){
- if (!rootBody->geometricalModel){
- FOREACH(const shared_ptr<Body>& b, *rootBody->bodies){
- if(!b || !b->physicalParameters->isDisplayed) continue;
- if (b->geometricalModel->shadowCaster && ( (b->getGroupMask() & Draw_mask) || b->getGroupMask()==0 ))
- shadowVolumeDispatcher(b->geometricalModel,b->physicalParameters,Light_position);
- }
- }
- else shadowVolumeDispatcher(rootBody->geometricalModel,rootBody->physicalParameters,Light_position);
-}
-#endif
-
void OpenGLRenderingEngine::renderDOF_ID(const shared_ptr<MetaBody>& rootBody){
const GLfloat ambientColorSelected[4]={10.0,0.0,0.0,1.0};
const GLfloat ambientColorUnselected[4]={0.5,0.5,0.5,1.0};
@@ -505,7 +320,7 @@
FOREACH(const shared_ptr<Interaction>& I, *rootBody->interactions){
if(!I->interactionGeometry) continue;
const shared_ptr<Body>& b1=Body::byId(I->getId1(),rootBody), b2=Body::byId(I->getId2(),rootBody);
- //FIXME: if(!(b1->physicalParameters->isDisplayed||b2->physicalParameters->isDisplayed)) continue;
+ if(!(bodyDisp[I->getId1()].isDisplayed||bodyDisp[I->getId2()].isDisplayed)) continue;
glPushMatrix(); interactionGeometryDispatcher(I->interactionGeometry,I,b1,b2,Interaction_wire); glPopMatrix();
}
}
@@ -518,9 +333,8 @@
FOREACH(const shared_ptr<Interaction>& I, *rootBody->interactions){
if(!I->interactionPhysics) continue;
const shared_ptr<Body>& b1=Body::byId(I->getId1(),rootBody), b2=Body::byId(I->getId2(),rootBody);
- // FIXME:
- // if(!b1->physicalParameters||!b2->physicalParameters) continue;
- // if(!(b1->physicalParameters->isDisplayed||b2->physicalParameters->isDisplayed)) continue;
+ body_id_t id1=I->getId1(), id2=I->getId2();
+ if(!(bodyDisp[id1].isDisplayed||bodyDisp[id2].isDisplayed)) continue;
glPushMatrix(); interactionPhysicsDispatcher(I->interactionPhysics,I,b1,b2,Interaction_wire); glPopMatrix();
}
}
@@ -541,9 +355,8 @@
void OpenGLRenderingEngine::renderBoundingVolume(const shared_ptr<MetaBody>& rootBody){
FOREACH(const shared_ptr<Body>& b, *rootBody->bodies){
- if(!b) continue;
- // FIXME:
- // if(b->physicalParameters && !b->physicalParameters->isDisplayed) continue;
+ if(!b || !b->boundingVolume) continue;
+ if(!bodyDisp[b->getId()].isDisplayed) continue;
if(b->boundingVolume && ((b->getGroupMask()&Draw_mask) || b->getGroupMask()==0)){
glPushMatrix(); boundingVolumeDispatcher(b->boundingVolume); glPopMatrix();
}
@@ -565,14 +378,14 @@
const GLfloat ambientColorUnselected[4]={0.5,0.5,0.5,1.0};
FOREACH(const shared_ptr<Body>& b, *rootBody->bodies){
- if(!b) continue;
- // FIXME:
- //if(b->physicalParameters && !b->physicalParameters->isDisplayed) continue;
- const Se3r& se3=b->state->se3; // dispSe3
+ if(!b || !b->interactingGeometry) continue;
+ if(!bodyDisp[b->getId()].isDisplayed) continue;
+ Vector3r pos=bodyDisp[b->getId()].pos;
+ Quaternionr ori=bodyDisp[b->getId()].ori;
if(b->interactingGeometry && ((b->getGroupMask()&Draw_mask) || b->getGroupMask()==0)){
glPushMatrix();
- Real angle; Vector3r axis; se3.orientation.ToAxisAngle(axis,angle);
- glTranslatef(se3.position[0],se3.position[1],se3.position[2]);
+ Real angle; Vector3r axis; ori.ToAxisAngle(axis,angle);
+ glTranslatef(pos[0],pos[1],pos[2]);
glRotatef(angle*Mathr::RAD_TO_DEG,axis[0],axis[1],axis[2]);
if(current_selection==b->getId() || b->interactingGeometry->highlight){
glLightModelfv(GL_LIGHT_MODEL_AMBIENT,ambientColorSelected);
@@ -592,46 +405,72 @@
}
glPopMatrix();
if(current_selection==b->getId() || b->interactingGeometry->highlight){
- if(!b->boundingVolume || Body_wire || b->interactingGeometry->wire) GLUtils::GLDrawInt(b->getId(),se3.position);
+ if(!b->boundingVolume || Body_wire || b->interactingGeometry->wire) GLUtils::GLDrawInt(b->getId(),pos);
else {
// move the label towards the camera by the bounding box so that it is not hidden inside the body
- const Vector3r& mn=b->boundingVolume->min; const Vector3r& mx=b->boundingVolume->max; const Vector3r& p=se3.position;
+ const Vector3r& mn=b->boundingVolume->min; const Vector3r& mx=b->boundingVolume->max; const Vector3r& p=pos;
Vector3r ext(viewDirection[0]>0?p[0]-mn[0]:p[0]-mx[0],viewDirection[1]>0?p[1]-mn[1]:p[1]-mx[1],viewDirection[2]>0?p[2]-mn[2]:p[2]-mx[2]); // signed extents towards the camera
Vector3r dr=-1.01*(viewDirection.Dot(ext)*viewDirection);
- GLUtils::GLDrawInt(b->getId(),se3.position+dr,Vector3r::ONE);
+ GLUtils::GLDrawInt(b->getId(),pos+dr,Vector3r::ONE);
}
}
// if the body goes over the cell margin, draw it in all other positions with wire
// this could be done in a nicer way perhaps...
if(b->boundingVolume && rootBody->isPeriodic){
const Vector3r& cellMin(rootBody->cellMin); const Vector3r& cellMax(rootBody->cellMax); Vector3r cellSize=cellMax-cellMin;
+ // traverse all periodic cells around the body, to see if any of them touches
Vector3<int> bodyPer,minPer,maxPer;
for(int i=0; i<3; i++){
- bodyPer[i]=(int)floor((b->state->pos[i]-cellMin[i])/cellSize[i]);
- minPer[i]=(int)floor((b->boundingVolume->min[i]-cellMin[i])/cellSize[i]);
- maxPer[i]=(int)floor((b->boundingVolume->max[i]-cellMin[i])/cellSize[i]);
- //assert(bodyPer[i]<=maxPer[i]); assert(bodyPer[i]>=minPer[i]);
- }
- /* m is bitmask from 3 couples (0…64=2^6) */
- for(int m=0; m<64; m++){
- // any mask containing 00 couple is invalid
- if((!(m&1) && (!(m&2))) || (!(m&4) && (!(m&8))) || (!(m&16) && (!(m&32)))) continue;
- Vector3r pt(se3.position);
- bool isInside=false;
- for(int j=0; j<3; j++){
- if(m&(1<<(2*j))) {
- if(m&(1<<(2*j+1))) { if(bodyPer[j]>=maxPer[j]) {isInside=true; break; } pt[j]-=cellSize[j]; }
- else { if(bodyPer[j]<=minPer[j]){ isInside=true; break; } pt[j]+=cellSize[j]; }
- }
- }
- if(isInside) continue;
- if(pt==se3.position) continue; // shouldn't happen, but it happens :-(
- glPushMatrix();
- glTranslatev(pt);
- glRotatef(angle*Mathr::RAD_TO_DEG,axis[0],axis[1],axis[2]);
- interactingGeometryDispatcher(b->interactingGeometry,b->state,/*Body_wire*/ true, viewInfo);
- glPopMatrix();
- }
+ bodyPer[i]=(int)floor((pos[i]-cellMin[i])/cellSize[i]); // period of the center
+ minPer[i]=(int)floor((b->boundingVolume->min[i]-cellMin[i])/cellSize[i]); // period of the minimum point
+ maxPer[i]=(int)floor((b->boundingVolume->max[i]-cellMin[i])/cellSize[i]); // period of the maximum point
+ }
+ if (minPer[0]!=maxPer[0] || minPer[1]!=maxPer[1] || minPer[2]!=maxPer[2]){ // crosses cell boundary?
+ Vector3<int> i;
+ for(i[0]=-1; i[0]<=1; i[0]++) for(i[1]=-1;i[1]<=1; i[1]++) for(i[2]=-1; i[2]<=1; i[2]++){
+ if(i[0]==0 && i[1]==0 && i[2]==0) continue; // middle; already rendered above
+ if((i[0]==0 && minPer[0]==maxPer[0]) || (i[1]==0 && minPer[0]==maxPer[0]) || (i[2]==0 && minPer[0]==maxPer[0])) continue;
+ if((i[0]<=0 && i[1]<=0 && i[2]<=0 && minPer[0]==bodyPer[0]+i[0] && minPer[1]==bodyPer[1]+i[1] && minPer[2]==bodyPer[2]+i[2]) ||
+ (i[0]>=0 && i[1]>=0 && i[2]>=0 && maxPer[0]==bodyPer[0]+i[0] && maxPer[1]==bodyPer[1]+i[1] && maxPer[2]==bodyPer[2]+i[2])){
+ Vector3r pt(pos[0]+cellSize[0]*i[0],pos[1]+cellSize[1]*i[1],pos[2]+cellSize[2]*i[2]);
+ glPushMatrix();
+ glTranslatev(pt);
+ glRotatef(angle*Mathr::RAD_TO_DEG,axis[0],axis[1],axis[2]);
+ interactingGeometryDispatcher(b->interactingGeometry,b->state,/*Body_wire*/ true, viewInfo);
+ glPopMatrix();
+ }
+ }
+ }
+ #if 0
+ Vector3<int> bodyPer,minPer,maxPer;
+ for(int i=0; i<3; i++){
+ bodyPer[i]=(int)floor((b->state->pos[i]-cellMin[i])/cellSize[i]); // period of the center
+ minPer[i]=(int)floor((b->boundingVolume->min[i]-cellMin[i])/cellSize[i]); // period of the minimum point
+ maxPer[i]=(int)floor((b->boundingVolume->max[i]-cellMin[i])/cellSize[i]); // period of the maximum point
+ //assert(bodyPer[i]<=maxPer[i]); assert(bodyPer[i]>=minPer[i]);
+ }
+
+ /* m is bitmask from 3 couples (0…64=2^6) */
+ for(int m=0; m<64; m++){
+ // any mask containing 00 couple is invalid
+ if((!(m&1) && (!(m&2))) || (!(m&4) && (!(m&8))) || (!(m&16) && (!(m&32)))) continue;
+ Vector3r pt(se3.position);
+ bool isInside=false;
+ for(int j=0; j<3; j++){
+ if(m&(1<<(2*j))) {
+ if(m&(1<<(2*j+1))) { if(bodyPer[j]>=maxPer[j]) {isInside=true; break; } pt[j]-=cellSize[j]; }
+ else { if(bodyPer[j]<=minPer[j]){ isInside=true; break; } pt[j]+=cellSize[j]; }
+ }
+ }
+ if(!isInside) continue;
+ if(pt==se3.position) continue; // shouldn't happen, but it happens :-(
+ glPushMatrix();
+ glTranslatev(pt);
+ glRotatef(angle*Mathr::RAD_TO_DEG,axis[0],axis[1],axis[2]);
+ interactingGeometryDispatcher(b->interactingGeometry,b->state,/*Body_wire*/ true, viewInfo);
+ glPopMatrix();
+ }
+ #endif
}
}
}
=== modified file 'pkg/common/RenderingEngine/OpenGLRenderingEngine.hpp'
--- pkg/common/RenderingEngine/OpenGLRenderingEngine.hpp 2009-11-21 16:46:58 +0000
+++ pkg/common/RenderingEngine/OpenGLRenderingEngine.hpp 2009-11-23 14:56:10 +0000
@@ -31,9 +31,7 @@
bool pointClipped(const Vector3r& p);
vector<Vector3r> clipPlaneNormals;
- void setBodiesRefSe3(const shared_ptr<MetaBody>& rootBody);
- void setBodiesDispSe3(const shared_ptr<MetaBody>& rootBody);
- long numBodiesWhenRefSe3LastSet,numIterWhenRefSe3LastSet;
+ void setBodiesDispInfo(const shared_ptr<MetaBody>& rootBody);
static bool glutInitDone;
static size_t selectBodyLimit;
Vector3r viewDirection; // updated from GLViewer regularly
@@ -50,6 +48,16 @@
Vector3r wrapCellPt(const Vector3r& pt, MetaBody* rb);
void drawPeriodicCell(MetaBody*);
+ void setBodiesRefSe3(const shared_ptr<MetaBody>& rootBody);
+
+ struct BodyDisp{
+ Vector3r pos;
+ Quaternionr ori;
+ bool isDisplayed;
+ };
+
+ vector<BodyDisp> bodyDisp;
+
private :
DynLibDispatcher< InteractionGeometry , GLDrawInteractionGeometryFunctor, void , TYPELIST_5(const shared_ptr<InteractionGeometry>&, const shared_ptr<Interaction>& , const shared_ptr<Body>&, const shared_ptr<Body>&, bool) > interactionGeometryDispatcher;
=== modified file 'pkg/dem/Engine/StandAloneEngine/VTKRecorder.cpp'
--- pkg/dem/Engine/StandAloneEngine/VTKRecorder.cpp 2009-11-21 16:46:58 +0000
+++ pkg/dem/Engine/StandAloneEngine/VTKRecorder.cpp 2009-11-23 14:56:10 +0000
@@ -150,9 +150,9 @@
spheresVelocity->InsertNextTupleValue((float*)(&b->state->vel));
}
if (recActive[REC_CPM]) {
- cpmDamage->InsertNextValue(YADE_PTR_CAST<CpmMat>(b->material)->normDmg);
- const Vector3r& ss=YADE_PTR_CAST<CpmMat>(b->material)->sigma;
- const Vector3r& tt=YADE_PTR_CAST<CpmMat>(b->material)->tau;
+ cpmDamage->InsertNextValue(YADE_PTR_CAST<CpmState>(b->material)->normDmg);
+ const Vector3r& ss=YADE_PTR_CAST<CpmState>(b->material)->sigma;
+ const Vector3r& tt=YADE_PTR_CAST<CpmState>(b->material)->tau;
float s[3]={ss[0],ss[1],ss[2]};
float t[3]={tt[0],tt[1],tt[2]};
cpmSigma->InsertNextTupleValue(s);
=== modified file 'pkg/dem/meta/ConcretePM.cpp'
--- pkg/dem/meta/ConcretePM.cpp 2009-11-21 16:46:58 +0000
+++ pkg/dem/meta/ConcretePM.cpp 2009-11-23 14:56:10 +0000
@@ -4,7 +4,7 @@
#include<yade/pkg-dem/DemXDofGeom.hpp>
#include<yade/pkg-dem/Shop.hpp>
-YADE_PLUGIN((CpmMat)(Ip2_CpmMat_CpmMat_CpmPhys)(CpmPhys)(Law2_Dem3DofGeom_CpmPhys_Cpm)(CpmGlobalCharacteristics)
+YADE_PLUGIN((CpmState)(CpmMat)(Ip2_CpmMat_CpmMat_CpmPhys)(CpmPhys)(Law2_Dem3DofGeom_CpmPhys_Cpm)(CpmGlobalCharacteristics)
#ifdef YADE_OPENGL
(GLDrawCpmPhys)
#endif
@@ -21,11 +21,11 @@
Dem3DofGeom* contGeom=YADE_CAST<Dem3DofGeom*>(interaction->interactionGeometry.get());
assert(contGeom);
- const shared_ptr<CpmMat>& elast1=YADE_PTR_CAST<CpmMat>(pp1);
- const shared_ptr<CpmMat>& elast2=YADE_PTR_CAST<CpmMat>(pp2);
+ const shared_ptr<CpmMat>& mat1=YADE_PTR_CAST<CpmMat>(pp1);
+ const shared_ptr<CpmMat>& mat2=YADE_PTR_CAST<CpmMat>(pp2);
- Real E12=2*elast1->young*elast2->young/(elast1->young+elast2->young); // harmonic Young's modulus average
- //Real nu12=2*elast1->poisson*elast2->poisson/(elast1->poisson+elast2->poisson); // dtto for Poisson ratio
+ Real E12=2*mat1->young*mat2->young/(mat1->young+mat2->young); // harmonic Young's modulus average
+ //Real nu12=2*mat1->poisson*mat2->poisson/(mat1->poisson+mat2->poisson); // dtto for Poisson ratio
Real minRad=(contGeom->refR1<=0?contGeom->refR2:(contGeom->refR2<=0?contGeom->refR1:min(contGeom->refR1,contGeom->refR2)));
Real S12=Mathr::PI*pow(minRad,2); // "surface" of interaction
//Real E=(E12 /* was here for Kn: *S12/d0 */)*((1+alpha)/(beta*(1+nu12)+gamma*(1-alpha*nu12)));
@@ -37,7 +37,7 @@
contPhys->E=E12;
contPhys->G=E12*G_over_E;
- contPhys->tanFrictionAngle=tan(.5*(elast1->frictionAngle+elast2->frictionAngle));
+ contPhys->tanFrictionAngle=tan(.5*(mat1->frictionAngle+mat2->frictionAngle));
contPhys->undamagedCohesion=sigmaT;
contPhys->crossSection=S12;
contPhys->epsCrackOnset=epsCrackOnset;
@@ -199,11 +199,10 @@
if(epsN>0. && ((isCohesive && omega>omegaThreshold) || !isCohesive)){
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<CpmMat>& rbp1=YADE_PTR_CAST<CpmMat>(body1->material), rbp2=YADE_PTR_CAST<CpmMat>(body2->material);
+ const shared_ptr<CpmState>& st1=YADE_PTR_CAST<CpmState>(body1->state), st2=YADE_PTR_CAST<CpmState>(body2->state);
// nice article about openMP::critical vs. scoped locks: http://www.thinkingparallel.com/2006/08/21/scoped-locking-vs-critical-in-openmp-a-personal-shootout/
- // FIXME: those state things must be moved outside CpmMat into something like CpmState; then the updateMutex will also work
- { /* boost::mutex::scoped_lock lock(rbp1->updateMutex); */ rbp1->numBrokenCohesive+=1; rbp1->epsPlBroken+=epsPlSum; }
- { /* boost::mutex::scoped_lock lock(rbp2->updateMutex); */ rbp2->numBrokenCohesive+=1; rbp2->epsPlBroken+=epsPlSum; }
+ { boost::mutex::scoped_lock lock(st1->updateMutex); st1->numBrokenCohesive+=1; st1->epsPlBroken+=epsPlSum; }
+ { boost::mutex::scoped_lock lock(st2->updateMutex); st2->numBrokenCohesive+=1; st2->epsPlBroken+=epsPlSum; }
}
rootBody->interactions->requestErase(I->getId1(),I->getId2());
return;
@@ -377,19 +376,19 @@
FOREACH(shared_ptr<Body> B, *rootBody->bodies){
const body_id_t& id=B->getId();
// add damaged contacts that have already been deleted
- CpmMat* bpp=dynamic_cast<CpmMat*>(B->material.get());
- if(!bpp) continue;
- bpp->sigma=bodyStats[id].sigma;
- bpp->tau=bodyStats[id].tau;
- int cohLinksWhenever=bodyStats[id].nCohLinks+bpp->numBrokenCohesive;
+ CpmState* state=dynamic_cast<CpmState*>(B->material.get());
+ if(!state) continue;
+ state->sigma=bodyStats[id].sigma;
+ state->tau=bodyStats[id].tau;
+ int cohLinksWhenever=bodyStats[id].nCohLinks+state->numBrokenCohesive;
if(cohLinksWhenever>0){
- bpp->normDmg=(bodyStats[id].dmgSum+bpp->numBrokenCohesive)/cohLinksWhenever;
- bpp->normEpsPl=(bodyStats[id].epsPlSum+bpp->epsPlBroken)/cohLinksWhenever;
- if(bpp->normDmg>1){
- LOG_WARN("#"<<id<<" normDmg="<<bpp->normDmg<<" nCohLinks="<<bodyStats[id].nCohLinks<<", numBrokenCohesive="<<bpp->numBrokenCohesive<<", dmgSum="<<bodyStats[id].dmgSum<<", numAllCohLinks"<<cohLinksWhenever);
+ state->normDmg=(bodyStats[id].dmgSum+state->numBrokenCohesive)/cohLinksWhenever;
+ state->normEpsPl=(bodyStats[id].epsPlSum+state->epsPlBroken)/cohLinksWhenever;
+ if(state->normDmg>1){
+ LOG_WARN("#"<<id<<" normDmg="<<state->normDmg<<" nCohLinks="<<bodyStats[id].nCohLinks<<", numBrokenCohesive="<<state->numBrokenCohesive<<", dmgSum="<<bodyStats[id].dmgSum<<", numAllCohLinks"<<cohLinksWhenever);
}
}
- else { bpp->normDmg=0; bpp->normEpsPl=0;}
- B->interactingGeometry->diffuseColor=Vector3r(bpp->normDmg,1-bpp->normDmg,B->isDynamic?0:1);
+ else { state->normDmg=0; state->normEpsPl=0;}
+ B->interactingGeometry->diffuseColor=Vector3r(state->normDmg,1-state->normDmg,B->isDynamic?0:1);
}
}
=== modified file 'pkg/dem/meta/ConcretePM.hpp'
--- pkg/dem/meta/ConcretePM.hpp 2009-11-21 23:50:06 +0000
+++ pkg/dem/meta/ConcretePM.hpp 2009-11-23 14:56:10 +0000
@@ -53,8 +53,12 @@
#include<yade/pkg-common/NormalShearInteractions.hpp>
#include<yade/pkg-common/ConstitutiveLaw.hpp>
-/* This class holds information associated with each body */
-class CpmMat: public GranularMat {
+/* Cpm state information about each body.
+
+None of that is used for computation (at least not now), only for post-processing.
+
+*/
+class CpmState: public State {
public:
//! volumetric strain around this body
Real epsVolumetric;
@@ -70,8 +74,18 @@
Real normEpsPl;
//! stresses on the particle
Vector3r sigma,tau;
- CpmMat(): epsVolumetric(0.), numBrokenCohesive(0), numContacts(0), normDmg(0.), epsPlBroken(0.), normEpsPl(0.), sigma(Vector3r::ZERO), tau(Vector3r::ZERO) { createIndex(); density=4800; };
- REGISTER_ATTRIBUTES(GranularMat, (epsVolumetric) (numBrokenCohesive) (numContacts) (normDmg) (epsPlBroken) (normEpsPl) (sigma) (tau));
+
+ CpmState(): epsVolumetric(0.), numBrokenCohesive(0), numContacts(0), normDmg(0.), epsPlBroken(0.), normEpsPl(0.), sigma(Vector3r::ZERO), tau(Vector3r::ZERO) { };
+ REGISTER_ATTRIBUTES(State, (epsVolumetric) (numBrokenCohesive) (numContacts) (normDmg) (epsPlBroken) (normEpsPl) (sigma) (tau));
+ REGISTER_CLASS_AND_BASE(CpmState,State);
+};
+REGISTER_SERIALIZABLE(CpmState);
+
+/* This class holds information associated with each body */
+class CpmMat: public GranularMat {
+ public:
+ CpmMat() { createIndex(); density=4800; };
+ REGISTER_ATTRIBUTES(GranularMat,);
REGISTER_CLASS_AND_BASE(CpmMat,GranularMat);
REGISTER_CLASS_INDEX(CpmMat,GranularMat);
};
=== modified file 'py/pack.py'
--- py/pack.py 2009-11-20 08:25:37 +0000
+++ py/pack.py 2009-11-23 14:56:10 +0000
@@ -309,7 +309,8 @@
#print x1,y1,z1,radius,rRelFuzz
num=sp.makeCloud(O.periodicCell[0],O.periodicCell[1],radius,rRelFuzz,spheresInCell,True)
O.engines=[BexResetter(),BoundingVolumeMetaEngine([InteractingSphere2AABB()]),PeriodicInsertionSortCollider(),InteractionDispatchers([ef2_Sphere_Sphere_Dem3DofGeom()],[SimpleElasticRelationships()],[Law2_Dem3Dof_Elastic_Elastic()]),PeriIsoCompressor(charLen=radius/5.,stresses=[100e9,1e8],maxUnbalanced=1e-2,doneHook='O.pause();',globalUpdateInt=5,keepProportions=True),NewtonsDampedLaw(damping=.6)]
- for s in sp: O.bodies.append(utils.sphere(s[0],s[1],density=1000))
+ O.materials.append(GranularMat(young=30e9,frictionAngle=.5,poisson=.3,density=1e3))
+ for s in sp: O.bodies.append(utils.sphere(s[0],s[1]))
O.dt=utils.PWaveTimeStep()
O.run(); O.wait()
sp=SpherePack(); sp.fromSimulation()
