From d21c2dd440e9672771bcf8ec9d79b064e6d681ce Mon Sep 17 00:00:00 2001 From: Denis Davydov Date: Thu, 17 Sep 2015 22:25:37 +0200 Subject: [PATCH] PArpack: extend unit tests to check that we get correct eigenvectors --- tests/arpack/step-36_parpack.cc | 20 ++++++++++++++++++++ tests/arpack/step-36_parpack_trilinos.cc | 20 ++++++++++++++++++++ 2 files changed, 40 insertions(+) diff --git a/tests/arpack/step-36_parpack.cc b/tests/arpack/step-36_parpack.cc index 98387375eb..cb8ec9ab59 100644 --- a/tests/arpack/step-36_parpack.cc +++ b/tests/arpack/step-36_parpack.cc @@ -296,6 +296,26 @@ void test () for (unsigned int i=0; i < eigenvalues.size(); i++) dealii::deallog << eigenvalues[i] << std::endl; + // make sure that we have eigenvectors and they are mass-orthonormal: + // a) (A*x_i-\lambda*B*x_i).L2() == 0 + // b) x_i*B*y_i=\delta_{ij} + { + const double precision = 1e-7; + PETScWrappers::MPI::Vector Ax(eigenfunctions[0]), Bx(eigenfunctions[0]); + for (unsigned int i=0; i < eigenfunctions.size(); ++i) + { + mass_matrix.vmult(Bx,eigenfunctions[i]); + + for (unsigned int j=0; j < eigenfunctions.size(); j++) + Assert( std::abs( eigenfunctions[j] * Bx - (i==j))< precision, + ExcMessage(std::to_string(eigenfunctions[j] * Bx))); + + stiffness_matrix.vmult(Ax,eigenfunctions[i]); + Ax.add(-1.0*std::real(lambda[i]),Bx); + Assert (Ax.l2_norm() < precision, + ExcMessage(std::to_string(Ax.l2_norm()))); + } + } } diff --git a/tests/arpack/step-36_parpack_trilinos.cc b/tests/arpack/step-36_parpack_trilinos.cc index 17d18057f2..31396a8b28 100644 --- a/tests/arpack/step-36_parpack_trilinos.cc +++ b/tests/arpack/step-36_parpack_trilinos.cc @@ -290,6 +290,26 @@ void test () for (unsigned int i=0; i < eigenvalues.size(); i++) deallog << eigenvalues[i] << std::endl; + // make sure that we have eigenvectors and they are mass-orthonormal: + // a) (A*x_i-\lambda*B*x_i).L2() == 0 + // b) x_i*B*y_i=\delta_{ij} + { + const double precision = 1e-7; + TrilinosWrappers::MPI::Vector Ax(eigenfunctions[0]), Bx(eigenfunctions[0]); + for (unsigned int i=0; i < eigenfunctions.size(); ++i) + { + mass_matrix.vmult(Bx,eigenfunctions[i]); + + for (unsigned int j=0; j < eigenfunctions.size(); j++) + Assert( std::abs( eigenfunctions[j] * Bx - (i==j))< precision, + ExcMessage(std::to_string(eigenfunctions[j] * Bx))); + + stiffness_matrix.vmult(Ax,eigenfunctions[i]); + Ax.add(-1.0*std::real(lambda[i]),Bx); + Assert (Ax.l2_norm() < precision, + ExcMessage(std::to_string(Ax.l2_norm()))); + } + } } -- 2.39.5