Author: eudoxos
Date: 2008-09-18 11:50:19 +0200 (Thu, 18 Sep 2008)
New Revision: 1517
Modified:
trunk/extra/Brefcom.cpp
trunk/extra/clump/Shop.cpp
trunk/extra/clump/Shop.hpp
trunk/lib/opengl/GLUtils.hpp
trunk/pkg/dem/DataClass/InteractionGeometry/SpheresContactGeometry.cpp
trunk/pkg/dem/DataClass/InteractionGeometry/SpheresContactGeometry.hpp
trunk/pkg/dem/Engine/EngineUnit/InteractingSphere2InteractingSphere4SpheresContactGeometry.cpp
trunk/pkg/dem/Engine/EngineUnit/InteractingSphere2InteractingSphere4SpheresContactGeometry.hpp
trunk/pkg/dem/RenderingEngine/GLDrawSpheresContactGeometry/GLDrawSpheresContactGeometry.cpp
Log:
1. Remove GL things from Shop, moved to lib/opengl/GLUtils.hpp
2. Cleanup of SpheresContactGeometry code, rename exactRot to hasShear
Modified: trunk/extra/Brefcom.cpp
===================================================================
--- trunk/extra/Brefcom.cpp 2008-09-17 21:48:17 UTC (rev 1516)
+++ trunk/extra/Brefcom.cpp 2008-09-18 09:50:19 UTC (rev 1517)
@@ -4,6 +4,7 @@
#include<yade/pkg-dem/BodyMacroParameters.hpp>
#include<yade/pkg-common/Sphere.hpp>
#include<yade/lib-QGLViewer/qglviewer.h>
+#include<yade/lib-opengl/GLUtils.hpp>
YADE_PLUGIN("BrefcomMakeContact","BrefcomContact","BrefcomLaw","GLDrawBrefcomContact","BrefcomDamageColorizer",
"BrefcomPhysParams", "BrefcomGlobalCharacteristics" /*
,"BrefcomStiffnessComputer"*/ );
@@ -314,9 +315,9 @@
min(1.,max(0.,BC->epsTrans/BC->epsCrackOnset)),
min(1.,max(0.,abs(BC->epsTrans)/BC->epsCrackOnset-1)));
- if(contactLine)
Shop::GLDrawLine(b1->physicalParameters->dispSe3.position,b2->physicalParameters->dispSe3.position,lineColor);
- if(dmgLabel){
Shop::GLDrawNum(BC->omega,0.5*(b1->physicalParameters->dispSe3.position+b2->physicalParameters->dispSe3.position),lineColor);
}
- else if(epsNLabel){
Shop::GLDrawNum(BC->epsN,0.5*(b1->physicalParameters->dispSe3.position+b2->physicalParameters->dispSe3.position),lineColor);
}
+ if(contactLine)
GLUtils::GLDrawLine(b1->physicalParameters->dispSe3.position,b2->physicalParameters->dispSe3.position,lineColor);
+ if(dmgLabel){
GLUtils::GLDrawNum(BC->omega,0.5*(b1->physicalParameters->dispSe3.position+b2->physicalParameters->dispSe3.position),lineColor);
}
+ else if(epsNLabel){
GLUtils::GLDrawNum(BC->epsN,0.5*(b1->physicalParameters->dispSe3.position+b2->physicalParameters->dispSe3.position),lineColor);
}
const Vector3r&
cp=static_pointer_cast<SpheresContactGeometry>(i->interactionGeometry)->contactPoint;
if(epsT){
@@ -325,16 +326,16 @@
Real scale=.5*BC->equilibriumDist;
Vector3r dirShear=BC->epsT; dirShear.Normalize();
if(epsTAxes){
-
Shop::GLDrawLine(cp-Vector3r(scale,0,0),cp+Vector3r(scale,0,0));
-
Shop::GLDrawLine(cp-Vector3r(0,scale,0),cp+Vector3r(0,scale,0));
-
Shop::GLDrawLine(cp-Vector3r(0,0,scale),cp+Vector3r(0,0,scale));
+
GLUtils::GLDrawLine(cp-Vector3r(scale,0,0),cp+Vector3r(scale,0,0));
+
GLUtils::GLDrawLine(cp-Vector3r(0,scale,0),cp+Vector3r(0,scale,0));
+
GLUtils::GLDrawLine(cp-Vector3r(0,0,scale),cp+Vector3r(0,0,scale));
}
-
Shop::GLDrawArrow(cp,cp+dirShear*relShear*scale,Vector3r(1.,0.,0.));
-
Shop::GLDrawLine(cp+dirShear*relShear*scale,cp+dirShear*scale,Vector3r(.3,.3,.3));
+
GLUtils::GLDrawArrow(cp,cp+dirShear*relShear*scale,Vector3r(1.,0.,0.));
+
GLUtils::GLDrawLine(cp+dirShear*relShear*scale,cp+dirShear*scale,Vector3r(.3,.3,.3));
- /* normal strain */
Shop::GLDrawArrow(cp,cp+BC->prevNormal*(BC->epsN/maxShear),Vector3r(0.,1.,0.));
+ /* normal strain */
GLUtils::GLDrawArrow(cp,cp+BC->prevNormal*(BC->epsN/maxShear),Vector3r(0.,1.,0.));
}
- //if(normal)
Shop::GLDrawArrow(cp,cp+BC->prevNormal*.5*BC->equilibriumDist,Vector3r(0.,1.,0.));
+ //if(normal)
GLUtils::GLDrawArrow(cp,cp+BC->prevNormal*.5*BC->equilibriumDist,Vector3r(0.,1.,0.));
}
/********************** BrefcomDamageColorizer ****************************/
Modified: trunk/extra/clump/Shop.cpp
===================================================================
--- trunk/extra/clump/Shop.cpp 2008-09-17 21:48:17 UTC (rev 1516)
+++ trunk/extra/clump/Shop.cpp 2008-09-18 09:50:19 UTC (rev 1517)
@@ -1071,56 +1071,3 @@
{0.370512,0.415453,0.055970,0.010546},
{-1.,-1.,-1.,-1. } /* sentinel: non-positive radius */
};
-#if 0
-Real Shop::ElasticWaveTimestepEstimate(shared_ptr<MetaBody> rootBody){
- Real minDt=std::numeric_limits<Real>::infinity();
- FOREACH(const shared_ptr<Body>& b, *rootBody->bodies){
- shared_ptr<Sphere>
sphere=dynamic_pointer_cast<Sphere>(b->geometricalModel);
- shared_ptr<ElasticBodyParameters>
elast=dynamic_pointer_cast<ElasticBodyParameters>(b->physicalParameters);
- if(!sphere || !elast) continue;
- Real
density=elast->mass/((4./3.)*pow(Mathr::PI*sphere->radius,3));
- Real thisDt=2*sphere->radius/sqrt(elast->young/density);
- minDt=std::min(minDt,thisDt);
- }
- return minDt;
-}
-#endif
-
-void Shop::GLDrawLine(Vector3r from, Vector3r to, Vector3r color){
- glEnable(GL_LIGHTING); glColor3v(color);
- glBegin(GL_LINES); glVertex3v(from); glVertex3v(to); glEnd();
-}
-
-void Shop::GLDrawArrow(Vector3r from, Vector3r to, Vector3r color){
- glEnable(GL_LIGHTING); glColor3v(color); qglviewer::Vec
a(from[0],from[1],from[2]),b(to[0],to[1],to[2]); QGLViewer::drawArrow(a,b);
-}
-
-void Shop::GLDrawNum(Real n, Vector3r pos, Vector3r color, unsigned precision){
- ostringstream oss; oss<<setprecision(precision)<< /* "w="<< */
(double)n;
- Shop::GLDrawText(oss.str(),pos,color);
-}
-
-void Shop::GLDrawText(std::string txt, Vector3r pos, Vector3r color){
- glPushMatrix();
- glTranslatev(pos);
- glColor3(color[0],color[1],color[2]);
- glRasterPos2i(0,0);
- for(unsigned int i=0;i<txt.length();i++)
- glutBitmapCharacter(GLUT_BITMAP_HELVETICA_12, txt[i]);
- glPopMatrix();
-
-}
-Real Shop::PWaveTimeStep(shared_ptr<MetaBody> _rb){
- shared_ptr<MetaBody> rb=_rb;
- if(!rb)rb=Omega::instance().getRootBody();
- Real dt=std::numeric_limits<Real>::infinity();
- FOREACH(const shared_ptr<Body>& b, *rb->bodies){
- if(!b->physicalParameters || !b->geometricalModel) continue;
- shared_ptr<ElasticBodyParameters>
ebp=dynamic_pointer_cast<ElasticBodyParameters>(b->physicalParameters);
- shared_ptr<Sphere>
s=dynamic_pointer_cast<Sphere>(b->geometricalModel);
- if(!ebp || !s) continue;
- Real density=ebp->mass/((4/3.)*Mathr::PI*pow(s->radius,3));
- dt=min(dt,s->radius/sqrt(ebp->young/density));
- }
- return dt;
-}
Modified: trunk/extra/clump/Shop.hpp
===================================================================
--- trunk/extra/clump/Shop.hpp 2008-09-17 21:48:17 UTC (rev 1516)
+++ trunk/extra/clump/Shop.hpp 2008-09-18 09:50:19 UTC (rev 1517)
@@ -79,11 +79,6 @@
//! Save spheres in the current simulation into a text file
static void saveSpheresToFile(string fileName);
- static void GLDrawLine(Vector3r from, Vector3r to, Vector3r
color=Vector3r(1,1,1));
- static void GLDrawArrow(Vector3r from, Vector3r to, Vector3r
color=Vector3r(1,1,1));
- static void GLDrawText(std::string text, Vector3r pos, Vector3r
color=Vector3r(1,1,1));
- static void GLDrawNum(Real n, Vector3r pos, Vector3r
color=Vector3r(1,1,1), unsigned precision=3);
-
/*! Cache for class indices for physical actions (body external
variables, Bex)
*
* It is necessary to populate the cache by calling
initCache(); then supported
Modified: trunk/lib/opengl/GLUtils.hpp
===================================================================
--- trunk/lib/opengl/GLUtils.hpp 2008-09-17 21:48:17 UTC (rev 1516)
+++ trunk/lib/opengl/GLUtils.hpp 2008-09-18 09:50:19 UTC (rev 1517)
@@ -12,20 +12,20 @@
struct GLUtils{
- static void GLDrawArrow(Vector3r from, Vector3r to, Vector3r color){
+ static void GLDrawArrow(const Vector3r& from, const Vector3r& to, const
Vector3r& color=Vector3r(1,1,1)){
glEnable(GL_LIGHTING); glColor3v(color); qglviewer::Vec
a(from[0],from[1],from[2]),b(to[0],to[1],to[2]); QGLViewer::drawArrow(a,b);
}
- static void GLDrawLine(Vector3r from, Vector3r to, Vector3r color){
+ static void GLDrawLine(const Vector3r& from, const Vector3r& to, const
Vector3r& color=Vector3r(1,1,1)){
glEnable(GL_LIGHTING); glColor3v(color);
glBegin(GL_LINES); glVertex3v(from); glVertex3v(to); glEnd();
}
- static void GLDrawNum(Real n, Vector3r pos, Vector3r color, unsigned
precision){
+ static void GLDrawNum(const Real& n, const Vector3r& pos, const
Vector3r& color=Vector3r(1,1,1), unsigned precision=3){
std::ostringstream oss; oss<<std::setprecision(precision)<< /*
"w="<< */ (double)n;
GLUtils::GLDrawText(oss.str(),pos,color);
}
- static void GLDrawText(std::string txt, Vector3r pos, Vector3r color){
+ static void GLDrawText(const std::string& txt, const Vector3r& pos,
const Vector3r& color=Vector3r(1,1,1)){
glPushMatrix();
glTranslatev(pos);
glColor3(color[0],color[1],color[2]);
Modified: trunk/pkg/dem/DataClass/InteractionGeometry/SpheresContactGeometry.cpp
===================================================================
--- trunk/pkg/dem/DataClass/InteractionGeometry/SpheresContactGeometry.cpp
2008-09-17 21:48:17 UTC (rev 1516)
+++ trunk/pkg/dem/DataClass/InteractionGeometry/SpheresContactGeometry.cpp
2008-09-18 09:50:19 UTC (rev 1517)
@@ -9,6 +9,7 @@
* (should be on the plane passing through origin and oriented with normal;
not checked!)
*/
void SpheresContactGeometry::setTgPlanePts(Vector3r p1new, Vector3r p2new){
+ assert(hasShear);
cp1rel=ori1.Conjugate()*rollPlanePtToSphere(p1new,d1,normal);
cp2rel=ori2.Conjugate()*rollPlanePtToSphere(p2new,d2,-normal);
}
@@ -19,6 +20,7 @@
* flip axis and the point would project on other side of the tangent plane
piece.
*/
void SpheresContactGeometry::relocateContactPoints(){
+ assert(hasShear);
Vector3r p1=contPtInTgPlane1(), p2=contPtInTgPlane2();
Vector3r midPt=.5*(p1+p2);
if((p1.SquaredLength()>4*d1 || p2.SquaredLength()>4*d2) &&
midPt.SquaredLength()>.5*min(d1,d2)){
@@ -29,15 +31,17 @@
}
/*! Perform slip of the projected contact points so that their distance
becomes equal (or remains smaller) than the given one.
- * TODO: not yet tested
+ * The slipped distance is returned.
*/
-void SpheresContactGeometry::slipToDistIfNeeded(Real dist){
+Real SpheresContactGeometry::slipToDisplacementTMax(Real displacementTMax){
+ assert(hasShear);
+ assert(displacementTMax>Mathr::ZERO_TOLERANCE);
Vector3r p1=contPtInTgPlane1(), p2=contPtInTgPlane2();
Real currDist=(p2-p1).Length();
- if(currDist<dist) return; // close enough, no slip needed
- assert(dist>Mathr::ZERO_TOLERANCE);
- Vector3r diff=.5*(currDist/dist-1)*(p2-p1);
+ if(currDist<displacementTMax) return 0; // close enough, no slip needed
+ Vector3r diff=.5*(currDist/displacementTMax-1)*(p2-p1);
setTgPlanePts(p1+diff,p2-diff);
+ return 2*diff.Length();
}
/*! Project point from sphere surface to tangent plane,
Modified: trunk/pkg/dem/DataClass/InteractionGeometry/SpheresContactGeometry.hpp
===================================================================
--- trunk/pkg/dem/DataClass/InteractionGeometry/SpheresContactGeometry.hpp
2008-09-17 21:48:17 UTC (rev 1516)
+++ trunk/pkg/dem/DataClass/InteractionGeometry/SpheresContactGeometry.hpp
2008-09-18 09:50:19 UTC (rev 1517)
@@ -30,7 +30,7 @@
contactPoint;
Real radius1,radius2,penetrationDepth;
- bool exactRot; // whether the exact rotation code is being used
-- following fields are needed for that
+ bool hasShear; // whether the exact rotation code is being used
-- following fields are needed for that
//! positions and orientations of both spheres -- must be
updated at every iteration
Vector3r pos1, pos2; Quaternionr ori1, ori2;
/*! Orientation of the contact point relative to each
sphere-local coordinates.
@@ -61,18 +61,21 @@
Vector3r displacementT(bool relocate=true){ return
contPtInTgPlane2()-contPtInTgPlane1();}
Vector3r epsT(){return displacementT()*(1./d0);}
- void slipToDistIfNeeded(Real slip);
+ Real slipToDisplacementTMax(Real displacementTMax);
+ //! slip to epsTMax if current epsT is greater; return the
relative slip magnitude
+ Real slipToEpsTMax(Real epsTMax){return
(1/d0)*slipToDisplacementTMax(epsTMax*d0);}
+
void relocateContactPoints();
-
SpheresContactGeometry():InteractionGeometry(),contactPoint(Vector3r::ZERO),radius1(0),radius2(0),exactRot(false){createIndex();}
+
SpheresContactGeometry():InteractionGeometry(),contactPoint(Vector3r::ZERO),radius1(0),radius2(0),hasShear(false){createIndex();}
virtual ~SpheresContactGeometry(){};
protected :
virtual void registerAttributes(){
REGISTER_ATTRIBUTE(radius1);
REGISTER_ATTRIBUTE(radius2);
REGISTER_ATTRIBUTE(contactPoint); // to allow access
from python
- // exactRot
- REGISTER_ATTRIBUTE(exactRot);
+ // hasShear
+ REGISTER_ATTRIBUTE(hasShear);
REGISTER_ATTRIBUTE(pos1);
REGISTER_ATTRIBUTE(pos2);
REGISTER_ATTRIBUTE(ori1);
Modified:
trunk/pkg/dem/Engine/EngineUnit/InteractingSphere2InteractingSphere4SpheresContactGeometry.cpp
===================================================================
---
trunk/pkg/dem/Engine/EngineUnit/InteractingSphere2InteractingSphere4SpheresContactGeometry.cpp
2008-09-17 21:48:17 UTC (rev 1516)
+++
trunk/pkg/dem/Engine/EngineUnit/InteractingSphere2InteractingSphere4SpheresContactGeometry.cpp
2008-09-18 09:50:19 UTC (rev 1517)
@@ -24,7 +24,7 @@
void
InteractingSphere2InteractingSphere4SpheresContactGeometry::registerAttributes()
{
REGISTER_ATTRIBUTE(interactionDetectionFactor);
- REGISTER_ATTRIBUTE(exactRot);
+ REGISTER_ATTRIBUTE(hasShear);
}
bool InteractingSphere2InteractingSphere4SpheresContactGeometry::go( const
shared_ptr<InteractingGeometry>& cm1,
@@ -51,8 +51,8 @@
scm->radius1 = s1->radius;
scm->radius2 = s2->radius;
if (!c->interactionGeometry) c->interactionGeometry = scm;
- if(exactRot){
- scm->exactRot=true;
+ if(hasShear){
+ scm->hasShear=true;
scm->pos1=se31.position; scm->pos2=se32.position;
scm->ori1=se31.orientation; scm->ori2=se32.orientation;
if(c->isNew){
@@ -64,8 +64,13 @@
scm->cp2rel.Align(Vector3r::UNIT_X,se32.orientation.Conjugate()*(-normal));
scm->cp1rel.Normalize();
scm->cp2rel.Normalize();
}
- // FIXME: temporary, testing only!
- if((Omega::instance().getCurrentIteration()%500)==0)
scm->relocateContactPoints();
+ // testing only
+ #if 0
+ if((Omega::instance().getCurrentIteration()%10000)==0)
scm->relocateContactPoints();
+ if((Omega::instance().getCurrentIteration()%10000)==0) {
+ Real slip=scm->slipToEpsNMax(1.); if(slip>0)
cerr<<"SLIP by Δε_N "<<slip<<" → ε_N="<<scm->epsN()<<endl;
+ }
+ #endif
}
return true;
} else return false;
Modified:
trunk/pkg/dem/Engine/EngineUnit/InteractingSphere2InteractingSphere4SpheresContactGeometry.hpp
===================================================================
---
trunk/pkg/dem/Engine/EngineUnit/InteractingSphere2InteractingSphere4SpheresContactGeometry.hpp
2008-09-17 21:48:17 UTC (rev 1516)
+++
trunk/pkg/dem/Engine/EngineUnit/InteractingSphere2InteractingSphere4SpheresContactGeometry.hpp
2008-09-18 09:50:19 UTC (rev 1517)
@@ -26,8 +26,8 @@
* @note This parameter is functionally coupled with
InteractinSphere2AABB::aabbEnlargeFactor,
* which will create larger bounding boxes and should be of the
same value. */
double interactionDetectionFactor;
- /*! Whether we create SpheresContactGeometry with data
necessary for exact rotation computation */
- bool exactRot;
+ /*! Whether we create SpheresContactGeometry with data
necessary for (exact) shear computation */
+ bool hasShear;
REGISTER_CLASS_NAME(InteractingSphere2InteractingSphere4SpheresContactGeometry);
REGISTER_BASE_CLASS_NAME(InteractionGeometryEngineUnit);
Modified:
trunk/pkg/dem/RenderingEngine/GLDrawSpheresContactGeometry/GLDrawSpheresContactGeometry.cpp
===================================================================
---
trunk/pkg/dem/RenderingEngine/GLDrawSpheresContactGeometry/GLDrawSpheresContactGeometry.cpp
2008-09-17 21:48:17 UTC (rev 1516)
+++
trunk/pkg/dem/RenderingEngine/GLDrawSpheresContactGeometry/GLDrawSpheresContactGeometry.cpp
2008-09-18 09:50:19 UTC (rev 1517)
@@ -74,7 +74,7 @@
#endif
}
- if(sc->exactRot){
+ if(sc->hasShear){
GLUtils::GLDrawLine(sc->pos1,sc->pos2,Vector3r(.5,.5,.5));
// sphere center to point on the sphere
//GLUtils::GLDrawText("[1]",sc->pos1,Vector3r(1,1,1));
GLUtils::GLDrawText("[2]",sc->pos2,Vector3r(1,1,1));
_______________________________________________
Mailing list: https://launchpad.net/~yade-dev
Post to : [email protected]
Unsubscribe : https://launchpad.net/~yade-dev
More help : https://help.launchpad.net/ListHelp
