div / 0 in BuildJacobianMatrix

Report any bugs here and we'll post fixes

Moderators: Sascha Willems, Thomas

Re: div / 0 in BuildJacobianMatrix

Postby arkeon » Sat Sep 20, 2014 6:35 am

In the dgWorldDynamicUpdate::CalculateForcesGameMode
When I add a conditionnal break point where "linearMomentum.m_y < (-10.0f)" on line 1047

I can make it stop when the huge acceleration happen.
what if strange is there is a condition "if (timestepRK != dgFloat32 (0.0f)) {"
but when the break point is triggered timestepRK is equal to 0.0;

how can it be ?!

.. never mind I'm not awake the 1047 line is in the else...
arkeon
 
Posts: 261
Joined: Sat Sep 13, 2014 5:25 pm

Re: div / 0 in BuildJacobianMatrix

Postby Julio Jerez » Sat Sep 20, 2014 6:43 am

yes that is why I can not debugged it, I renamed but I got this error
Untitled.png
Untitled.png (15.2 KiB) Viewed 4231 times


Si I thought that was the old car from newton 2.00, hence my assumption.
I any case I will not be able to check the code with a different DLL.

can you post the part that implement the call back in the joint you are using them.
The bug is there.
Julio Jerez
Moderator
Moderator
 
Posts: 11039
Joined: Sun Sep 14, 2003 2:18 pm
Location: Los Angeles

Re: div / 0 in BuildJacobianMatrix

Postby arkeon » Sat Sep 20, 2014 6:53 am

my slider joint code
Code: Select all
void CustomSlider::SubmitConstraints (dFloat timestep, int threadIndex)
{
   dFloat dist;
   dMatrix matrix0;
   dMatrix matrix1;
        dFloat invTimestep = (timestep > 0.0f) ? 1.0f / timestep: 1.0f;

   // calculate the position of the pivot point and the Jacobian direction vectors, in global space.
   CalculateGlobalMatrix (matrix0, matrix1);

   // Restrict the movement on the pivot point along all tree orthonormal direction
   dVector p0 (matrix0.m_posit);
   dVector p1 (matrix1.m_posit + matrix1.m_front.Scale ((p0 - matrix1.m_posit) % matrix1.m_front));
   NewtonUserJointAddLinearRow (m_joint, &p0[0], &p1[0], &matrix0.m_up[0]);
   NewtonUserJointAddLinearRow (m_joint, &p0[0], &p1[0], &matrix0.m_right[0]);
   
   // two constraints row perpendicular to the pin
   // get a point along the ping axis at some reasonable large distance from the pivot
   dVector q0 (p0 + matrix0.m_front.Scale(MIN_JOINT_PIN_LENGTH));
   dVector q1 (p1 + matrix1.m_front.Scale(MIN_JOINT_PIN_LENGTH));

   NewtonUserJointAddLinearRow (m_joint, &q0[0], &q1[0], &matrix0.m_up[0]);
   NewtonUserJointAddLinearRow (m_joint, &q0[0], &q1[0], &matrix0.m_right[0]);

   // one constraint row perpendicular to the pin
   // get a point along the ping axis at some reasonable large distance from the pivot
   dVector r0 (p0 + matrix0.m_up.Scale(MIN_JOINT_PIN_LENGTH));
   dVector r1 (p1 + matrix1.m_up.Scale(MIN_JOINT_PIN_LENGTH));
    NewtonUserJointAddLinearRow (m_joint, &r0[0], &r1[0], &matrix1.m_right[0]);

   // calculate position and speed   
   dVector veloc0(0.0f, 0.0f, 0.0f, 0.0f);
   dVector veloc1(0.0f, 0.0f, 0.0f, 0.0f); 
   if (m_body0) {
      NewtonBodyGetVelocity(m_body0, &veloc0[0]);
   }
   if (m_body1) {
      NewtonBodyGetVelocity(m_body1, &veloc1[0]);
   }
   m_posit = (matrix0.m_posit - matrix1.m_posit) % matrix1.m_front;
   m_speed = (veloc0 - veloc1) % matrix1.m_front;

   dist = (matrix0.m_posit - matrix1.m_posit) % matrix0.m_front;
   dFloat nextDist = dist + (m_motorVelocity * timestep);
   m_currentDistance = dist;

   // if limit are enable ...
   m_hitLimitOnLastUpdate = false;
   if (m_limitsOn)
  {
      if (dist < m_minDist && (!m_motorOn || (m_motorOn && (nextDist < m_minDist))))
    {
         // indicate that this row hit a limit
         m_hitLimitOnLastUpdate = true;

         // get a point along the up vector and set a constraint 
         const dVector& p0 = matrix0.m_posit;
         dVector p1 (p0 + matrix0.m_front.Scale (m_minDist - dist));
         NewtonUserJointAddLinearRow (m_joint, &p0[0], &p1[0], &matrix0.m_front[0]);
         // allow the object to return but not to kick going forward
         NewtonUserJointSetRowMinimumFriction (m_joint, 0.0f);         
      }
    else if (dist > m_maxDist && (!m_motorOn || (m_motorOn && (nextDist > m_maxDist))))
    {
         // indicate that this row hit a limit
         m_hitLimitOnLastUpdate = true;

         // get a point along the up vector and set a constraint 

         const dVector& p0 = matrix0.m_posit;
         dVector p1 (p0 + matrix0.m_front.Scale (m_maxDist - dist));
         NewtonUserJointAddLinearRow (m_joint, &p0[0], &p1[0], &matrix0.m_front[0]);
         // allow the object to return but not to kick going forward
         NewtonUserJointSetRowMaximumFriction (m_joint, 0.0f);

      }
    else if (m_friction > 0.0f)
    {
         // take any point on body0 (origin)
         const dVector& p0 = matrix0.m_posit;

         dVector veloc0(0.0f, 0.0f, 0.0f);
        dVector veloc1(0.0f, 0.0f, 0.0f);
        dVector omega1(0.0f, 0.0f, 0.0f);

        NewtonBodyGetVelocity(m_body0, &veloc0[0]);
      if (m_body1)
      {
          NewtonBodyGetVelocity(m_body1, &veloc1[0]);
          NewtonBodyGetOmega(m_body1, &omega1[0]);
      }

         // this assumes the origin of the bodies the matrix pivot are the same
         veloc1 += omega1 * (matrix1.m_posit - p0);
     
         dFloat relAccel;
      if (m_motorOn)
        relAccel = (m_motorVelocity - ((veloc0 - veloc1) % matrix0.m_front)) * m_motorStrength;
      else
        relAccel = ((veloc0 - veloc1) % matrix0.m_front);     

         NewtonUserJointAddLinearRow (m_joint, &p0[0], &p0[0], &matrix0.m_front[0]);
         NewtonUserJointSetRowAcceleration (m_joint, relAccel * invTimestep);
         NewtonUserJointSetRowMinimumFriction (m_joint, -m_friction);
         NewtonUserJointSetRowMaximumFriction(m_joint, m_friction);
      }
    else if (m_motorOn)
    {
      // take any point on body0 (origin)
        const dVector& p0 = matrix0.m_posit;

        dVector veloc0(0.0f, 0.0f, 0.0f);
        dVector veloc1(0.0f, 0.0f, 0.0f);
        dVector omega1(0.0f, 0.0f, 0.0f);

        NewtonBodyGetVelocity(m_body0, &veloc0[0]);
      if (m_body1)
      {
          NewtonBodyGetVelocity(m_body1, &veloc1[0]);
          NewtonBodyGetOmega(m_body1, &omega1[0]);
      }

        // this assumes the origin of the bodies the matrix pivot are the same
        veloc1 += omega1 * (matrix1.m_posit - p0);
     
        dFloat relAccel = (m_motorVelocity - ((veloc0 - veloc1) % matrix0.m_front)) * m_motorStrength;

        NewtonUserJointAddLinearRow (m_joint, &p0[0], &p0[0], &matrix0.m_front[0]);
        NewtonUserJointSetRowAcceleration (m_joint, relAccel * invTimestep);
    }
   }
  else if (m_friction > 0.0f)
  {
      // take any point on body0 (origin)
      const dVector& p0 = matrix0.m_posit;

      dVector veloc0(0.0f, 0.0f, 0.0f);
      dVector veloc1(0.0f, 0.0f, 0.0f);
      dVector omega1(0.0f, 0.0f, 0.0f);

      NewtonBodyGetVelocity(m_body0, &veloc0[0]);
    if (m_body1)
    {
        NewtonBodyGetVelocity(m_body1, &veloc1[0]);
        NewtonBodyGetOmega(m_body1, &omega1[0]);
    }

      // this assumes the origin of the bodies the matrix pivot are the same
      veloc1 += omega1 * (matrix1.m_posit - p0);
   
      dFloat relAccel;
    if (m_motorOn)
      relAccel = (m_motorVelocity - ((veloc0 - veloc1) % matrix0.m_front)) * m_motorStrength;
    else
      relAccel = ((veloc0 - veloc1) % matrix0.m_front); 

      NewtonUserJointAddLinearRow (m_joint, &p0[0], &p0[0], &matrix0.m_front[0]);
      NewtonUserJointSetRowAcceleration (m_joint, relAccel * invTimestep);
      NewtonUserJointSetRowMinimumFriction (m_joint, -m_friction);
      NewtonUserJointSetRowMaximumFriction(m_joint, m_friction);
   }
  else if (m_motorOn)
  {
    // take any point on body0 (origin)
      const dVector& p0 = matrix0.m_posit;

      dVector veloc0(0.0f, 0.0f, 0.0f);
      dVector veloc1(0.0f, 0.0f, 0.0f);
      dVector omega1(0.0f, 0.0f, 0.0f);

      NewtonBodyGetVelocity(m_body0, &veloc0[0]);
    if (m_body1)
    {
        NewtonBodyGetVelocity(m_body1, &veloc1[0]);
        NewtonBodyGetOmega(m_body1, &omega1[0]);
    }

      // this assumes the origin of the bodies the matrix pivot are the same
      veloc1 += omega1 * (matrix1.m_posit - p0);
   
      dFloat relAccel = (m_motorVelocity - ((veloc0 - veloc1) % matrix0.m_front)) * m_motorStrength;

      NewtonUserJointAddLinearRow (m_joint, &p0[0], &p0[0], &matrix0.m_front[0]);
      NewtonUserJointSetRowAcceleration (m_joint, relAccel * invTimestep);
  }

  if (m_useSpring)
  {
    NewtonUserJointSetRowStiffness(m_joint, 0.8f);
    NewtonUserJointSetRowSpringDamperAcceleration(m_joint, m_springStrength, m_springDamping);
  };
 }


specially on this joint I use high spring strength and damper value to get a smooth effect when I crouch.
arkeon
 
Posts: 261
Joined: Sat Sep 13, 2014 5:25 pm

Re: div / 0 in BuildJacobianMatrix

Postby Julio Jerez » Sat Sep 20, 2014 6:55 am

the bug happen on a ccd step in function
void dgWorldDynamicUpdate::CalculateIslandReactionForces (dgIsland* const island, dgFloat32 timestep, dgInt32 threadID) const

line line 204
Code: Select all
               CalculateIslandContacts (island, timeRemaining, lru, threadID);
               BuildJacobianMatrix (island, threadID, 0.0f);
               CalculateReactionsForces (island, threadID, 0.0f, DG_SOLVER_MAX_ERROR);

               bool islandResinding = true;
               for (dgInt32 k = 0; (k < DG_MAX_CONTINUE_COLLISON_STEPS) && islandResinding; k ++) {
                  dgFloat32 smallTimeStep = dgMin (timestep * dgFloat32 (1.0f / 8.0f), timeRemaining);
                  timeRemaining -= smallTimeStep;
                  for (dgInt32 j = 1; j < bodyCount; j ++) {
                     dgDynamicBody* const body = (dgDynamicBody*) bodyArray[j].m_body;
                     if (body->IsRTTIType (dgBody::m_dynamicBodyRTTI)) {
                        body->IntegrateVelocity (smallTimeStep);
                        body->UpdateWorlCollisionMatrix();
                     }
                  }


The call BuildJacobianMatrix (island, threadID, 0.0f);

is where the Joint callback are call from and in one of your joints the acceleration is growing exponentially, generation a huge velocity.
Is no an error is jut huge acceleration whi is wrong.
Julio Jerez
Moderator
Moderator
 
Posts: 11039
Joined: Sun Sep 14, 2003 2:18 pm
Location: Los Angeles

Re: div / 0 in BuildJacobianMatrix

Postby arkeon » Sat Sep 20, 2014 6:58 am

I think the acceleration is set by NewtonUserJointSetRowSpringDamperAcceleration ?
arkeon
 
Posts: 261
Joined: Sat Sep 13, 2014 5:25 pm

Re: div / 0 in BuildJacobianMatrix

Postby Julio Jerez » Sat Sep 20, 2014 7:00 am

even when I set a break point on this function
dFloat NewtonCalculateSpringDamperAcceleration(dFloat dt, dFloat ks, dFloat x, dFloat kd, dFloat s)

the function is not called but the joint.
There is a bug on that joint, and with out being able to step in the code I can not help.

My suggestion Is that you use the new library.
Julio Jerez
Moderator
Moderator
 
Posts: 11039
Joined: Sun Sep 14, 2003 2:18 pm
Location: Los Angeles

Re: div / 0 in BuildJacobianMatrix

Postby Julio Jerez » Sat Sep 20, 2014 7:17 am

ha, ok
Is set a break point there, and I can clear see the bug.
the bug was here in the newton.dll side

this function was not handlin the impulse called correctly.
Code: Select all
dgFloat32 dgBilateralConstraint::CalculateSpringDamperAcceleration (
   dgInt32 index,
   const dgContraintDescritor& desc,
   dgFloat32 jointAngle,
   const dgVector& p0Global,
   const dgVector& p1Global,
   dgFloat32 springK,
   dgFloat32 springD)
{


I fixed now and it seems to work, please sync and try again.

after that, I still recommend you use the new dJlointLibrary.dll
there vehicle is better, and there si also a player controller that resolve all those issue that you are having.
Newton 3.00 has a feature called listened that is use to implmenet lot of functionality.
Julio Jerez
Moderator
Moderator
 
Posts: 11039
Joined: Sun Sep 14, 2003 2:18 pm
Location: Los Angeles

Re: div / 0 in BuildJacobianMatrix

Postby arkeon » Sat Sep 20, 2014 7:51 am

cool :) well done this one was tricky to produce and debug ^^
Everything seems to works well now.

My first goal is to make all working the same as before, I'll take the time to clean and update everything later. Maybe when soft bodies will be finished ;)

and my weird configuration permit to find hidden bugs :)
Thanks a lot ! (again).
arkeon
 
Posts: 261
Joined: Sat Sep 13, 2014 5:25 pm

Re: div / 0 in BuildJacobianMatrix

Postby arkeon » Sat Sep 20, 2014 7:58 am

ho I forgot the weird collision shape on the hand object using NewtonCreateCompoundCollisionFromMesh.
the vertices created are out of the origin mesh and the shape is weird.

You remember the one fro the F8 demo on http://arkeon.dyndns.org/scol/physic_demo.rar
arkeon
 
Posts: 261
Joined: Sat Sep 13, 2014 5:25 pm

Re: div / 0 in BuildJacobianMatrix

Postby Julio Jerez » Sat Sep 20, 2014 8:21 am

can you make the demo so that it just do one hand, not tow, also only one demo not multiple demos.
Is the hand model as single mesh?


The one thong I noticed is that you use then engine with multi thread and asynchronous, very cool!
I assume you do have frame interpolation.
Julio Jerez
Moderator
Moderator
 
Posts: 11039
Joined: Sun Sep 14, 2003 2:18 pm
Location: Los Angeles

Re: div / 0 in BuildJacobianMatrix

Postby Julio Jerez » Sat Sep 20, 2014 9:05 am

when I set a break point on this function

Code: Select all
NewtonMesh* NewtonMeshApproximateConvexDecomposition (const NewtonMesh* const mesh, dFloat maxConcavity, dFloat backFaceDistanceFactor, int maxCount, int maxVertexPerHull, NewtonReportProgress progressReportCallback, void* const reportProgressUserData)


It is never called from program: physic_demo.exe
what function are you calling to create the convex decomposition?
Julio Jerez
Moderator
Moderator
 
Posts: 11039
Joined: Sun Sep 14, 2003 2:18 pm
Location: Los Angeles

Re: div / 0 in BuildJacobianMatrix

Postby arkeon » Sat Sep 20, 2014 10:52 am

you can download the sample with only one hand here : http://arkeon.dyndns.org/scol/handshape.rar

(the memory assert on program close was because I forgot to destroy the mesh ^^)
the code to create the collision look like that :
Code: Select all
          //first create a tree collision
          TreeCollision treecol(world, obj, false, id, FW_DEFAULT, scale);

          //get the mesh from the collision
          NewtonMesh* const mesh = NewtonMeshCreateFromCollision(treecol.getNewtonCollision());
          if (m_col)
            NewtonDestroyCollision(m_col);
          m_col = NewtonCreateCompoundCollisionFromMesh(m_world->getNewtonWorld(), mesh, tolerance, id, id);
          if (!m_col)
            m_col = NewtonCreateNull(m_world->getNewtonWorld());

          if (mesh)
            NewtonMeshDestroy(mesh);


and yes multi threaded solver is really cool :), I had strange behavior with frame interpolation. I'll see if I can enable it again from the sample in dnewton with GetInterpolationParam
arkeon
 
Posts: 261
Joined: Sat Sep 13, 2014 5:25 pm

Re: div / 0 in BuildJacobianMatrix

Postby Julio Jerez » Sat Sep 20, 2014 11:28 am

well the are lot of problem there. firs off all the input mesh is wrong I save the mesh as an
.off file using this function
NewtonMeshSaveOFF (convexAproximation, "xxxxx.off");
this is how it looks like, you are getting is what is being input
Untitled.png
Untitled.png (15.77 KiB) Viewed 4220 times


The second problem is that the mesh I tiny is if the order of centimeters is size and you using the process all backward. You are making a collision tree and converting it to a newtonmesh, however the process that builds a collision tree, condition the mesh for better collision, so I am guessing the all the missing polygons are filtered out by the collision tree builder because the as too small.

you should do it the other way around, use the NewtonMesh to create the mesh. the once you Mesh you can do few thing to the mesh, you can bake a scale on it, or you can optimize it if you want.

the to make the Compound collision you nee to create a convex decomposition by calling this function first
NEWTON_API NewtonMesh* NewtonMeshApproximateConvexDecomposition (const NewtonMesh* const mesh, dFloat maxConcavity, dFloat backFaceDistanceFactor, int maxCount, int maxVertexPerHull, NewtonReportProgress reportProgressCallback, void* const reportProgressUserData);

then after your mesh has being processed, the you can now call
NewtonCreateCompoundCollisionFromMesh to make a compound.

you are doing the worse possible way, there is no need to make a collision tree as intermediate step
Julio Jerez
Moderator
Moderator
 
Posts: 11039
Joined: Sun Sep 14, 2003 2:18 pm
Location: Los Angeles

Re: div / 0 in BuildJacobianMatrix

Postby arkeon » Sat Sep 20, 2014 11:38 am

doh ! it was working in most cases, but sure the way you explain must be better :)
The collision tree is build as it without optimization, so I thought it could work that way.
arkeon
 
Posts: 261
Joined: Sat Sep 13, 2014 5:25 pm

Re: div / 0 in BuildJacobianMatrix

Postby arkeon » Sat Sep 20, 2014 12:07 pm

wow it's a lot better :)
thanks again ^^

Image
arkeon
 
Posts: 261
Joined: Sat Sep 13, 2014 5:25 pm

PreviousNext

Return to Bugs and Fixes

Who is online

Users browsing this forum: No registered users and 2 guests

cron