| 37 | |
| 38 | template <class T> |
| 39 | void GivensR(T* S_, const size_t dim[2], size_t m, T a, T b) { |
| 40 | T r = sqrt(a * a + b * b); |
| 41 | if (fabs(r) < 1e-7) |
| 42 | return; |
| 43 | T c = a / r; |
| 44 | T s = -b / r; |
| 45 | |
| 46 | for (size_t i = 0; i < dim[0]; i++) { |
| 47 | T S0 = S(i, m + 0); |
| 48 | T S1 = S(i, m + 1); |
| 49 | S(i, m) += S0 * (c - 1); |
| 50 | S(i, m) += S1 * (-s); |
| 51 | |
| 52 | S(i, m + 1) += S0 * (s); |
| 53 | S(i, m + 1) += S1 * (c - 1); |
| 54 | } |
| 55 | } |
| 56 | |
| 57 | template <class T> |
| 58 | void SVD(const size_t dim[2], T* U_, T* S_, T* V_, T eps = -1) { |