From f1f9065aa3aee80ed59e8bb1c5eb40364cac1084 Mon Sep 17 00:00:00 2001 From: Daniel Arndt Date: Sun, 7 Oct 2018 04:37:41 +0200 Subject: [PATCH] Adress Wolfgang's comments a.k.a add more documentation and switches --- include/deal.II/base/utilities.h | 8 ++++ .../deal.II/lac/read_write_vector.templates.h | 45 ++++++++++--------- include/deal.II/lac/trilinos_sparse_matrix.h | 16 +++++++ source/base/utilities.cc | 2 + 4 files changed, 50 insertions(+), 21 deletions(-) diff --git a/include/deal.II/base/utilities.h b/include/deal.II/base/utilities.h index 71d05faf24..dc7414587e 100644 --- a/include/deal.II/base/utilities.h +++ b/include/deal.II/base/utilities.h @@ -888,6 +888,14 @@ namespace Utilities const Epetra_Comm & comm_self(); + /** + * Return a Teuchos::Comm object needed for creation of Tpetra::Maps. + * + * If deal.II has been configured to use a compiler that does not support + * MPI then the resulting communicator will be a serial one. Otherwise, + * the communicator will correspond to MPI_COMM_SELF, i.e. a communicator + * that comprises only this one processor. + */ const Teuchos::RCP> & tpetra_comm_self(); diff --git a/include/deal.II/lac/read_write_vector.templates.h b/include/deal.II/lac/read_write_vector.templates.h index 87139aedd7..9da3c32e68 100644 --- a/include/deal.II/lac/read_write_vector.templates.h +++ b/include/deal.II/lac/read_write_vector.templates.h @@ -565,30 +565,33 @@ namespace LinearAlgebra Assert(size == 0 || values != nullptr, ExcInternalError("Export failed.")); AssertDimension(size, stored_elements.n_elements()); - if (operation == VectorOperation::insert) - { - for (size_type i = 0; i < size; ++i) - values[i] = new_values[i]; - } - else if (operation == VectorOperation::add) - { - for (size_type i = 0; i < size; ++i) - values[i] += new_values[i]; - } - else if (operation == VectorOperation::min) + switch (operation) { - for (size_type i = 0; i < size; ++i) - if (std::real(new_values[i]) - std::real(values[i]) < 0.0) - values[i] = new_values[i]; - } - else if (operation == VectorOperation::max) - { - for (size_type i = 0; i < size; ++i) - if (std::real(new_values[i]) - std::real(values[i]) > 0.0) + case VectorOperation::insert: + for (size_type i = 0; i < size; ++i) values[i] = new_values[i]; + break; + + case VectorOperation::add: + for (size_type i = 0; i < size; ++i) + values[i] += new_values[i]; + break; + + case VectorOperation::min: + for (size_type i = 0; i < size; ++i) + if (std::real(new_values[i]) - std::real(values[i]) < 0.0) + values[i] = new_values[i]; + break; + + case VectorOperation::max: + for (size_type i = 0; i < size; ++i) + if (std::real(new_values[i]) - std::real(values[i]) > 0.0) + values[i] = new_values[i]; + break; + + default: + AssertThrow(false, ExcNotImplemented()); } - else - AssertThrow(false, ExcNotImplemented()); } diff --git a/include/deal.II/lac/trilinos_sparse_matrix.h b/include/deal.II/lac/trilinos_sparse_matrix.h index 2ae42c5ef1..368fc70a12 100644 --- a/include/deal.II/lac/trilinos_sparse_matrix.h +++ b/include/deal.II/lac/trilinos_sparse_matrix.h @@ -1576,6 +1576,9 @@ namespace TrilinosWrappers * initialized with the same IndexSet that was used for the column indices * of the matrix. * + * This function will be called when the underlying number type for the + * matrix object and the one one for the vector object are the same. + * * In case of a serial vector, this function will only work when * running on one processor, since the matrix object is inherently * distributed. Otherwise, an exception will be thrown. @@ -1585,6 +1588,11 @@ namespace TrilinosWrappers TrilinosScalar>::value>::type vmult(VectorType &dst, const VectorType &src) const; + /** + * Same as the function above for the case that the underlying number type + * for the matrix object and the one for the vector object do not coincide. + * This case is not implemented. Calling it will result in an runtime error. + */ template typename std::enable_if::value>::type @@ -1599,12 +1607,20 @@ namespace TrilinosWrappers * * This function can be called with several types of vector objects, * see the discussion about @p VectorType in vmult(). + * + * This function will be called when the underlying number type for the + * matrix object and the one one for the vector object are the same. */ template typename std::enable_if::value>::type Tvmult(VectorType &dst, const VectorType &src) const; + /** + * Same as the function above for the case that the underlying number type + * for the matrix object and the one for the vector object do not coincide. + * This case is not implemented. Calling it will result in an runtime error. + */ template typename std::enable_if::value>::type diff --git a/source/base/utilities.cc b/source/base/utilities.cc index ff4d0d9573..d5f120192b 100644 --- a/source/base/utilities.cc +++ b/source/base/utilities.cc @@ -993,6 +993,8 @@ namespace Utilities return communicator; } + + const Epetra_Comm & comm_self() { -- 2.39.5