1c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamathnamespace Eigen { 2c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 3c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamathnamespace internal { 4c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 5c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath// TODO : move this to GivensQR once there's such a thing in Eigen 6c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 7c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamathtemplate <typename Scalar> 8c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamathvoid r1mpyq(DenseIndex m, DenseIndex n, Scalar *a, const std::vector<JacobiRotation<Scalar> > &v_givens, const std::vector<JacobiRotation<Scalar> > &w_givens) 9c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath{ 10c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typedef DenseIndex Index; 11c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 12c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /* apply the first set of givens rotations to a. */ 13c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath for (Index j = n-2; j>=0; --j) 14c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath for (Index i = 0; i<m; ++i) { 15c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath Scalar temp = v_givens[j].c() * a[i+m*j] - v_givens[j].s() * a[i+m*(n-1)]; 16c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath a[i+m*(n-1)] = v_givens[j].s() * a[i+m*j] + v_givens[j].c() * a[i+m*(n-1)]; 17c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath a[i+m*j] = temp; 18c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 19c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /* apply the second set of givens rotations to a. */ 20c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath for (Index j = 0; j<n-1; ++j) 21c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath for (Index i = 0; i<m; ++i) { 22c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath Scalar temp = w_givens[j].c() * a[i+m*j] + w_givens[j].s() * a[i+m*(n-1)]; 23c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath a[i+m*(n-1)] = -w_givens[j].s() * a[i+m*j] + w_givens[j].c() * a[i+m*(n-1)]; 24c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath a[i+m*j] = temp; 25c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 26c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath} 27c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 28c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath} // end namespace internal 29c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 30c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath} // end namespace Eigen 31