From 7b34b5f6f98115efbda9d127378aea5d403a4905 Mon Sep 17 00:00:00 2001
From: Gael Guennebaud <g.gael@free.fr>
Date: Tue, 24 Jul 2012 18:19:56 +0200
Subject: do not apply plane rotation when it is exactly the identity

---
 Eigen/src/Jacobi/Jacobi.h | 26 +++++++++++++++-----------
 1 file changed, 15 insertions(+), 11 deletions(-)

(limited to 'Eigen/src/Jacobi')

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;
     }
-- 
cgit v1.2.3