53 template<
class mat_t,
class vec_t>
54 void QR_decomp(
size_t M,
size_t N, mat_t &A, vec_t &tau) {
59 for(
size_t i=0;i<imax;i++) {
70 template<
class mat_t,
class vec_t,
class vec2_t>
72 const mat_t &QR,
const vec_t &tau, vec2_t &v) {
78 for(
size_t i=0;i<imax;i++) {
86 template<
class mat1_t,
class mat2_t,
class mat3_t,
class vec_t>
88 const mat1_t &QR,
const vec_t &tau, mat2_t &Q,
97 if (i==j) O2SCL_IX2(Q,i,j)=1.0;
98 else O2SCL_IX2(Q,i,j)=0.0;
106 for (i=istart; i-- > 0;) {
112 for (i=0; i < M; i++) {
113 for (j=0; j < i && j < N; j++) {
114 O2SCL_IX2(R,i,j)=0.0;
116 for (j=i; j < N; j++) {
117 O2SCL_IX2(R,i,j)=O2SCL_IX2(QR,i,j);
126 template<
class mat_t,
class vec_t,
class vec2_t>
127 void QR_svx(
size_t M,
size_t N,
const mat_t &QR,
const vec_t &tau,
131 o2scl_cblas::o2cblas_Upper,
132 o2scl_cblas::o2cblas_NoTrans,
133 o2scl_cblas::o2cblas_NonUnit,
140 template<
class mat_t,
class vec_t,
class vec2_t,
class vec3_t>
141 void QR_solve(
size_t N,
const mat_t &QR,
const vec_t &tau,
142 const vec2_t &b, vec3_t &x) {
164 template<
class mat1_t,
class mat2_t,
class vec1_t,
class vec2_t>
165 void QR_update(
size_t M,
size_t N, mat1_t &Q, mat2_t &R,
166 vec1_t &w, vec2_t &v) {
182 for (k=(
int)(M-1); k > 0; k--) {
185 double wk=O2SCL_IX(w,k);
186 double wkm1=O2SCL_IX(w,k-1);
197 for (j=0; j < N; j++) {
198 double r0j=O2SCL_IX2(R,0,j);
199 double vj=O2SCL_IX(v,j);
200 O2SCL_IX2(R,0,j)=r0j+w0*vj;
208 for (k=1;k<kmax;k++) {
211 double diag=O2SCL_IX2(R,k-1,k-1);
212 double offdiag=O2SCL_IX2(R,k,k-1);
217 O2SCL_IX2(R,k,k-1)=0.0;
230 template<
class mat_t,
class mat2_t,
class mat3_t>
232 mat_t &A, mat2_t &Q, mat3_t &R) {
void QR_solve(size_t N, const mat_t &QR, const vec_t &tau, const vec2_t &b, vec3_t &x)
Solve the system A x = b using the QR factorization.
void dtrsv(const enum o2cblas_order order, const enum o2cblas_uplo Uplo, const enum o2cblas_transpose TransA, const enum o2cblas_diag Diag, const size_t M, const size_t N, const mat_t &A, vec_t &X)
Compute .
void QR_unpack(const size_t M, const size_t N, const mat1_t &QR, const vec_t &tau, mat2_t &Q, mat3_t &R)
Unpack the QR matrix to the individual Q and R components.
void apply_givens_vec(vec_t &v, size_t i, size_t j, double c, double s)
Apply a rotation to a vector, .
void QR_decomp(size_t M, size_t N, mat_t &A, vec_t &tau)
Compute the QR decomposition of matrix A.
void QR_decomp_unpack(const size_t M, const size_t N, mat_t &A, mat2_t &Q, mat3_t &R)
Compute the unpacked QR decomposition of matrix A.
double householder_transform_subcol(mat_t &A, const size_t ir, const size_t ic, const size_t M)
Compute the Householder transform of a vector formed with n rows of a column of a matrix...
void create_givens(const double a, const double b, double &c, double &s)
Create a Givens rotation matrix.
void vector_copy(const vec_t &src, vec2_t &dest)
Simple vector copy.
void householder_hm_subcol(mat_t &M, const size_t ir, const size_t ic, const size_t nr, const size_t nc, const mat_t &M2, const size_t ir2, const size_t ic2, double tau)
Apply a Householder transformation to the lower-right part of M when the transformation is stored in ...
The namespace for linear algebra classes and functions.
void QR_update(size_t M, size_t N, mat1_t &Q, mat2_t &R, vec1_t &w, vec2_t &v)
Update a QR factorisation for A= Q R, A' = A + u v^T,.
void apply_givens_qr(size_t M, size_t N, mat1_t &Q, mat2_t &R, size_t i, size_t j, double c, double s)
Apply a rotation to matrices from the QR decomposition.
void QR_QTvec(const size_t M, const size_t N, const mat_t &QR, const vec_t &tau, vec2_t &v)
Form the product Q^T v from a QR factorized matrix.
void QR_svx(size_t M, size_t N, const mat_t &QR, const vec_t &tau, vec2_t &x)
Solve the system A x = b in place using the QR factorization.
void householder_hv_subcol(const mat_t &A, vec_t &w, double tau, const size_t ie, const size_t N)
Apply a Householder transformation v to vector w where v is stored as a column in a matrix A...