diff --git a/include/matrix.hpp b/include/matrix.hpp index 3ed9290..a44b20b 100644 --- a/include/matrix.hpp +++ b/include/matrix.hpp @@ -13,4 +13,6 @@ arma::mat create_tridiagonal(int n, double a, double d, double e); arma::mat create_symmetric_tridiagonal(int n, double a, double d); +double max_offdiag_symmetric(arma::mat& A, int& k, int& l); + #endif diff --git a/src/matrix.cpp b/src/matrix.cpp index f1099b7..53734d8 100644 --- a/src/matrix.cpp +++ b/src/matrix.cpp @@ -36,3 +36,20 @@ arma::mat create_symmetric_tridiagonal(int n, double a, double d) { return create_tridiagonal(n, a, d, a); } + +double max_offdiag_symmetric(arma::mat& A, int& k, int& l) +{ + int m = A.n_rows; + int n = A.n_cols; + double item, res = 0.; + for (int i=0; i < m; i++) { + for (int j=i+1; j < n; j++) { + if ((item = std::abs(A(i,j))) > res) { + res = item; + k = i; + l = j; + } + } + } + return res; +} diff --git a/src/test_suite.cpp b/src/test_suite.cpp index dc13845..becc2f3 100644 --- a/src/test_suite.cpp +++ b/src/test_suite.cpp @@ -48,12 +48,23 @@ void test_create_symmetric_tridiagonal() { } void test_max_off_diag_symmetric() { + arma::mat A = arma::mat(4,4,arma::fill::eye); + double tol = 10e-8; + A(0, 3) = .5; + A(3, 0) = .5; + A(1, 2) = -.7; + A(2, 1) = -.7; + + int k,l; + assert(max_offdiag_symmetric(A, k, l) - .7 < tol); + assert(k == 1 && l == 2); } int main() { test_create_symmetric_tridiagonal(); + test_max_off_diag_symmetric(); return 0; }