yade-dev team mailing list archive
-
yade-dev team
-
Mailing list archive
-
Message #08268
[Branch ~yade-dev/yade/trunk] Rev 3009: - dispatching: don't make bboxes larger for static bodies
------------------------------------------------------------
revno: 3009
committer: Bruno Chareyre <bruno.chareyre@xxxxxxxxxxx>
branch nick: trunk
timestamp: Mon 2012-01-30 21:53:49 +0100
message:
- dispatching: don't make bboxes larger for static bodies
- gravity : deprec. warning
- collider: try and fix Anton's crash (https://lists.launchpad.net/yade-dev/msg08256.html and next in the thread)
- Newton: workaround the resetForce bug in OMP build (https://bugs.launchpad.net/yade/+bug/923929)
modified:
pkg/common/Dispatching.cpp
pkg/common/GravityEngines.hpp
pkg/common/InsertionSortCollider.cpp
pkg/dem/NewtonIntegrator.cpp
pkg/dem/NewtonIntegrator.hpp
--
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 'pkg/common/Dispatching.cpp'
--- pkg/common/Dispatching.cpp 2012-01-23 14:43:54 +0000
+++ pkg/common/Dispatching.cpp 2012-01-30 20:53:49 +0000
@@ -28,9 +28,11 @@
if (targetInterv>=0) {
Vector3r disp = b->state->pos-b->bound->refPos;
Real dist = max(abs(disp[0]),max(abs(disp[1]),abs(disp[2])));
- Real newLength = dist*targetInterv/(scene->iter-b->bound->lastUpdateIter);
- newLength = max(0.9*sweepLength,newLength);//don't decrease size too fast to prevent time consuming oscillations
- sweepLength=max(minSweepDistFactor*sweepDist,min(newLength,sweepDist));
+ if (dist){
+ Real newLength = dist*targetInterv/(scene->iter-b->bound->lastUpdateIter);
+ newLength = max(0.9*sweepLength,newLength);//don't decrease size too fast to prevent time consuming oscillations
+ sweepLength=max(minSweepDistFactor*sweepDist,min(newLength,sweepDist));}
+ else sweepLength=0;
} else sweepLength=sweepDist;
}
#ifdef BV_FUNCTOR_CACHE
@@ -40,11 +42,10 @@
operator()(shape,b->bound,b->state->se3,b.get());
#endif
if(!b->bound) continue; // the functor did not create new bound
-
- if(sweepDist>0){
- b->bound->refPos=b->state->pos;
- b->bound->lastUpdateIter=scene->iter;
- const Real& sweepLength = b->bound->sweepLength;
+ b->bound->refPos=b->state->pos;
+ b->bound->lastUpdateIter=scene->iter;
+ const Real& sweepLength = b->bound->sweepLength;
+ if(sweepLength>0){
Aabb* aabb=YADE_CAST<Aabb*>(b->bound.get());
aabb->min-=Vector3r(sweepLength,sweepLength,sweepLength);
aabb->max+=Vector3r(sweepLength,sweepLength,sweepLength);
=== modified file 'pkg/common/GravityEngines.hpp'
--- pkg/common/GravityEngines.hpp 2012-01-30 10:42:13 +0000
+++ pkg/common/GravityEngines.hpp 2012-01-30 20:53:49 +0000
@@ -9,7 +9,7 @@
class GravityEngine: public FieldApplier{
public:
virtual void action();
- YADE_CLASS_BASE_DOC_ATTRS_CTOR_PY(GravityEngine,FieldApplier,"Engine applying constant acceleration to all bodies.",
+ YADE_CLASS_BASE_DOC_ATTRS_CTOR_PY(GravityEngine,FieldApplier,"Engine applying constant acceleration to all bodies. DEPRECATED, use :yref:`Newton::gravity` unless you need energy tracking or selective gravity application using groupMask).",
((Vector3r,gravity,Vector3r::Zero(),,"Acceleration [kgms⁻²]"))
((int,gravPotIx,-1,(Attr::noSave|Attr::hidden),"Index for gravPot energy"))
((int,mask,0,,"If mask defined, only bodies with corresponding groupMask will be affected by this engine. If 0, all bodies will be affected."))
=== modified file 'pkg/common/InsertionSortCollider.cpp'
--- pkg/common/InsertionSortCollider.cpp 2012-01-30 10:42:13 +0000
+++ pkg/common/InsertionSortCollider.cpp 2012-01-30 20:53:49 +0000
@@ -64,7 +64,7 @@
if (!it->flags.isMin || !it->flags.hasBB) continue;
int offset = 3*it->id;
const shared_ptr<Body>& b=Body::byId(it->id,scene);
- if(unlikely(!b)) continue;
+ if(unlikely(!b || !b->bound)) continue;
const Real& sweepLength = b->bound->sweepLength;
Vector3r disp = b->state->pos - b->bound->refPos;
if (!(maxima[offset]-sweepLength+disp[0] < bv.min[0] ||
=== modified file 'pkg/dem/NewtonIntegrator.cpp'
--- pkg/dem/NewtonIntegrator.cpp 2012-01-30 10:42:13 +0000
+++ pkg/dem/NewtonIntegrator.cpp 2012-01-30 20:53:49 +0000
@@ -81,8 +81,24 @@
#endif
}
+#ifdef YADE_OPENMP
+void NewtonIntegrator::ensureSync()
+{
+ if (syncEnsured) return;
+ YADE_PARALLEL_FOREACH_BODY_BEGIN(const shared_ptr<Body>& b, scene->bodies){
+ if(b->isClump()) continue;
+ scene->forces.addForce(b->getId(),Vector3r(0,0,0));
+ } YADE_PARALLEL_FOREACH_BODY_END();
+ syncEnsured=true;
+}
+#endif
+
void NewtonIntegrator::action()
{
+ #ifdef YADE_OPENMP
+ //prevent https://bugs.launchpad.net/yade/+bug/923929
+ ensureSync();
+ #endif
scene->forces.sync();
bodySelected=(scene->selectedBody>=0);
if(warnNoForceReset && scene->forces.lastReset<scene->iter) LOG_WARN("O.forces last reset in step "<<scene->forces.lastReset<<", while the current step is "<<scene->iter<<". Did you forget to include ForceResetter in O.engines?");
@@ -157,7 +173,7 @@
Vector3r linAccel=computeAccel(f,state->mass,state->blockedDOFs);
if(state->isDamped) cundallDamp2nd(dt,fluctVel,linAccel);
//This is the convective term, appearing in the time derivation of Cundall/Thornton expression (dx/dt=velGrad*pos -> d²x/dt²=dvelGrad/dt*pos+velGrad*vel), negligible in many cases but not for high speed large deformations (gaz or turbulent flow).
- linAccel+=prevVelGrad*state->vel;
+ if (isPeriodic) linAccel+=prevVelGrad*state->vel;
//finally update velocity
state->vel+=dt*linAccel;
// angular acceleration
@@ -222,11 +238,11 @@
void NewtonIntegrator::leapfrogAsphericalRotate(State* state, const Body::id_t& id, const Real& dt, const Vector3r& M){
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_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
const Vector3r angVel_b_n = l_b_n.cwise()/state->inertia; // local angular velocity at time n
const Quaternionr dotQ_n=DotQ(angVel_b_n,state->ori); // dQ/dt at time n
- const Quaternionr Q_half = state->ori + dt/2 * dotQ_n; // Q at time n+1/2
+ const Quaternionr Q_half = state->ori + dt/2. * dotQ_n; // Q at time n+1/2
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
=== modified file 'pkg/dem/NewtonIntegrator.hpp'
--- pkg/dem/NewtonIntegrator.hpp 2012-01-30 10:42:13 +0000
+++ pkg/dem/NewtonIntegrator.hpp 2012-01-30 20:53:49 +0000
@@ -33,7 +33,9 @@
Vector3r computeAngAccel(const Vector3r& torque, const Vector3r& inertia, int blockedDOFs);
void updateEnergy(const shared_ptr<Body>&b, const State* state, const Vector3r& fluctVel, const Vector3r& f, const Vector3r& m);
-
+ #ifdef YADE_OPENMP
+ void ensureSync(); bool syncEnsured;
+ #endif
// whether the cell has changed from the previous step
bool cellChanged;
@@ -70,7 +72,7 @@
,
/*ctor*/
#ifdef YADE_OPENMP
- threadMaxVelocitySq.resize(omp_get_max_threads());
+ threadMaxVelocitySq.resize(omp_get_max_threads()); syncEnsured=false;
#endif
);
DECLARE_LOGGER;
