From: Doug Shi-Dong Date: Fri, 24 Apr 2020 17:12:52 +0000 (-0400) Subject: Allow LinOp(TrilinosMatrix). X-Git-Tag: v9.2.0-rc1~179^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=de1035509b3dc2b7266af6c4b7423ff969370145;p=dealii.git Allow LinOp(TrilinosMatrix). --- diff --git a/include/deal.II/lac/la_parallel_vector.h b/include/deal.II/lac/la_parallel_vector.h index 015f184fc4..7f4374099a 100644 --- a/include/deal.II/lac/la_parallel_vector.h +++ b/include/deal.II/lac/la_parallel_vector.h @@ -1759,26 +1759,137 @@ namespace internal class ReinitHelper> { public: - template + // A helper type-trait that leverage SFINAE to figure out if type T has + // void T::get_mpi_communicator() + template + struct has_get_mpi_communicator + { + private: + static bool + detect(...); + + template + static decltype(std::declval().get_mpi_communicator()) + detect(const U &); + + public: + static const bool value = + !std::is_same()))>::value; + }; + + // A helper type-trait that leverage SFINAE to figure out if type T has + // void T::locally_owned_domain_indices() + template + struct has_locally_owned_domain_indices + { + private: + static bool + detect(...); + + template + static decltype(std::declval().locally_owned_domain_indices()) + detect(const U &); + + public: + static const bool value = + !std::is_same()))>::value; + }; + + // A helper type-trait that leverage SFINAE to figure out if type T has + // void T::locally_owned_range_indices() + template + struct has_locally_owned_range_indices + { + private: + static bool + detect(...); + + template + static decltype(std::declval().locally_owned_range_indices()) + detect(const U &); + + public: + static const bool value = + !std::is_same()))>::value; + }; + + // A helper type-trait that leverage SFINAE to figure out if type T has + // void T::initialize_dof_vector(VectorType v) + template + struct has_initialize_dof_vector + { + private: + static bool + detect(...); + + template + static decltype(std::declval().initialize_dof_vector( + LinearAlgebra::distributed::Vector())) + detect(const U &); + + public: + static const bool value = + !std::is_same()))>::value; + }; + + // Used for (Trilinos/PETSc)Wrappers::SparseMatrix + template ::value && + has_locally_owned_domain_indices::value, + MatrixType>::type * = nullptr> static void - reinit_range_vector(const Matrix & matrix, - LinearAlgebra::distributed::Vector &v, - bool omit_zeroing_entries) + reinit_domain_vector(MatrixType & mat, + LinearAlgebra::distributed::Vector &vec, + bool /*omit_zeroing_entries*/) { - matrix.initialize_dof_vector(v); - if (!omit_zeroing_entries) - v = Number(); + vec.reinit(mat.locally_owned_domain_indices(), + mat.get_mpi_communicator()); } - template + // Used for MatrixFree and DiagonalMatrix + template < + typename MatrixType, + typename std::enable_if::value, + MatrixType>::type * = nullptr> static void - reinit_domain_vector(const Matrix & matrix, - LinearAlgebra::distributed::Vector &v, + reinit_domain_vector(MatrixType & mat, + LinearAlgebra::distributed::Vector &vec, bool omit_zeroing_entries) { - matrix.initialize_dof_vector(v); + mat.initialize_dof_vector(vec); + if (!omit_zeroing_entries) + vec = Number(); + } + + // Used for (Trilinos/PETSc)Wrappers::SparseMatrix + template ::value && + has_locally_owned_range_indices::value, + MatrixType>::type * = nullptr> + static void + reinit_range_vector(MatrixType & mat, + LinearAlgebra::distributed::Vector &vec, + bool /*omit_zeroing_entries*/) + { + vec.reinit(mat.locally_owned_range_indices(), + mat.get_mpi_communicator()); + } + + // Used for MatrixFree and DiagonalMatrix + template < + typename MatrixType, + typename std::enable_if::value, + MatrixType>::type * = nullptr> + static void + reinit_range_vector(MatrixType & mat, + LinearAlgebra::distributed::Vector &vec, + bool omit_zeroing_entries) + { + mat.initialize_dof_vector(vec); if (!omit_zeroing_entries) - v = Number(); + vec = Number(); } }; diff --git a/tests/lac/linear_operator_18.cc b/tests/lac/linear_operator_18.cc new file mode 100644 index 0000000000..00a962feb0 --- /dev/null +++ b/tests/lac/linear_operator_18.cc @@ -0,0 +1,53 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2016 - 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. +// +// --------------------------------------------------------------------- + +// Test that it is possible to use a PackagedOperation created by +// operator*() of a LinearOperator object for Trilinos matrices and vectors + +#include +#include +#include +#include + +#include "../tests.h" + + +int +main(int argc, char *argv[]) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + + initlog(); + deallog << std::setprecision(10); + + typedef LinearAlgebra::distributed::Vector vector_t; + typedef TrilinosWrappers::SparseMatrix matrix_t; + + matrix_t a(5U, 4U, 3U); + a.compress(VectorOperation::add); + + auto op_a = linear_operator(a); + vector_t u, res; + op_a.reinit_domain_vector(u, false); + res = op_a * u; + + auto op_a_transpose = transpose_operator(op_a); + op_a.reinit_range_vector(u, false); + res = op_a_transpose * u; + + deallog << "OK" << std::endl; + + return 0; +} diff --git a/tests/lac/linear_operator_18.with_trilinos=on.output b/tests/lac/linear_operator_18.with_trilinos=on.output new file mode 100644 index 0000000000..0fd8fc12f0 --- /dev/null +++ b/tests/lac/linear_operator_18.with_trilinos=on.output @@ -0,0 +1,2 @@ + +DEAL::OK