From 80fc788ef61141be89828e61508668f290e27151 Mon Sep 17 00:00:00 2001 From: rainmaker6 Date: Wed, 9 Jun 2021 12:45:10 +0530 Subject: [PATCH 1/3] Added partial qr solver --- include/solvers/qr.hpp | 48 ++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 48 insertions(+) create mode 100644 include/solvers/qr.hpp diff --git a/include/solvers/qr.hpp b/include/solvers/qr.hpp new file mode 100644 index 0000000..420e512 --- /dev/null +++ b/include/solvers/qr.hpp @@ -0,0 +1,48 @@ +#pragma once +#include +#include + +namespace sw +{ + namespace hprblas + { + template + auto qr(const universal::blas::matrix &A, + universal::blas::matrix &Q, + universal::blas::matrix &R) + { + using size_type = typename universal::blas::matrix::size_type; + size_type n = num_rows(A), m = num_cols(A); + for (size_type i = 0; i < m; ++i) + { + universal::blas::matrix v; // Store i'th column of A (A[(rows 0->n), i'th column]) + /* + Need a function which return's a matrix consisting + of i'th column of A (similar to slice in matlab) + */ + if (i > 0) + { + /* + Need a function similar to slice of a matrix in matlab + + R[(rows 0->i), i'th column] = Q[(rows 0->n, columns 0->i)]' * A[rows 0->n, + i'th column] v= A[(rows 0->n), ith column] - Q[(rows 0->n), (column 0->k)] + * R[(rows 0->i), i't column] + */ + } + R(i, i) = frobenius_norm(v); + // i'th column of Q = v/ R(i,i) + } + return std::make_pair(Q, R); + } + template + auto qr(const universal::blas::matrix &A) + { + size_type m = num_cols(A), n = num_rows(A); + universal::blas::matrix Q(n, m), R(n, n); + qr(A, Q, R); + return std::make_pair(Q, R); + } + + } // namespace hprblas +} // namespace sw From 890967f6868a0d1c9ae4d0a2c44c12c40a15b865 Mon Sep 17 00:00:00 2001 From: rainmaker6 Date: Thu, 10 Jun 2021 19:15:05 +0530 Subject: [PATCH 2/3] Added slice function --- include/solvers/qr.hpp | 36 +++++++++++++++++++++++------------- 1 file changed, 23 insertions(+), 13 deletions(-) diff --git a/include/solvers/qr.hpp b/include/solvers/qr.hpp index 420e512..fa45878 100644 --- a/include/solvers/qr.hpp +++ b/include/solvers/qr.hpp @@ -15,34 +15,44 @@ namespace sw size_type n = num_rows(A), m = num_cols(A); for (size_type i = 0; i < m; ++i) { - universal::blas::matrix v; // Store i'th column of A (A[(rows 0->n), i'th column]) - /* - Need a function which return's a matrix consisting - of i'th column of A (similar to slice in matlab) - */ + universal::blas::matrix v = help(A, 0, n, q); + if (i > 0) { - /* - Need a function similar to slice of a matrix in matlab - - R[(rows 0->i), i'th column] = Q[(rows 0->n, columns 0->i)]' * A[rows 0->n, - i'th column] v= A[(rows 0->n), ith column] - Q[(rows 0->n), (column 0->k)] - * R[(rows 0->i), i't column] - */ + R(0, i, i) = transpose(Q(0, n, 0, i)) * A(0, n, i, i);//To solve problem + v = A(0, n, k, k) - Q(0, n, 0, i) * R(0, i, i, i); } R(i, i) = frobenius_norm(v); - // i'th column of Q = v/ R(i,i) + Q(0, n, i, i) = v / R(i, i); } return std::make_pair(Q, R); } template auto qr(const universal::blas::matrix &A) { + using size_type = typename universal::blas::matrix::size_type; size_type m = num_cols(A), n = num_rows(A); universal::blas::matrix Q(n, m), R(n, n); qr(A, Q, R); return std::make_pair(Q, R); } + template + universal::blas::matrix help(const universal::blas::matrix &A, Scalar rb, Scalar re, Scalar cb, Scalar ce) + { + using size_type = typename universal::blas::matrix::size_type; + universal::blas::matrix B(re - rb + 1, ce - cb + 1); + size_type p = 1, q = 1; + for (size_type i = rb; i <= re; ++i) + { + for (size_type j = cb; j <= ce; ++j) + { + B(p, q) = A(i, j); + ++q; + } + ++p; + } + return B; + } } // namespace hprblas } // namespace sw From df0d3b7922bff2584b6a8745a70b246b6af33c74 Mon Sep 17 00:00:00 2001 From: rainmaker6 Date: Thu, 17 Jun 2021 09:11:17 +0530 Subject: [PATCH 3/3] fix --- include/solvers/qr.hpp | 7 ++++--- 1 file changed, 4 insertions(+), 3 deletions(-) diff --git a/include/solvers/qr.hpp b/include/solvers/qr.hpp index fa45878..e83f493 100644 --- a/include/solvers/qr.hpp +++ b/include/solvers/qr.hpp @@ -19,11 +19,12 @@ namespace sw if (i > 0) { - R(0, i, i) = transpose(Q(0, n, 0, i)) * A(0, n, i, i);//To solve problem - v = A(0, n, k, k) - Q(0, n, 0, i) * R(0, i, i, i); + help(R, 0, i, i, i) = transpose(help(R, 0, n, 0, i)) * help(A, 0, n, i, i); + //Problem => function should change value of the given matrix + v = help(A, 0, n, i, i) - help(Q, 0, n, 0, i) * help(R, 0, i, i, i); } R(i, i) = frobenius_norm(v); - Q(0, n, i, i) = v / R(i, i); + help(Q, 0, n, i, i) = v / R(i, i); } return std::make_pair(Q, R); }