]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Avoid the use of ArrayView for Tensors in MPI::sum(). 12685/head
authorWolfgang Bangerth <bangerth@colostate.edu>
Fri, 20 Aug 2021 02:19:55 +0000 (20:19 -0600)
committerWolfgang Bangerth <bangerth@colostate.edu>
Fri, 20 Aug 2021 19:27:13 +0000 (13:27 -0600)
include/deal.II/base/mpi.templates.h

index 6c9e2eb178caff261d0e0cfa1e5e3471f3d0961f..fe5d640253e3ac01519baf8794792398155cffb4 100644 (file)
@@ -323,12 +323,20 @@ namespace Utilities
 
     template <int rank, int dim, typename Number>
     Tensor<rank, dim, Number>
-    sum(const Tensor<rank, dim, Number> &local,
-        const MPI_Comm &                 mpi_communicator)
+    sum(const Tensor<rank, dim, Number> &t, const MPI_Comm &mpi_communicator)
     {
-      Tensor<rank, dim, Number> sums;
-      sum(local, mpi_communicator, sums);
-      return sums;
+      // Copy the tensor into a C-style array with which we can then
+      // call the other sum() function.
+      Number array[Tensor<rank, dim, Number>::n_independent_components];
+      for (unsigned int i = 0;
+           i < Tensor<rank, dim, Number>::n_independent_components;
+           ++i)
+        array[i] =
+          t[Tensor<rank, dim, Number>::unrolled_to_component_indices(i)];
+
+      sum(array, mpi_communicator, array);
+
+      return Tensor<rank, dim, Number>(make_array_view(array));
     }
 
 
@@ -338,6 +346,8 @@ namespace Utilities
     sum(const SymmetricTensor<rank, dim, Number> &local,
         const MPI_Comm &                          mpi_communicator)
     {
+      // Copy the tensor into a C-style array with which we can then
+      // call the other sum() function.
       const unsigned int n_entries =
         SymmetricTensor<rank, dim, Number>::n_independent_components;
       Number

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.