aboutsummaryrefslogtreecommitdiffhomepage
path: root/Eigen/src/Jacobi
diff options
context:
space:
mode:
authorGravatar Gael Guennebaud <g.gael@free.fr>2012-07-24 18:19:56 +0200
committerGravatar Gael Guennebaud <g.gael@free.fr>2012-07-24 18:19:56 +0200
commit7b34b5f6f98115efbda9d127378aea5d403a4905 (patch)
tree2f19605b569a48b5e233723201e2f5eb148fe242 /Eigen/src/Jacobi
parente7c07de549ace45f11c000c15f6a9fcb1c9b8170 (diff)
do not apply plane rotation when it is exactly the identity
Diffstat (limited to 'Eigen/src/Jacobi')
-rw-r--r--Eigen/src/Jacobi/Jacobi.h26
1 files changed, 15 insertions, 11 deletions
diff --git a/Eigen/src/Jacobi/Jacobi.h b/Eigen/src/Jacobi/Jacobi.h
index a9c17dcdf..7eb19a023 100644
--- a/Eigen/src/Jacobi/Jacobi.h
+++ b/Eigen/src/Jacobi/Jacobi.h
@@ -207,7 +207,6 @@ void JacobiRotation<Scalar>::makeGivens(const Scalar& p, const Scalar& q, Scalar
template<typename Scalar>
void JacobiRotation<Scalar>::makeGivens(const Scalar& p, const Scalar& q, Scalar* r, internal::false_type)
{
-
if(q==Scalar(0))
{
m_c = p<Scalar(0) ? Scalar(-1) : Scalar(1);
@@ -303,6 +302,11 @@ void /*EIGEN_DONT_INLINE*/ apply_rotation_in_the_plane(VectorX& _x, VectorY& _y,
Scalar* EIGEN_RESTRICT x = &_x.coeffRef(0);
Scalar* EIGEN_RESTRICT y = &_y.coeffRef(0);
+
+ OtherScalar c = j.c();
+ OtherScalar s = j.s();
+ if (c==OtherScalar(1) && s==OtherScalar(0))
+ return;
/*** dynamic-size vectorized paths ***/
@@ -316,16 +320,16 @@ void /*EIGEN_DONT_INLINE*/ apply_rotation_in_the_plane(VectorX& _x, VectorY& _y,
Index alignedStart = internal::first_aligned(y, size);
Index alignedEnd = alignedStart + ((size-alignedStart)/PacketSize)*PacketSize;
- const Packet pc = pset1<Packet>(j.c());
- const Packet ps = pset1<Packet>(j.s());
+ const Packet pc = pset1<Packet>(c);
+ const Packet ps = pset1<Packet>(s);
conj_helper<Packet,Packet,NumTraits<Scalar>::IsComplex,false> pcj;
for(Index i=0; i<alignedStart; ++i)
{
Scalar xi = x[i];
Scalar yi = y[i];
- x[i] = j.c() * xi + conj(j.s()) * yi;
- y[i] = -j.s() * xi + conj(j.c()) * yi;
+ x[i] = c * xi + conj(s) * yi;
+ y[i] = -s * xi + conj(c) * yi;
}
Scalar* EIGEN_RESTRICT px = x + alignedStart;
@@ -372,8 +376,8 @@ void /*EIGEN_DONT_INLINE*/ apply_rotation_in_the_plane(VectorX& _x, VectorY& _y,
{
Scalar xi = x[i];
Scalar yi = y[i];
- x[i] = j.c() * xi + conj(j.s()) * yi;
- y[i] = -j.s() * xi + conj(j.c()) * yi;
+ x[i] = c * xi + conj(s) * yi;
+ y[i] = -s * xi + conj(c) * yi;
}
}
@@ -382,8 +386,8 @@ void /*EIGEN_DONT_INLINE*/ apply_rotation_in_the_plane(VectorX& _x, VectorY& _y,
(VectorX::Flags & VectorY::Flags & PacketAccessBit) &&
(VectorX::Flags & VectorY::Flags & AlignedBit))
{
- const Packet pc = pset1<Packet>(j.c());
- const Packet ps = pset1<Packet>(j.s());
+ const Packet pc = pset1<Packet>(c);
+ const Packet ps = pset1<Packet>(s);
conj_helper<Packet,Packet,NumTraits<Scalar>::IsComplex,false> pcj;
Scalar* EIGEN_RESTRICT px = x;
Scalar* EIGEN_RESTRICT py = y;
@@ -405,8 +409,8 @@ void /*EIGEN_DONT_INLINE*/ apply_rotation_in_the_plane(VectorX& _x, VectorY& _y,
{
Scalar xi = *x;
Scalar yi = *y;
- *x = j.c() * xi + conj(j.s()) * yi;
- *y = -j.s() * xi + conj(j.c()) * yi;
+ *x = c * xi + conj(s) * yi;
+ *y = -s * xi + conj(c) * yi;
x += incrx;
y += incry;
}