From 9f1edfe36399c1389a047e7038000c4195e823f7 Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Mon, 19 Feb 2024 16:11:47 -0700 Subject: [PATCH] Add tests. --- .../sparse_matrix_vector_05.cc | 119 +++++++++++++++++ .../sparse_matrix_vector_05.output | 2 + .../sparse_matrix_vector_06.cc | 109 ++++++++++++++++ .../sparse_matrix_vector_06.output | 2 + .../sparse_matrix_vector_07.cc | 123 ++++++++++++++++++ .../sparse_matrix_vector_07.output | 2 + 6 files changed, 357 insertions(+) create mode 100644 tests/trilinos_tpetra/sparse_matrix_vector_05.cc create mode 100644 tests/trilinos_tpetra/sparse_matrix_vector_05.output create mode 100644 tests/trilinos_tpetra/sparse_matrix_vector_06.cc create mode 100644 tests/trilinos_tpetra/sparse_matrix_vector_06.output create mode 100644 tests/trilinos_tpetra/sparse_matrix_vector_07.cc create mode 100644 tests/trilinos_tpetra/sparse_matrix_vector_07.output diff --git a/tests/trilinos_tpetra/sparse_matrix_vector_05.cc b/tests/trilinos_tpetra/sparse_matrix_vector_05.cc new file mode 100644 index 0000000000..af2185e2da --- /dev/null +++ b/tests/trilinos_tpetra/sparse_matrix_vector_05.cc @@ -0,0 +1,119 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2004 - 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.md at +// the top level directory of deal.II. +// +// --------------------------------------------------------------------- + + + +// check SparseMatrix::matrix_scalar_product + +#include + +#include +#include + +#include +#include + +#include "../tests.h" + + +void +test(LinearAlgebra::TpetraWrappers::Vector &v, + LinearAlgebra::TpetraWrappers::Vector &w) +{ + LinearAlgebra::TpetraWrappers::SparseMatrix m(v.size(), + v.size(), + v.size()); + for (unsigned int i = 0; i < m.m(); ++i) + for (unsigned int j = 0; j < m.m(); ++j) + m.set(i, j, i + 2 * j); + + for (unsigned int i = 0; i < v.size(); ++i) + { + v(i) = i; + w(i) = i + 1; + } + + m.compress(VectorOperation::insert); + v.compress(VectorOperation::insert); + w.compress(VectorOperation::insert); + + // + const TrilinosScalar s = m.matrix_scalar_product(w, v); + + // make sure we get the expected result + for (unsigned int i = 0; i < v.size(); ++i) + { + AssertThrow(v(i) == i, ExcInternalError()); + AssertThrow(w(i) == i + 1, ExcInternalError()); + } + + TrilinosScalar result = 0; + for (unsigned int i = 0; i < m.m(); ++i) + for (unsigned int j = 0; j < m.m(); ++j) + result += (i + 2 * j) * j * (i + 1); + + AssertThrow(s == result, ExcInternalError()); + + deallog << "OK" << std::endl; +} + + + +int +main(int argc, char **argv) +{ + initlog(); + + Utilities::MPI::MPI_InitFinalize mpi_initialization( + argc, argv, testing_max_num_threads()); + + + try + { + { + LinearAlgebra::TpetraWrappers::Vector v; + v.reinit(complete_index_set(100), MPI_COMM_WORLD); + LinearAlgebra::TpetraWrappers::Vector w; + w.reinit(complete_index_set(100), MPI_COMM_WORLD); + test(v, w); + } + } + catch (const std::exception &exc) + { + std::cerr << std::endl + << std::endl + << "----------------------------------------------------" + << std::endl; + std::cerr << "Exception on processing: " << std::endl + << exc.what() << std::endl + << "Aborting!" << std::endl + << "----------------------------------------------------" + << std::endl; + + return 1; + } + catch (...) + { + std::cerr << std::endl + << std::endl + << "----------------------------------------------------" + << std::endl; + std::cerr << "Unknown exception!" << std::endl + << "Aborting!" << std::endl + << "----------------------------------------------------" + << std::endl; + return 1; + }; +} diff --git a/tests/trilinos_tpetra/sparse_matrix_vector_05.output b/tests/trilinos_tpetra/sparse_matrix_vector_05.output new file mode 100644 index 0000000000..0fd8fc12f0 --- /dev/null +++ b/tests/trilinos_tpetra/sparse_matrix_vector_05.output @@ -0,0 +1,2 @@ + +DEAL::OK diff --git a/tests/trilinos_tpetra/sparse_matrix_vector_06.cc b/tests/trilinos_tpetra/sparse_matrix_vector_06.cc new file mode 100644 index 0000000000..d6dc7519ee --- /dev/null +++ b/tests/trilinos_tpetra/sparse_matrix_vector_06.cc @@ -0,0 +1,109 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2004 - 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.md at +// the top level directory of deal.II. +// +// --------------------------------------------------------------------- + + + +// check SparseMatrix::matrix_norm_square + +#include + +#include +#include + +#include +#include + +#include "../tests.h" + + +void +test(LinearAlgebra::TpetraWrappers::Vector &v) +{ + LinearAlgebra::TpetraWrappers::SparseMatrix m(v.size(), + v.size(), + v.size()); + for (unsigned int i = 0; i < m.m(); ++i) + for (unsigned int j = 0; j < m.m(); ++j) + m.set(i, j, i + 2 * j); + + for (unsigned int i = 0; i < v.size(); ++i) + v(i) = i; + + m.compress(VectorOperation::insert); + v.compress(VectorOperation::insert); + + // + const TrilinosScalar s = m.matrix_norm_square(v); + + // make sure we get the expected result + for (unsigned int i = 0; i < v.size(); ++i) + AssertThrow(v(i) == i, ExcInternalError()); + + TrilinosScalar result = 0; + for (unsigned int i = 0; i < m.m(); ++i) + for (unsigned int j = 0; j < m.m(); ++j) + result += (i + 2 * j) * j * i; + + AssertThrow(s == result, ExcInternalError()); + + deallog << "OK" << std::endl; +} + + + +int +main(int argc, char **argv) +{ + initlog(); + + Utilities::MPI::MPI_InitFinalize mpi_initialization( + argc, argv, testing_max_num_threads()); + + + try + { + { + LinearAlgebra::TpetraWrappers::Vector v; + v.reinit(complete_index_set(30), MPI_COMM_WORLD); + test(v); + } + } + catch (const std::exception &exc) + { + std::cerr << std::endl + << std::endl + << "----------------------------------------------------" + << std::endl; + std::cerr << "Exception on processing: " << std::endl + << exc.what() << std::endl + << "Aborting!" << std::endl + << "----------------------------------------------------" + << std::endl; + + return 1; + } + catch (...) + { + std::cerr << std::endl + << std::endl + << "----------------------------------------------------" + << std::endl; + std::cerr << "Unknown exception!" << std::endl + << "Aborting!" << std::endl + << "----------------------------------------------------" + << std::endl; + return 1; + }; +} diff --git a/tests/trilinos_tpetra/sparse_matrix_vector_06.output b/tests/trilinos_tpetra/sparse_matrix_vector_06.output new file mode 100644 index 0000000000..0fd8fc12f0 --- /dev/null +++ b/tests/trilinos_tpetra/sparse_matrix_vector_06.output @@ -0,0 +1,2 @@ + +DEAL::OK diff --git a/tests/trilinos_tpetra/sparse_matrix_vector_07.cc b/tests/trilinos_tpetra/sparse_matrix_vector_07.cc new file mode 100644 index 0000000000..93cd24099f --- /dev/null +++ b/tests/trilinos_tpetra/sparse_matrix_vector_07.cc @@ -0,0 +1,123 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2004 - 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.md at +// the top level directory of deal.II. +// +// --------------------------------------------------------------------- + + + +// check SparseMatrix::matrix_norm_square + +#include + +#include +#include + +#include +#include + +#include "../tests.h" + + +void +test(LinearAlgebra::TpetraWrappers::Vector &v, + LinearAlgebra::TpetraWrappers::Vector &w, + LinearAlgebra::TpetraWrappers::Vector &x) +{ + LinearAlgebra::TpetraWrappers::SparseMatrix m(v.size(), + v.size(), + v.size()); + for (unsigned int i = 0; i < m.m(); ++i) + for (unsigned int j = 0; j < m.m(); ++j) + m.set(i, j, i + 2 * j); + + for (unsigned int i = 0; i < v.size(); ++i) + { + v(i) = i; + w(i) = i + 1; + } + + m.compress(VectorOperation::insert); + v.compress(VectorOperation::insert); + w.compress(VectorOperation::insert); + + // x=w-Mv + const double s = m.residual(x, v, w); + + // make sure we get the expected result + for (unsigned int i = 0; i < v.size(); ++i) + { + AssertThrow(v(i) == i, ExcInternalError()); + AssertThrow(w(i) == i + 1, ExcInternalError()); + + double result = i + 1; + for (unsigned int j = 0; j < m.m(); ++j) + result -= (i + 2 * j) * j; + + AssertThrow(x(i) == result, ExcInternalError()); + } + + AssertThrow(s == x.l2_norm(), ExcInternalError()); + + deallog << "OK" << std::endl; +} + + + +int +main(int argc, char **argv) +{ + initlog(); + + Utilities::MPI::MPI_InitFinalize mpi_initialization( + argc, argv, testing_max_num_threads()); + + + try + { + { + LinearAlgebra::TpetraWrappers::Vector v; + v.reinit(complete_index_set(100), MPI_COMM_WORLD); + LinearAlgebra::TpetraWrappers::Vector w; + w.reinit(complete_index_set(100), MPI_COMM_WORLD); + LinearAlgebra::TpetraWrappers::Vector x; + x.reinit(complete_index_set(100), MPI_COMM_WORLD); + test(v, w, x); + } + } + catch (const std::exception &exc) + { + std::cerr << std::endl + << std::endl + << "----------------------------------------------------" + << std::endl; + std::cerr << "Exception on processing: " << std::endl + << exc.what() << std::endl + << "Aborting!" << std::endl + << "----------------------------------------------------" + << std::endl; + + return 1; + } + catch (...) + { + std::cerr << std::endl + << std::endl + << "----------------------------------------------------" + << std::endl; + std::cerr << "Unknown exception!" << std::endl + << "Aborting!" << std::endl + << "----------------------------------------------------" + << std::endl; + return 1; + }; +} diff --git a/tests/trilinos_tpetra/sparse_matrix_vector_07.output b/tests/trilinos_tpetra/sparse_matrix_vector_07.output new file mode 100644 index 0000000000..0fd8fc12f0 --- /dev/null +++ b/tests/trilinos_tpetra/sparse_matrix_vector_07.output @@ -0,0 +1,2 @@ + +DEAL::OK -- 2.39.5