yade-dev team mailing list archive
-
yade-dev team
-
Mailing list archive
-
Message #06614
[Branch ~yade-dev/yade/trunk] Rev 2625: 1. Fix NewtonIntegrator for clumps, blockedDOFs &c; msg will be sent to yade-dev
------------------------------------------------------------
revno: 2625
committer: Václav Šmilauer <eu@xxxxxxxx>
branch nick: yade
timestamp: Tue 2010-12-21 23:50:34 +0100
message:
1. Fix NewtonIntegrator for clumps, blockedDOFs &c; msg will be sent to yade-dev
2. Make static attributes properly hyperlinked and documented in the ui
3. Fix yref in L3Geom
modified:
gui/qt4/Inspector.py
gui/qt4/SerializableEditor.py
pkg/dem/L3Geom.hpp
pkg/dem/NewtonIntegrator.cpp
--
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 'gui/qt4/Inspector.py'
--- gui/qt4/Inspector.py 2010-12-21 15:31:40 +0000
+++ gui/qt4/Inspector.py 2010-12-21 22:50:34 +0000
@@ -67,15 +67,17 @@
return -1
class BodyInspector(QWidget):
- def __init__(self,bodyId=0,parent=None,bodyLinkCallback=None,intrLinkCallback=None):
+ def __init__(self,bodyId=-1,parent=None,bodyLinkCallback=None,intrLinkCallback=None):
QWidget.__init__(self,parent)
+ v=yade.qt.views()
self.bodyId=bodyId
- self.idGlSync=-1
+ if bodyId<0 and len(v)>0 and v[0].selection>0: self.bodyId=v[0].selection
+ self.idGlSync=self.bodyId
self.bodyLinkCallback,self.intrLinkCallback=bodyLinkCallback,intrLinkCallback
self.bodyIdBox=QSpinBox(self)
self.bodyIdBox.setMinimum(0)
self.bodyIdBox.setMaximum(100000000)
- self.bodyIdBox.setValue(0)
+ self.bodyIdBox.setValue(self.bodyId)
self.intrWithCombo=QComboBox(self);
self.gotoBodyButton=QPushButton(u'â #',self)
self.gotoIntrButton=QPushButton(u'â #+#',self)
@@ -166,6 +168,7 @@
# no body shown yet, try to get the first one...
if self.bodyId<0 and len(O.bodies)>0:
try:
+ print 'SET ZERO'
b=O.bodies[0]; self.bodyIdBox.setValue(0)
except IndexError: pass
v=yade.qt.views()
=== modified file 'gui/qt4/SerializableEditor.py'
--- gui/qt4/SerializableEditor.py 2010-12-21 15:31:40 +0000
+++ gui/qt4/SerializableEditor.py 2010-12-21 22:50:34 +0000
@@ -15,14 +15,15 @@
seqSerializableShowType=True # show type headings in serializable sequences (takes vertical space, but makes the type hyperlinked)
-def makeWrapperHref(text,className,attr=None):
+def makeWrapperHref(text,className,attr=None,static=False):
"""Create clickable HTML hyperlink to a Yade class or its attribute.
:param className: name of the class to link to.
:param attr: attribute to link to. If given, must exist directly in given *className*; if not given or empty, link to the class itself is created and *attr* is ignored.
:return: HTML with the hyperref.
"""
- return '<a href="%s#yade.wrapper.%s%s">%s</a>'%(yade.qt.sphinxDocWrapperPage,className,(('.'+attr) if attr else ''),text)
+ if not static: return '<a href="%s#yade.wrapper.%s%s">%s</a>'%(yade.qt.sphinxDocWrapperPage,className,(('.'+attr) if attr else ''),text)
+ else: return '<a href="%s#ystaticattr-%s.%s">%s</a>'%(yade.qt.sphinxDocWrapperPage,className,attr,text)
def serializableHref(ser,attr=None,text=None):
"""Return HTML href to a *ser* optionally to the attribute *attr*.
@@ -48,7 +49,7 @@
else:
klass=ser.__class__
if not text: text=klass.__name__
- return makeWrapperHref(text,klass.__name__,attr)
+ return makeWrapperHref(text,klass.__name__,attr,static=(attr and getattr(klass,attr)==getattr(ser,attr)))
class AttrEditor():
"""Abstract base class handing some aspects common to all attribute editors.
@@ -346,8 +347,7 @@
t=None
doc=getattr(self.ser.__class__,attr).__doc__;
if '|yhidden|' in doc: continue
- if attr in self.ignoredAttrs or attr=='blockedDOFs': # HACK here
- continue
+ if attr in self.ignoredAttrs: continue
if isinstance(val,list):
t=self.getListTypeFromDocstring(attr)
if not t and len(val)==0: t=(val[0].__class__,) # 1-tuple is list of the contained type
@@ -364,10 +364,28 @@
doc=(getattr(self.ser.__class__,attr).__doc__ if attr else self.ser.__class__.__doc__)
if not doc: return ''
doc=re.sub(':y(attrtype|default|attrflags):`[^`]*`','',doc)
+ statAttr=re.compile('^.. ystaticattr::.*$',re.MULTILINE|re.DOTALL)
+ doc=re.sub(statAttr,'',doc) # static classes have their proper docs at the beginning, discard static memeber docs
+ # static: attribute of the type is the same object as attribute of the instance
+ # in that case, get docstring from the class documentation by parsing it
+ if attr and getattr(self.ser.__class__,attr)==getattr(self.ser,attr): doc=self.getStaticAttrDocstring(attr)
doc=re.sub(':yref:`([^`]*)`','\\1',doc)
import textwrap
wrapper=textwrap.TextWrapper(replace_whitespace=False)
return wrapper.fill(textwrap.dedent(doc))
+ def getStaticAttrDocstring(self,attr):
+ ret=''; c=self.ser.__class__
+ while hasattr(c,attr) and hasattr(c.__base__,attr): c=c.__base__
+ start='.. ystaticattr:: %s.%s('%(c.__name__,attr)
+ if start in c.__doc__:
+ ll=c.__doc__.split('\n')
+ for i in range(len(ll)):
+ if ll[i].startswith(start): break
+ for i in range(i+1,len(ll)):
+ if len(ll[i])>0 and ll[i][0] not in ' \t': break
+ ret+=ll[i]
+ return ret
+ else: return '[no documentation found]'
def mkWidget(self,entry):
if not entry.T: return None
# single fundamental object
@@ -407,8 +425,6 @@
grid.setVerticalSpacing(0)
grid.setLabelAlignment(Qt.AlignRight)
if self.showType:
- #lab=QLabel(makeSerializableLabel(self.ser,addr=True,href=True))
- #lab.setToolTip(('<b>'+self.path+'</b><br>' if self.path else '')+self.getDocstring())
lab=SerQLabel(self,makeSerializableLabel(self.ser,addr=True,href=True),tooltip=self.getDocstring(),path=self.path)
lab.setFrameShape(QFrame.Box); lab.setFrameShadow(QFrame.Sunken); lab.setLineWidth(2); lab.setAlignment(Qt.AlignHCenter); lab.linkActivated.connect(yade.qt.openUrl)
grid.setWidget(0,QFormLayout.SpanningRole,lab)
@@ -416,8 +432,6 @@
entry.widget=self.mkWidget(entry)
objPath=(self.path+'.'+entry.name) if self.path else None
label=SerQLabel(self,serializableHref(self.ser,entry.name),tooltip=self.getDocstring(entry.name),path=objPath)
- #label=SerAttrQLabel(self,serializableHref(self.ser,entry.name),tooltip=('<b>'+self.path+'.'+entry.name+'</b><br>' if self.path else '')+self.getDocstring(entry.name),
- #label=QLabel(self); label.setText(); label.setToolTip(; label.linkActivated.connect(yade.qt.openUrl)
grid.addRow(label,entry.widget if entry.widget else QLabel('<i>unhandled type</i>'))
self.setLayout(grid)
self.refreshEvent()
=== modified file 'pkg/dem/L3Geom.hpp'
--- pkg/dem/L3Geom.hpp 2010-12-13 12:11:43 +0000
+++ pkg/dem/L3Geom.hpp 2010-12-21 22:50:34 +0000
@@ -74,7 +74,7 @@
FUNCTOR1D(L6Geom);
void go(const shared_ptr<IGeom>&, const shared_ptr<Interaction>&, const shared_ptr<Body>&, const shared_ptr<Body>&, bool);
YADE_CLASS_BASE_DOC_STATICATTRS(Gl1_L6Geom,Gl1_L3Geom,"Render :yref:`L6Geom` geometry.",
- ((Real,phiScale,1.,,"Scale local rotations (:yref:`phi<L3Geom.phi>` - :yref:`phi0<L3Geom.phi0>`). The default scale is to draw $\\pi$ rotation with length equal to minimum radius."))
+ ((Real,phiScale,1.,,"Scale local rotations (:yref:`phi<L6Geom.phi>` - :yref:`phi0<L6Geom.phi0>`). The default scale is to draw $\\pi$ rotation with length equal to minimum radius."))
);
};
REGISTER_SERIALIZABLE(Gl1_L6Geom);
=== modified file 'pkg/dem/NewtonIntegrator.cpp'
--- pkg/dem/NewtonIntegrator.cpp 2010-12-21 18:31:15 +0000
+++ pkg/dem/NewtonIntegrator.cpp 2010-12-21 22:50:34 +0000
@@ -98,11 +98,8 @@
if(unlikely(!b)) continue;
State* state=b->state.get();
const Body::id_t& id=b->getId();
- // clump members are non-dynamic; we only get their velocities here
- if(unlikely(b->isClumpMember())){
- saveMaximaVelocity(scene,id,state);
- continue;
- }
+ // clump members are handled inside clumps
+ if(unlikely(b->isClumpMember())) continue;
const Vector3r& f=scene->forces.getForce(id);
const Vector3r& m=scene->forces.getTorque(id);
#ifdef YADE_DEBUG
@@ -137,30 +134,27 @@
if (likely(b->isStandalone())){
// translate equation
- if(b->isDynamic()){
- state->accel=computeAccel(f,state->mass,state->blockedDOFs);
- cundallDamp(dt,f,fluctVel,state->accel);}
+ state->accel=computeAccel(f,state->mass,state->blockedDOFs);
+ cundallDamp(dt,f,fluctVel,state->accel);
leapfrogTranslate(scene,state,id,dt);
// rotate equation
// exactAsphericalRot is disabled or the body is spherical
- if (likely(!exactAsphericalRot || !b->isAspherical() || !b->isDynamic())){
- if(b->isDynamic()) {
- state->angAccel=computeAngAccel(m,state->inertia,state->blockedDOFs);
- cundallDamp(dt,m,state->angVel,state->angAccel);}
+ if (likely(!exactAsphericalRot || !b->isAspherical() || b->state->inertia.minCoeff()<=0.)){
+ state->angAccel=computeAngAccel(m,state->inertia,state->blockedDOFs);
+ cundallDamp(dt,m,state->angVel,state->angAccel);
leapfrogSphericalRotate(scene,state,id,dt);
} else { // exactAsphericalRot enabled & aspherical body
// no damping in this case
leapfrogAsphericalRotate(scene,state,id,dt,m);
}
} else if (b->isClump()){
- if(b->isDynamic()){
- // reset acceleration of the clump itself; computed from accels on constituents
- state->accel=state->angAccel=Vector3r::Zero();
- // clump mass forces
- Vector3r dLinAccel=computeAccel(f,state->mass,state->blockedDOFs);
- cundallDamp(dt,f,fluctVel,dLinAccel);
- state->accel+=dLinAccel;}
- Vector3r M(m);
+ // reset acceleration of the clump itself; computed from accels on constituents
+ state->accel=state->angAccel=Vector3r::Zero();
+ // clump mass forces
+ Vector3r dLinAccel=computeAccel(f,state->mass,state->blockedDOFs);
+ cundallDamp(dt,f,fluctVel,dLinAccel);
+ state->accel+=dLinAccel;
+ Vector3r M(m); // torque on the clump itself, will accumulate from members
// sum force on clump memebrs
// exactAsphericalRot enabled and clump is aspherical
if (exactAsphericalRot && b->isAspherical() && b->state->inertia.maxCoeff()>0){
@@ -168,38 +162,35 @@
const Body::id_t& memberId=mm.first;
const shared_ptr<Body>& member=Body::byId(memberId,scene); assert(member->isClumpMember());
State* memberState=member->state.get();
- if(b->isDynamic()){
- handleClumpMemberAccel(scene,memberId,memberState,state);
- handleClumpMemberTorque(scene,memberId,memberState,state,M);}
- //FIXME : we are saving max velocity of previous timestep here, and actually it's already saved above
- saveMaximaVelocity(scene,memberId,memberState);
+ handleClumpMemberAccel(scene,memberId,memberState,state);
+ handleClumpMemberTorque(scene,memberId,memberState,state,M);
}
// motion
leapfrogTranslate(scene,state,id,dt);
leapfrogAsphericalRotate(scene,state,id,dt,M);
} else { // exactAsphericalRot disabled or clump is spherical
- if(b->isDynamic()){
- Vector3r dAngAccel=computeAngAccel(M,state->inertia,state->blockedDOFs);
- cundallDamp(dt,M,state->angVel,dAngAccel);
- state->angAccel+=dAngAccel;}
+ Vector3r dAngAccel=computeAngAccel(M,state->inertia,state->blockedDOFs);
+ cundallDamp(dt,M,state->angVel,dAngAccel);
+ state->angAccel+=dAngAccel;
FOREACH(Clump::MemberMap::value_type mm, static_cast<Clump*>(b->shape.get())->members){
const Body::id_t& memberId=mm.first;
const shared_ptr<Body>& member=Body::byId(memberId,scene); assert(member->isClumpMember());
State* memberState=member->state.get();
- if(b->isDynamic()){
- handleClumpMemberAccel(scene,memberId,memberState,state);
- handleClumpMemberAngAccel(scene,memberId,memberState,state);}
- //FIXME : we are saving max velocity of previous timestep here, and actually it's already saved above
- saveMaximaVelocity(scene,memberId,memberState);
+ handleClumpMemberAccel(scene,memberId,memberState,state);
+ handleClumpMemberAngAccel(scene,memberId,memberState,state);
}
// motion
leapfrogTranslate(scene,state,id,dt);
leapfrogSphericalRotate(scene,state,id,dt);
}
Clump::moveMembers(b,scene);
+ // save max velocity for all members of the clump
+ FOREACH(Clump::MemberMap::value_type mm, static_cast<Clump*>(b->shape.get())->members){
+ const shared_ptr<Body>& member=Body::byId(mm.first,scene); assert(member->isClumpMember());
+ saveMaximaVelocity(scene,mm.first,member->state.get());
+ }
}
saveMaximaVelocity(scene,id,state);
-
// process callbacks
for(size_t i=0; i<callbacksSize; i++){
cerr<<"<"<<b->id<<",cb="<<callbacks[i]<<",scene="<<callbacks[i]->scene<<">"; // <<",force="<<callbacks[i]->scene->forces.getForce(b->id)<<">";
@@ -228,13 +219,13 @@
} else if (homoDeform==Cell::HOMO_POS){
state->pos+=scene->cell->velGrad*state->pos*dt;
}
- if(state->blockedDOFs!=State::DOF_ALL) state->vel+=dt*state->accel;
- state->pos += state->vel*dt;
+ state->vel+=dt*state->accel;
+ state->pos+=state->vel*dt;
}
inline void NewtonIntegrator::leapfrogSphericalRotate(Scene* scene, State* state, const Body::id_t& id, const Real& dt )
{
- if(state->blockedDOFs!=State::DOF_ALL) state->angVel+=dt*state->angAccel;
+ state->angVel+=dt*state->angAccel;
Vector3r axis = state->angVel;
if (axis!=Vector3r::Zero()) { //If we have an angular velocity, we make a rotation
@@ -252,6 +243,14 @@
}
void NewtonIntegrator::leapfrogAsphericalRotate(Scene* scene, State* state, const Body::id_t& id, const Real& dt, const Vector3r& M){
+ // if all DoFs are blocked then return; warn if there is angular acceleration, since it will be ignored
+ if((state->blockedDOFs & State::DOF_RXRYRZ)==State::DOF_RXRYRZ){
+ /* necessarily angular acceleration is zero, we can call the spherical integrator safely */
+ if(state->angVel!=Vector3r::Zero()) leapfrogSphericalRotate(scene,state,id,dt); /* in most cases, angVel will be zero; no need to integrate anything, just return */
+ return;
+ }
+ // if only some rotations are allowed, warn, since it is not handled (yet?)
+ if((state->blockedDOFs & State::DOF_RXRYRZ)!=0) LOG_WARN("Aspherical body #"<<id<<" with some rotational DoFs blocked is not supported (blocked DoFs ignored).");
Matrix3r A=state->ori.conjugate().toRotationMatrix(); // rotation matrix from global to local r.f.
const Vector3r l_n = state->angMom + dt/2 * M; // global angular momentum at time n
const Vector3r l_b_n = A*l_n; // local angular momentum at time n
@@ -261,7 +260,6 @@
state->angMom+=dt*M; // global angular momentum at time n+1/2
const Vector3r l_b_half = A*state->angMom; // local angular momentum at time n+1/2
Vector3r angVel_b_half = l_b_half.cwise()/state->inertia; // local angular velocity at time n+1/2
- //blockRotateDOFs( state->blockedDOFs, angVel_b_half ); // FIXME: angVel_b_half is the local velocity, but blockedDOFs need for the _global_ velocity.
const Quaternionr dotQ_half=DotQ(angVel_b_half,Q_half); // dQ/dt at time n+1/2
state->ori=state->ori+dt*dotQ_half; // Q at time n+1
state->angVel=state->ori*angVel_b_half; // global angular velocity at time n+1/2