mirror/eigen
2
0
Fork 0
mirror of https://gitlab.com/libeigen/eigen.git synced 2024-12-27 07:29:52 +08:00
Code Issues Packages Projects Releases Wiki Activity

do not apply plane rotation when it is exactly the identity

Browse Source
This commit is contained in:
Gael Guennebaud 2012-07-24 18:19:56 +02:00
parent e7c07de549
commit 7b34b5f6f9
1 changed files with 15 additions and 11 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

26
Eigen/src/Jacobi/Jacobi.h
View File

@ -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;
}