diff --git a/src/test_suite.cpp b/src/test_suite.cpp index becc2f3..b083671 100644 --- a/src/test_suite.cpp +++ b/src/test_suite.cpp @@ -1,15 +1,28 @@ +/** @file test_suite.cpp + * @brief Test suite for project 2 + * + * This is where all the tests for the project are implemented. + * The tests test the creation of symmetrical tridiagonals, + * finding the largest off-diagonal absolute value of a symmetric matrix, + * and the Jacobi rotation method. + * + * @author Cory Alexander Balaton (coryab) + * @author Janita Ovidie Sandtrøen Willumsen (janitaws) + * @bug No known bugs + */ #include #include #include #include "utils.hpp" #include "matrix.hpp" +#include "jacobi.hpp" void test_create_symmetric_tridiagonal() { double tol = 10e-8; double diff; int N = 6; - double h = 1./(double) (N+1); + double h = 1. / (double) (N+1); double a = -1. / (h*h), d = 2. / (h*h); arma::mat A = create_symmetric_tridiagonal(N, a, d); @@ -19,11 +32,13 @@ void test_create_symmetric_tridiagonal() { arma::eig_sym(eigval, eigvec, A); + // Test eigenvalues for (int i=0; i < N; i++) { diff = eigval(i) - (d + 2.*a*std::cos(((i+1.)*M_PI)/(N+1.))); assert(std::abs(diff) < tol); } + // Test eigenvectors arma::vec v, res, analytic_vec = arma::vec(N); for (int i=0; i < N; i++) { v = eigvec.col(i); @@ -35,7 +50,7 @@ void test_create_symmetric_tridiagonal() { analytic_vec = arma::normalise(analytic_vec); - // flip the sign of the analytic vector if they are different + // Flip the sign of the analytic vector if they are different if (analytic_vec(0)*v(0) < 0.) { analytic_vec *= -1; } @@ -61,10 +76,58 @@ void test_max_off_diag_symmetric() { assert(k == 1 && l == 2); } +void test_jacobi_eigensolver() +{ + double tol = 10e-8; + double diff; + int N = 6; + double h = 1. / (double) (N+1); + double a = -1. / (h*h), d = 2. / (h*h); + + arma::mat A = create_symmetric_tridiagonal(N, a, d); + + arma::vec eigval; + arma::mat eigvec; + + int iters; + bool converged; + jacobi_eigensolver(A, 10e-14, eigval, eigvec, 10000, iters, converged); + + // Test eigenvalues + for (int i=0; i < N; i++) { + diff = eigval(i) - (d + 2.*a*std::cos(((i+1.)*M_PI)/(N+1.))); + assert(std::abs(diff) < tol); + } + + // Test eigenvectors + arma::vec v, res, analytic_vec = arma::vec(N); + for (int i=0; i < N; i++) { + v = eigvec.col(i); + + // Build the analytic vector + for (int j=0; j < N; j++) { + analytic_vec(j) = std::sin(((j+1.)*(i+1.)*M_PI) / (N+1.)); + } + + analytic_vec = arma::normalise(analytic_vec); + + // Flip the sign of the analytic vector if they are different + if (analytic_vec(0)*v(0) < 0.) { + analytic_vec *= -1; + } + + res = v - analytic_vec; + for (int j=0; j < N; j++) { + assert(std::abs(res(j)) < tol); + } + } +} + int main() { test_create_symmetric_tridiagonal(); test_max_off_diag_symmetric(); + test_jacobi_eigensolver(); return 0; }