[kernel] add a default allocation of _jacobianWrenchTwist since _jacobianMGyrtwist is always present. Correct the size of _jacobianMGyrtwist
@@ -400,8 +400,14 @@ void NewtonEulerDS::init()
_wrench.reset(new SiconosVector(_n));
_mGyr.reset(new SiconosVector(3,0.0));
+ /** The follwing jacobian are always allocated since we have always
+ * Gyroscopical forces that has non linear forces
+ * This should be remove if the integration is explicit or _nullifyMGyr(false) is set to true ?
+ */
+
+ _jacobianMGyrtwist.reset(new SimpleMatrix(3, _n));
+ _jacobianWrenchTwist.reset(new SimpleMatrix(_n, _n));
- _jacobianMGyrtwist.reset(new SimpleMatrix(_n, _n));
//We initialize _z with a null vector of size 1, since z is required in plug-in functions call.
_z.reset(new SiconosVector(1));
@@ -963,7 +969,7 @@ void NewtonEulerDS::computeJacobianvForces(double time)
{
//computeJacobianMGyrtwistByFD(time,_q,_twist);
computeJacobianMGyrtwist(time);
- *_jacobianWrenchTwist -= *_jacobianMGyrtwist;
+ _jacobianWrenchTwist->setBlock(3,0,-1.0 * *_jacobianMGyrtwist);
}
}
// std::cout << "_jacobianWrenchTwist : "<< std::endl;
@@ -1000,7 +1006,7 @@ void NewtonEulerDS::computeJacobianMGyrtwist(double time)
cross_product(omega, Iei, omega_Iei);
cross_product(ei, Iomega, ei_Iomega);
for (int j = 0; j < 3; j++)
- _jacobianMGyrtwist->setValue(3 + j, 3 + i, ei_Iomega.getValue(j) + omega_Iei.getValue(j));
+ _jacobianMGyrtwist->setValue(j, 3 + i, ei_Iomega.getValue(j) + omega_Iei.getValue(j));
}
// Check if Jacobian is valid. Warning to the transpose operation in
// _jacobianMGyrtwist->setValue(3 + j, 3 + i, ei_Iomega.getValue(j) + omega_Iei.getValue(j));
