Consider denormals as zero in makeJacobi and 2x2 SVD.

This also fix serious issues with x387 for which values can be much smaller than the smallest denormal!
This commit is contained in:
Gael Guennebaud 2016-07-12 17:21:03 +02:00
parent 263993a7b6
commit e2aa58b631
2 changed files with 5 additions and 3 deletions

View File

@ -85,7 +85,8 @@ bool JacobiRotation<Scalar>::makeJacobi(const RealScalar& x, const Scalar& y, co
using std::sqrt;
using std::abs;
typedef typename NumTraits<Scalar>::Real RealScalar;
if(y == Scalar(0))
RealScalar deno = RealScalar(2)*abs(y);
if(deno < (std::numeric_limits<RealScalar>::min)())
{
m_c = Scalar(1);
m_s = Scalar(0);
@ -93,7 +94,7 @@ bool JacobiRotation<Scalar>::makeJacobi(const RealScalar& x, const Scalar& y, co
}
else
{
RealScalar tau = (x-z)/(RealScalar(2)*abs(y));
RealScalar tau = (x-z)/deno;
RealScalar w = sqrt(numext::abs2(tau) + RealScalar(1));
RealScalar t;
if(tau>RealScalar(0))

View File

@ -28,7 +28,8 @@ void real_2x2_jacobi_svd(const MatrixType& matrix, Index p, Index q,
JacobiRotation<RealScalar> rot1;
RealScalar t = m.coeff(0,0) + m.coeff(1,1);
RealScalar d = m.coeff(1,0) - m.coeff(0,1);
if(d == RealScalar(0))
if(abs(d) < (std::numeric_limits<RealScalar>::min)())
{
rot1.s() = RealScalar(0);
rot1.c() = RealScalar(1);