src/NonLinearOptimization/rwupdt.h

0001 namespace Eigen {
0002
0003 namespace internal {
0004
0005 template <typename Scalar>
0006 void rwupdt(
0007         Matrix< Scalar, Dynamic, Dynamic >  &r,
0008         const Matrix< Scalar, Dynamic, 1>  &w,
0009         Matrix< Scalar, Dynamic, 1>  &b,
0010         Scalar alpha)
0011 {
0012     typedef DenseIndex Index;
0013
0014     const Index n = r.cols();
0015     eigen_assert(r.rows()>=n);
0016     std::vector<JacobiRotation<Scalar> > givens(n);
0017
0018     /* Local variables */
0019     Scalar temp, rowj;
0020
0021     /* Function Body */
0022     for (Index j = 0; j < n; ++j) {
0023         rowj = w[j];
0024
0025         /* apply the previous transformations to */
0026         /* r(i,j), i=0,1,...,j-1, and to w(j). */
0027         for (Index i = 0; i < j; ++i) {
0028             temp = givens[i].c() * r(i,j) + givens[i].s() * rowj;
0029             rowj = -givens[i].s() * r(i,j) + givens[i].c() * rowj;
0030             r(i,j) = temp;
0031         }
0032
0033         /* determine a givens rotation which eliminates w(j). */
0034         givens[j].makeGivens(-r(j,j), rowj);
0035
0036         if (rowj == 0.)
0037             continue; // givens[j] is identity
0038
0039         /* apply the current transformation to r(j,j), b(j), and alpha. */
0040         r(j,j) = givens[j].c() * r(j,j) + givens[j].s() * rowj;
0041         temp = givens[j].c() * b[j] + givens[j].s() * alpha;
0042         alpha = -givens[j].s() * b[j] + givens[j].c() * alpha;
0043         b[j] = temp;
0044     }
0045 }
0046
0047 } // end namespace internal
0048
0049 } // end namespace Eigen