![]() |
Eigen-unsupported
3.3.3
|
00001 namespace Eigen { 00002 00003 namespace internal { 00004 00005 // TODO : move this to GivensQR once there's such a thing in Eigen 00006 00007 template <typename Scalar> 00008 void r1mpyq(DenseIndex m, DenseIndex n, Scalar *a, const std::vector<JacobiRotation<Scalar> > &v_givens, const std::vector<JacobiRotation<Scalar> > &w_givens) 00009 { 00010 typedef DenseIndex Index; 00011 00012 /* apply the first set of givens rotations to a. */ 00013 for (Index j = n-2; j>=0; --j) 00014 for (Index i = 0; i<m; ++i) { 00015 Scalar temp = v_givens[j].c() * a[i+m*j] - v_givens[j].s() * a[i+m*(n-1)]; 00016 a[i+m*(n-1)] = v_givens[j].s() * a[i+m*j] + v_givens[j].c() * a[i+m*(n-1)]; 00017 a[i+m*j] = temp; 00018 } 00019 /* apply the second set of givens rotations to a. */ 00020 for (Index j = 0; j<n-1; ++j) 00021 for (Index i = 0; i<m; ++i) { 00022 Scalar temp = w_givens[j].c() * a[i+m*j] + w_givens[j].s() * a[i+m*(n-1)]; 00023 a[i+m*(n-1)] = -w_givens[j].s() * a[i+m*j] + w_givens[j].c() * a[i+m*(n-1)]; 00024 a[i+m*j] = temp; 00025 } 00026 } 00027 00028 } // end namespace internal 00029 00030 } // end namespace Eigen