]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Make deal.II compile with Trilinos but without MPI. 4296/head
authorMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Sat, 22 Apr 2017 13:02:28 +0000 (15:02 +0200)
committerMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Sat, 22 Apr 2017 13:02:28 +0000 (15:02 +0200)
cmake/configure/configure_2_trilinos.cmake
include/deal.II/dofs/dof_accessor.templates.h
include/deal.II/fe/fe_tools_extrapolate.templates.h
include/deal.II/lac/read_write_vector.h
include/deal.II/lac/trilinos_epetra_vector.h
include/deal.II/lac/vector_element_access.h
include/deal.II/multigrid/mg_transfer.h
source/lac/trilinos_sparse_matrix.cc
source/numerics/vector_tools_mean_value.cc

index e00ed235341cdd963d6414afbfe0e1397b443bb5..5665d907f2a2aa56d5e03bbf4d40741e076f916b 100644 (file)
@@ -55,7 +55,7 @@ MACRO(FEATURE_TRILINOS_FIND_EXTERNAL var)
     ENDFOREACH()
 
     IF((TRILINOS_VERSION_MAJOR EQUAL 11 AND
-        NOT (TRILINOS_VERSION_MINOR LESS 14)) 
+        NOT (TRILINOS_VERSION_MINOR LESS 14))
        OR
        (NOT (TRILINOS_VERSION_MAJOR LESS 12)))
         ITEM_MATCHES(_module_found MueLu ${Trilinos_PACKAGE_LIST})
@@ -194,7 +194,9 @@ MACRO(FEATURE_TRILINOS_CONFIGURE_EXTERNAL)
   SET(DEAL_II_EXPAND_TRILINOS_BLOCK_SPARSITY_PATTERN "TrilinosWrappers::BlockSparsityPattern")
   SET(DEAL_II_EXPAND_TRILINOS_MPI_BLOCKVECTOR "TrilinosWrappers::MPI::BlockVector")
   SET(DEAL_II_EXPAND_TRILINOS_MPI_VECTOR "TrilinosWrappers::MPI::Vector")
-  SET(DEAL_II_EXPAND_EPETRA_VECTOR "LinearAlgebra::EpetraWrappers::Vector")
+  IF (TRILINOS_WITH_MPI)
+    SET(DEAL_II_EXPAND_EPETRA_VECTOR "LinearAlgebra::EpetraWrappers::Vector")
+  ENDIF()
 ENDMACRO()
 
 
index fa19f0c392df1d6bf3698d7f5d3cb0d14bb0efb4..5a082e200cfc0631bc0253dd67b2b2b08248f0f3 100644 (file)
@@ -1494,7 +1494,7 @@ namespace internal
 
 
 
-#ifdef DEAL_II_WITH_TRILINOS
+#if defined(DEAL_II_WITH_TRILINOS) && defined(DEAL_II_WITH_MPI)
       static
       std::vector<unsigned int> sort_indices(const types::global_dof_index *v_begin,
                                              const types::global_dof_index *v_end)
index eee5051ef43f36fcf877905a09fcd753d0405ef3..23c79145755ffc2e68b55d48ce894c5050d2cd0c 100644 (file)
@@ -1454,6 +1454,7 @@ namespace FETools
 
 
 
+#ifdef DEAL_II_WITH_MPI
       template <int dim, int spacedim>
       void reinit_distributed(const DoFHandler<dim, spacedim> &dh,
                               LinearAlgebra::EpetraWrappers::Vector &vector)
@@ -1465,6 +1466,7 @@ namespace FETools
         const IndexSet locally_owned_dofs = dh.locally_owned_dofs();
         vector.reinit(locally_owned_dofs, parallel_tria->get_communicator());
       }
+#endif
 #endif //DEAL_II_WITH_TRILINOS
 
       template <int dim, int spacedim, typename Number>
index 4ae7b29fd35b5150a7f444043c218e4b1802c227..631125cfb756e237083e4f08c267b077ac87daba 100644 (file)
@@ -295,6 +295,7 @@ namespace LinearAlgebra
                 std::shared_ptr<const CommunicationPatternBase> communication_pattern =
                   std::shared_ptr<const CommunicationPatternBase> ());
 
+#ifdef DEAL_II_WITH_MPI
     /**
      * Imports all the elements present in the vector's IndexSet from the input
      * vector @p epetra_vec. VectorOperation::values @p operation is used to
@@ -308,6 +309,7 @@ namespace LinearAlgebra
                 std::shared_ptr<const CommunicationPatternBase> communication_pattern =
                   std::shared_ptr<const CommunicationPatternBase> ());
 #endif
+#endif
 
 #ifdef DEAL_II_WITH_CUDA
     /**
index 803cc2c2301c483a28b74e9890bb840ab73d8eea..2be2a661b0d9e18f815e45f4e55dde0c4974a42b 100644 (file)
@@ -19,9 +19,7 @@
 
 #include <deal.II/base/config.h>
 
-#ifdef DEAL_II_WITH_TRILINOS
-
-#ifdef DEAL_II_WITH_MPI
+#if defined(DEAL_II_WITH_TRILINOS) && defined(DEAL_II_WITH_MPI)
 
 #include <deal.II/base/index_set.h>
 #include <deal.II/base/subscriptor.h>
@@ -351,5 +349,3 @@ DEAL_II_NAMESPACE_CLOSE
 #endif
 
 #endif
-
-#endif
index 1ad1054ad2bf85060eed9337a9cff0b8908417ef..2d79a16a45db8bdf7983d0840c3d13528560adda 100644 (file)
@@ -69,7 +69,7 @@ namespace internal
 
 
 
-#ifdef DEAL_II_WITH_TRILINOS
+#if defined(DEAL_II_WITH_TRILINOS) && defined(DEAL_II_WITH_MPI)
   template <>
   inline
   void ElementAccess<LinearAlgebra::EpetraWrappers::Vector>::
index b585b1756078203b51284c63e454fd4a01ac5cc3..211c0e6ea3cf1026f458499b6806c2612264c3b3 100644 (file)
@@ -89,6 +89,7 @@ namespace internal
 
   };
 
+#ifdef DEAL_II_WITH_MPI
   template <>
   struct MatrixSelector<dealii::LinearAlgebra::EpetraWrappers::Vector>
   {
@@ -102,8 +103,8 @@ namespace internal
                     dh.locally_owned_mg_dofs(level),
                     sp, MPI_COMM_WORLD, true);
     }
-
   };
+#endif
 
   template <>
   struct MatrixSelector<dealii::TrilinosWrappers::Vector>
index f97722f38fd48884f96d826b8140e3212360ec5c..9081b440bef26b3d26cc48f3613a7281b5eb18e0 100644 (file)
@@ -121,40 +121,41 @@ namespace TrilinosWrappers
       return V.begin();
     }
 
-    template <>
-    double *begin(LinearAlgebra::EpetraWrappers::Vector &V)
-    {
-      return V.trilinos_vector()[0];
-    }
-
     template <typename VectorType>
     const double *begin(const VectorType &V)
     {
       return V.begin();
     }
 
-    template <>
-    const double *begin(const LinearAlgebra::EpetraWrappers::Vector &V)
+    template <typename VectorType>
+    double *end(VectorType &V)
     {
-      return V.trilinos_vector()[0];
+      return V.end();
     }
 
     template <typename VectorType>
-    double *end(VectorType &V)
+    const double *end(const VectorType &V)
     {
       return V.end();
     }
 
+#ifdef DEAL_II_WITH_MPI
     template <>
-    double *end(LinearAlgebra::EpetraWrappers::Vector &V)
+    double *begin(LinearAlgebra::EpetraWrappers::Vector &V)
     {
-      return V.trilinos_vector()[0]+V.trilinos_vector().MyLength();
+      return V.trilinos_vector()[0];
     }
 
-    template <typename VectorType>
-    const double *end(const VectorType &V)
+    template <>
+    const double *begin(const LinearAlgebra::EpetraWrappers::Vector &V)
     {
-      return V.end();
+      return V.trilinos_vector()[0];
+    }
+
+    template <>
+    double *end(LinearAlgebra::EpetraWrappers::Vector &V)
+    {
+      return V.trilinos_vector()[0]+V.trilinos_vector().MyLength();
     }
 
     template <>
@@ -162,6 +163,7 @@ namespace TrilinosWrappers
     {
       return V.trilinos_vector()[0]+V.trilinos_vector().MyLength();
     }
+#endif
   }
 
 
@@ -3460,9 +3462,11 @@ namespace TrilinosWrappers
   template void
   SparseMatrix::vmult (dealii::LinearAlgebra::distributed::Vector<double> &,
                        const dealii::LinearAlgebra::distributed::Vector<double> &) const;
+#ifdef DEAL_II_WITH_MPI
   template void
   SparseMatrix::vmult (dealii::LinearAlgebra::EpetraWrappers::Vector &,
                        const dealii::LinearAlgebra::EpetraWrappers::Vector &) const;
+#endif
   template void
   SparseMatrix::Tvmult (VectorBase &,
                         const VectorBase &) const;
@@ -3478,9 +3482,11 @@ namespace TrilinosWrappers
   template void
   SparseMatrix::Tvmult (dealii::LinearAlgebra::distributed::Vector<double> &,
                         const dealii::LinearAlgebra::distributed::Vector<double> &) const;
+#ifdef DEAL_II_WITH_MPI
   template void
   SparseMatrix::Tvmult (dealii::LinearAlgebra::EpetraWrappers::Vector &,
                         const dealii::LinearAlgebra::EpetraWrappers::Vector &) const;
+#endif
   template void
   SparseMatrix::vmult_add (VectorBase &,
                            const VectorBase &) const;
@@ -3496,9 +3502,11 @@ namespace TrilinosWrappers
   template void
   SparseMatrix::vmult_add (dealii::LinearAlgebra::distributed::Vector<double> &,
                            const dealii::LinearAlgebra::distributed::Vector<double> &) const;
+#ifdef DEAL_II_WITH_MPI
   template void
   SparseMatrix::vmult_add (dealii::LinearAlgebra::EpetraWrappers::Vector &,
                            const dealii::LinearAlgebra::EpetraWrappers::Vector &) const;
+#endif
   template void
   SparseMatrix::Tvmult_add (VectorBase &,
                             const VectorBase &) const;
@@ -3514,10 +3522,11 @@ namespace TrilinosWrappers
   template void
   SparseMatrix::Tvmult_add (dealii::LinearAlgebra::distributed::Vector<double> &,
                             const dealii::LinearAlgebra::distributed::Vector<double> &) const;
+#ifdef DEAL_II_WITH_MPI
   template void
   SparseMatrix::Tvmult_add (dealii::LinearAlgebra::EpetraWrappers::Vector &,
                             const dealii::LinearAlgebra::EpetraWrappers::Vector &) const;
-
+#endif
 }
 
 DEAL_II_NAMESPACE_CLOSE
index 890ff4d2ce32bb356565854b7de9775d0e9e552e..2861b76a98ef8f1867a6b47148c18869a39c3870 100644 (file)
@@ -18,7 +18,7 @@
 
 DEAL_II_NAMESPACE_OPEN
 
-#ifdef DEAL_II_WITH_TRILINOS
+#if defined(DEAL_II_WITH_TRILINOS) && defined(DEAL_II_WITH_MPI)
 namespace VectorTools
 {
   template <>

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.