From: Benjamin Brands Date: Mon, 19 Mar 2018 15:07:34 +0000 (+0100) Subject: tests for the MRRR eigensolvers X-Git-Tag: v9.0.0-rc1~202^2~1 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=0e1270fe685a37fffa68f8ad4b7fc4613699c8d2;p=dealii.git tests for the MRRR eigensolvers --- diff --git a/tests/scalapack/scalapack_15.cc b/tests/scalapack/scalapack_15.cc new file mode 100644 index 0000000000..3b1ef5ceb5 --- /dev/null +++ b/tests/scalapack/scalapack_15.cc @@ -0,0 +1,158 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2017 - 2018 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + +#include "../tests.h" +#include "../lapack/create_matrix.h" + +// test eigenpairs_symmetric_by_index_MRRR(const std::pair &, const bool) +// for all eigenvalues with eigenvectors + +#include +#include +#include +#include +#include +#include + +#include +#include +#include +#include + +#include +#include +#include +#include + + +template +void test(const unsigned int size, const unsigned int block_size, const NumberType tol) +{ + MPI_Comm mpi_communicator(MPI_COMM_WORLD); + const unsigned int n_mpi_processes(Utilities::MPI::n_mpi_processes(mpi_communicator)); + const unsigned int this_mpi_process(Utilities::MPI::this_mpi_process(mpi_communicator)); + + ConditionalOStream pcout (std::cout, (this_mpi_process ==0)); + + std::shared_ptr grid = std::make_shared(mpi_communicator,size,size,block_size,block_size); + + pcout << size << " " << block_size << " " << grid->get_process_grid_rows() << " " << grid->get_process_grid_columns() << std::endl; + + const unsigned int n_eigenvalues = size; + const unsigned int max_n_eigenvalues = 5; + + // Create SPD matrices of requested size: + FullMatrix full_A(size); + std::vector eigenvalues_Lapack(size); + std::vector> s_eigenvectors_ (max_n_eigenvalues,Vector(size)); + std::vector> p_eigenvectors_ (max_n_eigenvalues,Vector(size)); + FullMatrix p_eigenvectors (size,size); + + ScaLAPACKMatrix scalapack_syevr (size, grid, block_size); + scalapack_syevr.set_property(LAPACKSupport::Property::symmetric); + + create_spd (full_A); + scalapack_syevr = full_A; + + //Lapack as reference + { + std::vector lapack_A(size*size); + for (unsigned int i = 0; i < size; ++i) + for (unsigned int j = 0; j < size; ++j) + lapack_A[i*size+j] = full_A(i,j); + + int info; //Variable containing information about the successful exit of the lapack routine + char jobz = 'V'; //'V': all eigenpairs of A are computed + char uplo = 'U'; //storage format of the matrix A; not so important as matrix is symmetric + char range = 'I'; //the il-th through iu-th eigenvalues will be found + int LDA = size; //leading dimension of the matrix A + NumberType vl=0, vu=0; + int il=1, iu=size; + char sign = 'S'; + NumberType abstol; + lamch(&sign,abstol); + abstol *= 2; + int m=0; + std::vector eigenvectors(size*size); + int LDZ=size; + std::vector isuppz(2*size); + int lwork=-1; //length of vector/array work + std::vector work(1); + int liwork=-1; //length of vector/array iwork + std::vector iwork(1); + //by setting lwork to -1 a workspace query for work is done + //as matrix is symmetric: LDA == size of matrix + syevr(&jobz, &range, &uplo, &LDA, lapack_A.data(), &LDA, &vl, &vu, &il, &iu, + &abstol, &m, eigenvalues_Lapack.data(), eigenvectors.data(), &LDZ, isuppz.data(), + work.data(), &lwork, iwork.data(), &liwork, &info); + + lwork=work[0]; + work.resize (lwork); + liwork=iwork[0]; + iwork.resize(liwork); + + syevr(&jobz, &range, &uplo, &LDA, lapack_A.data(), &LDA, &vl, &vu, &il, &iu, + &abstol, &m, eigenvalues_Lapack.data(), eigenvectors.data(), &LDZ, isuppz.data(), + work.data(), &lwork, iwork.data(), &liwork, &info); + + AssertThrow (info==0, LAPACKSupport::ExcErrorCode("syevr", info)); + for (int i=0; i eigenvalues_psyer = scalapack_syevr.eigenpairs_symmetric_by_index_MRRR(std::make_pair(0,size-1),true); + scalapack_syevr.copy_to(p_eigenvectors); + for (unsigned int i=0; i sizes = {{200,400,600}}; + const std::vector blocks = {{32,64}}; + + const double tol = 1e-10; + + for (const auto &s : sizes) + for (const auto &b : blocks) + if (b <= s) + test(s,b,tol); +} diff --git a/tests/scalapack/scalapack_15.mpirun=1.output b/tests/scalapack/scalapack_15.mpirun=1.output new file mode 100644 index 0000000000..b55bfbdf68 --- /dev/null +++ b/tests/scalapack/scalapack_15.mpirun=1.output @@ -0,0 +1,30 @@ +200 32 1 1 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +200 64 1 1 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +400 32 1 1 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +400 64 1 1 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +600 32 1 1 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +600 64 1 1 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + diff --git a/tests/scalapack/scalapack_15.mpirun=10.output b/tests/scalapack/scalapack_15.mpirun=10.output new file mode 100644 index 0000000000..b3ea1071a5 --- /dev/null +++ b/tests/scalapack/scalapack_15.mpirun=10.output @@ -0,0 +1,30 @@ +200 32 3 3 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +200 64 3 3 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +400 32 3 3 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +400 64 3 3 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +600 32 3 3 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +600 64 3 3 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + diff --git a/tests/scalapack/scalapack_15.mpirun=4.output b/tests/scalapack/scalapack_15.mpirun=4.output new file mode 100644 index 0000000000..10ea4f4a8c --- /dev/null +++ b/tests/scalapack/scalapack_15.mpirun=4.output @@ -0,0 +1,30 @@ +200 32 2 2 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +200 64 2 2 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +400 32 2 2 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +400 64 2 2 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +600 32 2 2 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +600 64 2 2 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + diff --git a/tests/scalapack/scalapack_15.mpirun=5.output b/tests/scalapack/scalapack_15.mpirun=5.output new file mode 100644 index 0000000000..10ea4f4a8c --- /dev/null +++ b/tests/scalapack/scalapack_15.mpirun=5.output @@ -0,0 +1,30 @@ +200 32 2 2 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +200 64 2 2 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +400 32 2 2 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +400 64 2 2 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +600 32 2 2 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +600 64 2 2 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + diff --git a/tests/scalapack/scalapack_15_a.cc b/tests/scalapack/scalapack_15_a.cc new file mode 100644 index 0000000000..d37d30af51 --- /dev/null +++ b/tests/scalapack/scalapack_15_a.cc @@ -0,0 +1,139 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2017 - 2018 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + +#include "../tests.h" +#include "../lapack/create_matrix.h" + +// test eigenpairs_symmetric_by_index_MRRR(const std::pair &, const bool) +// for all eigenvalues without eigenvectors + +#include +#include +#include +#include +#include +#include + +#include +#include +#include +#include + +#include +#include +#include +#include + + +template +void test(const unsigned int size, const unsigned int block_size, const NumberType tol) +{ + MPI_Comm mpi_communicator(MPI_COMM_WORLD); + const unsigned int n_mpi_processes(Utilities::MPI::n_mpi_processes(mpi_communicator)); + const unsigned int this_mpi_process(Utilities::MPI::this_mpi_process(mpi_communicator)); + + ConditionalOStream pcout (std::cout, (this_mpi_process ==0)); + + std::shared_ptr grid = std::make_shared(mpi_communicator,size,size,block_size,block_size); + + pcout << size << " " << block_size << " " << grid->get_process_grid_rows() << " " << grid->get_process_grid_columns() << std::endl; + + const unsigned int n_eigenvalues = size; + const unsigned int max_n_eigenvalues = 5; + + // Create SPD matrices of requested size: + FullMatrix full_A(size); + std::vector eigenvalues_Lapack(size); + std::vector> s_eigenvectors_ (max_n_eigenvalues,Vector(size)); + std::vector> p_eigenvectors_ (max_n_eigenvalues,Vector(size)); + FullMatrix p_eigenvectors (size,size); + + ScaLAPACKMatrix scalapack_syevr (size, grid, block_size); + scalapack_syevr.set_property(LAPACKSupport::Property::symmetric); + + create_spd (full_A); + scalapack_syevr = full_A; + + //Lapack as reference + { + std::vector lapack_A(size*size); + for (unsigned int i = 0; i < size; ++i) + for (unsigned int j = 0; j < size; ++j) + lapack_A[i*size+j] = full_A(i,j); + + int info; //Variable containing information about the successful exit of the lapack routine + char jobz = 'N'; //'N': all eigenvalues of A are computed + char uplo = 'U'; //storage format of the matrix A; not so important as matrix is symmetric + char range = 'I'; //the il-th through iu-th eigenvalues will be found + int LDA = size; //leading dimension of the matrix A + NumberType vl=0, vu=0; + int il=1, iu=size; + char sign = 'S'; + NumberType abstol; + lamch(&sign,abstol); + abstol *= 2; + int m=0; + std::vector eigenvectors(size*size); //will not be referenced + int LDZ=size; + std::vector isuppz(2*size); + int lwork=-1; //length of vector/array work + std::vector work(1); + int liwork=-1; //length of vector/array iwork + std::vector iwork(1); + //by setting lwork to -1 a workspace query for work is done + //as matrix is symmetric: LDA == size of matrix + syevr(&jobz, &range, &uplo, &LDA, lapack_A.data(), &LDA, &vl, &vu, &il, &iu, + &abstol, &m, eigenvalues_Lapack.data(), eigenvectors.data(), &LDZ, isuppz.data(), + work.data(), &lwork, iwork.data(), &liwork, &info); + + lwork=work[0]; + work.resize (lwork); + liwork=iwork[0]; + iwork.resize(liwork); + + syevr(&jobz, &range, &uplo, &LDA, lapack_A.data(), &LDA, &vl, &vu, &il, &iu, + &abstol, &m, eigenvalues_Lapack.data(), eigenvectors.data(), &LDZ, isuppz.data(), + work.data(), &lwork, iwork.data(), &liwork, &info); + } + + // the actual test: + + pcout << "comparing " << max_n_eigenvalues << " eigenvalues computed using LAPACK and ScaLAPACK p_syevr:" << std::endl; + const std::vector eigenvalues_psyer = scalapack_syevr.eigenpairs_symmetric_by_index_MRRR(std::make_pair(0,size-1),false); + scalapack_syevr.copy_to(p_eigenvectors); + for (unsigned int i=0; i sizes = {{200,400,600}}; + const std::vector blocks = {{32,64}}; + + const double tol = 1e-10; + + for (const auto &s : sizes) + for (const auto &b : blocks) + if (b <= s) + test(s,b,tol); +} diff --git a/tests/scalapack/scalapack_15_a.mpirun=1.output b/tests/scalapack/scalapack_15_a.mpirun=1.output new file mode 100644 index 0000000000..e2dec6161d --- /dev/null +++ b/tests/scalapack/scalapack_15_a.mpirun=1.output @@ -0,0 +1,18 @@ +200 32 1 1 +comparing 5 eigenvalues computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide +200 64 1 1 +comparing 5 eigenvalues computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide +400 32 1 1 +comparing 5 eigenvalues computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide +400 64 1 1 +comparing 5 eigenvalues computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide +600 32 1 1 +comparing 5 eigenvalues computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide +600 64 1 1 +comparing 5 eigenvalues computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide diff --git a/tests/scalapack/scalapack_15_a.mpirun=10.output b/tests/scalapack/scalapack_15_a.mpirun=10.output new file mode 100644 index 0000000000..624b588864 --- /dev/null +++ b/tests/scalapack/scalapack_15_a.mpirun=10.output @@ -0,0 +1,18 @@ +200 32 3 3 +comparing 5 eigenvalues computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide +200 64 3 3 +comparing 5 eigenvalues computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide +400 32 3 3 +comparing 5 eigenvalues computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide +400 64 3 3 +comparing 5 eigenvalues computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide +600 32 3 3 +comparing 5 eigenvalues computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide +600 64 3 3 +comparing 5 eigenvalues computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide diff --git a/tests/scalapack/scalapack_15_a.mpirun=4.output b/tests/scalapack/scalapack_15_a.mpirun=4.output new file mode 100644 index 0000000000..4881bde55c --- /dev/null +++ b/tests/scalapack/scalapack_15_a.mpirun=4.output @@ -0,0 +1,18 @@ +200 32 2 2 +comparing 5 eigenvalues computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide +200 64 2 2 +comparing 5 eigenvalues computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide +400 32 2 2 +comparing 5 eigenvalues computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide +400 64 2 2 +comparing 5 eigenvalues computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide +600 32 2 2 +comparing 5 eigenvalues computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide +600 64 2 2 +comparing 5 eigenvalues computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide diff --git a/tests/scalapack/scalapack_15_a.mpirun=5.output b/tests/scalapack/scalapack_15_a.mpirun=5.output new file mode 100644 index 0000000000..4881bde55c --- /dev/null +++ b/tests/scalapack/scalapack_15_a.mpirun=5.output @@ -0,0 +1,18 @@ +200 32 2 2 +comparing 5 eigenvalues computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide +200 64 2 2 +comparing 5 eigenvalues computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide +400 32 2 2 +comparing 5 eigenvalues computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide +400 64 2 2 +comparing 5 eigenvalues computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide +600 32 2 2 +comparing 5 eigenvalues computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide +600 64 2 2 +comparing 5 eigenvalues computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide diff --git a/tests/scalapack/scalapack_15_b.cc b/tests/scalapack/scalapack_15_b.cc new file mode 100644 index 0000000000..c8b675db82 --- /dev/null +++ b/tests/scalapack/scalapack_15_b.cc @@ -0,0 +1,158 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2017 - 2018 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + +#include "../tests.h" +#include "../lapack/create_matrix.h" + +// test eigenpairs_symmetric_by_index_MRRR(const std::pair &, const bool) +// for the largest eigenvalues with eigenvectors + +#include +#include +#include +#include +#include +#include + +#include +#include +#include +#include + +#include +#include +#include +#include + + +template +void test(const unsigned int size, const unsigned int block_size, const NumberType tol) +{ + MPI_Comm mpi_communicator(MPI_COMM_WORLD); + const unsigned int n_mpi_processes(Utilities::MPI::n_mpi_processes(mpi_communicator)); + const unsigned int this_mpi_process(Utilities::MPI::this_mpi_process(mpi_communicator)); + + ConditionalOStream pcout (std::cout, (this_mpi_process ==0)); + + std::shared_ptr grid = std::make_shared(mpi_communicator,size,size,block_size,block_size); + + pcout << size << " " << block_size << " " << grid->get_process_grid_rows() << " " << grid->get_process_grid_columns() << std::endl; + + const unsigned int n_eigenvalues = size; + const unsigned int max_n_eigenvalues = 5; + + // Create SPD matrices of requested size: + FullMatrix full_A(size); + std::vector eigenvalues_Lapack(size); + std::vector> s_eigenvectors_ (max_n_eigenvalues,Vector(size)); + std::vector> p_eigenvectors_ (max_n_eigenvalues,Vector(size)); + FullMatrix p_eigenvectors (size,size); + + ScaLAPACKMatrix scalapack_syevr (size, grid, block_size); + scalapack_syevr.set_property(LAPACKSupport::Property::symmetric); + + create_spd (full_A); + scalapack_syevr = full_A; + + //Lapack as reference + { + std::vector lapack_A(size*size); + for (unsigned int i = 0; i < size; ++i) + for (unsigned int j = 0; j < size; ++j) + lapack_A[i*size+j] = full_A(i,j); + + int info; //Variable containing information about the successful exit of the lapack routine + char jobz = 'V'; //'V': eigenpairs of A are computed + char uplo = 'U'; //storage format of the matrix A; not so important as matrix is symmetric + char range = 'I'; //the il-th through iu-th eigenvalues will be found + int LDA = size; //leading dimension of the matrix A + NumberType vl=0, vu=0; + int il=size-max_n_eigenvalues+1, iu=size; + char sign = 'S'; + NumberType abstol; + lamch(&sign,abstol); + abstol *= 2; + int m=0; + std::vector eigenvectors(size*size); + int LDZ=size; + std::vector isuppz(2*size); + int lwork=-1; //length of vector/array work + std::vector work(1); + int liwork=-1; //length of vector/array iwork + std::vector iwork(1); + //by setting lwork to -1 a workspace query for work is done + //as matrix is symmetric: LDA == size of matrix + syevr(&jobz, &range, &uplo, &LDA, lapack_A.data(), &LDA, &vl, &vu, &il, &iu, + &abstol, &m, eigenvalues_Lapack.data(), eigenvectors.data(), &LDZ, isuppz.data(), + work.data(), &lwork, iwork.data(), &liwork, &info); + + lwork=work[0]; + work.resize (lwork); + liwork=iwork[0]; + iwork.resize(liwork); + + syevr(&jobz, &range, &uplo, &LDA, lapack_A.data(), &LDA, &vl, &vu, &il, &iu, + &abstol, &m, eigenvalues_Lapack.data(), eigenvectors.data(), &LDZ, isuppz.data(), + work.data(), &lwork, iwork.data(), &liwork, &info); + + AssertThrow (info==0, LAPACKSupport::ExcErrorCode("syevr", info)); + for (int i=0; i eigenvalues_psyer = scalapack_syevr.eigenpairs_symmetric_by_index_MRRR(std::make_pair(size-max_n_eigenvalues,size-1),true); + scalapack_syevr.copy_to(p_eigenvectors); + for (unsigned int i=0; i sizes = {{200,400,600}}; + const std::vector blocks = {{32,64}}; + + const double tol = 1e-10; + + for (const auto &s : sizes) + for (const auto &b : blocks) + if (b <= s) + test(s,b,tol); +} diff --git a/tests/scalapack/scalapack_15_b.mpirun=1.output b/tests/scalapack/scalapack_15_b.mpirun=1.output new file mode 100644 index 0000000000..b55bfbdf68 --- /dev/null +++ b/tests/scalapack/scalapack_15_b.mpirun=1.output @@ -0,0 +1,30 @@ +200 32 1 1 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +200 64 1 1 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +400 32 1 1 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +400 64 1 1 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +600 32 1 1 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +600 64 1 1 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + diff --git a/tests/scalapack/scalapack_15_b.mpirun=10.output b/tests/scalapack/scalapack_15_b.mpirun=10.output new file mode 100644 index 0000000000..b3ea1071a5 --- /dev/null +++ b/tests/scalapack/scalapack_15_b.mpirun=10.output @@ -0,0 +1,30 @@ +200 32 3 3 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +200 64 3 3 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +400 32 3 3 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +400 64 3 3 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +600 32 3 3 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +600 64 3 3 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + diff --git a/tests/scalapack/scalapack_15_b.mpirun=4.output b/tests/scalapack/scalapack_15_b.mpirun=4.output new file mode 100644 index 0000000000..10ea4f4a8c --- /dev/null +++ b/tests/scalapack/scalapack_15_b.mpirun=4.output @@ -0,0 +1,30 @@ +200 32 2 2 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +200 64 2 2 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +400 32 2 2 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +400 64 2 2 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +600 32 2 2 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +600 64 2 2 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + diff --git a/tests/scalapack/scalapack_15_b.mpirun=5.output b/tests/scalapack/scalapack_15_b.mpirun=5.output new file mode 100644 index 0000000000..10ea4f4a8c --- /dev/null +++ b/tests/scalapack/scalapack_15_b.mpirun=5.output @@ -0,0 +1,30 @@ +200 32 2 2 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +200 64 2 2 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +400 32 2 2 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +400 64 2 2 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +600 32 2 2 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide + +600 64 2 2 +comparing 5 eigenvalues and eigenvectors computed using LAPACK and ScaLAPACK p_syevr: + with respect to the given tolerance the eigenvalues coincide + with respect to the given tolerance also the eigenvectors coincide +