From b8bf3d93af2c83224ddd698b1228b1d415050406 Mon Sep 17 00:00:00 2001
From: Daniel Arndt <arndtd@ornl.gov>
Date: Tue, 6 Aug 2024 14:50:35 -0400
Subject: [PATCH] Tpetra: Test with MemorySpace::Default

---
 cmake/configure/configure_20_trilinos.cmake   |  16 +-
 include/deal.II/dofs/dof_accessor.templates.h |  11 +-
 .../fe/fe_tools_extrapolate.templates.h       |   7 +-
 .../fe/fe_tools_interpolate.templates.h       |  14 +-
 .../lac/affine_constraints.templates.h        |  10 +-
 .../trilinos_tpetra_sparse_matrix.templates.h |  44 ++--
 .../lac/trilinos_tpetra_vector.templates.h    |   2 +-
 include/deal.II/multigrid/mg_transfer.h       |   5 +-
 .../numerics/data_out_dof_data.templates.h    |   6 +-
 .../vector_tools_interpolate.templates.h      |   6 +-
 .../vector_tools_mean_value.templates.h       |   4 +-
 source/lac/affine_constraints.cc              |  10 +-
 source/lac/affine_constraints.inst.in         |  89 +++++--
 source/lac/read_write_vector.cc               |  32 ++-
 source/lac/trilinos_sparse_matrix.cc          | 120 ++++++++--
 .../trilinos_tpetra_block_sparse_matrix.cc    |  28 +--
 source/lac/trilinos_tpetra_block_vector.cc    |  34 ++-
 source/lac/trilinos_tpetra_precondition.cc    | 217 +++++++++++++-----
 source/lac/trilinos_tpetra_solver_direct.cc   |  40 +++-
 source/lac/trilinos_tpetra_sparse_matrix.cc   | 110 +++++++--
 source/lac/trilinos_tpetra_vector.cc          |  85 ++++---
 source/lac/vector.cc                          |  48 +++-
 tests/trilinos_tpetra/01.cc                   |   7 +-
 tests/trilinos_tpetra/02.cc                   |   7 +-
 tests/trilinos_tpetra/03a.cc                  |   7 +-
 tests/trilinos_tpetra/04.cc                   |   7 +-
 tests/trilinos_tpetra/05.cc                   |   7 +-
 tests/trilinos_tpetra/08.cc                   |   7 +-
 tests/trilinos_tpetra/09.cc                   |   7 +-
 tests/trilinos_tpetra/10.cc                   |   7 +-
 tests/trilinos_tpetra/11.cc                   |   4 +-
 tests/trilinos_tpetra/12.cc                   |   4 +-
 tests/trilinos_tpetra/13.cc                   |   4 +-
 tests/trilinos_tpetra/14.cc                   |   4 +-
 tests/trilinos_tpetra/15.cc                   |   4 +-
 tests/trilinos_tpetra/16.cc                   |   4 +-
 tests/trilinos_tpetra/17.cc                   |   4 +-
 tests/trilinos_tpetra/18.cc                   |   4 +-
 tests/trilinos_tpetra/19.cc                   |   4 +-
 tests/trilinos_tpetra/20.cc                   |   4 +-
 tests/trilinos_tpetra/21.cc                   |   4 +-
 tests/trilinos_tpetra/22.cc                   |   8 +-
 tests/trilinos_tpetra/23.cc                   |   8 +-
 tests/trilinos_tpetra/24.cc                   |   4 +-
 tests/trilinos_tpetra/27.cc                   |   6 +-
 tests/trilinos_tpetra/28.cc                   |   6 +-
 tests/trilinos_tpetra/29.cc                   |   4 +-
 tests/trilinos_tpetra/30.cc                   |   4 +-
 tests/trilinos_tpetra/31.cc                   |   4 +-
 tests/trilinos_tpetra/32.cc                   |   4 +-
 tests/trilinos_tpetra/34.cc                   |   4 +-
 tests/trilinos_tpetra/35.cc                   |   8 +-
 tests/trilinos_tpetra/36.cc                   |   8 +-
 tests/trilinos_tpetra/37.cc                   |   4 +-
 tests/trilinos_tpetra/39.cc                   |   8 +-
 tests/trilinos_tpetra/40.cc                   |  12 +-
 tests/trilinos_tpetra/42.cc                   |   8 +-
 tests/trilinos_tpetra/45.cc                   |   8 +-
 tests/trilinos_tpetra/46.cc                   |   8 +-
 tests/trilinos_tpetra/49.cc                   |   4 +-
 tests/trilinos_tpetra/51.cc                   |   6 +-
 tests/trilinos_tpetra/52.cc                   |  12 +-
 tests/trilinos_tpetra/53.cc                   |   4 +-
 tests/trilinos_tpetra/54.cc                   |   4 +-
 tests/trilinos_tpetra/55.cc                   |   4 +-
 tests/trilinos_tpetra/56.cc                   |   4 +-
 tests/trilinos_tpetra/57.cc                   |   4 +-
 tests/trilinos_tpetra/62.cc                   |   7 +-
 tests/trilinos_tpetra/63.cc                   |   7 +-
 tests/trilinos_tpetra/64.cc                   |  15 +-
 tests/trilinos_tpetra/66.cc                   |   7 +-
 tests/trilinos_tpetra/67.cc                   |   7 +-
 tests/trilinos_tpetra/68.cc                   |   7 +-
 tests/trilinos_tpetra/69.cc                   |   7 +-
 tests/trilinos_tpetra/add_matrices_06.cc      |  14 +-
 tests/trilinos_tpetra/deal_solver_01.cc       |  17 +-
 tests/trilinos_tpetra/deal_solver_01.output   |   2 +-
 tests/trilinos_tpetra/deal_solver_02.cc       |  15 +-
 tests/trilinos_tpetra/deal_solver_02.output   |   2 +-
 tests/trilinos_tpetra/deal_solver_03.cc       |  23 +-
 tests/trilinos_tpetra/deal_solver_03.output   |   2 +-
 tests/trilinos_tpetra/deal_solver_04.cc       |  17 +-
 tests/trilinos_tpetra/deal_solver_04.output   |   2 +-
 tests/trilinos_tpetra/deal_solver_05.cc       |  17 +-
 tests/trilinos_tpetra/deal_solver_05.output   |   2 +-
 tests/trilinos_tpetra/deal_solver_06.cc       |  18 +-
 tests/trilinos_tpetra/deal_solver_06.output   |   2 +-
 tests/trilinos_tpetra/direct_solver_01.cc     |  46 ++--
 tests/trilinos_tpetra/direct_solver_02.cc     |  36 +--
 tests/trilinos_tpetra/index_set_01.cc         |   4 +-
 tests/trilinos_tpetra/n_nonzeros_01.cc        |   6 +-
 tests/trilinos_tpetra/n_nonzeros_02.cc        |   4 +-
 .../parallel_sparse_matrix_01.cc              |   2 +-
 tests/trilinos_tpetra/precondition.cc         | 158 ++++++++-----
 tests/trilinos_tpetra/slowness_02.cc          |   7 +-
 tests/trilinos_tpetra/slowness_03.cc          |  11 +-
 tests/trilinos_tpetra/slowness_04.cc          |  11 +-
 tests/trilinos_tpetra/sparse_matrix_02.cc     |   3 +-
 tests/trilinos_tpetra/sparse_matrix_03.cc     |   3 +-
 tests/trilinos_tpetra/sparse_matrix_04.cc     |   3 +-
 .../sparse_matrix_copy_from_02.cc             |  12 +-
 .../sparse_matrix_iterator_01.cc              |   7 +-
 tests/trilinos_tpetra/sparse_matrix_set_01.cc |  12 +-
 tests/trilinos_tpetra/sparse_matrix_set_02.cc |  10 +-
 .../sparse_matrix_vector_01.cc                |  13 +-
 .../sparse_matrix_vector_02.cc                |  13 +-
 .../sparse_matrix_vector_03.cc                |  13 +-
 .../sparse_matrix_vector_04.cc                |  13 +-
 .../sparse_matrix_vector_05.cc                |  13 +-
 .../sparse_matrix_vector_06.cc                |   9 +-
 .../sparse_matrix_vector_07.cc                |  17 +-
 .../sparse_matrix_vector_08.cc                |   5 +-
 .../sparse_matrix_vector_14.cc                |  17 +-
 .../sparse_matrix_vector_15.cc                |  17 +-
 tests/trilinos_tpetra/tpetra_vector_01.cc     |   8 +-
 tests/trilinos_tpetra/tpetra_vector_02.cc     |  14 +-
 tests/trilinos_tpetra/tpetra_vector_03.cc     |   4 +-
 .../trilinos_sparse_matrix_iterator_01.cc     |  11 +-
 .../trilinos_sparse_matrix_iterator_02.cc     |  11 +-
 .../trilinos_sparse_matrix_iterator_03.cc     |  11 +-
 .../trilinos_sparse_matrix_iterator_04.cc     |  11 +-
 .../trilinos_sparse_matrix_iterator_05.cc     |  11 +-
 .../trilinos_sparse_matrix_iterator_06.cc     |  11 +-
 .../trilinos_sparse_matrix_iterator_07.cc     |  11 +-
 .../trilinos_sparse_matrix_iterator_08.cc     |  11 +-
 .../trilinos_sparse_matrix_iterator_09.cc     |  11 +-
 .../trilinos_sparse_matrix_iterator_10.cc     |  11 +-
 .../trilinos_sparse_matrix_iterator_11.cc     |  13 +-
 .../trilinos_sparse_matrix_iterator_12.cc     |  12 +-
 tests/trilinos_tpetra/vector_assign_01.cc     |   8 +-
 tests/trilinos_tpetra/vector_assign_02.cc     |   8 +-
 tests/trilinos_tpetra/vector_equality_1.cc    |   8 +-
 tests/trilinos_tpetra/vector_equality_2.cc    |   8 +-
 tests/trilinos_tpetra/vector_equality_3.cc    |   8 +-
 tests/trilinos_tpetra/vector_equality_4.cc    |   8 +-
 tests/trilinos_tpetra/vector_reference_01.cc  |   2 +-
 tests/trilinos_tpetra/vector_reinit.cc        |   2 +-
 tests/trilinos_tpetra/vector_swap_01.cc       |   6 +-
 138 files changed, 1346 insertions(+), 724 deletions(-)

diff --git a/cmake/configure/configure_20_trilinos.cmake b/cmake/configure/configure_20_trilinos.cmake
index 76ef7f1008..068dd8c6ad 100644
--- a/cmake/configure/configure_20_trilinos.cmake
+++ b/cmake/configure/configure_20_trilinos.cmake
@@ -446,20 +446,28 @@ macro(feature_trilinos_configure_external)
   if(${DEAL_II_TRILINOS_WITH_TPETRA})
     if(_tpetra_inst_double)
       set(DEAL_II_EXPAND_TPETRA_TYPES "double")
-      set(DEAL_II_EXPAND_TPETRA_VECTOR_DOUBLE "LinearAlgebra::TpetraWrappers::Vector<double>")
+      set(DEAL_II_EXPAND_TPETRA_VECTOR_DOUBLE
+        "LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Host>"
+        "LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>")
     endif()
 
     if(_tpetra_inst_float)
-      set(DEAL_II_EXPAND_TPETRA_VECTOR_FLOAT "LinearAlgebra::TpetraWrappers::Vector<float>")
+      set(DEAL_II_EXPAND_TPETRA_VECTOR_FLOAT
+        "LinearAlgebra::TpetraWrappers::Vector<float, MemorySpace::Host>"
+        "LinearAlgebra::TpetraWrappers::Vector<float, MemorySpace::Default>")
     endif()
 
     if(${DEAL_II_WITH_COMPLEX_NUMBERS})
       if(_tpetra_inst_complex_double)
-        set(DEAL_II_EXPAND_TPETRA_VECTOR_COMPLEX_DOUBLE "LinearAlgebra::TpetraWrappers::Vector<std::complex<double>>")
+        set(DEAL_II_EXPAND_TPETRA_VECTOR_COMPLEX_DOUBLE
+          "LinearAlgebra::TpetraWrappers::Vector<std::complex<double>, MemorySpace::Host>"
+          "LinearAlgebra::TpetraWrappers::Vector<std::complex<double>, MemorySpace::Default>")
       endif()
 
       if(_tpetra_inst_complex_float)
-        set(DEAL_II_EXPAND_TPETRA_VECTOR_COMPLEX_FLOAT "LinearAlgebra::TpetraWrappers::Vector<std::complex<float>>")
+        set(DEAL_II_EXPAND_TPETRA_VECTOR_COMPLEX_FLOAT
+          "LinearAlgebra::TpetraWrappers::Vector<std::complex<float>, MemorySpace::Host>"
+          "LinearAlgebra::TpetraWrappers::Vector<std::complex<float>, MemorySpace::Default>")
       endif()
     endif()
   endif()
diff --git a/include/deal.II/dofs/dof_accessor.templates.h b/include/deal.II/dofs/dof_accessor.templates.h
index c3afa34c61..851b3e1bd4 100644
--- a/include/deal.II/dofs/dof_accessor.templates.h
+++ b/include/deal.II/dofs/dof_accessor.templates.h
@@ -821,13 +821,14 @@ namespace internal
 
 
 #  ifdef DEAL_II_TRILINOS_WITH_TPETRA
-      template <typename ForwardIterator, typename Number>
+      template <typename ForwardIterator, typename Number, typename MemorySpace>
       static void
       extract_subvector_to(
-        const LinearAlgebra::TpetraWrappers::Vector<Number> &values,
-        const types::global_dof_index                       *cache_begin,
-        const types::global_dof_index                       *cache_end,
-        ForwardIterator                                      local_values_begin)
+        const LinearAlgebra::TpetraWrappers::Vector<Number, MemorySpace>
+                                      &values,
+        const types::global_dof_index *cache_begin,
+        const types::global_dof_index *cache_end,
+        ForwardIterator                local_values_begin)
       {
         std::vector<unsigned int> sorted_indices_pos =
           sort_indices(cache_begin, cache_end);
diff --git a/include/deal.II/fe/fe_tools_extrapolate.templates.h b/include/deal.II/fe/fe_tools_extrapolate.templates.h
index 6ebcb4dd46..54fd178c79 100644
--- a/include/deal.II/fe/fe_tools_extrapolate.templates.h
+++ b/include/deal.II/fe/fe_tools_extrapolate.templates.h
@@ -1516,10 +1516,11 @@ namespace FETools
 
 #  ifdef DEAL_II_WITH_MPI
 #    ifdef DEAL_II_TRILINOS_WITH_TPETRA
-    template <int dim, int spacedim, typename Number>
+    template <int dim, int spacedim, typename Number, typename MemorySpace>
     void
-    reinit_distributed(const DoFHandler<dim, spacedim>               &dh,
-                       LinearAlgebra::TpetraWrappers::Vector<Number> &vector)
+    reinit_distributed(
+      const DoFHandler<dim, spacedim>                            &dh,
+      LinearAlgebra::TpetraWrappers::Vector<Number, MemorySpace> &vector)
     {
       const parallel::distributed::Triangulation<dim, spacedim> *parallel_tria =
         dynamic_cast<
diff --git a/include/deal.II/fe/fe_tools_interpolate.templates.h b/include/deal.II/fe/fe_tools_interpolate.templates.h
index 2bfc84f004..bfb50b7b1d 100644
--- a/include/deal.II/fe/fe_tools_interpolate.templates.h
+++ b/include/deal.II/fe/fe_tools_interpolate.templates.h
@@ -492,17 +492,17 @@ namespace FETools
     }
 
 #  ifdef DEAL_II_TRILINOS_WITH_TPETRA
-    template <int dim, int spacedim, typename Number>
+    template <int dim, int spacedim, typename Number, typename MemorySpace>
     void
     back_interpolate(
       const DoFHandler<dim, spacedim> &,
-      const AffineConstraints<
-        typename LinearAlgebra::TpetraWrappers::Vector<Number>::value_type> &,
-      const LinearAlgebra::TpetraWrappers::Vector<Number> &,
+      const AffineConstraints<typename LinearAlgebra::TpetraWrappers::
+                                Vector<Number, MemorySpace>::value_type> &,
+      const LinearAlgebra::TpetraWrappers::Vector<Number, MemorySpace> &,
       const DoFHandler<dim, spacedim> &,
-      const AffineConstraints<
-        typename LinearAlgebra::TpetraWrappers::Vector<Number>::value_type> &,
-      LinearAlgebra::TpetraWrappers::Vector<Number> &)
+      const AffineConstraints<typename LinearAlgebra::TpetraWrappers::
+                                Vector<Number, MemorySpace>::value_type> &,
+      LinearAlgebra::TpetraWrappers::Vector<Number, MemorySpace> &)
     {
       AssertThrow(false, ExcNotImplemented());
     }
diff --git a/include/deal.II/lac/affine_constraints.templates.h b/include/deal.II/lac/affine_constraints.templates.h
index f0edc5f0da..05880bf333 100644
--- a/include/deal.II/lac/affine_constraints.templates.h
+++ b/include/deal.II/lac/affine_constraints.templates.h
@@ -2819,13 +2819,13 @@ namespace internal
 
 
 #ifdef DEAL_II_TRILINOS_WITH_TPETRA
-  template <typename Number>
+  template <typename Number, typename MemorySpace>
   inline void
   import_vector_with_ghost_elements(
-    const LinearAlgebra::TpetraWrappers::Vector<Number> &vec,
-    const IndexSet                                      &locally_owned_elements,
-    const IndexSet                                      &needed_elements,
-    LinearAlgebra::TpetraWrappers::Vector<Number>       &output,
+    const LinearAlgebra::TpetraWrappers::Vector<Number, MemorySpace> &vec,
+    const IndexSet &locally_owned_elements,
+    const IndexSet &needed_elements,
+    LinearAlgebra::TpetraWrappers::Vector<Number, MemorySpace> &output,
     const std::bool_constant<false> /*is_block_vector*/)
   {
     Assert(!vec.has_ghost_elements(), ExcGhostsPresent());
diff --git a/include/deal.II/lac/trilinos_tpetra_sparse_matrix.templates.h b/include/deal.II/lac/trilinos_tpetra_sparse_matrix.templates.h
index b12c3f394d..23cae8d196 100644
--- a/include/deal.II/lac/trilinos_tpetra_sparse_matrix.templates.h
+++ b/include/deal.II/lac/trilinos_tpetra_sparse_matrix.templates.h
@@ -46,28 +46,37 @@ namespace LinearAlgebra
             Number alpha = Teuchos::ScalarTraits<Number>::one(),
             Number beta  = Teuchos::ScalarTraits<Number>::zero())
       {
-        Assert(&src != &dst,
-               SparseMatrix<double>::ExcSourceEqualsDestination());
+        Assert(
+          &src != &dst,
+          (typename SparseMatrix<double,
+                                 MemorySpace>::ExcSourceEqualsDestination()));
         Assert(M.trilinos_matrix().isFillComplete(),
-               SparseMatrix<double>::ExcMatrixNotCompressed());
+               (typename SparseMatrix<double,
+                                      MemorySpace>::ExcMatrixNotCompressed()));
 
         if (mode == Teuchos::NO_TRANS)
           {
             Assert(src.trilinos_vector().getMap()->isSameAs(
                      *M.trilinos_matrix().getDomainMap()),
-                   SparseMatrix<double>::ExcColMapMismatch());
-            Assert(dst.trilinos_vector().getMap()->isSameAs(
-                     *M.trilinos_matrix().getRangeMap()),
-                   SparseMatrix<double>::ExcDomainMapMismatch());
+                   (typename SparseMatrix<double,
+                                          MemorySpace>::ExcColMapMismatch()));
+            Assert(
+              dst.trilinos_vector().getMap()->isSameAs(
+                *M.trilinos_matrix().getRangeMap()),
+              (typename SparseMatrix<double,
+                                     MemorySpace>::ExcDomainMapMismatch()));
           }
         else
           {
             Assert(dst.trilinos_vector().getMap()->isSameAs(
                      *M.trilinos_matrix().getDomainMap()),
-                   SparseMatrix<double>::ExcColMapMismatch());
-            Assert(src.trilinos_vector().getMap()->isSameAs(
-                     *M.trilinos_matrix().getRangeMap()),
-                   SparseMatrix<double>::ExcDomainMapMismatch());
+                   (typename SparseMatrix<double,
+                                          MemorySpace>::ExcColMapMismatch()));
+            Assert(
+              src.trilinos_vector().getMap()->isSameAs(
+                *M.trilinos_matrix().getRangeMap()),
+              (typename SparseMatrix<double,
+                                     MemorySpace>::ExcDomainMapMismatch()));
           }
 
         M.trilinos_matrix().apply(
@@ -85,10 +94,13 @@ namespace LinearAlgebra
             Number alpha = Teuchos::ScalarTraits<Number>::one(),
             Number beta  = Teuchos::ScalarTraits<Number>::zero())
       {
-        Assert(&src != &dst,
-               SparseMatrix<double>::ExcSourceEqualsDestination());
+        Assert(
+          &src != &dst,
+          (typename SparseMatrix<double,
+                                 MemorySpace>::ExcSourceEqualsDestination()));
         Assert(M.trilinos_matrix().isFillComplete(),
-               SparseMatrix<double>::ExcMatrixNotCompressed());
+               (typename SparseMatrix<double,
+                                      MemorySpace>::ExcMatrixNotCompressed()));
 
         // get the size of the input vectors:
         const size_type dst_local_size = dst.end() - dst.begin();
@@ -189,7 +201,7 @@ namespace LinearAlgebra
         // correct number of indices right from the start
         if (exchange_data)
           {
-            SparsityPattern trilinos_sparsity;
+            SparsityPattern<MemorySpace> trilinos_sparsity;
             trilinos_sparsity.reinit(row_parallel_partitioning,
                                      column_parallel_partitioning,
                                      sparsity_pattern,
@@ -307,7 +319,7 @@ namespace LinearAlgebra
         // correct number of indices right from the start
         if (exchange_data)
           {
-            SparsityPattern trilinos_sparsity;
+            SparsityPattern<MemorySpace> trilinos_sparsity;
             trilinos_sparsity.reinit(row_parallel_partitioning,
                                      column_parallel_partitioning,
                                      sparsity_pattern,
diff --git a/include/deal.II/lac/trilinos_tpetra_vector.templates.h b/include/deal.II/lac/trilinos_tpetra_vector.templates.h
index 98f3a318b8..acb36e1e97 100644
--- a/include/deal.II/lac/trilinos_tpetra_vector.templates.h
+++ b/include/deal.II/lac/trilinos_tpetra_vector.templates.h
@@ -1279,7 +1279,7 @@ namespace LinearAlgebra
       source_stored_elements = source_index_set;
 
       tpetra_comm_pattern =
-        Teuchos::rcp(new TpetraWrappers::CommunicationPattern(
+        Teuchos::rcp(new TpetraWrappers::CommunicationPattern<MemorySpace>(
           locally_owned_elements(), source_index_set, mpi_comm));
     }
   } // namespace TpetraWrappers
diff --git a/include/deal.II/multigrid/mg_transfer.h b/include/deal.II/multigrid/mg_transfer.h
index d020c8cc61..7117c97371 100644
--- a/include/deal.II/multigrid/mg_transfer.h
+++ b/include/deal.II/multigrid/mg_transfer.h
@@ -121,8 +121,9 @@ namespace internal
 
 #  ifdef DEAL_II_WITH_MPI
 #    ifdef DEAL_II_TRILINOS_WITH_TPETRA
-  template <typename Number>
-  struct MatrixSelector<dealii::LinearAlgebra::TpetraWrappers::Vector<Number>>
+  template <typename Number, typename MemorySpace>
+  struct MatrixSelector<
+    dealii::LinearAlgebra::TpetraWrappers::Vector<Number, MemorySpace>>
   {
     using Sparsity = ::dealii::TrilinosWrappers::SparsityPattern;
     using Matrix   = ::dealii::TrilinosWrappers::SparseMatrix;
diff --git a/include/deal.II/numerics/data_out_dof_data.templates.h b/include/deal.II/numerics/data_out_dof_data.templates.h
index 8e22d75b0e..4265f49d1c 100644
--- a/include/deal.II/numerics/data_out_dof_data.templates.h
+++ b/include/deal.II/numerics/data_out_dof_data.templates.h
@@ -847,11 +847,11 @@ namespace internal
 #endif
 
 #ifdef DEAL_II_TRILINOS_WITH_TPETRA
-      template <typename Number>
+      template <typename Number, typename MemorySpace>
       void
       copy_locally_owned_data_from(
-        const LinearAlgebra::TpetraWrappers::Vector<Number> &src,
-        LinearAlgebra::distributed::Vector<Number>          &dst)
+        const LinearAlgebra::TpetraWrappers::Vector<Number, MemorySpace> &src,
+        LinearAlgebra::distributed::Vector<Number>                       &dst)
       {
         // ReadWriteVector does not work for ghosted
         // TrilinosWrappers::MPI::Vector objects. Fall back to copy the
diff --git a/include/deal.II/numerics/vector_tools_interpolate.templates.h b/include/deal.II/numerics/vector_tools_interpolate.templates.h
index 83a93b42eb..21b1dea516 100644
--- a/include/deal.II/numerics/vector_tools_interpolate.templates.h
+++ b/include/deal.II/numerics/vector_tools_interpolate.templates.h
@@ -1147,11 +1147,11 @@ namespace VectorTools
 #endif
 
 #ifdef DEAL_II_TRILINOS_WITH_TPETRA
-    template <typename Number>
+    template <typename Number, typename MemorySpace>
     void
     copy_locally_owned_data_from(
-      const LinearAlgebra::TpetraWrappers::Vector<Number> &src,
-      LinearAlgebra::distributed::Vector<Number>          &dst)
+      const LinearAlgebra::TpetraWrappers::Vector<Number, MemorySpace> &src,
+      LinearAlgebra::distributed::Vector<Number>                       &dst)
     {
       // ReadWriteVector does not work for ghosted
       // TrilinosWrappers::MPI::Vector objects. Fall back to copy the
diff --git a/include/deal.II/numerics/vector_tools_mean_value.templates.h b/include/deal.II/numerics/vector_tools_mean_value.templates.h
index a87c1c3d52..2005fd05c8 100644
--- a/include/deal.II/numerics/vector_tools_mean_value.templates.h
+++ b/include/deal.II/numerics/vector_tools_mean_value.templates.h
@@ -273,9 +273,9 @@ namespace VectorTools
 
 
 #  ifdef DEAL_II_TRILINOS_WITH_TPETRA
-  template <int dim, int spacedim, typename ValueType>
+  template <int dim, int spacedim, typename ValueType, typename MemorySpace>
   void
-  add_constant(LinearAlgebra::TpetraWrappers::Vector<ValueType> &,
+  add_constant(LinearAlgebra::TpetraWrappers::Vector<ValueType, MemorySpace> &,
                const DoFHandler<dim, spacedim> &,
                const unsigned int,
                const ValueType)
diff --git a/source/lac/affine_constraints.cc b/source/lac/affine_constraints.cc
index 99483b5cc7..e051e03b83 100644
--- a/source/lac/affine_constraints.cc
+++ b/source/lac/affine_constraints.cc
@@ -143,8 +143,14 @@ INSTANTIATE_DLTG_MATRIX(TrilinosWrappers::BlockSparseMatrix);
 // FIXME: This mixed variant is needed for multigrid and matrix free.
 template void
 dealii::AffineConstraints<double>::distribute<
-  dealii::LinearAlgebra::TpetraWrappers::Vector<float>>(
-  dealii::LinearAlgebra::TpetraWrappers::Vector<float> &) const;
+  dealii::LinearAlgebra::TpetraWrappers::Vector<float, MemorySpace::Host>>(
+  dealii::LinearAlgebra::TpetraWrappers::Vector<float, MemorySpace::Host> &)
+  const;
+template void
+dealii::AffineConstraints<double>::distribute<
+  dealii::LinearAlgebra::TpetraWrappers::Vector<float, MemorySpace::Default>>(
+  dealii::LinearAlgebra::TpetraWrappers::Vector<float, MemorySpace::Default> &)
+  const;
 #    endif
 #  endif
 #endif
diff --git a/source/lac/affine_constraints.inst.in b/source/lac/affine_constraints.inst.in
index a40947ebeb..bffe8f17bc 100644
--- a/source/lac/affine_constraints.inst.in
+++ b/source/lac/affine_constraints.inst.in
@@ -246,53 +246,110 @@ for (S : REAL_AND_COMPLEX_SCALARS)
 for (S : TRILINOS_SCALARS)
   {
     template void AffineConstraints<S>::distribute_local_to_global<
-      LinearAlgebra::TpetraWrappers::SparseMatrix<S>,
-      LinearAlgebra::TpetraWrappers::Vector<S>>(
+      LinearAlgebra::TpetraWrappers::SparseMatrix<S, MemorySpace::Host>,
+      LinearAlgebra::TpetraWrappers::Vector<S, MemorySpace::Host>>(
       const FullMatrix<S> &,
       const Vector<S> &,
       const std::vector<AffineConstraints<S>::size_type> &,
-      LinearAlgebra::TpetraWrappers::SparseMatrix<S> &,
-      LinearAlgebra::TpetraWrappers::Vector<S> &,
+      LinearAlgebra::TpetraWrappers::SparseMatrix<S, MemorySpace::Host> &,
+      LinearAlgebra::TpetraWrappers::Vector<S, MemorySpace::Host> &,
       bool,
       std::integral_constant<bool, false>) const;
 
     // BlockSparseMatrix
     template void AffineConstraints<S>::distribute_local_to_global<
-      LinearAlgebra::TpetraWrappers::BlockSparseMatrix<S>,
-      LinearAlgebra::TpetraWrappers::Vector<S>>(
+      LinearAlgebra::TpetraWrappers::BlockSparseMatrix<S, MemorySpace::Host>,
+      LinearAlgebra::TpetraWrappers::Vector<S, MemorySpace::Host>>(
       const FullMatrix<S> &,
       const Vector<S> &,
       const std::vector<AffineConstraints<S>::size_type> &,
-      LinearAlgebra::TpetraWrappers::BlockSparseMatrix<S> &,
-      LinearAlgebra::TpetraWrappers::Vector<S> &,
+      LinearAlgebra::TpetraWrappers::BlockSparseMatrix<S, MemorySpace::Host> &,
+      LinearAlgebra::TpetraWrappers::Vector<S, MemorySpace::Host> &,
       bool,
       std::bool_constant<true>) const;
 
     template void AffineConstraints<S>::distribute_local_to_global<
-      LinearAlgebra::TpetraWrappers::BlockSparseMatrix<S>,
-      LinearAlgebra::TpetraWrappers::BlockVector<S>>(
+      LinearAlgebra::TpetraWrappers::BlockSparseMatrix<S, MemorySpace::Host>,
+      LinearAlgebra::TpetraWrappers::BlockVector<S, MemorySpace::Host>>(
       const FullMatrix<S> &,
       const Vector<S> &,
       const std::vector<AffineConstraints<S>::size_type> &,
-      LinearAlgebra::TpetraWrappers::BlockSparseMatrix<S> &,
-      LinearAlgebra::TpetraWrappers::BlockVector<S> &,
+      LinearAlgebra::TpetraWrappers::BlockSparseMatrix<S, MemorySpace::Host> &,
+      LinearAlgebra::TpetraWrappers::BlockVector<S, MemorySpace::Host> &,
       bool,
       std::bool_constant<true>) const;
 
     template void AffineConstraints<S>::distribute_local_to_global<
-      LinearAlgebra::TpetraWrappers::BlockSparseMatrix<S>>(
+      LinearAlgebra::TpetraWrappers::BlockSparseMatrix<S, MemorySpace::Host>>(
       const FullMatrix<S> &,
       const std::vector<AffineConstraints<S>::size_type> &,
       const std::vector<AffineConstraints<S>::size_type> &,
-      LinearAlgebra::TpetraWrappers::BlockSparseMatrix<S> &) const;
+      LinearAlgebra::TpetraWrappers::BlockSparseMatrix<S, MemorySpace::Host> &)
+      const;
+
+    template void AffineConstraints<S>::distribute_local_to_global<
+      LinearAlgebra::TpetraWrappers::BlockSparseMatrix<S, MemorySpace::Host>>(
+      const FullMatrix<S> &,
+      const std::vector<AffineConstraints<S>::size_type> &,
+      const AffineConstraints<S> &,
+      const std::vector<AffineConstraints<S>::size_type> &,
+      LinearAlgebra::TpetraWrappers::BlockSparseMatrix<S, MemorySpace::Host> &)
+      const;
+
+    template void AffineConstraints<S>::distribute_local_to_global<
+      LinearAlgebra::TpetraWrappers::SparseMatrix<S, MemorySpace::Default>,
+      LinearAlgebra::TpetraWrappers::Vector<S, MemorySpace::Default>>(
+      const FullMatrix<S> &,
+      const Vector<S> &,
+      const std::vector<AffineConstraints<S>::size_type> &,
+      LinearAlgebra::TpetraWrappers::SparseMatrix<S, MemorySpace::Default> &,
+      LinearAlgebra::TpetraWrappers::Vector<S, MemorySpace::Default> &,
+      bool,
+      std::integral_constant<bool, false>) const;
+
+    // BlockSparseMatrix
+    template void AffineConstraints<S>::distribute_local_to_global<
+      LinearAlgebra::TpetraWrappers::BlockSparseMatrix<S, MemorySpace::Default>,
+      LinearAlgebra::TpetraWrappers::Vector<S, MemorySpace::Default>>(
+      const FullMatrix<S> &,
+      const Vector<S> &,
+      const std::vector<AffineConstraints<S>::size_type> &,
+      LinearAlgebra::TpetraWrappers::BlockSparseMatrix<S, MemorySpace::Default>
+        &,
+      LinearAlgebra::TpetraWrappers::Vector<S, MemorySpace::Default> &,
+      bool,
+      std::bool_constant<true>) const;
+
+    template void AffineConstraints<S>::distribute_local_to_global<
+      LinearAlgebra::TpetraWrappers::BlockSparseMatrix<S, MemorySpace::Default>,
+      LinearAlgebra::TpetraWrappers::BlockVector<S, MemorySpace::Default>>(
+      const FullMatrix<S> &,
+      const Vector<S> &,
+      const std::vector<AffineConstraints<S>::size_type> &,
+      LinearAlgebra::TpetraWrappers::BlockSparseMatrix<S, MemorySpace::Default>
+        &,
+      LinearAlgebra::TpetraWrappers::BlockVector<S, MemorySpace::Default> &,
+      bool,
+      std::bool_constant<true>) const;
+
+    template void AffineConstraints<S>::distribute_local_to_global<
+      LinearAlgebra::TpetraWrappers::BlockSparseMatrix<S,
+                                                       MemorySpace::Default>>(
+      const FullMatrix<S> &,
+      const std::vector<AffineConstraints<S>::size_type> &,
+      const std::vector<AffineConstraints<S>::size_type> &,
+      LinearAlgebra::TpetraWrappers::BlockSparseMatrix<S, MemorySpace::Default>
+        &) const;
 
     template void AffineConstraints<S>::distribute_local_to_global<
-      LinearAlgebra::TpetraWrappers::BlockSparseMatrix<S>>(
+      LinearAlgebra::TpetraWrappers::BlockSparseMatrix<S,
+                                                       MemorySpace::Default>>(
       const FullMatrix<S> &,
       const std::vector<AffineConstraints<S>::size_type> &,
       const AffineConstraints<S> &,
       const std::vector<AffineConstraints<S>::size_type> &,
-      LinearAlgebra::TpetraWrappers::BlockSparseMatrix<S> &) const;
+      LinearAlgebra::TpetraWrappers::BlockSparseMatrix<S, MemorySpace::Default>
+        &) const;
   }
 
 
diff --git a/source/lac/read_write_vector.cc b/source/lac/read_write_vector.cc
index 508789b343..a46256e693 100644
--- a/source/lac/read_write_vector.cc
+++ b/source/lac/read_write_vector.cc
@@ -82,14 +82,24 @@ namespace LinearAlgebra
 #  ifdef HAVE_TPETRA_INST_FLOAT
   template void
   ReadWriteVector<float>::import_elements(
-    const LinearAlgebra::TpetraWrappers::Vector<float> &,
+    const LinearAlgebra::TpetraWrappers::Vector<float, MemorySpace::Host> &,
+    VectorOperation::values,
+    const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase> &);
+  template void
+  ReadWriteVector<float>::import_elements(
+    const LinearAlgebra::TpetraWrappers::Vector<float, MemorySpace::Default> &,
     VectorOperation::values,
     const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase> &);
 #  endif
 #  ifdef HAVE_TPETRA_INST_DOUBLE
   template void
   ReadWriteVector<double>::import_elements(
-    const LinearAlgebra::TpetraWrappers::Vector<double> &,
+    const LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Host> &,
+    VectorOperation::values,
+    const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase> &);
+  template void
+  ReadWriteVector<double>::import_elements(
+    const LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &,
     VectorOperation::values,
     const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase> &);
 #  endif
@@ -97,14 +107,28 @@ namespace LinearAlgebra
 #    ifdef HAVE_TPETRA_INST_COMPLEX_FLOAT
   template void
   ReadWriteVector<std::complex<float>>::import_elements(
-    const LinearAlgebra::TpetraWrappers::Vector<std::complex<float>> &,
+    const LinearAlgebra::TpetraWrappers::Vector<std::complex<float>,
+                                                MemorySpace::Host> &,
+    VectorOperation::values,
+    const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase> &);
+  template void
+  ReadWriteVector<std::complex<float>>::import_elements(
+    const LinearAlgebra::TpetraWrappers::Vector<std::complex<float>,
+                                                MemorySpace::Default> &,
     VectorOperation::values,
     const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase> &);
 #    endif
 #    ifdef HAVE_TPETRA_INST_COMPLEX_DOUBLE
   template void
   ReadWriteVector<std::complex<double>>::import_elements(
-    const LinearAlgebra::TpetraWrappers::Vector<std::complex<double>> &,
+    const LinearAlgebra::TpetraWrappers::Vector<std::complex<double>,
+                                                MemorySpace::Host> &,
+    VectorOperation::values,
+    const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase> &);
+  template void
+  ReadWriteVector<std::complex<double>>::import_elements(
+    const LinearAlgebra::TpetraWrappers::Vector<std::complex<double>,
+                                                MemorySpace::Default> &,
     VectorOperation::values,
     const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase> &);
 #    endif
diff --git a/source/lac/trilinos_sparse_matrix.cc b/source/lac/trilinos_sparse_matrix.cc
index 99bd58cb91..20f6029a40 100644
--- a/source/lac/trilinos_sparse_matrix.cc
+++ b/source/lac/trilinos_sparse_matrix.cc
@@ -102,31 +102,31 @@ namespace TrilinosWrappers
     }
 
 #  ifdef DEAL_II_TRILINOS_WITH_TPETRA
-    template <typename Number>
+    template <typename Number, typename MemorySpace>
     Number *
-    begin(LinearAlgebra::TpetraWrappers::Vector<Number> &V)
+    begin(LinearAlgebra::TpetraWrappers::Vector<Number, MemorySpace> &V)
     {
       return V.trilinos_vector().getDataNonConst().get();
     }
 
-    template <typename Number>
+    template <typename Number, typename MemorySpace>
     const Number *
-    begin(const LinearAlgebra::TpetraWrappers::Vector<Number> &V)
+    begin(const LinearAlgebra::TpetraWrappers::Vector<Number, MemorySpace> &V)
     {
       return V.trilinos_vector().getData().get();
     }
 
-    template <typename Number>
+    template <typename Number, typename MemorySpace>
     Number *
-    end(LinearAlgebra::TpetraWrappers::Vector<Number> &V)
+    end(LinearAlgebra::TpetraWrappers::Vector<Number, MemorySpace> &V)
     {
       return V.trilinos_vector().getDataNonConst().get() +
              V.trilinos_vector().getLocalLength();
     }
 
-    template <typename Number>
+    template <typename Number, typename MemorySpace>
     const Number *
-    end(const LinearAlgebra::TpetraWrappers::Vector<Number> &V)
+    end(const LinearAlgebra::TpetraWrappers::Vector<Number, MemorySpace> &V)
     {
       return V.trilinos_vector().getData().get() +
              V.trilinos_vector().getLocalLength();
@@ -3049,15 +3049,33 @@ namespace TrilinosWrappers
 #      if defined(HAVE_TPETRA_INST_DOUBLE)
   template void
   SparseMatrix::vmult(
-    dealii::LinearAlgebra::TpetraWrappers::Vector<double> &,
-    const dealii::LinearAlgebra::TpetraWrappers::Vector<double> &) const;
+    dealii::LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Host> &,
+    const dealii::LinearAlgebra::TpetraWrappers::Vector<double,
+                                                        MemorySpace::Host> &)
+    const;
+  template void
+  SparseMatrix::vmult(
+    dealii::LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>
+      &,
+    const dealii::LinearAlgebra::TpetraWrappers::Vector<double,
+                                                        MemorySpace::Default> &)
+    const;
 #      endif
 
 #      if defined(HAVE_TPETRA_INST_FLOAT)
   template void
   SparseMatrix::vmult(
-    dealii::LinearAlgebra::TpetraWrappers::Vector<float> &,
-    const dealii::LinearAlgebra::TpetraWrappers::Vector<float> &) const;
+    dealii::LinearAlgebra::TpetraWrappers::Vector<float, MemorySpace::Host> &,
+    const dealii::LinearAlgebra::TpetraWrappers::Vector<float,
+                                                        MemorySpace::Host> &)
+    const;
+  template void
+  SparseMatrix::vmult(
+    dealii::LinearAlgebra::TpetraWrappers::Vector<float, MemorySpace::Default>
+      &,
+    const dealii::LinearAlgebra::TpetraWrappers::Vector<float,
+                                                        MemorySpace::Default> &)
+    const;
 #      endif
 #    endif
 
@@ -3082,15 +3100,33 @@ namespace TrilinosWrappers
 #      if defined(HAVE_TPETRA_INST_DOUBLE)
   template void
   SparseMatrix::Tvmult(
-    dealii::LinearAlgebra::TpetraWrappers::Vector<double> &,
-    const dealii::LinearAlgebra::TpetraWrappers::Vector<double> &) const;
+    dealii::LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Host> &,
+    const dealii::LinearAlgebra::TpetraWrappers::Vector<double,
+                                                        MemorySpace::Host> &)
+    const;
+  template void
+  SparseMatrix::Tvmult(
+    dealii::LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>
+      &,
+    const dealii::LinearAlgebra::TpetraWrappers::Vector<double,
+                                                        MemorySpace::Default> &)
+    const;
 #      endif
 
 #      if defined(HAVE_TPETRA_INST_FLOAT)
   template void
   SparseMatrix::Tvmult(
-    dealii::LinearAlgebra::TpetraWrappers::Vector<float> &,
-    const dealii::LinearAlgebra::TpetraWrappers::Vector<float> &) const;
+    dealii::LinearAlgebra::TpetraWrappers::Vector<float, MemorySpace::Host> &,
+    const dealii::LinearAlgebra::TpetraWrappers::Vector<float,
+                                                        MemorySpace::Host> &)
+    const;
+  template void
+  SparseMatrix::Tvmult(
+    dealii::LinearAlgebra::TpetraWrappers::Vector<float, MemorySpace::Default>
+      &,
+    const dealii::LinearAlgebra::TpetraWrappers::Vector<float,
+                                                        MemorySpace::Default> &)
+    const;
 #      endif
 #    endif
 
@@ -3115,15 +3151,33 @@ namespace TrilinosWrappers
 #      if defined(HAVE_TPETRA_INST_DOUBLE)
   template void
   SparseMatrix::vmult_add(
-    dealii::LinearAlgebra::TpetraWrappers::Vector<double> &,
-    const dealii::LinearAlgebra::TpetraWrappers::Vector<double> &) const;
+    dealii::LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Host> &,
+    const dealii::LinearAlgebra::TpetraWrappers::Vector<double,
+                                                        MemorySpace::Host> &)
+    const;
+  template void
+  SparseMatrix::vmult_add(
+    dealii::LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>
+      &,
+    const dealii::LinearAlgebra::TpetraWrappers::Vector<double,
+                                                        MemorySpace::Default> &)
+    const;
 #      endif
 
 #      if defined(HAVE_TPETRA_INST_FLOAT)
   template void
   SparseMatrix::vmult_add(
-    dealii::LinearAlgebra::TpetraWrappers::Vector<float> &,
-    const dealii::LinearAlgebra::TpetraWrappers::Vector<float> &) const;
+    dealii::LinearAlgebra::TpetraWrappers::Vector<float, MemorySpace::Host> &,
+    const dealii::LinearAlgebra::TpetraWrappers::Vector<float,
+                                                        MemorySpace::Host> &)
+    const;
+  template void
+  SparseMatrix::vmult_add(
+    dealii::LinearAlgebra::TpetraWrappers::Vector<float, MemorySpace::Default>
+      &,
+    const dealii::LinearAlgebra::TpetraWrappers::Vector<float,
+                                                        MemorySpace::Default> &)
+    const;
 #      endif
 #    endif
 
@@ -3148,15 +3202,33 @@ namespace TrilinosWrappers
 #      if defined(HAVE_TPETRA_INST_DOUBLE)
   template void
   SparseMatrix::Tvmult_add(
-    dealii::LinearAlgebra::TpetraWrappers::Vector<double> &,
-    const dealii::LinearAlgebra::TpetraWrappers::Vector<double> &) const;
+    dealii::LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Host> &,
+    const dealii::LinearAlgebra::TpetraWrappers::Vector<double,
+                                                        MemorySpace::Host> &)
+    const;
+  template void
+  SparseMatrix::Tvmult_add(
+    dealii::LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>
+      &,
+    const dealii::LinearAlgebra::TpetraWrappers::Vector<double,
+                                                        MemorySpace::Default> &)
+    const;
 #      endif
 
 #      if defined(HAVE_TPETRA_INST_FLOAT)
   template void
   SparseMatrix::Tvmult_add(
-    dealii::LinearAlgebra::TpetraWrappers::Vector<float> &,
-    const dealii::LinearAlgebra::TpetraWrappers::Vector<float> &) const;
+    dealii::LinearAlgebra::TpetraWrappers::Vector<float, MemorySpace::Host> &,
+    const dealii::LinearAlgebra::TpetraWrappers::Vector<float,
+                                                        MemorySpace::Host> &)
+    const;
+  template void
+  SparseMatrix::Tvmult_add(
+    dealii::LinearAlgebra::TpetraWrappers::Vector<float, MemorySpace::Default>
+      &,
+    const dealii::LinearAlgebra::TpetraWrappers::Vector<float,
+                                                        MemorySpace::Default> &)
+    const;
 #      endif
 #    endif
 
diff --git a/source/lac/trilinos_tpetra_block_sparse_matrix.cc b/source/lac/trilinos_tpetra_block_sparse_matrix.cc
index e41ad48a61..a51c197ce9 100644
--- a/source/lac/trilinos_tpetra_block_sparse_matrix.cc
+++ b/source/lac/trilinos_tpetra_block_sparse_matrix.cc
@@ -26,31 +26,31 @@ namespace LinearAlgebra
 {
   namespace TpetraWrappers
   {
-    template class BlockSparseMatrix<double>;
+    template class BlockSparseMatrix<double, MemorySpace::Host>;
 
     template void
-    BlockSparseMatrix<double>::reinit(
+    BlockSparseMatrix<double, MemorySpace::Host>::reinit(
       const ::dealii::BlockDynamicSparsityPattern &);
 
     template void
-    BlockSparseMatrix<double>::vmult(
-      TpetraWrappers::Vector<double> &,
-      const TpetraWrappers::Vector<double> &) const;
+    BlockSparseMatrix<double, MemorySpace::Host>::vmult(
+      TpetraWrappers::Vector<double, MemorySpace::Host> &,
+      const TpetraWrappers::Vector<double, MemorySpace::Host> &) const;
 
     template void
-    BlockSparseMatrix<double>::Tvmult(
-      TpetraWrappers::Vector<double> &,
-      const TpetraWrappers::Vector<double> &) const;
+    BlockSparseMatrix<double, MemorySpace::Host>::Tvmult(
+      TpetraWrappers::Vector<double, MemorySpace::Host> &,
+      const TpetraWrappers::Vector<double, MemorySpace::Host> &) const;
 
     template void
-    BlockSparseMatrix<double>::vmult(
-      TpetraWrappers::BlockVector<double> &,
-      const TpetraWrappers::BlockVector<double> &) const;
+    BlockSparseMatrix<double, MemorySpace::Host>::vmult(
+      TpetraWrappers::BlockVector<double, MemorySpace::Host> &,
+      const TpetraWrappers::BlockVector<double, MemorySpace::Host> &) const;
 
     template void
-    BlockSparseMatrix<double>::Tvmult(
-      TpetraWrappers::BlockVector<double> &,
-      const TpetraWrappers::BlockVector<double> &) const;
+    BlockSparseMatrix<double, MemorySpace::Host>::Tvmult(
+      TpetraWrappers::BlockVector<double, MemorySpace::Host> &,
+      const TpetraWrappers::BlockVector<double, MemorySpace::Host> &) const;
   } // namespace TpetraWrappers
 } // namespace LinearAlgebra
 #  endif // DOXYGEN
diff --git a/source/lac/trilinos_tpetra_block_vector.cc b/source/lac/trilinos_tpetra_block_vector.cc
index 673c61b590..d1a114203a 100644
--- a/source/lac/trilinos_tpetra_block_vector.cc
+++ b/source/lac/trilinos_tpetra_block_vector.cc
@@ -32,30 +32,44 @@ namespace LinearAlgebra
 
 #  ifdef HAVE_TPETRA_INST_FLOAT
 #    ifdef DEAL_II_HAVE_CXX20
-    static_assert(concepts::is_vector_space_vector<BlockVector<float>>);
+    static_assert(
+      concepts::is_vector_space_vector<BlockVector<float, MemorySpace::Host>>);
+    static_assert(concepts::is_vector_space_vector<
+                  BlockVector<float, MemorySpace::Default>>);
 #    endif
-    template class BlockVector<float>;
+    template class BlockVector<float, MemorySpace::Host>;
+    template class BlockVector<float, MemorySpace::Default>;
 #  endif
 #  ifdef HAVE_TPETRA_INST_DOUBLE
 #    ifdef DEAL_II_HAVE_CXX20
-    static_assert(concepts::is_vector_space_vector<BlockVector<double>>);
+    static_assert(
+      concepts::is_vector_space_vector<BlockVector<double, MemorySpace::Host>>);
+    static_assert(concepts::is_vector_space_vector<
+                  BlockVector<double, MemorySpace::Default>>);
 #    endif
-    template class BlockVector<double>;
+    template class BlockVector<double, MemorySpace::Host>;
+    template class BlockVector<double, MemorySpace::Default>;
 #  endif
 #  ifdef DEAL_II_WITH_COMPLEX_VALUES
 #    ifdef HAVE_TPETRA_INST_COMPLEX_FLOAT
 #      ifdef DEAL_II_HAVE_CXX20
-    static_assert(
-      concepts::is_vector_space_vector<BlockVector<std::complex<float>>>);
+    static_assert(concepts::is_vector_space_vector<
+                  BlockVector<std::complex<float>, MemorySpace::Host>>);
+    static_assert(concepts::is_vector_space_vector<
+                  BlockVector<std::complex<float>, MemorySpace::Default>>);
 #      endif
-    template class BlockVector<std::complex<float>>;
+    template class BlockVector<std::complex<float>, MemorySpace::Host>;
+    template class BlockVector<std::complex<float>, MemorySpace::Default>;
 #    endif
 #    ifdef HAVE_TPETRA_INST_COMPLEX_DOUBLE
 #      ifdef DEAL_II_HAVE_CXX20
-    static_assert(
-      concepts::is_vector_space_vector<BlockVector<std::complex<double>>>);
+    static_assert(concepts::is_vector_space_vector<
+                  BlockVector<std::complex<double>, MemorySpace::Host>>);
+    static_assert(concepts::is_vector_space_vector<
+                  BlockVector<std::complex<double>, MemorySpace::Default>>);
 #      endif
-    template class BlockVector<std::complex<double>>;
+    template class BlockVector<std::complex<double>, MemorySpace::Host>;
+    template class BlockVector<std::complex<double>, MemorySpace::Default>;
 #    endif
 #  endif
   } // namespace TpetraWrappers
diff --git a/source/lac/trilinos_tpetra_precondition.cc b/source/lac/trilinos_tpetra_precondition.cc
index 21f7ee2562..0f8d3f77ee 100644
--- a/source/lac/trilinos_tpetra_precondition.cc
+++ b/source/lac/trilinos_tpetra_precondition.cc
@@ -27,74 +27,171 @@ namespace LinearAlgebra
   namespace TpetraWrappers
   {
 #      ifdef HAVE_TPETRA_INST_FLOAT
-    template class PreconditionBase<float>;
-    template class PreconditionIdentity<float>;
-    template class PreconditionIfpackBase<float>;
-    template class PreconditionIfpack<float>;
-    template class PreconditionJacobi<float>;
-    template class PreconditionL1Jacobi<float>;
-    template class PreconditionL1GaussSeidel<float>;
-    template class PreconditionSOR<float>;
-    template class PreconditionSSOR<float>;
-    template class PreconditionBlockJacobi<float>;
-    template class PreconditionBlockSOR<float>;
-    template class PreconditionBlockSSOR<float>;
-    template class PreconditionChebyshev<float>;
-    template class PreconditionILU<float>;
-    template class PreconditionILUT<float>;
+    template class PreconditionBase<float, MemorySpace::Host>;
+    template class PreconditionIdentity<float, MemorySpace::Host>;
+    template class PreconditionIfpackBase<float, MemorySpace::Host>;
+    template class PreconditionIfpack<float, MemorySpace::Host>;
+    template class PreconditionJacobi<float, MemorySpace::Host>;
+    template class PreconditionL1Jacobi<float, MemorySpace::Host>;
+    template class PreconditionL1GaussSeidel<float, MemorySpace::Host>;
+    template class PreconditionSOR<float, MemorySpace::Host>;
+    template class PreconditionSSOR<float, MemorySpace::Host>;
+    template class PreconditionBlockJacobi<float, MemorySpace::Host>;
+    template class PreconditionBlockSOR<float, MemorySpace::Host>;
+    template class PreconditionBlockSSOR<float, MemorySpace::Host>;
+    template class PreconditionChebyshev<float, MemorySpace::Host>;
+    template class PreconditionILU<float, MemorySpace::Host>;
+    template class PreconditionILUT<float, MemorySpace::Host>;
+
+    template class PreconditionBase<float, MemorySpace::Default>;
+    template class PreconditionIdentity<float, MemorySpace::Default>;
+    template class PreconditionIfpackBase<float, MemorySpace::Default>;
+    template class PreconditionIfpack<float, MemorySpace::Default>;
+    template class PreconditionJacobi<float, MemorySpace::Default>;
+    template class PreconditionL1Jacobi<float, MemorySpace::Default>;
+    template class PreconditionL1GaussSeidel<float, MemorySpace::Default>;
+    template class PreconditionSOR<float, MemorySpace::Default>;
+    template class PreconditionSSOR<float, MemorySpace::Default>;
+    template class PreconditionBlockJacobi<float, MemorySpace::Default>;
+    template class PreconditionBlockSOR<float, MemorySpace::Default>;
+    template class PreconditionBlockSSOR<float, MemorySpace::Default>;
+    template class PreconditionChebyshev<float, MemorySpace::Default>;
+    template class PreconditionILU<float, MemorySpace::Default>;
+    template class PreconditionILUT<float, MemorySpace::Default>;
 #      endif
 
 #      ifdef HAVE_TPETRA_INST_DOUBLE
-    template class PreconditionBase<double>;
-    template class PreconditionIdentity<double>;
-    template class PreconditionIfpackBase<double>;
-    template class PreconditionIfpack<double>;
-    template class PreconditionJacobi<double>;
-    template class PreconditionL1Jacobi<double>;
-    template class PreconditionL1GaussSeidel<double>;
-    template class PreconditionSOR<double>;
-    template class PreconditionSSOR<double>;
-    template class PreconditionBlockJacobi<double>;
-    template class PreconditionBlockSOR<double>;
-    template class PreconditionBlockSSOR<double>;
-    template class PreconditionChebyshev<double>;
-    template class PreconditionILU<double>;
-    template class PreconditionILUT<double>;
+    template class PreconditionBase<double, MemorySpace::Host>;
+    template class PreconditionIdentity<double, MemorySpace::Host>;
+    template class PreconditionIfpackBase<double, MemorySpace::Host>;
+    template class PreconditionIfpack<double, MemorySpace::Host>;
+    template class PreconditionJacobi<double, MemorySpace::Host>;
+    template class PreconditionL1Jacobi<double, MemorySpace::Host>;
+    template class PreconditionL1GaussSeidel<double, MemorySpace::Host>;
+    template class PreconditionSOR<double, MemorySpace::Host>;
+    template class PreconditionSSOR<double, MemorySpace::Host>;
+    template class PreconditionBlockJacobi<double, MemorySpace::Host>;
+    template class PreconditionBlockSOR<double, MemorySpace::Host>;
+    template class PreconditionBlockSSOR<double, MemorySpace::Host>;
+    template class PreconditionChebyshev<double, MemorySpace::Host>;
+    template class PreconditionILU<double, MemorySpace::Host>;
+    template class PreconditionILUT<double, MemorySpace::Host>;
+
+    template class PreconditionBase<double, MemorySpace::Default>;
+    template class PreconditionIdentity<double, MemorySpace::Default>;
+    template class PreconditionIfpackBase<double, MemorySpace::Default>;
+    template class PreconditionIfpack<double, MemorySpace::Default>;
+    template class PreconditionJacobi<double, MemorySpace::Default>;
+    template class PreconditionL1Jacobi<double, MemorySpace::Default>;
+    template class PreconditionL1GaussSeidel<double, MemorySpace::Default>;
+    template class PreconditionSOR<double, MemorySpace::Default>;
+    template class PreconditionSSOR<double, MemorySpace::Default>;
+    template class PreconditionBlockJacobi<double, MemorySpace::Default>;
+    template class PreconditionBlockSOR<double, MemorySpace::Default>;
+    template class PreconditionBlockSSOR<double, MemorySpace::Default>;
+    template class PreconditionChebyshev<double, MemorySpace::Default>;
+    template class PreconditionILU<double, MemorySpace::Default>;
+    template class PreconditionILUT<double, MemorySpace::Default>;
 #      endif
 #      ifdef DEAL_II_WITH_COMPLEX_VALUES
 #        ifdef HAVE_TPETRA_INST_COMPLEX_FLOAT
-    template class PreconditionBase<std::complex<float>>;
-    template class PreconditionIdentity<std::complex<float>>;
-    template class PreconditionIfpackBase<std::complex<float>>;
-    template class PreconditionIfpack<std::complex<float>>;
-    template class PreconditionJacobi<std::complex<float>>;
-    template class PreconditionL1Jacobi<std::complex<float>>;
-    template class PreconditionL1GaussSeidel<std::complex<float>>;
-    template class PreconditionSOR<std::complex<float>>;
-    template class PreconditionSSOR<std::complex<float>>;
-    template class PreconditionBlockJacobi<std::complex<float>>;
-    template class PreconditionBlockSOR<std::complex<float>>;
-    template class PreconditionBlockSSOR<std::complex<float>>;
-    template class PreconditionChebyshev<std::complex<float>>;
-    template class PreconditionILU<std::complex<float>>;
-    template class PreconditionILUT<std::complex<float>>;
+    template class PreconditionBase<std::complex<float>, MemorySpace::Host>;
+    template class PreconditionIdentity<std::complex<float>, MemorySpace::Host>;
+    template class PreconditionIfpackBase<std::complex<float>,
+                                          MemorySpace::Host>;
+    template class PreconditionIfpack<std::complex<float>, MemorySpace::Host>;
+    template class PreconditionJacobi<std::complex<float>, MemorySpace::Host>;
+    template class PreconditionL1Jacobi<std::complex<float>, MemorySpace::Host>;
+    template class PreconditionL1GaussSeidel<std::complex<float>,
+                                             MemorySpace::Host>;
+    template class PreconditionSOR<std::complex<float>, MemorySpace::Host>;
+    template class PreconditionSSOR<std::complex<float>, MemorySpace::Host>;
+    template class PreconditionBlockJacobi<std::complex<float>,
+                                           MemorySpace::Host>;
+    template class PreconditionBlockSOR<std::complex<float>, MemorySpace::Host>;
+    template class PreconditionBlockSSOR<std::complex<float>,
+                                         MemorySpace::Host>;
+    template class PreconditionChebyshev<std::complex<float>,
+                                         MemorySpace::Host>;
+    template class PreconditionILU<std::complex<float>, MemorySpace::Host>;
+    template class PreconditionILUT<std::complex<float>, MemorySpace::Host>;
+
+    template class PreconditionBase<std::complex<float>, MemorySpace::Default>;
+    template class PreconditionIdentity<std::complex<float>,
+                                        MemorySpace::Default>;
+    template class PreconditionIfpackBase<std::complex<float>,
+                                          MemorySpace::Default>;
+    template class PreconditionIfpack<std::complex<float>,
+                                      MemorySpace::Default>;
+    template class PreconditionJacobi<std::complex<float>,
+                                      MemorySpace::Default>;
+    template class PreconditionL1Jacobi<std::complex<float>,
+                                        MemorySpace::Default>;
+    template class PreconditionL1GaussSeidel<std::complex<float>,
+                                             MemorySpace::Default>;
+    template class PreconditionSOR<std::complex<float>, MemorySpace::Default>;
+    template class PreconditionSSOR<std::complex<float>, MemorySpace::Default>;
+    template class PreconditionBlockJacobi<std::complex<float>,
+                                           MemorySpace::Default>;
+    template class PreconditionBlockSOR<std::complex<float>,
+                                        MemorySpace::Default>;
+    template class PreconditionBlockSSOR<std::complex<float>,
+                                         MemorySpace::Default>;
+    template class PreconditionChebyshev<std::complex<float>,
+                                         MemorySpace::Default>;
+    template class PreconditionILU<std::complex<float>, MemorySpace::Default>;
+    template class PreconditionILUT<std::complex<float>, MemorySpace::Default>;
 #        endif
 #        ifdef HAVE_TPETRA_INST_COMPLEX_DOUBLE
-    template class PreconditionBase<std::complex<double>>;
-    template class PreconditionIdentity<std::complex<double>>;
-    template class PreconditionIfpackBase<std::complex<double>>;
-    template class PreconditionIfpack<std::complex<double>>;
-    template class PreconditionJacobi<std::complex<double>>;
-    template class PreconditionL1Jacobi<std::complex<double>>;
-    template class PreconditionL1GaussSeidel<std::complex<double>>;
-    template class PreconditionSOR<std::complex<double>>;
-    template class PreconditionSSOR<std::complex<double>>;
-    template class PreconditionBlockJacobi<std::complex<double>>;
-    template class PreconditionBlockSOR<std::complex<double>>;
-    template class PreconditionBlockSSOR<std::complex<double>>;
-    template class PreconditionChebyshev<std::complex<double>>;
-    template class PreconditionILU<std::complex<double>>;
-    template class PreconditionILUT<std::complex<double>>;
+    template class PreconditionBase<std::complex<double>, MemorySpace::Host>;
+    template class PreconditionIdentity<std::complex<double>,
+                                        MemorySpace::Host>;
+    template class PreconditionIfpackBase<std::complex<double>,
+                                          MemorySpace::Host>;
+    template class PreconditionIfpack<std::complex<double>, MemorySpace::Host>;
+    template class PreconditionJacobi<std::complex<double>, MemorySpace::Host>;
+    template class PreconditionL1Jacobi<std::complex<double>,
+                                        MemorySpace::Host>;
+    template class PreconditionL1GaussSeidel<std::complex<double>,
+                                             MemorySpace::Host>;
+    template class PreconditionSOR<std::complex<double>, MemorySpace::Host>;
+    template class PreconditionSSOR<std::complex<double>, MemorySpace::Host>;
+    template class PreconditionBlockJacobi<std::complex<double>,
+                                           MemorySpace::Host>;
+    template class PreconditionBlockSOR<std::complex<double>,
+                                        MemorySpace::Host>;
+    template class PreconditionBlockSSOR<std::complex<double>,
+                                         MemorySpace::Host>;
+    template class PreconditionChebyshev<std::complex<double>,
+                                         MemorySpace::Host>;
+    template class PreconditionILU<std::complex<double>, MemorySpace::Host>;
+    template class PreconditionILUT<std::complex<double>, MemorySpace::Host>;
+
+    template class PreconditionBase<std::complex<double>, MemorySpace::Default>;
+    template class PreconditionIdentity<std::complex<double>,
+                                        MemorySpace::Default>;
+    template class PreconditionIfpackBase<std::complex<double>,
+                                          MemorySpace::Default>;
+    template class PreconditionIfpack<std::complex<double>,
+                                      MemorySpace::Default>;
+    template class PreconditionJacobi<std::complex<double>,
+                                      MemorySpace::Default>;
+    template class PreconditionL1Jacobi<std::complex<double>,
+                                        MemorySpace::Default>;
+    template class PreconditionL1GaussSeidel<std::complex<double>,
+                                             MemorySpace::Default>;
+    template class PreconditionSOR<std::complex<double>, MemorySpace::Default>;
+    template class PreconditionSSOR<std::complex<double>, MemorySpace::Default>;
+    template class PreconditionBlockJacobi<std::complex<double>,
+                                           MemorySpace::Default>;
+    template class PreconditionBlockSOR<std::complex<double>,
+                                        MemorySpace::Default>;
+    template class PreconditionBlockSSOR<std::complex<double>,
+                                         MemorySpace::Default>;
+    template class PreconditionChebyshev<std::complex<double>,
+                                         MemorySpace::Default>;
+    template class PreconditionILU<std::complex<double>, MemorySpace::Default>;
+    template class PreconditionILUT<std::complex<double>, MemorySpace::Default>;
 #        endif
 #      endif
 
diff --git a/source/lac/trilinos_tpetra_solver_direct.cc b/source/lac/trilinos_tpetra_solver_direct.cc
index 095ef21a60..97a201a496 100644
--- a/source/lac/trilinos_tpetra_solver_direct.cc
+++ b/source/lac/trilinos_tpetra_solver_direct.cc
@@ -26,26 +26,42 @@ namespace LinearAlgebra
   namespace TpetraWrappers
   {
 #      ifdef HAVE_TPETRA_INST_FLOAT
-    template class SolverDirectBase<float>;
-    template class SolverDirect<float>;
-    template class SolverDirectKLU2<float>;
+    template class SolverDirectBase<float, MemorySpace::Host>;
+    template class SolverDirect<float, MemorySpace::Host>;
+    template class SolverDirectKLU2<float, MemorySpace::Host>;
+
+    template class SolverDirectBase<float, MemorySpace::Default>;
+    template class SolverDirect<float, MemorySpace::Default>;
+    template class SolverDirectKLU2<float, MemorySpace::Default>;
 #      endif
 
 #      ifdef HAVE_TPETRA_INST_DOUBLE
-    template class SolverDirectBase<double>;
-    template class SolverDirect<double>;
-    template class SolverDirectKLU2<double>;
+    template class SolverDirectBase<double, MemorySpace::Host>;
+    template class SolverDirect<double, MemorySpace::Host>;
+    template class SolverDirectKLU2<double, MemorySpace::Host>;
+
+    template class SolverDirectBase<double, MemorySpace::Default>;
+    template class SolverDirect<double, MemorySpace::Default>;
+    template class SolverDirectKLU2<double, MemorySpace::Default>;
 #      endif
 #      ifdef DEAL_II_WITH_COMPLEX_VALUES
 #        ifdef HAVE_TPETRA_INST_COMPLEX_FLOAT
-    template class SolverDirectBase<std::complex<float>>;
-    template class SolverDirect<std::complex<float>>;
-    template class SolverDirectKLU2<std::complex<float>>;
+    template class SolverDirectBase<std::complex<float>, MemorySpace::Host>;
+    template class SolverDirect<std::complex<float>, MemorySpace::Host>;
+    template class SolverDirectKLU2<std::complex<float>, MemorySpace::Host>;
+
+    template class SolverDirectBase<std::complex<float>, MemorySpace::Default>;
+    template class SolverDirect<std::complex<float>, MemorySpace::Default>;
+    template class SolverDirectKLU2<std::complex<float>, MemorySpace::Default>;
 #        endif
 #        ifdef HAVE_TPETRA_INST_COMPLEX_DOUBLE
-    template class SolverDirectBase<std::complex<double>>;
-    template class SolverDirect<std::complex<double>>;
-    template class SolverDirectKLU2<std::complex<double>>;
+    template class SolverDirectBase<std::complex<double>, MemorySpace::Host>;
+    template class SolverDirect<std::complex<double>, MemorySpace::Host>;
+    template class SolverDirectKLU2<std::complex<double>, MemorySpace::Host>;
+
+    template class SolverDirectBase<std::complex<double>, MemorySpace::Default>;
+    template class SolverDirect<std::complex<double>, MemorySpace::Default>;
+    template class SolverDirectKLU2<std::complex<double>, MemorySpace::Default>;
 #        endif
 #      endif
 
diff --git a/source/lac/trilinos_tpetra_sparse_matrix.cc b/source/lac/trilinos_tpetra_sparse_matrix.cc
index 3cc43d6226..61f113372d 100644
--- a/source/lac/trilinos_tpetra_sparse_matrix.cc
+++ b/source/lac/trilinos_tpetra_sparse_matrix.cc
@@ -26,17 +26,17 @@ namespace LinearAlgebra
 {
   namespace TpetraWrappers
   {
-    template class SparseMatrix<double>;
+    template class SparseMatrix<double, MemorySpace::Host>;
 
     template void
-    SparseMatrix<double>::reinit(
+    SparseMatrix<double, MemorySpace::Host>::reinit(
       const IndexSet                       &parallel_partitioning,
       const dealii::DynamicSparsityPattern &sparsity_pattern,
       const MPI_Comm                        communicator,
       const bool                            exchange_data);
 
     template void
-    SparseMatrix<double>::reinit(
+    SparseMatrix<double, MemorySpace::Host>::reinit(
       const IndexSet                       &row_parallel_partitioning,
       const IndexSet                       &col_parallel_partitioning,
       const dealii::DynamicSparsityPattern &sparsity_pattern,
@@ -44,39 +44,109 @@ namespace LinearAlgebra
       const bool                            exchange_data);
 
     template void
-    SparseMatrix<double>::reinit(const dealii::DynamicSparsityPattern &);
+    SparseMatrix<double, MemorySpace::Host>::reinit(
+      const dealii::DynamicSparsityPattern &);
 
     template void
-    SparseMatrix<double>::vmult(Vector<double>       &dst,
-                                const Vector<double> &src) const;
+    SparseMatrix<double, MemorySpace::Host>::vmult(
+      Vector<double, MemorySpace::Host>       &dst,
+      const Vector<double, MemorySpace::Host> &src) const;
 
     template void
-    SparseMatrix<double>::Tvmult(Vector<double>       &dst,
-                                 const Vector<double> &src) const;
+    SparseMatrix<double, MemorySpace::Host>::Tvmult(
+      Vector<double, MemorySpace::Host>       &dst,
+      const Vector<double, MemorySpace::Host> &src) const;
 
     template void
-    SparseMatrix<double>::vmult_add(Vector<double>       &dst,
-                                    const Vector<double> &src) const;
+    SparseMatrix<double, MemorySpace::Host>::vmult_add(
+      Vector<double, MemorySpace::Host>       &dst,
+      const Vector<double, MemorySpace::Host> &src) const;
 
     template void
-    SparseMatrix<double>::Tvmult_add(Vector<double>       &dst,
-                                     const Vector<double> &src) const;
+    SparseMatrix<double, MemorySpace::Host>::Tvmult_add(
+      Vector<double, MemorySpace::Host>       &dst,
+      const Vector<double, MemorySpace::Host> &src) const;
 
     template void
-    SparseMatrix<double>::vmult(::dealii::Vector<double>       &dst,
-                                const ::dealii::Vector<double> &src) const;
+    SparseMatrix<double, MemorySpace::Host>::vmult(
+      ::dealii::Vector<double>       &dst,
+      const ::dealii::Vector<double> &src) const;
 
     template void
-    SparseMatrix<double>::Tvmult(::dealii::Vector<double>       &dst,
-                                 const ::dealii::Vector<double> &src) const;
+    SparseMatrix<double, MemorySpace::Host>::Tvmult(
+      ::dealii::Vector<double>       &dst,
+      const ::dealii::Vector<double> &src) const;
 
     template void
-    SparseMatrix<double>::vmult_add(::dealii::Vector<double>       &dst,
-                                    const ::dealii::Vector<double> &src) const;
+    SparseMatrix<double, MemorySpace::Host>::vmult_add(
+      ::dealii::Vector<double>       &dst,
+      const ::dealii::Vector<double> &src) const;
 
     template void
-    SparseMatrix<double>::Tvmult_add(::dealii::Vector<double>       &dst,
-                                     const ::dealii::Vector<double> &src) const;
+    SparseMatrix<double, MemorySpace::Host>::Tvmult_add(
+      ::dealii::Vector<double>       &dst,
+      const ::dealii::Vector<double> &src) const;
+
+    template class SparseMatrix<double, MemorySpace::Default>;
+
+    template void
+    SparseMatrix<double, MemorySpace::Default>::reinit(
+      const IndexSet                       &parallel_partitioning,
+      const dealii::DynamicSparsityPattern &sparsity_pattern,
+      const MPI_Comm                        communicator,
+      const bool                            exchange_data);
+
+    template void
+    SparseMatrix<double, MemorySpace::Default>::reinit(
+      const IndexSet                       &row_parallel_partitioning,
+      const IndexSet                       &col_parallel_partitioning,
+      const dealii::DynamicSparsityPattern &sparsity_pattern,
+      const MPI_Comm                        communicator,
+      const bool                            exchange_data);
+
+    template void
+    SparseMatrix<double, MemorySpace::Default>::reinit(
+      const dealii::DynamicSparsityPattern &);
+
+    template void
+    SparseMatrix<double, MemorySpace::Default>::vmult(
+      Vector<double, MemorySpace::Default>       &dst,
+      const Vector<double, MemorySpace::Default> &src) const;
+
+    template void
+    SparseMatrix<double, MemorySpace::Default>::Tvmult(
+      Vector<double, MemorySpace::Default>       &dst,
+      const Vector<double, MemorySpace::Default> &src) const;
+
+    template void
+    SparseMatrix<double, MemorySpace::Default>::vmult_add(
+      Vector<double, MemorySpace::Default>       &dst,
+      const Vector<double, MemorySpace::Default> &src) const;
+
+    template void
+    SparseMatrix<double, MemorySpace::Default>::Tvmult_add(
+      Vector<double, MemorySpace::Default>       &dst,
+      const Vector<double, MemorySpace::Default> &src) const;
+
+    template void
+    SparseMatrix<double, MemorySpace::Default>::vmult(
+      ::dealii::Vector<double>       &dst,
+      const ::dealii::Vector<double> &src) const;
+
+    template void
+    SparseMatrix<double, MemorySpace::Default>::Tvmult(
+      ::dealii::Vector<double>       &dst,
+      const ::dealii::Vector<double> &src) const;
+
+    template void
+    SparseMatrix<double, MemorySpace::Default>::vmult_add(
+      ::dealii::Vector<double>       &dst,
+      const ::dealii::Vector<double> &src) const;
+
+    template void
+    SparseMatrix<double, MemorySpace::Default>::Tvmult_add(
+      ::dealii::Vector<double>       &dst,
+      const ::dealii::Vector<double> &src) const;
   } // namespace TpetraWrappers
 } // namespace LinearAlgebra
 #  endif // DOXYGEN
diff --git a/source/lac/trilinos_tpetra_vector.cc b/source/lac/trilinos_tpetra_vector.cc
index f7f2755692..5b1e4b113b 100644
--- a/source/lac/trilinos_tpetra_vector.cc
+++ b/source/lac/trilinos_tpetra_vector.cc
@@ -32,59 +32,92 @@ namespace LinearAlgebra
 
 #  ifdef HAVE_TPETRA_INST_FLOAT
 #    ifdef DEAL_II_HAVE_CXX20
-    static_assert(concepts::is_vector_space_vector<Vector<float>>);
+    static_assert(
+      concepts::is_vector_space_vector<Vector<float, MemorySpace::Host>>);
+    static_assert(
+      concepts::is_vector_space_vector<Vector<float, MemorySpace::Default>>);
 #    endif
-    template class Vector<float>;
-    template Vector<float> &
-    Vector<float>::operator=<float>(const dealii::Vector<float> &);
+    template class Vector<float, MemorySpace::Host>;
+    template class Vector<float, MemorySpace::Default>;
+    template Vector<float, MemorySpace::Host> &
+    Vector<float, MemorySpace::Host>::operator=
+      <float>(const dealii::Vector<float> &);
+    template Vector<float, MemorySpace::Default> &
+    Vector<float, MemorySpace::Default>::operator=
+      <float>(const dealii::Vector<float> &);
     namespace internal
     {
-      template class VectorReference<float>;
-    }
+      template class VectorReference<float, MemorySpace::Host>;
+      template class VectorReference<float, MemorySpace::Default>;
+    } // namespace internal
 #  endif
 
 #  ifdef HAVE_TPETRA_INST_DOUBLE
 #    ifdef DEAL_II_HAVE_CXX20
-    static_assert(concepts::is_vector_space_vector<Vector<double>>);
+    static_assert(
+      concepts::is_vector_space_vector<Vector<double, MemorySpace::Host>>);
+    static_assert(
+      concepts::is_vector_space_vector<Vector<double, MemorySpace::Default>>);
 #    endif
-    template class Vector<double>;
-    template Vector<double> &
-    Vector<double>::operator=<double>(const dealii::Vector<double> &);
+    template class Vector<double, MemorySpace::Host>;
+    template class Vector<double, MemorySpace::Default>;
+    template Vector<double, MemorySpace::Host> &
+    Vector<double, MemorySpace::Host>::operator=
+      <double>(const dealii::Vector<double> &);
+    template Vector<double, MemorySpace::Default> &
+    Vector<double, MemorySpace::Default>::operator=
+      <double>(const dealii::Vector<double> &);
     namespace internal
     {
-      template class VectorReference<double>;
-    }
+      template class VectorReference<double, MemorySpace::Host>;
+      template class VectorReference<double, MemorySpace::Default>;
+    } // namespace internal
 #  endif
 
 #  ifdef DEAL_II_WITH_COMPLEX_VALUES
 #    ifdef HAVE_TPETRA_INST_COMPLEX_FLOAT
 #      ifdef DEAL_II_HAVE_CXX20
-    static_assert(
-      concepts::is_vector_space_vector<Vector<std::complex<float>>>);
+    static_assert(concepts::is_vector_space_vector<
+                  Vector<std::complex<float>, MemorySpace::Host>>);
+    static_assert(concepts::is_vector_space_vector<
+                  Vector<std::complex<float>, MemorySpace::Default>>);
 #      endif
-    template class Vector<std::complex<float>>;
-    template Vector<std::complex<float>> &
-    Vector<std::complex<float>>::operator=
+    template class Vector<std::complex<float>, MemorySpace::Host>;
+    template class Vector<std::complex<float>, MemorySpace::Default>;
+    template Vector<std::complex<float>, MemorySpace::Host> &
+    Vector<std::complex<float>, MemorySpace::Host>::operator=
+      <std::complex<float>>(const dealii::Vector<std::complex<float>> &);
+    template Vector<std::complex<float>, MemorySpace::Default> &
+    Vector<std::complex<float>, MemorySpace::Default>::operator=
       <std::complex<float>>(const dealii::Vector<std::complex<float>> &);
     namespace internal
     {
-      template class VectorReference<std::complex<float>>;
-    }
+      template class VectorReference<std::complex<float>, MemorySpace::Host>;
+      template class VectorReference<std::complex<float>, MemorySpace::Default>;
+    } // namespace internal
 #    endif
 
 #    ifdef HAVE_TPETRA_INST_COMPLEX_DOUBLE
 #      ifdef DEAL_II_HAVE_CXX20
-    static_assert(
-      concepts::is_vector_space_vector<Vector<std::complex<double>>>);
+    static_assert(concepts::is_vector_space_vector<
+                  Vector<std::complex<double>, MemorySpace::Host>>);
+    static_assert(concepts::is_vector_space_vector<
+                  Vector<std::complex<double>, MemorySpace::Default>>);
 #      endif
-    template class Vector<std::complex<double>>;
-    template Vector<std::complex<double>> &
-    Vector<std::complex<double>>::operator=
+    template class Vector<std::complex<double>, MemorySpace::Host>;
+    template class Vector<std::complex<double>, MemorySpace::Default>;
+    template Vector<std::complex<double>, MemorySpace::Host> &
+    Vector<std::complex<double>, MemorySpace::Host>::operator=
+      <std::complex<double>>(const dealii::Vector<std::complex<double>> &);
+    template Vector<std::complex<double>, MemorySpace::Default> &
+    Vector<std::complex<double>, MemorySpace::Default>::operator=
       <std::complex<double>>(const dealii::Vector<std::complex<double>> &);
     namespace internal
     {
-      template class VectorReference<std::complex<double>>;
-    }
+      template class VectorReference<std::complex<double>, MemorySpace::Host>;
+      template class VectorReference<std::complex<double>,
+                                     MemorySpace::Default>;
+    } // namespace internal
 #    endif
 #  endif
   } // namespace TpetraWrappers
diff --git a/source/lac/vector.cc b/source/lac/vector.cc
index d50a270c0e..7bbed57d59 100644
--- a/source/lac/vector.cc
+++ b/source/lac/vector.cc
@@ -85,35 +85,63 @@ Vector<int>::lp_norm(const real_type) const
 #ifdef DEAL_II_TRILINOS_WITH_TPETRA
 #  ifdef HAVE_TPETRA_INST_FLOAT
 template Vector<float>::Vector(
-  const LinearAlgebra::TpetraWrappers::Vector<float> &);
+  const LinearAlgebra::TpetraWrappers::Vector<float, MemorySpace::Host> &);
+template Vector<float>::Vector(
+  const LinearAlgebra::TpetraWrappers::Vector<float, MemorySpace::Default> &);
+template Vector<float> &
+Vector<float>::operator=<float>(
+  const LinearAlgebra::TpetraWrappers::Vector<float, MemorySpace::Host> &);
 template Vector<float> &
-Vector<float>::operator=
-  <float>(const LinearAlgebra::TpetraWrappers::Vector<float> &);
+Vector<float>::operator=<float>(
+  const LinearAlgebra::TpetraWrappers::Vector<float, MemorySpace::Default> &);
 #  endif
 
 #  ifdef HAVE_TPETRA_INST_DOUBLE
 template Vector<double>::Vector(
-  const LinearAlgebra::TpetraWrappers::Vector<double> &);
+  const LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Host> &);
+template Vector<double>::Vector(
+  const LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &);
 template Vector<double> &
-Vector<double>::operator=
-  <double>(const LinearAlgebra::TpetraWrappers::Vector<double> &);
+Vector<double>::operator=<double>(
+  const LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Host> &);
+template Vector<double> &
+Vector<double>::operator=<double>(
+  const LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &);
 #  endif
 
 #  ifdef DEAL_II_WITH_COMPLEX_VALUES
 #    ifdef HAVE_TPETRA_INST_COMPLEX_FLOAT
 template Vector<std::complex<float>>::Vector(
-  const LinearAlgebra::TpetraWrappers::Vector<std::complex<float>> &);
+  const LinearAlgebra::TpetraWrappers::Vector<std::complex<float>,
+                                              MemorySpace::Host> &);
+template Vector<std::complex<float>>::Vector(
+  const LinearAlgebra::TpetraWrappers::Vector<std::complex<float>,
+                                              MemorySpace::Default> &);
 template Vector<std::complex<float>> &
 Vector<std::complex<float>>::operator=<std::complex<float>>(
-  const LinearAlgebra::TpetraWrappers::Vector<std::complex<float>> &);
+  const LinearAlgebra::TpetraWrappers::Vector<std::complex<float>,
+                                              MemorySpace::Host> &);
+template Vector<std::complex<float>> &
+Vector<std::complex<float>>::operator=<std::complex<float>>(
+  const LinearAlgebra::TpetraWrappers::Vector<std::complex<float>,
+                                              MemorySpace::Default> &);
 #    endif
 
 #    ifdef HAVE_TPETRA_INST_COMPLEX_DOUBLE
 template Vector<std::complex<double>>::Vector(
-  const LinearAlgebra::TpetraWrappers::Vector<std::complex<double>> &);
+  const LinearAlgebra::TpetraWrappers::Vector<std::complex<double>,
+                                              MemorySpace::Host> &);
+template Vector<std::complex<double>>::Vector(
+  const LinearAlgebra::TpetraWrappers::Vector<std::complex<double>,
+                                              MemorySpace::Default> &);
+template Vector<std::complex<double>> &
+Vector<std::complex<double>>::operator=<std::complex<double>>(
+  const LinearAlgebra::TpetraWrappers::Vector<std::complex<double>,
+                                              MemorySpace::Host> &);
 template Vector<std::complex<double>> &
 Vector<std::complex<double>>::operator=<std::complex<double>>(
-  const LinearAlgebra::TpetraWrappers::Vector<std::complex<double>> &);
+  const LinearAlgebra::TpetraWrappers::Vector<std::complex<double>,
+                                              MemorySpace::Default> &);
 #    endif
 #  endif
 #endif
diff --git a/tests/trilinos_tpetra/01.cc b/tests/trilinos_tpetra/01.cc
index e8ca662ded..a59a0cf1b2 100644
--- a/tests/trilinos_tpetra/01.cc
+++ b/tests/trilinos_tpetra/01.cc
@@ -27,7 +27,8 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::SparseMatrix<double> &m)
+test(
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> &m)
 {
   // first set a few entries
   for (unsigned int i = 0; i < m.m(); ++i)
@@ -66,7 +67,9 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::SparseMatrix<double> m(5U, 5U, 3U);
+        LinearAlgebra::TpetraWrappers::SparseMatrix<double,
+                                                    MemorySpace::Default>
+          m(5U, 5U, 3U);
         test(m);
       }
     }
diff --git a/tests/trilinos_tpetra/02.cc b/tests/trilinos_tpetra/02.cc
index cafa94b0d7..c3d6529a62 100644
--- a/tests/trilinos_tpetra/02.cc
+++ b/tests/trilinos_tpetra/02.cc
@@ -27,7 +27,8 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::SparseMatrix<double> &m)
+test(
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> &m)
 {
   // first set a few entries
   for (unsigned int i = 0; i < m.m(); ++i)
@@ -72,7 +73,9 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::SparseMatrix<double> m(5U, 5U, 3U);
+        LinearAlgebra::TpetraWrappers::SparseMatrix<double,
+                                                    MemorySpace::Default>
+          m(5U, 5U, 3U);
         test(m);
       }
     }
diff --git a/tests/trilinos_tpetra/03a.cc b/tests/trilinos_tpetra/03a.cc
index 104de096a2..c0716ea928 100644
--- a/tests/trilinos_tpetra/03a.cc
+++ b/tests/trilinos_tpetra/03a.cc
@@ -31,7 +31,8 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::SparseMatrix<double> &m)
+test(
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> &m)
 {
   // first set a few entries
   for (unsigned int i = 0; i < m.m(); ++i)
@@ -76,7 +77,9 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::SparseMatrix<double> m(5U, 5U, 3U);
+        LinearAlgebra::TpetraWrappers::SparseMatrix<double,
+                                                    MemorySpace::Default>
+          m(5U, 5U, 3U);
         test(m);
       }
     }
diff --git a/tests/trilinos_tpetra/04.cc b/tests/trilinos_tpetra/04.cc
index d488188e95..11c8f4e6ac 100644
--- a/tests/trilinos_tpetra/04.cc
+++ b/tests/trilinos_tpetra/04.cc
@@ -26,7 +26,8 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::SparseMatrix<double> &m)
+test(
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> &m)
 {
   AssertThrow(m.m() == 5, ExcInternalError());
   AssertThrow(m.n() == 5, ExcInternalError());
@@ -48,7 +49,9 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::SparseMatrix<double> m(5U, 5U, 3U);
+        LinearAlgebra::TpetraWrappers::SparseMatrix<double,
+                                                    MemorySpace::Default>
+          m(5U, 5U, 3U);
         test(m);
       }
     }
diff --git a/tests/trilinos_tpetra/05.cc b/tests/trilinos_tpetra/05.cc
index 42067d5840..def86ef8cd 100644
--- a/tests/trilinos_tpetra/05.cc
+++ b/tests/trilinos_tpetra/05.cc
@@ -27,7 +27,8 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::SparseMatrix<double> &m)
+test(
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> &m)
 {
   // first set a few entries. count how many
   // entries we have
@@ -62,7 +63,9 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::SparseMatrix<double> m(5U, 5U, 3U);
+        LinearAlgebra::TpetraWrappers::SparseMatrix<double,
+                                                    MemorySpace::Default>
+          m(5U, 5U, 3U);
         test(m);
       }
     }
diff --git a/tests/trilinos_tpetra/08.cc b/tests/trilinos_tpetra/08.cc
index 6b08722953..40bd1f3fb9 100644
--- a/tests/trilinos_tpetra/08.cc
+++ b/tests/trilinos_tpetra/08.cc
@@ -25,7 +25,8 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::SparseMatrix<double> &m)
+test(
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> &m)
 {
   // first set a few entries. count how many
   // entries we have
@@ -63,7 +64,9 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::SparseMatrix<double> m(5U, 5U, 3U);
+        LinearAlgebra::TpetraWrappers::SparseMatrix<double,
+                                                    MemorySpace::Default>
+          m(5U, 5U, 3U);
         test(m);
       }
     }
diff --git a/tests/trilinos_tpetra/09.cc b/tests/trilinos_tpetra/09.cc
index acef4d7c63..2e8dab2246 100644
--- a/tests/trilinos_tpetra/09.cc
+++ b/tests/trilinos_tpetra/09.cc
@@ -26,7 +26,8 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::SparseMatrix<double> &m)
+test(
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> &m)
 {
   // first set a few entries
   for (unsigned int i = 0; i < m.m(); ++i)
@@ -71,7 +72,9 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::SparseMatrix<double> m(5U, 5U, 3U);
+        LinearAlgebra::TpetraWrappers::SparseMatrix<double,
+                                                    MemorySpace::Default>
+          m(5U, 5U, 3U);
         test(m);
       }
     }
diff --git a/tests/trilinos_tpetra/10.cc b/tests/trilinos_tpetra/10.cc
index 37bcbe50ca..7779511fd8 100644
--- a/tests/trilinos_tpetra/10.cc
+++ b/tests/trilinos_tpetra/10.cc
@@ -26,7 +26,8 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::SparseMatrix<double> &m)
+test(
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> &m)
 {
   // first set a few entries
   for (unsigned int i = 0; i < m.m(); ++i)
@@ -71,7 +72,9 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::SparseMatrix<double> m(5U, 5U, 3U);
+        LinearAlgebra::TpetraWrappers::SparseMatrix<double,
+                                                    MemorySpace::Default>
+          m(5U, 5U, 3U);
         test(m);
       }
     }
diff --git a/tests/trilinos_tpetra/11.cc b/tests/trilinos_tpetra/11.cc
index 8cddd57f5a..b233008443 100644
--- a/tests/trilinos_tpetra/11.cc
+++ b/tests/trilinos_tpetra/11.cc
@@ -26,7 +26,7 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v)
 {
   // set only certain elements of the vector
   for (unsigned int i = 0; i < v.size(); i += 1 + i)
@@ -53,7 +53,7 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v);
       }
diff --git a/tests/trilinos_tpetra/12.cc b/tests/trilinos_tpetra/12.cc
index c23bc5eddc..c74507af58 100644
--- a/tests/trilinos_tpetra/12.cc
+++ b/tests/trilinos_tpetra/12.cc
@@ -27,7 +27,7 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v)
 {
   // set only certain elements of the
   // vector. have a bit pattern of where we
@@ -65,7 +65,7 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v);
       }
diff --git a/tests/trilinos_tpetra/13.cc b/tests/trilinos_tpetra/13.cc
index 6bbfbdb607..abe64983db 100644
--- a/tests/trilinos_tpetra/13.cc
+++ b/tests/trilinos_tpetra/13.cc
@@ -27,7 +27,7 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v)
 {
   // set only certain elements of the
   // vector. have a bit pattern of where we
@@ -65,7 +65,7 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v);
       }
diff --git a/tests/trilinos_tpetra/14.cc b/tests/trilinos_tpetra/14.cc
index 17508e1274..1c98a9fa16 100644
--- a/tests/trilinos_tpetra/14.cc
+++ b/tests/trilinos_tpetra/14.cc
@@ -41,7 +41,7 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v)
 {
   // set only certain elements of the
   // vector. have a bit pattern of where we
@@ -83,7 +83,7 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v);
       }
diff --git a/tests/trilinos_tpetra/15.cc b/tests/trilinos_tpetra/15.cc
index b1d6a676d9..8fd451ec5b 100644
--- a/tests/trilinos_tpetra/15.cc
+++ b/tests/trilinos_tpetra/15.cc
@@ -44,7 +44,7 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v)
 {
   // set only certain elements of the
   // vector. have a bit pattern of where we
@@ -83,7 +83,7 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v);
       }
diff --git a/tests/trilinos_tpetra/16.cc b/tests/trilinos_tpetra/16.cc
index ea9e199435..9d8ca8108f 100644
--- a/tests/trilinos_tpetra/16.cc
+++ b/tests/trilinos_tpetra/16.cc
@@ -44,7 +44,7 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v)
 {
   // set only certain elements of the
   // vector. have a bit pattern of where we
@@ -84,7 +84,7 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v);
       }
diff --git a/tests/trilinos_tpetra/17.cc b/tests/trilinos_tpetra/17.cc
index d511d00a97..ea5e6639c4 100644
--- a/tests/trilinos_tpetra/17.cc
+++ b/tests/trilinos_tpetra/17.cc
@@ -27,7 +27,7 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v)
 {
   // set some elements of the vector
   double norm = 0;
@@ -58,7 +58,7 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v);
       }
diff --git a/tests/trilinos_tpetra/18.cc b/tests/trilinos_tpetra/18.cc
index bafa7bbb71..4c131958ea 100644
--- a/tests/trilinos_tpetra/18.cc
+++ b/tests/trilinos_tpetra/18.cc
@@ -27,7 +27,7 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v)
 {
   // set some elements of the vector
   TrilinosScalar norm = 0;
@@ -59,7 +59,7 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v);
       }
diff --git a/tests/trilinos_tpetra/19.cc b/tests/trilinos_tpetra/19.cc
index fac8224eff..1ba2ce6e8c 100644
--- a/tests/trilinos_tpetra/19.cc
+++ b/tests/trilinos_tpetra/19.cc
@@ -27,7 +27,7 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v)
 {
   // set some elements of the vector
   double norm = 0;
@@ -58,7 +58,7 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v);
       }
diff --git a/tests/trilinos_tpetra/20.cc b/tests/trilinos_tpetra/20.cc
index 95d12048f9..54f399c6fb 100644
--- a/tests/trilinos_tpetra/20.cc
+++ b/tests/trilinos_tpetra/20.cc
@@ -27,7 +27,7 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v)
 {
   // set only certain elements of the
   // vector. have a bit pattern of where we
@@ -66,7 +66,7 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v);
       }
diff --git a/tests/trilinos_tpetra/21.cc b/tests/trilinos_tpetra/21.cc
index c619e7f137..683028a9b9 100644
--- a/tests/trilinos_tpetra/21.cc
+++ b/tests/trilinos_tpetra/21.cc
@@ -27,7 +27,7 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v)
 {
   // set only certain elements of the
   // vector. have a bit pattern of where we
@@ -66,7 +66,7 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v);
       }
diff --git a/tests/trilinos_tpetra/22.cc b/tests/trilinos_tpetra/22.cc
index 941f89099d..2cedfe4773 100644
--- a/tests/trilinos_tpetra/22.cc
+++ b/tests/trilinos_tpetra/22.cc
@@ -28,8 +28,8 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v,
-     LinearAlgebra::TpetraWrappers::Vector<double> &w)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v,
+     LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &w)
 {
   // set only certain elements of each
   // vector, but disjoint sets of elements
@@ -61,9 +61,9 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
-        LinearAlgebra::TpetraWrappers::Vector<double> w;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> w;
         w.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v, w);
       }
diff --git a/tests/trilinos_tpetra/23.cc b/tests/trilinos_tpetra/23.cc
index 20681b2584..07b98841cb 100644
--- a/tests/trilinos_tpetra/23.cc
+++ b/tests/trilinos_tpetra/23.cc
@@ -28,8 +28,8 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v,
-     LinearAlgebra::TpetraWrappers::Vector<double> &w)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v,
+     LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &w)
 {
   // set only certain elements of each
   // vector, and record the expected scalar
@@ -68,9 +68,9 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
-        LinearAlgebra::TpetraWrappers::Vector<double> w;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> w;
         w.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v, w);
       }
diff --git a/tests/trilinos_tpetra/24.cc b/tests/trilinos_tpetra/24.cc
index bf0dd6c537..41d919cdd1 100644
--- a/tests/trilinos_tpetra/24.cc
+++ b/tests/trilinos_tpetra/24.cc
@@ -30,7 +30,7 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v)
 {
   // set some entries of the vector
   for (unsigned int i = 0; i < v.size(); ++i)
@@ -62,7 +62,7 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v);
       }
diff --git a/tests/trilinos_tpetra/27.cc b/tests/trilinos_tpetra/27.cc
index 4844bb9108..b2f29c4f71 100644
--- a/tests/trilinos_tpetra/27.cc
+++ b/tests/trilinos_tpetra/27.cc
@@ -27,7 +27,7 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v)
 {
   // set some entries of the vector
   for (unsigned int i = 0; i < v.size(); ++i)
@@ -36,7 +36,7 @@ test(LinearAlgebra::TpetraWrappers::Vector<double> &v)
   v.compress(VectorOperation::insert);
 
   // then copy it
-  LinearAlgebra::TpetraWrappers::Vector<double> w;
+  LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> w;
   w.reinit(complete_index_set(v.size()), MPI_COMM_WORLD);
   w = v;
 
@@ -64,7 +64,7 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v);
       }
diff --git a/tests/trilinos_tpetra/28.cc b/tests/trilinos_tpetra/28.cc
index c15f6307e9..ecc9bb5ac3 100644
--- a/tests/trilinos_tpetra/28.cc
+++ b/tests/trilinos_tpetra/28.cc
@@ -28,7 +28,7 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v)
 {
   // set some entries of the vector
   for (unsigned int i = 0; i < v.size(); ++i)
@@ -38,7 +38,7 @@ test(LinearAlgebra::TpetraWrappers::Vector<double> &v)
 
   // then copy it to a vector of different
   // size
-  LinearAlgebra::TpetraWrappers::Vector<double> w;
+  LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> w;
   w.reinit(complete_index_set(1), MPI_COMM_WORLD);
   w = v;
 
@@ -66,7 +66,7 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v);
       }
diff --git a/tests/trilinos_tpetra/29.cc b/tests/trilinos_tpetra/29.cc
index 04ab5f280f..9063d439c2 100644
--- a/tests/trilinos_tpetra/29.cc
+++ b/tests/trilinos_tpetra/29.cc
@@ -27,7 +27,7 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v)
 {
   v.reinit(complete_index_set(13), MPI_COMM_WORLD);
 
@@ -50,7 +50,7 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v);
       }
diff --git a/tests/trilinos_tpetra/30.cc b/tests/trilinos_tpetra/30.cc
index 5af7dd695c..703606c115 100644
--- a/tests/trilinos_tpetra/30.cc
+++ b/tests/trilinos_tpetra/30.cc
@@ -27,7 +27,7 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v)
 {
   // set some entries of the vector
   for (unsigned int i = 0; i < v.size(); ++i)
@@ -58,7 +58,7 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v);
       }
diff --git a/tests/trilinos_tpetra/31.cc b/tests/trilinos_tpetra/31.cc
index 8035958cc8..6e2d35ec16 100644
--- a/tests/trilinos_tpetra/31.cc
+++ b/tests/trilinos_tpetra/31.cc
@@ -27,7 +27,7 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v)
 {
   // set some elements of the vector
   TrilinosScalar norm = 0;
@@ -59,7 +59,7 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v);
       }
diff --git a/tests/trilinos_tpetra/32.cc b/tests/trilinos_tpetra/32.cc
index 608052c97b..56436456e6 100644
--- a/tests/trilinos_tpetra/32.cc
+++ b/tests/trilinos_tpetra/32.cc
@@ -27,7 +27,7 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v)
 {
   // set some elements of the vector
   TrilinosScalar sum = 0;
@@ -60,7 +60,7 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v);
       }
diff --git a/tests/trilinos_tpetra/34.cc b/tests/trilinos_tpetra/34.cc
index d1c5b853f4..3cf6dcf963 100644
--- a/tests/trilinos_tpetra/34.cc
+++ b/tests/trilinos_tpetra/34.cc
@@ -27,7 +27,7 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v)
 {
   // set some elements of the vector
   double sum = 0;
@@ -61,7 +61,7 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v);
       }
diff --git a/tests/trilinos_tpetra/35.cc b/tests/trilinos_tpetra/35.cc
index 0901b9a13b..57b3e4aa53 100644
--- a/tests/trilinos_tpetra/35.cc
+++ b/tests/trilinos_tpetra/35.cc
@@ -27,8 +27,8 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v,
-     LinearAlgebra::TpetraWrappers::Vector<double> &w)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v,
+     LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &w)
 {
   // set only certain elements of each
   // vector
@@ -77,9 +77,9 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
-        LinearAlgebra::TpetraWrappers::Vector<double> w;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> w;
         w.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v, w);
       }
diff --git a/tests/trilinos_tpetra/36.cc b/tests/trilinos_tpetra/36.cc
index 8a90ae46fe..a309dd99dc 100644
--- a/tests/trilinos_tpetra/36.cc
+++ b/tests/trilinos_tpetra/36.cc
@@ -27,8 +27,8 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v,
-     LinearAlgebra::TpetraWrappers::Vector<double> &w)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v,
+     LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &w)
 {
   // set only certain elements of each
   // vector
@@ -77,9 +77,9 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
-        LinearAlgebra::TpetraWrappers::Vector<double> w;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> w;
         w.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v, w);
       }
diff --git a/tests/trilinos_tpetra/37.cc b/tests/trilinos_tpetra/37.cc
index cefd5b163a..e2a7574b20 100644
--- a/tests/trilinos_tpetra/37.cc
+++ b/tests/trilinos_tpetra/37.cc
@@ -27,7 +27,7 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v)
 {
   for (unsigned int i = 0; i < v.size(); ++i)
     v(i) = i;
@@ -57,7 +57,7 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v);
       }
diff --git a/tests/trilinos_tpetra/39.cc b/tests/trilinos_tpetra/39.cc
index 6e388e90d6..b660b5bbe6 100644
--- a/tests/trilinos_tpetra/39.cc
+++ b/tests/trilinos_tpetra/39.cc
@@ -27,8 +27,8 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v,
-     LinearAlgebra::TpetraWrappers::Vector<double> &w)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v,
+     LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &w)
 {
   for (unsigned int i = 0; i < v.size(); ++i)
     {
@@ -65,9 +65,9 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
-        LinearAlgebra::TpetraWrappers::Vector<double> w;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> w;
         w.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v, w);
       }
diff --git a/tests/trilinos_tpetra/40.cc b/tests/trilinos_tpetra/40.cc
index 0a1253e86a..37d65edd31 100644
--- a/tests/trilinos_tpetra/40.cc
+++ b/tests/trilinos_tpetra/40.cc
@@ -27,9 +27,9 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v,
-     LinearAlgebra::TpetraWrappers::Vector<double> &w,
-     LinearAlgebra::TpetraWrappers::Vector<double> &x)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v,
+     LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &w,
+     LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &x)
 {
   for (unsigned int i = 0; i < v.size(); ++i)
     {
@@ -69,11 +69,11 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
-        LinearAlgebra::TpetraWrappers::Vector<double> w;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> w;
         w.reinit(complete_index_set(100), MPI_COMM_WORLD);
-        LinearAlgebra::TpetraWrappers::Vector<double> x;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> x;
         x.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v, w, x);
       }
diff --git a/tests/trilinos_tpetra/42.cc b/tests/trilinos_tpetra/42.cc
index 0afc4ce54d..15f2fb0fa5 100644
--- a/tests/trilinos_tpetra/42.cc
+++ b/tests/trilinos_tpetra/42.cc
@@ -28,8 +28,8 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v,
-     LinearAlgebra::TpetraWrappers::Vector<double> &w)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v,
+     LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &w)
 {
   for (unsigned int i = 0; i < v.size(); ++i)
     {
@@ -66,9 +66,9 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
-        LinearAlgebra::TpetraWrappers::Vector<double> w;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> w;
         w.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v, w);
       }
diff --git a/tests/trilinos_tpetra/45.cc b/tests/trilinos_tpetra/45.cc
index 30e741349c..8d0b684c82 100644
--- a/tests/trilinos_tpetra/45.cc
+++ b/tests/trilinos_tpetra/45.cc
@@ -27,8 +27,8 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v,
-     LinearAlgebra::TpetraWrappers::Vector<double> &w)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v,
+     LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &w)
 {
   for (unsigned int i = 0; i < v.size(); ++i)
     {
@@ -65,9 +65,9 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
-        LinearAlgebra::TpetraWrappers::Vector<double> w;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> w;
         w.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v, w);
       }
diff --git a/tests/trilinos_tpetra/46.cc b/tests/trilinos_tpetra/46.cc
index 80c500fdf5..3f500cff23 100644
--- a/tests/trilinos_tpetra/46.cc
+++ b/tests/trilinos_tpetra/46.cc
@@ -27,8 +27,8 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v,
-     LinearAlgebra::TpetraWrappers::Vector<double> &w)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v,
+     LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &w)
 {
   for (unsigned int i = 0; i < v.size(); ++i)
     {
@@ -65,9 +65,9 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
-        LinearAlgebra::TpetraWrappers::Vector<double> w;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> w;
         w.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v, w);
       }
diff --git a/tests/trilinos_tpetra/49.cc b/tests/trilinos_tpetra/49.cc
index 56e31989b2..8f50f0ac58 100644
--- a/tests/trilinos_tpetra/49.cc
+++ b/tests/trilinos_tpetra/49.cc
@@ -29,7 +29,7 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v)
 {
   dealii::Vector<TrilinosScalar> w(v.size());
 
@@ -63,7 +63,7 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v);
       }
diff --git a/tests/trilinos_tpetra/51.cc b/tests/trilinos_tpetra/51.cc
index d1419e9157..382369a41b 100644
--- a/tests/trilinos_tpetra/51.cc
+++ b/tests/trilinos_tpetra/51.cc
@@ -28,7 +28,7 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v)
 {
   // set some entries of the vector
   for (unsigned int i = 0; i < v.size(); ++i)
@@ -37,7 +37,7 @@ test(LinearAlgebra::TpetraWrappers::Vector<double> &v)
   v.compress(VectorOperation::insert);
 
   // then copy it
-  LinearAlgebra::TpetraWrappers::Vector<double> w(v);
+  LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> w(v);
   w.compress(VectorOperation::insert);
 
   // make sure they're equal
@@ -64,7 +64,7 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v);
       }
diff --git a/tests/trilinos_tpetra/52.cc b/tests/trilinos_tpetra/52.cc
index 593d4e4b78..21fa7182e6 100644
--- a/tests/trilinos_tpetra/52.cc
+++ b/tests/trilinos_tpetra/52.cc
@@ -28,7 +28,8 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::SparseMatrix<double> &m)
+test(
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> &m)
 {
   // first set a few entries
   for (unsigned int i = 0; i < m.m(); ++i)
@@ -70,10 +71,11 @@ main(int argc, char **argv)
       {
         std::vector<unsigned int> row_lengths(5, 3U);
         row_lengths.back() = 2;
-        LinearAlgebra::TpetraWrappers::SparseMatrix<double> m(
-          static_cast<types::global_dof_index>(5U),
-          static_cast<types::global_dof_index>(5U),
-          row_lengths);
+        LinearAlgebra::TpetraWrappers::SparseMatrix<double,
+                                                    MemorySpace::Default>
+          m(static_cast<types::global_dof_index>(5U),
+            static_cast<types::global_dof_index>(5U),
+            row_lengths);
         test(m);
       }
     }
diff --git a/tests/trilinos_tpetra/53.cc b/tests/trilinos_tpetra/53.cc
index cad10618b4..07e21aa2f4 100644
--- a/tests/trilinos_tpetra/53.cc
+++ b/tests/trilinos_tpetra/53.cc
@@ -28,7 +28,7 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v)
 {
   // set only certain elements of the
   // vector. have a bit pattern of where we
@@ -69,7 +69,7 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v);
       }
diff --git a/tests/trilinos_tpetra/54.cc b/tests/trilinos_tpetra/54.cc
index dd0d9b88ef..e1a46fa3ad 100644
--- a/tests/trilinos_tpetra/54.cc
+++ b/tests/trilinos_tpetra/54.cc
@@ -28,7 +28,7 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v)
 {
   // set only certain elements of the
   // vector. have a bit pattern of where we
@@ -69,7 +69,7 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v);
       }
diff --git a/tests/trilinos_tpetra/55.cc b/tests/trilinos_tpetra/55.cc
index 2eec167758..2a6e01690c 100644
--- a/tests/trilinos_tpetra/55.cc
+++ b/tests/trilinos_tpetra/55.cc
@@ -28,7 +28,7 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v)
 {
   // set only certain elements of the
   // vector. have a bit pattern of where we
@@ -68,7 +68,7 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v);
       }
diff --git a/tests/trilinos_tpetra/56.cc b/tests/trilinos_tpetra/56.cc
index 1a02312904..754c5fbb11 100644
--- a/tests/trilinos_tpetra/56.cc
+++ b/tests/trilinos_tpetra/56.cc
@@ -28,7 +28,7 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v)
 {
   // set only certain elements of the
   // vector. have a bit pattern of where we
@@ -68,7 +68,7 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v);
       }
diff --git a/tests/trilinos_tpetra/57.cc b/tests/trilinos_tpetra/57.cc
index 9612b8cb33..b929352d69 100644
--- a/tests/trilinos_tpetra/57.cc
+++ b/tests/trilinos_tpetra/57.cc
@@ -27,7 +27,7 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v)
 {
   // set only certain elements of the
   // vector. they are all positive
@@ -66,7 +66,7 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v);
       }
diff --git a/tests/trilinos_tpetra/62.cc b/tests/trilinos_tpetra/62.cc
index 870543ae89..c67f5c164d 100644
--- a/tests/trilinos_tpetra/62.cc
+++ b/tests/trilinos_tpetra/62.cc
@@ -28,7 +28,8 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::SparseMatrix<double> &m)
+test(
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> &m)
 {
   AssertThrow(m.m() != 0, ExcInternalError());
   AssertThrow(m.n() != 0, ExcInternalError());
@@ -55,7 +56,9 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::SparseMatrix<double> v(100U, 100U, 5U);
+        LinearAlgebra::TpetraWrappers::SparseMatrix<double,
+                                                    MemorySpace::Default>
+          v(100U, 100U, 5U);
         test(v);
       }
     }
diff --git a/tests/trilinos_tpetra/63.cc b/tests/trilinos_tpetra/63.cc
index f8d8688642..cfe4031d0e 100644
--- a/tests/trilinos_tpetra/63.cc
+++ b/tests/trilinos_tpetra/63.cc
@@ -26,7 +26,8 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::SparseMatrix<double> &m)
+test(
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> &m)
 {
   AssertThrow(m.m() == 100, ExcInternalError());
   AssertThrow(m.n() == 100, ExcInternalError());
@@ -53,7 +54,9 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::SparseMatrix<double> v(100U, 100U, 5U);
+        LinearAlgebra::TpetraWrappers::SparseMatrix<double,
+                                                    MemorySpace::Default>
+          v(100U, 100U, 5U);
         test(v);
       }
     }
diff --git a/tests/trilinos_tpetra/64.cc b/tests/trilinos_tpetra/64.cc
index a294ae4a45..2a93c37e10 100644
--- a/tests/trilinos_tpetra/64.cc
+++ b/tests/trilinos_tpetra/64.cc
@@ -61,10 +61,11 @@ main(int argc, char **argv)
       {
         const unsigned int n_dofs = 420;
         // check
-        // LinearAlgebra::TpetraWrappers::SparseMatrix<double>
-        LinearAlgebra::TpetraWrappers::SparseMatrix<double> v1(n_dofs,
-                                                               n_dofs,
-                                                               5U);
+        // LinearAlgebra::TpetraWrappers::SparseMatrix<double,
+        // MemorySpace::Default>
+        LinearAlgebra::TpetraWrappers::SparseMatrix<double,
+                                                    MemorySpace::Default>
+          v1(n_dofs, n_dofs, 5U);
         test(v1);
 
         // check
@@ -77,9 +78,9 @@ main(int argc, char **argv)
         const unsigned int n_local_dofs = n_dofs / n_jobs;
         IndexSet           local_rows(n_dofs);
         local_rows.add_range(n_local_dofs * my_id, n_local_dofs * (my_id + 1));
-        LinearAlgebra::TpetraWrappers::SparseMatrix<double> v2(local_rows,
-                                                               MPI_COMM_WORLD,
-                                                               5);
+        LinearAlgebra::TpetraWrappers::SparseMatrix<double,
+                                                    MemorySpace::Default>
+          v2(local_rows, MPI_COMM_WORLD, 5);
         test(v2);
       }
     }
diff --git a/tests/trilinos_tpetra/66.cc b/tests/trilinos_tpetra/66.cc
index e28488e43f..ad15dec3c6 100644
--- a/tests/trilinos_tpetra/66.cc
+++ b/tests/trilinos_tpetra/66.cc
@@ -28,7 +28,8 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::SparseMatrix<double> &m)
+test(
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> &m)
 {
   Assert(m.m() != 0, ExcInternalError());
   Assert(m.n() != 0, ExcInternalError());
@@ -107,7 +108,9 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::SparseMatrix<double> v(14U, 14U, 3U);
+        LinearAlgebra::TpetraWrappers::SparseMatrix<double,
+                                                    MemorySpace::Default>
+          v(14U, 14U, 3U);
         test(v);
       }
     }
diff --git a/tests/trilinos_tpetra/67.cc b/tests/trilinos_tpetra/67.cc
index d59543e27d..9c630615b0 100644
--- a/tests/trilinos_tpetra/67.cc
+++ b/tests/trilinos_tpetra/67.cc
@@ -28,7 +28,8 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::SparseMatrix<double> &m)
+test(
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> &m)
 {
   Assert(m.m() != 0, ExcInternalError());
   Assert(m.n() != 0, ExcInternalError());
@@ -113,7 +114,9 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::SparseMatrix<double> v(14U, 14U, 3U);
+        LinearAlgebra::TpetraWrappers::SparseMatrix<double,
+                                                    MemorySpace::Default>
+          v(14U, 14U, 3U);
         test(v);
       }
     }
diff --git a/tests/trilinos_tpetra/68.cc b/tests/trilinos_tpetra/68.cc
index ffdccc1a8f..7c18509c1c 100644
--- a/tests/trilinos_tpetra/68.cc
+++ b/tests/trilinos_tpetra/68.cc
@@ -29,7 +29,8 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::SparseMatrix<double> &m)
+test(
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> &m)
 {
   Assert(m.m() != 0, ExcInternalError());
   Assert(m.n() != 0, ExcInternalError());
@@ -114,7 +115,9 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::SparseMatrix<double> v(14U, 14U, 3U);
+        LinearAlgebra::TpetraWrappers::SparseMatrix<double,
+                                                    MemorySpace::Default>
+          v(14U, 14U, 3U);
         test(v);
       }
     }
diff --git a/tests/trilinos_tpetra/69.cc b/tests/trilinos_tpetra/69.cc
index 45d93f09a4..876a981db4 100644
--- a/tests/trilinos_tpetra/69.cc
+++ b/tests/trilinos_tpetra/69.cc
@@ -29,7 +29,8 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::SparseMatrix<double> &m)
+test(
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> &m)
 {
   Assert(m.m() != 0, ExcInternalError());
   Assert(m.n() != 0, ExcInternalError());
@@ -122,7 +123,9 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::SparseMatrix<double> v(14U, 14U, 3U);
+        LinearAlgebra::TpetraWrappers::SparseMatrix<double,
+                                                    MemorySpace::Default>
+          v(14U, 14U, 3U);
         test(v);
       }
     }
diff --git a/tests/trilinos_tpetra/add_matrices_06.cc b/tests/trilinos_tpetra/add_matrices_06.cc
index 180281b000..a431be2c85 100644
--- a/tests/trilinos_tpetra/add_matrices_06.cc
+++ b/tests/trilinos_tpetra/add_matrices_06.cc
@@ -30,9 +30,11 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::SparseMatrix<double> &m)
+test(
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> &m)
 {
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double> m2(m.m(), m.n(), 3U),
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> m2(
+    m.m(), m.n(), 3U),
     m3(m.m(), m.n(), 3U);
 
   // first set a few entries one-by-one
@@ -59,7 +61,9 @@ test(LinearAlgebra::TpetraWrappers::SparseMatrix<double> &m)
   Vector<double> src(m.m()), dst(m.m());
   src = 1.;
 
-  LinearAlgebra::TpetraWrappers::PreconditionJacobi<double> prec;
+  LinearAlgebra::TpetraWrappers::PreconditionJacobi<double,
+                                                    MemorySpace::Default>
+    prec;
   prec.initialize(m);
 
   m.vmult(dst, src);
@@ -93,7 +97,9 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::SparseMatrix<double> m(5U, 5U, 3U);
+        LinearAlgebra::TpetraWrappers::SparseMatrix<double,
+                                                    MemorySpace::Default>
+          m(5U, 5U, 3U);
 
         test(m);
       }
diff --git a/tests/trilinos_tpetra/deal_solver_01.cc b/tests/trilinos_tpetra/deal_solver_01.cc
index 0ba5045435..9d209b39c1 100644
--- a/tests/trilinos_tpetra/deal_solver_01.cc
+++ b/tests/trilinos_tpetra/deal_solver_01.cc
@@ -58,23 +58,26 @@ main(int argc, char **argv)
     FDMatrix               testproblem(size, size);
     DynamicSparsityPattern csp(dim, dim);
     testproblem.five_point_structure(csp);
-    LinearAlgebra::TpetraWrappers::SparseMatrix<double> A;
+    LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> A;
     A.reinit(csp);
     testproblem.five_point(A);
 
-    LinearAlgebra::TpetraWrappers::Vector<double> f;
+    LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> f;
     f.reinit(complete_index_set(dim), MPI_COMM_WORLD);
-    LinearAlgebra::TpetraWrappers::Vector<double> u;
+    LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> u;
     u.reinit(complete_index_set(dim), MPI_COMM_WORLD);
     f.trilinos_rcp()->putScalar(1.0);
     A.compress(VectorOperation::insert);
     f.compress(VectorOperation::insert);
     u.compress(VectorOperation::insert);
 
-    GrowingVectorMemory<LinearAlgebra::TpetraWrappers::Vector<double>> mem;
-    SolverCG<LinearAlgebra::TpetraWrappers::Vector<double>> solver(control,
-                                                                   mem);
-    PreconditionIdentity                                    preconditioner;
+    GrowingVectorMemory<
+      LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>>
+      mem;
+    SolverCG<
+      LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>>
+                         solver(control, mem);
+    PreconditionIdentity preconditioner;
 
     deallog << "Solver type: " << typeid(solver).name() << std::endl;
 
diff --git a/tests/trilinos_tpetra/deal_solver_01.output b/tests/trilinos_tpetra/deal_solver_01.output
index 9ba3a2dae4..065c0a5ae2 100644
--- a/tests/trilinos_tpetra/deal_solver_01.output
+++ b/tests/trilinos_tpetra/deal_solver_01.output
@@ -1,4 +1,4 @@
 
 DEAL::Size 32 Unknowns 961
-DEAL::Solver type: N6dealii8SolverCGINS_13LinearAlgebra14TpetraWrappers6VectorIdNS_11MemorySpace4HostEEEEE
+DEAL::Solver type: N6dealii8SolverCGINS_13LinearAlgebra14TpetraWrappers6VectorIdNS_11MemorySpace7DefaultEEEEE
 DEAL::Solver stopped within 42 - 44 iterations
diff --git a/tests/trilinos_tpetra/deal_solver_02.cc b/tests/trilinos_tpetra/deal_solver_02.cc
index 01350ad3c3..b403fbac88 100644
--- a/tests/trilinos_tpetra/deal_solver_02.cc
+++ b/tests/trilinos_tpetra/deal_solver_02.cc
@@ -59,22 +59,25 @@ main(int argc, char **argv)
     FDMatrix               testproblem(size, size);
     DynamicSparsityPattern csp(dim, dim);
     testproblem.five_point_structure(csp);
-    LinearAlgebra::TpetraWrappers::SparseMatrix<double> A;
+    LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> A;
     A.reinit(csp);
     testproblem.five_point(A);
 
-    LinearAlgebra::TpetraWrappers::Vector<double> f;
+    LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> f;
     f.reinit(complete_index_set(dim), MPI_COMM_WORLD);
-    LinearAlgebra::TpetraWrappers::Vector<double> u;
+    LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> u;
     u.reinit(complete_index_set(dim), MPI_COMM_WORLD);
     f.trilinos_rcp()->putScalar(1.0);
     A.compress(VectorOperation::insert);
     f.compress(VectorOperation::insert);
     u.compress(VectorOperation::insert);
 
-    GrowingVectorMemory<LinearAlgebra::TpetraWrappers::Vector<double>> mem;
-    SolverBicgstab<LinearAlgebra::TpetraWrappers::Vector<double>>      solver(
-      control, mem);
+    GrowingVectorMemory<
+      LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>>
+      mem;
+    SolverBicgstab<
+      LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>>
+                         solver(control, mem);
     PreconditionIdentity preconditioner;
 
     deallog << "Solver type: " << typeid(solver).name() << std::endl;
diff --git a/tests/trilinos_tpetra/deal_solver_02.output b/tests/trilinos_tpetra/deal_solver_02.output
index ceafaf1d09..7ce9609c32 100644
--- a/tests/trilinos_tpetra/deal_solver_02.output
+++ b/tests/trilinos_tpetra/deal_solver_02.output
@@ -1,4 +1,4 @@
 
 DEAL::Size 32 Unknowns 961
-DEAL::Solver type: N6dealii14SolverBicgstabINS_13LinearAlgebra14TpetraWrappers6VectorIdNS_11MemorySpace4HostEEEEE
+DEAL::Solver type: N6dealii14SolverBicgstabINS_13LinearAlgebra14TpetraWrappers6VectorIdNS_11MemorySpace7DefaultEEEEE
 DEAL::Solver stopped within 49 - 51 iterations
diff --git a/tests/trilinos_tpetra/deal_solver_03.cc b/tests/trilinos_tpetra/deal_solver_03.cc
index 1b926dae5c..57fa49f96d 100644
--- a/tests/trilinos_tpetra/deal_solver_03.cc
+++ b/tests/trilinos_tpetra/deal_solver_03.cc
@@ -61,26 +61,29 @@ main(int argc, char **argv)
     FDMatrix               testproblem(size, size);
     DynamicSparsityPattern csp(dim, dim);
     testproblem.five_point_structure(csp);
-    LinearAlgebra::TpetraWrappers::SparseMatrix<double> A;
+    LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> A;
     A.reinit(csp);
     testproblem.five_point(A);
 
-    LinearAlgebra::TpetraWrappers::Vector<double> f;
+    LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> f;
     f.reinit(complete_index_set(dim), MPI_COMM_WORLD);
-    LinearAlgebra::TpetraWrappers::Vector<double> u;
+    LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> u;
     u.reinit(complete_index_set(dim), MPI_COMM_WORLD);
     f.trilinos_rcp()->putScalar(1.0);
     A.compress(VectorOperation::insert);
     f.compress(VectorOperation::insert);
     u.compress(VectorOperation::insert);
 
-    GrowingVectorMemory<LinearAlgebra::TpetraWrappers::Vector<double>> mem;
-    SolverGMRES<LinearAlgebra::TpetraWrappers::Vector<double>>::AdditionalData
-                                                               sol_data(28);
-    SolverGMRES<LinearAlgebra::TpetraWrappers::Vector<double>> solver(control,
-                                                                      mem,
-                                                                      sol_data);
-    PreconditionIdentity                                       preconditioner;
+    GrowingVectorMemory<
+      LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>>
+      mem;
+    SolverGMRES<
+      LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>>::
+      AdditionalData sol_data(28);
+    SolverGMRES<
+      LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>>
+                         solver(control, mem, sol_data);
+    PreconditionIdentity preconditioner;
 
     deallog << "Solver type: " << typeid(solver).name() << std::endl;
 
diff --git a/tests/trilinos_tpetra/deal_solver_03.output b/tests/trilinos_tpetra/deal_solver_03.output
index 9d72b0e9c5..9e919b6b2f 100644
--- a/tests/trilinos_tpetra/deal_solver_03.output
+++ b/tests/trilinos_tpetra/deal_solver_03.output
@@ -1,4 +1,4 @@
 
 DEAL::Size 32 Unknowns 961
-DEAL::Solver type: N6dealii11SolverGMRESINS_13LinearAlgebra14TpetraWrappers6VectorIdNS_11MemorySpace4HostEEEEE
+DEAL::Solver type: N6dealii11SolverGMRESINS_13LinearAlgebra14TpetraWrappers6VectorIdNS_11MemorySpace7DefaultEEEEE
 DEAL::Solver stopped within 74 - 76 iterations
diff --git a/tests/trilinos_tpetra/deal_solver_04.cc b/tests/trilinos_tpetra/deal_solver_04.cc
index e4c7ae7c73..8487d66d30 100644
--- a/tests/trilinos_tpetra/deal_solver_04.cc
+++ b/tests/trilinos_tpetra/deal_solver_04.cc
@@ -60,23 +60,26 @@ main(int argc, char **argv)
     FDMatrix               testproblem(size, size);
     DynamicSparsityPattern csp(dim, dim);
     testproblem.five_point_structure(csp);
-    LinearAlgebra::TpetraWrappers::SparseMatrix<double> A;
+    LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> A;
     A.reinit(csp);
     testproblem.five_point(A);
 
-    LinearAlgebra::TpetraWrappers::Vector<double> f;
+    LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> f;
     f.reinit(complete_index_set(dim), MPI_COMM_WORLD);
-    LinearAlgebra::TpetraWrappers::Vector<double> u;
+    LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> u;
     u.reinit(complete_index_set(dim), MPI_COMM_WORLD);
     f.trilinos_rcp()->putScalar(1.0);
     A.compress(VectorOperation::insert);
     f.compress(VectorOperation::insert);
     u.compress(VectorOperation::insert);
 
-    GrowingVectorMemory<LinearAlgebra::TpetraWrappers::Vector<double>> mem;
-    SolverMinRes<LinearAlgebra::TpetraWrappers::Vector<double>> solver(control,
-                                                                       mem);
-    PreconditionIdentity                                        preconditioner;
+    GrowingVectorMemory<
+      LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>>
+      mem;
+    SolverMinRes<
+      LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>>
+                         solver(control, mem);
+    PreconditionIdentity preconditioner;
 
     deallog << "Solver type: " << typeid(solver).name() << std::endl;
 
diff --git a/tests/trilinos_tpetra/deal_solver_04.output b/tests/trilinos_tpetra/deal_solver_04.output
index 662bfbeceb..59257f2894 100644
--- a/tests/trilinos_tpetra/deal_solver_04.output
+++ b/tests/trilinos_tpetra/deal_solver_04.output
@@ -1,4 +1,4 @@
 
 DEAL::Size 32 Unknowns 961
-DEAL::Solver type: N6dealii12SolverMinResINS_13LinearAlgebra14TpetraWrappers6VectorIdNS_11MemorySpace4HostEEEEE
+DEAL::Solver type: N6dealii12SolverMinResINS_13LinearAlgebra14TpetraWrappers6VectorIdNS_11MemorySpace7DefaultEEEEE
 DEAL::Solver stopped after 43 iterations
diff --git a/tests/trilinos_tpetra/deal_solver_05.cc b/tests/trilinos_tpetra/deal_solver_05.cc
index 2cd5b919ab..02253e3c74 100644
--- a/tests/trilinos_tpetra/deal_solver_05.cc
+++ b/tests/trilinos_tpetra/deal_solver_05.cc
@@ -59,23 +59,26 @@ main(int argc, char **argv)
     FDMatrix               testproblem(size, size);
     DynamicSparsityPattern csp(dim, dim);
     testproblem.five_point_structure(csp);
-    LinearAlgebra::TpetraWrappers::SparseMatrix<double> A;
+    LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> A;
     A.reinit(csp);
     testproblem.five_point(A);
 
-    LinearAlgebra::TpetraWrappers::Vector<double> f;
+    LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> f;
     f.reinit(complete_index_set(dim), MPI_COMM_WORLD);
-    LinearAlgebra::TpetraWrappers::Vector<double> u;
+    LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> u;
     u.reinit(complete_index_set(dim), MPI_COMM_WORLD);
     f.trilinos_rcp()->putScalar(1.0);
     A.compress(VectorOperation::insert);
     f.compress(VectorOperation::insert);
     u.compress(VectorOperation::insert);
 
-    GrowingVectorMemory<LinearAlgebra::TpetraWrappers::Vector<double>> mem;
-    SolverQMRS<LinearAlgebra::TpetraWrappers::Vector<double>> solver(control,
-                                                                     mem);
-    PreconditionIdentity                                      preconditioner;
+    GrowingVectorMemory<
+      LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>>
+      mem;
+    SolverQMRS<
+      LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>>
+                         solver(control, mem);
+    PreconditionIdentity preconditioner;
     deallog << "Solver type: " << typeid(solver).name() << std::endl;
 
     check_solver_within_range(solver.solve(A, u, f, preconditioner),
diff --git a/tests/trilinos_tpetra/deal_solver_05.output b/tests/trilinos_tpetra/deal_solver_05.output
index d15b2da96f..9ffde89206 100644
--- a/tests/trilinos_tpetra/deal_solver_05.output
+++ b/tests/trilinos_tpetra/deal_solver_05.output
@@ -1,4 +1,4 @@
 
 DEAL::Size 32 Unknowns 961
-DEAL::Solver type: N6dealii10SolverQMRSINS_13LinearAlgebra14TpetraWrappers6VectorIdNS_11MemorySpace4HostEEEEE
+DEAL::Solver type: N6dealii10SolverQMRSINS_13LinearAlgebra14TpetraWrappers6VectorIdNS_11MemorySpace7DefaultEEEEE
 DEAL::Solver stopped within 45 - 47 iterations
diff --git a/tests/trilinos_tpetra/deal_solver_06.cc b/tests/trilinos_tpetra/deal_solver_06.cc
index a5a104fcac..6581a6fee2 100644
--- a/tests/trilinos_tpetra/deal_solver_06.cc
+++ b/tests/trilinos_tpetra/deal_solver_06.cc
@@ -58,25 +58,29 @@ main(int argc, char **argv)
     FDMatrix               testproblem(size, size);
     DynamicSparsityPattern csp(dim, dim);
     testproblem.five_point_structure(csp);
-    LinearAlgebra::TpetraWrappers::SparseMatrix<double> A;
+    LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> A;
     A.reinit(csp);
     testproblem.five_point(A);
 
-    LinearAlgebra::TpetraWrappers::Vector<double> f;
+    LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> f;
     f.reinit(complete_index_set(dim), MPI_COMM_WORLD);
-    LinearAlgebra::TpetraWrappers::Vector<double> u;
+    LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> u;
     u.reinit(complete_index_set(dim), MPI_COMM_WORLD);
     f.trilinos_rcp()->putScalar(1.0);
     A.compress(VectorOperation::insert);
     f.compress(VectorOperation::insert);
     u.compress(VectorOperation::insert);
 
-    GrowingVectorMemory<LinearAlgebra::TpetraWrappers::Vector<double>> mem;
-    SolverRichardson<LinearAlgebra::TpetraWrappers::Vector<double>>::
+    GrowingVectorMemory<
+      LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>>
+      mem;
+    SolverRichardson<
+      LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>>::
       AdditionalData data(
         /*omega=*/0.1);
-    SolverRichardson<LinearAlgebra::TpetraWrappers::Vector<double>> solver(
-      control, mem, data);
+    SolverRichardson<
+      LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>>
+                         solver(control, mem, data);
     PreconditionIdentity preconditioner;
     deallog << "Solver type: " << typeid(solver).name() << std::endl;
 
diff --git a/tests/trilinos_tpetra/deal_solver_06.output b/tests/trilinos_tpetra/deal_solver_06.output
index 039bff27a1..c82a0d9a71 100644
--- a/tests/trilinos_tpetra/deal_solver_06.output
+++ b/tests/trilinos_tpetra/deal_solver_06.output
@@ -1,4 +1,4 @@
 
 DEAL::Size 32 Unknowns 961
-DEAL::Solver type: N6dealii16SolverRichardsonINS_13LinearAlgebra14TpetraWrappers6VectorIdNS_11MemorySpace4HostEEEEE
+DEAL::Solver type: N6dealii16SolverRichardsonINS_13LinearAlgebra14TpetraWrappers6VectorIdNS_11MemorySpace7DefaultEEEEE
 DEAL::Solver stopped within 5269 - 5271 iterations
diff --git a/tests/trilinos_tpetra/direct_solver_01.cc b/tests/trilinos_tpetra/direct_solver_01.cc
index 24f71ee43f..bc965f3607 100644
--- a/tests/trilinos_tpetra/direct_solver_01.cc
+++ b/tests/trilinos_tpetra/direct_solver_01.cc
@@ -71,11 +71,14 @@ private:
   AffineConstraints<double> constraints;
   SparsityPattern           sparsity_pattern;
 
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double> system_matrix;
-
-  LinearAlgebra::TpetraWrappers::Vector<double> solution;
-  LinearAlgebra::TpetraWrappers::Vector<double> system_rhs;
-  LinearAlgebra::TpetraWrappers::Vector<double> system_rhs_two;
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default>
+    system_matrix;
+
+  LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> solution;
+  LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>
+    system_rhs;
+  LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>
+    system_rhs_two;
 };
 
 
@@ -300,20 +303,23 @@ void
 Step4<dim>::solve()
 {
   // Compute 'reference' solution with CG solver and SSOR preconditioner
-  LinearAlgebra::TpetraWrappers::PreconditionSSOR<double> preconditioner;
+  LinearAlgebra::TpetraWrappers::PreconditionSSOR<double, MemorySpace::Default>
+    preconditioner;
   preconditioner.initialize(system_matrix);
-  LinearAlgebra::TpetraWrappers::Vector<double> temp_solution(system_rhs);
+  LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>
+    temp_solution(system_rhs);
   temp_solution = 0;
   SolverControl solver_control(1000, 1e-12);
-  SolverCG<LinearAlgebra::TpetraWrappers::Vector<double>> solver(
-    solver_control);
+  SolverCG<LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>>
+    solver(solver_control);
 
   solver.solve(system_matrix, temp_solution, system_rhs, preconditioner);
 
   constraints.distribute(temp_solution);
   solution = temp_solution;
 
-  LinearAlgebra::TpetraWrappers::Vector<double> output(temp_solution);
+  LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> output(
+    temp_solution);
 
   // do CG solve for new rhs
   temp_solution = 0;
@@ -323,7 +329,8 @@ Step4<dim>::solve()
   constraints.distribute(temp_solution);
   solution = temp_solution;
 
-  LinearAlgebra::TpetraWrappers::Vector<double> output_two(temp_solution);
+  LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>
+    output_two(temp_solution);
 
   // factorize matrix for direct solver
   temp_solution = 0;
@@ -331,10 +338,10 @@ Step4<dim>::solve()
 
   {
     deallog.push("DirectKLU2");
-    LinearAlgebra::TpetraWrappers::SolverDirect<double>::AdditionalData data(
-      "KLU2");
-    LinearAlgebra::TpetraWrappers::SolverDirect<double> direct_solver(
-      solver_control, data);
+    LinearAlgebra::TpetraWrappers::SolverDirect<double, MemorySpace::Default>::
+      AdditionalData data("KLU2");
+    LinearAlgebra::TpetraWrappers::SolverDirect<double, MemorySpace::Default>
+      direct_solver(solver_control, data);
 
     Teuchos::ParameterList amesos2_params("Amesos2");
     Teuchos::ParameterList klu2_params = amesos2_params.sublist("KLU2");
@@ -384,12 +391,13 @@ Step4<dim>::solve()
     // MEASLES is a non-existent solver so it throws the exception
     // independent of the Amesos2 configuration.
     deallog.push("DirectMEASLES");
-    LinearAlgebra::TpetraWrappers::SolverDirect<double>::AdditionalData data(
-      "MEASLES");
+    LinearAlgebra::TpetraWrappers::SolverDirect<double, MemorySpace::Default>::
+      AdditionalData data("MEASLES");
     try
       {
-        LinearAlgebra::TpetraWrappers::SolverDirect<double> direct_solver(
-          solver_control, data);
+        LinearAlgebra::TpetraWrappers::SolverDirect<double,
+                                                    MemorySpace::Default>
+          direct_solver(solver_control, data);
       }
     catch (dealii::ExceptionBase &exc)
       {
diff --git a/tests/trilinos_tpetra/direct_solver_02.cc b/tests/trilinos_tpetra/direct_solver_02.cc
index e384b2fe0b..cf940b4a59 100644
--- a/tests/trilinos_tpetra/direct_solver_02.cc
+++ b/tests/trilinos_tpetra/direct_solver_02.cc
@@ -70,11 +70,14 @@ private:
   AffineConstraints<double> constraints;
   SparsityPattern           sparsity_pattern;
 
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double> system_matrix;
-
-  LinearAlgebra::TpetraWrappers::Vector<double> solution;
-  LinearAlgebra::TpetraWrappers::Vector<double> system_rhs;
-  LinearAlgebra::TpetraWrappers::Vector<double> system_rhs_two;
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default>
+    system_matrix;
+
+  LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> solution;
+  LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>
+    system_rhs;
+  LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>
+    system_rhs_two;
 };
 
 
@@ -299,20 +302,23 @@ void
 Step4<dim>::solve()
 {
   // Compute 'reference' solution with CG solver and SSOR preconditioner
-  LinearAlgebra::TpetraWrappers::PreconditionSSOR<double> preconditioner;
+  LinearAlgebra::TpetraWrappers::PreconditionSSOR<double, MemorySpace::Default>
+    preconditioner;
   preconditioner.initialize(system_matrix);
-  LinearAlgebra::TpetraWrappers::Vector<double> temp_solution(system_rhs);
+  LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>
+    temp_solution(system_rhs);
   temp_solution = 0;
   SolverControl solver_control(1000, 1e-12);
-  SolverCG<LinearAlgebra::TpetraWrappers::Vector<double>> solver(
-    solver_control);
+  SolverCG<LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>>
+    solver(solver_control);
 
   solver.solve(system_matrix, temp_solution, system_rhs, preconditioner);
 
   constraints.distribute(temp_solution);
   solution = temp_solution;
 
-  LinearAlgebra::TpetraWrappers::Vector<double> output(temp_solution);
+  LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> output(
+    temp_solution);
 
   // do CG solve for new rhs
   temp_solution = 0;
@@ -322,16 +328,18 @@ Step4<dim>::solve()
   constraints.distribute(temp_solution);
   solution = temp_solution;
 
-  LinearAlgebra::TpetraWrappers::Vector<double> output_two(temp_solution);
+  LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>
+    output_two(temp_solution);
 
   // factorize matrix for direct solver
   temp_solution = 0;
   solution      = 0;
 
   deallog.push("DirectKLU2");
-  LinearAlgebra::TpetraWrappers::SolverDirectKLU2<double>::AdditionalData data;
-  LinearAlgebra::TpetraWrappers::SolverDirectKLU2<double> direct_solver(
-    solver_control, data);
+  LinearAlgebra::TpetraWrappers::
+    SolverDirectKLU2<double, MemorySpace::Default>::AdditionalData data;
+  LinearAlgebra::TpetraWrappers::SolverDirectKLU2<double, MemorySpace::Default>
+    direct_solver(solver_control, data);
   direct_solver.initialize(system_matrix);
 
   // do solve 1
diff --git a/tests/trilinos_tpetra/index_set_01.cc b/tests/trilinos_tpetra/index_set_01.cc
index 3c2d965178..94452a969e 100644
--- a/tests/trilinos_tpetra/index_set_01.cc
+++ b/tests/trilinos_tpetra/index_set_01.cc
@@ -49,7 +49,9 @@ test()
       my_set.add_range(400, 529);
     }
   my_set.compress();
-  my_set.make_tpetra_map(MPI_COMM_WORLD, false);
+  my_set.template make_tpetra_map<
+    LinearAlgebra::TpetraWrappers::TpetraTypes::NodeType<MemorySpace::Default>>(
+    MPI_COMM_WORLD, false);
 }
 
 
diff --git a/tests/trilinos_tpetra/n_nonzeros_01.cc b/tests/trilinos_tpetra/n_nonzeros_01.cc
index cca9edf331..fb17aae461 100644
--- a/tests/trilinos_tpetra/n_nonzeros_01.cc
+++ b/tests/trilinos_tpetra/n_nonzeros_01.cc
@@ -26,11 +26,13 @@ void
 test()
 {
   // create an empty sparsity pattern
-  LinearAlgebra::TpetraWrappers::SparsityPattern sparsity(4, 5, 1);
+  LinearAlgebra::TpetraWrappers::SparsityPattern<MemorySpace::Default> sparsity(
+    4, 5, 1);
   sparsity.compress();
 
   // attach a sparse matrix to it
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double> A(sparsity);
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> A(
+    sparsity);
 
   // see how many nonzero elements it reports
   deallog << A.n_nonzero_elements() << std::endl;
diff --git a/tests/trilinos_tpetra/n_nonzeros_02.cc b/tests/trilinos_tpetra/n_nonzeros_02.cc
index 9089ed4729..0d8d8dd579 100644
--- a/tests/trilinos_tpetra/n_nonzeros_02.cc
+++ b/tests/trilinos_tpetra/n_nonzeros_02.cc
@@ -44,7 +44,7 @@ test()
   coupling.fill(DoFTools::none);
 
   // create an empty sparsity pattern
-  LinearAlgebra::TpetraWrappers::SparsityPattern sparsity;
+  LinearAlgebra::TpetraWrappers::SparsityPattern<MemorySpace::Default> sparsity;
   sparsity.reinit(dof_handler.n_dofs(), dof_handler.n_dofs());
   DoFTools::make_sparsity_pattern(dof_handler,
                                   coupling,
@@ -56,7 +56,7 @@ test()
   sparsity.compress();
 
   // attach a sparse matrix to it
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double> A;
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> A;
   A.reinit(sparsity);
 
   // see how many nonzero elements it reports
diff --git a/tests/trilinos_tpetra/parallel_sparse_matrix_01.cc b/tests/trilinos_tpetra/parallel_sparse_matrix_01.cc
index 5136775452..eacea604a3 100644
--- a/tests/trilinos_tpetra/parallel_sparse_matrix_01.cc
+++ b/tests/trilinos_tpetra/parallel_sparse_matrix_01.cc
@@ -87,7 +87,7 @@ test()
   IndexSet           locally_owned_dofs(N);
   locally_owned_dofs.add_range(start_row[my_id], local_rows_per_process[my_id]);
 
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double> m;
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> m;
   m.reinit(locally_owned_dofs, locally_owned_dofs, csp, MPI_COMM_WORLD);
   // now write into the exact same matrix
   // entries as have been created by the
diff --git a/tests/trilinos_tpetra/precondition.cc b/tests/trilinos_tpetra/precondition.cc
index 7ecb8536bf..0d7cf4a032 100644
--- a/tests/trilinos_tpetra/precondition.cc
+++ b/tests/trilinos_tpetra/precondition.cc
@@ -70,10 +70,12 @@ private:
 
   AffineConstraints<double> constraints;
 
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double> system_matrix;
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default>
+    system_matrix;
 
-  LinearAlgebra::TpetraWrappers::Vector<double> solution;
-  LinearAlgebra::TpetraWrappers::Vector<double> system_rhs;
+  LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> solution;
+  LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>
+    system_rhs;
 };
 
 
@@ -233,11 +235,14 @@ Step4<dim>::solve(int cycle)
   {
     deallog.push("Identity");
     static constexpr std::array<unsigned int, 2> lower{{49, 100}};
-    LinearAlgebra::TpetraWrappers::PreconditionIdentity<double> preconditioner;
+    LinearAlgebra::TpetraWrappers::PreconditionIdentity<double,
+                                                        MemorySpace::Default>
+      preconditioner;
     solution = 0;
     SolverControl solver_control(1000, 1e-10);
-    SolverCG<LinearAlgebra::TpetraWrappers::Vector<double>> solver(
-      solver_control);
+    SolverCG<
+      LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>>
+                           solver(solver_control);
     Teuchos::ParameterList precon_params;
     preconditioner.initialize(system_matrix);
     check_solver_within_range(
@@ -253,13 +258,15 @@ Step4<dim>::solve(int cycle)
     // For this we use the point relaxation (SSOR)
     // with default parameters, resulting in symmetric Gauss-Seidel (SGS)
     deallog.push("CustomSGS");
-    static constexpr std::array<unsigned int, 2>              lower{{49, 100}};
-    LinearAlgebra::TpetraWrappers::PreconditionIfpack<double> preconditioner(
-      "RELAXATION");
+    static constexpr std::array<unsigned int, 2> lower{{49, 100}};
+    LinearAlgebra::TpetraWrappers::PreconditionIfpack<double,
+                                                      MemorySpace::Default>
+      preconditioner("RELAXATION");
     solution = 0;
     SolverControl solver_control(1000, 1e-10);
-    SolverCG<LinearAlgebra::TpetraWrappers::Vector<double>> solver(
-      solver_control);
+    SolverCG<
+      LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>>
+                           solver(solver_control);
     Teuchos::ParameterList precon_params;
     precon_params.set("relaxation: type", "Symmetric Gauss-Seidel");
     preconditioner.set_parameter_list(precon_params);
@@ -274,13 +281,16 @@ Step4<dim>::solve(int cycle)
 
   {
     deallog.push("Jacobi");
-    static constexpr std::array<unsigned int, 2>              lower{{49, 100}};
-    LinearAlgebra::TpetraWrappers::PreconditionJacobi<double> preconditioner;
+    static constexpr std::array<unsigned int, 2> lower{{49, 100}};
+    LinearAlgebra::TpetraWrappers::PreconditionJacobi<double,
+                                                      MemorySpace::Default>
+      preconditioner;
     solution = 0;
 
     SolverControl solver_control(1000, 1e-10);
-    SolverCG<LinearAlgebra::TpetraWrappers::Vector<double>> solver(
-      solver_control);
+    SolverCG<
+      LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>>
+      solver(solver_control);
     preconditioner.initialize(system_matrix);
     check_solver_within_range(
       solver.solve(system_matrix, solution, system_rhs, preconditioner),
@@ -293,12 +303,15 @@ Step4<dim>::solve(int cycle)
   {
     deallog.push("l1Jacobi");
     static constexpr std::array<unsigned int, 2> lower{{49, 100}};
-    LinearAlgebra::TpetraWrappers::PreconditionL1Jacobi<double> preconditioner;
+    LinearAlgebra::TpetraWrappers::PreconditionL1Jacobi<double,
+                                                        MemorySpace::Default>
+      preconditioner;
     solution = 0;
 
     SolverControl solver_control(1000, 1e-10);
-    SolverCG<LinearAlgebra::TpetraWrappers::Vector<double>> solver(
-      solver_control);
+    SolverCG<
+      LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>>
+      solver(solver_control);
     preconditioner.initialize(system_matrix);
     check_solver_within_range(
       solver.solve(system_matrix, solution, system_rhs, preconditioner),
@@ -311,13 +324,15 @@ Step4<dim>::solve(int cycle)
   {
     deallog.push("l1GaussSeidel");
     static constexpr std::array<unsigned int, 2> lower{{31, 62}};
-    LinearAlgebra::TpetraWrappers::PreconditionL1GaussSeidel<double>
-      preconditioner;
+    LinearAlgebra::TpetraWrappers::
+      PreconditionL1GaussSeidel<double, MemorySpace::Default>
+        preconditioner;
     solution = 0;
 
     SolverControl solver_control(1000, 1e-10);
-    SolverBicgstab<LinearAlgebra::TpetraWrappers::Vector<double>> solver(
-      solver_control);
+    SolverBicgstab<
+      LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>>
+      solver(solver_control);
     preconditioner.initialize(system_matrix);
     check_solver_within_range(
       solver.solve(system_matrix, solution, system_rhs, preconditioner),
@@ -329,13 +344,15 @@ Step4<dim>::solve(int cycle)
 
   {
     deallog.push("SOR");
-    static constexpr std::array<unsigned int, 2>           lower{{31, 62}};
-    LinearAlgebra::TpetraWrappers::PreconditionSOR<double> preconditioner;
+    static constexpr std::array<unsigned int, 2> lower{{31, 62}};
+    LinearAlgebra::TpetraWrappers::PreconditionSOR<double, MemorySpace::Default>
+      preconditioner;
     solution = 0;
 
     SolverControl solver_control(1000, 1e-5);
-    SolverBicgstab<LinearAlgebra::TpetraWrappers::Vector<double>> solver(
-      solver_control);
+    SolverBicgstab<
+      LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>>
+      solver(solver_control);
     preconditioner.initialize(system_matrix);
     check_solver_within_range(
       solver.solve(system_matrix, solution, system_rhs, preconditioner),
@@ -347,13 +364,16 @@ Step4<dim>::solve(int cycle)
 
   {
     deallog.push("SSOR");
-    static constexpr std::array<unsigned int, 2>            lower{{40, 77}};
-    LinearAlgebra::TpetraWrappers::PreconditionSSOR<double> preconditioner;
+    static constexpr std::array<unsigned int, 2> lower{{40, 77}};
+    LinearAlgebra::TpetraWrappers::PreconditionSSOR<double,
+                                                    MemorySpace::Default>
+      preconditioner;
     solution = 0;
 
     SolverControl solver_control(1000, 1e-10);
-    SolverCG<LinearAlgebra::TpetraWrappers::Vector<double>> solver(
-      solver_control);
+    SolverCG<
+      LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>>
+      solver(solver_control);
     preconditioner.initialize(system_matrix);
     check_solver_within_range(
       solver.solve(system_matrix, solution, system_rhs, preconditioner),
@@ -366,15 +386,18 @@ Step4<dim>::solve(int cycle)
   {
     deallog.push("BlockJacobi");
     static constexpr std::array<unsigned int, 2> lower{{73, 145}};
-    LinearAlgebra::TpetraWrappers::PreconditionBlockJacobi<double>
+    LinearAlgebra::TpetraWrappers::PreconditionBlockJacobi<double,
+                                                           MemorySpace::Default>
       preconditioner;
     LinearAlgebra::TpetraWrappers::PreconditionBlockJacobi<
-      double>::AdditionalData data;
+      double,
+      MemorySpace::Default>::AdditionalData data;
     data.n_local_parts = 16;
     solution           = 0;
     SolverControl solver_control(1000, 1e-10);
-    SolverCG<LinearAlgebra::TpetraWrappers::Vector<double>> solver(
-      solver_control);
+    SolverCG<
+      LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>>
+      solver(solver_control);
     preconditioner.initialize(system_matrix, data);
     check_solver_within_range(
       solver.solve(system_matrix, solution, system_rhs, preconditioner),
@@ -386,16 +409,19 @@ Step4<dim>::solve(int cycle)
 
   {
     deallog.push("BlockSOR");
-    static constexpr std::array<unsigned int, 2>                lower{{18, 37}};
-    LinearAlgebra::TpetraWrappers::PreconditionBlockSOR<double> preconditioner;
-    LinearAlgebra::TpetraWrappers::PreconditionBlockSOR<double>::AdditionalData
-      data;
+    static constexpr std::array<unsigned int, 2> lower{{18, 37}};
+    LinearAlgebra::TpetraWrappers::PreconditionBlockSOR<double,
+                                                        MemorySpace::Default>
+      preconditioner;
+    LinearAlgebra::TpetraWrappers::
+      PreconditionBlockSOR<double, MemorySpace::Default>::AdditionalData data;
     data.n_local_parts = 16;
     data.omega         = 0.8;
     solution           = 0;
     SolverControl solver_control(1000, 1e-5);
-    SolverBicgstab<LinearAlgebra::TpetraWrappers::Vector<double>> solver(
-      solver_control);
+    SolverBicgstab<
+      LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>>
+      solver(solver_control);
     preconditioner.initialize(system_matrix, data);
     check_solver_within_range(
       solver.solve(system_matrix, solution, system_rhs, preconditioner),
@@ -408,15 +434,18 @@ Step4<dim>::solve(int cycle)
   {
     deallog.push("BlockSSOR");
     static constexpr std::array<unsigned int, 2> lower{{30, 59}};
-    LinearAlgebra::TpetraWrappers::PreconditionBlockSSOR<double> preconditioner;
-    LinearAlgebra::TpetraWrappers::PreconditionBlockSSOR<double>::AdditionalData
-      data;
+    LinearAlgebra::TpetraWrappers::PreconditionBlockSSOR<double,
+                                                         MemorySpace::Default>
+      preconditioner;
+    LinearAlgebra::TpetraWrappers::
+      PreconditionBlockSSOR<double, MemorySpace::Default>::AdditionalData data;
     data.n_local_parts = 16;
     data.omega         = 1.2;
     solution           = 0;
     SolverControl solver_control(1000, 1e-10);
-    SolverCG<LinearAlgebra::TpetraWrappers::Vector<double>> solver(
-      solver_control);
+    SolverCG<
+      LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>>
+      solver(solver_control);
     preconditioner.initialize(system_matrix, data);
     check_solver_within_range(
       solver.solve(system_matrix, solution, system_rhs, preconditioner),
@@ -429,11 +458,13 @@ Step4<dim>::solve(int cycle)
   {
     static constexpr std::array<unsigned int, 2> lower{{30, 56}};
     deallog.push("ILU");
-    LinearAlgebra::TpetraWrappers::PreconditionILU<double> preconditioner;
+    LinearAlgebra::TpetraWrappers::PreconditionILU<double, MemorySpace::Default>
+      preconditioner;
     solution = 0;
     SolverControl solver_control(1000, 1e-10);
-    SolverCG<LinearAlgebra::TpetraWrappers::Vector<double>> solver(
-      solver_control);
+    SolverCG<
+      LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>>
+      solver(solver_control);
     preconditioner.initialize(system_matrix);
     check_solver_within_range(
       solver.solve(system_matrix, solution, system_rhs, preconditioner),
@@ -445,16 +476,19 @@ Step4<dim>::solve(int cycle)
 
   {
     deallog.push("ILUT");
-    static constexpr std::array<unsigned int, 2>            lower{{11, 19}};
-    LinearAlgebra::TpetraWrappers::PreconditionILUT<double> preconditioner;
-    LinearAlgebra::TpetraWrappers::PreconditionILUT<double>::AdditionalData
-      data;
+    static constexpr std::array<unsigned int, 2> lower{{11, 19}};
+    LinearAlgebra::TpetraWrappers::PreconditionILUT<double,
+                                                    MemorySpace::Default>
+      preconditioner;
+    LinearAlgebra::TpetraWrappers::
+      PreconditionILUT<double, MemorySpace::Default>::AdditionalData data;
     data.ilut_drop = 1e-6;
     data.ilut_fill = 3;
     solution       = 0;
     SolverControl solver_control(1000, 1e-5);
-    SolverBicgstab<LinearAlgebra::TpetraWrappers::Vector<double>> solver(
-      solver_control);
+    SolverBicgstab<
+      LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>>
+      solver(solver_control);
     preconditioner.initialize(system_matrix, data);
     check_solver_within_range(
       solver.solve(system_matrix, solution, system_rhs, preconditioner),
@@ -467,15 +501,18 @@ Step4<dim>::solve(int cycle)
   {
     deallog.push("Chebyshev");
     static constexpr std::array<unsigned int, 2> lower{{23, 46}};
-    LinearAlgebra::TpetraWrappers::PreconditionChebyshev<double> preconditioner;
-    LinearAlgebra::TpetraWrappers::PreconditionChebyshev<double>::AdditionalData
-      data;
+    LinearAlgebra::TpetraWrappers::PreconditionChebyshev<double,
+                                                         MemorySpace::Default>
+      preconditioner;
+    LinearAlgebra::TpetraWrappers::
+      PreconditionChebyshev<double, MemorySpace::Default>::AdditionalData data;
     data.max_eigenvalue = 2.5;
     data.degree         = 3;
     solution            = 0;
     SolverControl solver_control(1000, 1e-10);
-    SolverCG<LinearAlgebra::TpetraWrappers::Vector<double>> solver(
-      solver_control);
+    SolverCG<
+      LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default>>
+      solver(solver_control);
     preconditioner.initialize(system_matrix, data);
     check_solver_within_range(
       solver.solve(system_matrix, solution, system_rhs, preconditioner),
@@ -510,7 +547,8 @@ Step4<dim>::solve(int cycle)
     deallog.push("CustomIGLOO");
     try
       {
-        LinearAlgebra::TpetraWrappers::PreconditionIfpack<double>
+        LinearAlgebra::TpetraWrappers::PreconditionIfpack<double,
+                                                          MemorySpace::Default>
           preconditioner("IGLOO");
       }
     catch (dealii::ExceptionBase &exc)
diff --git a/tests/trilinos_tpetra/slowness_02.cc b/tests/trilinos_tpetra/slowness_02.cc
index 0055228577..95efd95335 100644
--- a/tests/trilinos_tpetra/slowness_02.cc
+++ b/tests/trilinos_tpetra/slowness_02.cc
@@ -38,7 +38,8 @@ test()
   const unsigned int N = 200;
 
   // build the sparse matrix
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double> matrix(N * N, N * N, 5U);
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default>
+    matrix(N * N, N * N, 5U);
   for (unsigned int i = 0; i < N; ++i)
     for (unsigned int j = 0; j < N; ++j)
       {
@@ -71,9 +72,9 @@ test()
   // then do a single matrix-vector
   // multiplication with subsequent formation
   // of the matrix norm
-  LinearAlgebra::TpetraWrappers::Vector<double> v1;
+  LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v1;
   v1.reinit(complete_index_set(N * N), MPI_COMM_WORLD);
-  LinearAlgebra::TpetraWrappers::Vector<double> v2;
+  LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v2;
   v2.reinit(complete_index_set(N * N), MPI_COMM_WORLD);
   for (unsigned int i = 0; i < N * N; ++i)
     v1(i) = i;
diff --git a/tests/trilinos_tpetra/slowness_03.cc b/tests/trilinos_tpetra/slowness_03.cc
index 75a7b9a16d..ebe48b6efb 100644
--- a/tests/trilinos_tpetra/slowness_03.cc
+++ b/tests/trilinos_tpetra/slowness_03.cc
@@ -39,9 +39,8 @@ test()
 
   // build the sparse matrix
   IndexSet indices = complete_index_set(N * N);
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double> matrix(indices,
-                                                             MPI_COMM_WORLD,
-                                                             5);
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default>
+    matrix(indices, MPI_COMM_WORLD, 5);
   for (unsigned int i = 0; i < N; ++i)
     for (unsigned int j = 0; j < N; ++j)
       {
@@ -73,8 +72,10 @@ test()
   // then do a single matrix-vector
   // multiplication with subsequent formation
   // of the matrix norm
-  LinearAlgebra::TpetraWrappers::Vector<double> v1(indices, MPI_COMM_WORLD);
-  LinearAlgebra::TpetraWrappers::Vector<double> v2(indices, MPI_COMM_WORLD);
+  LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v1(
+    indices, MPI_COMM_WORLD);
+  LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v2(
+    indices, MPI_COMM_WORLD);
   for (unsigned int i = 0; i < N * N; ++i)
     v1(i) = i;
   matrix.vmult(v2, v1);
diff --git a/tests/trilinos_tpetra/slowness_04.cc b/tests/trilinos_tpetra/slowness_04.cc
index 79ca15c824..93c2fe7bba 100644
--- a/tests/trilinos_tpetra/slowness_04.cc
+++ b/tests/trilinos_tpetra/slowness_04.cc
@@ -67,9 +67,8 @@ test()
 
   // build the sparse matrix
   IndexSet indices = complete_index_set(N * N);
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double> matrix(indices,
-                                                             MPI_COMM_WORLD,
-                                                             5);
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default>
+    matrix(indices, MPI_COMM_WORLD, 5);
   for (unsigned int i_ = 0; i_ < N; ++i_)
     for (unsigned int j_ = 0; j_ < N; ++j_)
       {
@@ -104,8 +103,10 @@ test()
   // then do a single matrix-vector
   // multiplication with subsequent formation
   // of the matrix norm
-  LinearAlgebra::TpetraWrappers::Vector<double> v1(indices, MPI_COMM_WORLD);
-  LinearAlgebra::TpetraWrappers::Vector<double> v2(indices, MPI_COMM_WORLD);
+  LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v1(
+    indices, MPI_COMM_WORLD);
+  LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v2(
+    indices, MPI_COMM_WORLD);
   for (unsigned int i = 0; i < N * N; ++i)
     v1(i) = i;
   matrix.vmult(v2, v1);
diff --git a/tests/trilinos_tpetra/sparse_matrix_02.cc b/tests/trilinos_tpetra/sparse_matrix_02.cc
index d01a631773..8df7a82c48 100644
--- a/tests/trilinos_tpetra/sparse_matrix_02.cc
+++ b/tests/trilinos_tpetra/sparse_matrix_02.cc
@@ -54,7 +54,8 @@ main(int argc, char **argv)
 
   // now copy everything into a Trilinos matrix
   const auto local_rows = complete_index_set(5);
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double> tmatrix;
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default>
+    tmatrix;
   tmatrix.reinit(local_rows, local_rows, matrix, MPI_COMM_SELF);
 
   deallog << "Copy with all values:" << std::endl;
diff --git a/tests/trilinos_tpetra/sparse_matrix_03.cc b/tests/trilinos_tpetra/sparse_matrix_03.cc
index d7bf2e2292..c6585b17ea 100644
--- a/tests/trilinos_tpetra/sparse_matrix_03.cc
+++ b/tests/trilinos_tpetra/sparse_matrix_03.cc
@@ -54,7 +54,8 @@ main(int argc, char **argv)
 
   // now copy everything into a Trilinos matrix
   const auto local_rows = complete_index_set(5);
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double> tmatrix;
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default>
+    tmatrix;
   tmatrix.reinit(local_rows, local_rows, matrix, MPI_COMM_SELF, 4.);
 
   deallog << "Copy with a drop tolerance of 4:" << std::endl;
diff --git a/tests/trilinos_tpetra/sparse_matrix_04.cc b/tests/trilinos_tpetra/sparse_matrix_04.cc
index 304762a628..3b600247aa 100644
--- a/tests/trilinos_tpetra/sparse_matrix_04.cc
+++ b/tests/trilinos_tpetra/sparse_matrix_04.cc
@@ -54,7 +54,8 @@ main(int argc, char **argv)
 
   // now copy everything into a Trilinos matrix
   const auto local_rows = complete_index_set(5);
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double> tmatrix;
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default>
+    tmatrix;
   tmatrix.reinit(local_rows, local_rows, matrix, MPI_COMM_SELF, 0, false);
 
   deallog << "Copy structure only:" << std::endl;
diff --git a/tests/trilinos_tpetra/sparse_matrix_copy_from_02.cc b/tests/trilinos_tpetra/sparse_matrix_copy_from_02.cc
index d64e1a9658..86755ccba3 100644
--- a/tests/trilinos_tpetra/sparse_matrix_copy_from_02.cc
+++ b/tests/trilinos_tpetra/sparse_matrix_copy_from_02.cc
@@ -54,8 +54,8 @@ test()
   // with reordered dofs by its components, such that the rows in the
   // final matrix are locally not in a contiguous set.
 
-  dealii::LinearAlgebra::TpetraWrappers::SparsityPattern sp_M(
-    parallel_partitioning, MPI_COMM_WORLD, 2);
+  dealii::LinearAlgebra::TpetraWrappers::SparsityPattern<MemorySpace::Default>
+    sp_M(parallel_partitioning, MPI_COMM_WORLD, 2);
 
   sp_M.add(MyPID, MyPID);
   sp_M.add(MyPID, NumProc + MyPID);
@@ -65,7 +65,9 @@ test()
   sp_M.compress();
 
   // create matrix with dummy entries on the diagonal
-  dealii::LinearAlgebra::TpetraWrappers::SparseMatrix<double> M0;
+  dealii::LinearAlgebra::TpetraWrappers::SparseMatrix<double,
+                                                      MemorySpace::Default>
+    M0;
   M0.reinit(sp_M);
   M0 = 0;
 
@@ -78,7 +80,9 @@ test()
   // test ::add(TrilinosScalar, SparseMatrix)
   //
 
-  dealii::LinearAlgebra::TpetraWrappers::SparseMatrix<double> M1;
+  dealii::LinearAlgebra::TpetraWrappers::SparseMatrix<double,
+                                                      MemorySpace::Default>
+    M1;
   M1.reinit(sp_M); // avoid deep copy
   M1 = 0;
   M1.copy_from(M0);
diff --git a/tests/trilinos_tpetra/sparse_matrix_iterator_01.cc b/tests/trilinos_tpetra/sparse_matrix_iterator_01.cc
index 0657d1c22e..5bfa1ab68d 100644
--- a/tests/trilinos_tpetra/sparse_matrix_iterator_01.cc
+++ b/tests/trilinos_tpetra/sparse_matrix_iterator_01.cc
@@ -28,13 +28,14 @@
 void
 test()
 {
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double> m(5U, 5U, 5U);
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> m(
+    5U, 5U, 5U);
   m.set(0, 0, 1);
   m.set(1, 1, 2);
   m.set(1, 2, 3);
   m.compress(VectorOperation::insert);
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double>::const_iterator i =
-    m.begin();
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default>::
+    const_iterator i = m.begin();
   deallog << i->row() << ' ' << i->column() << ' ' << i->value() << std::endl;
   ++i;
   deallog << i->row() << ' ' << i->column() << ' ' << i->value() << std::endl;
diff --git a/tests/trilinos_tpetra/sparse_matrix_set_01.cc b/tests/trilinos_tpetra/sparse_matrix_set_01.cc
index 6f8ac7249a..00aab2cb0d 100644
--- a/tests/trilinos_tpetra/sparse_matrix_set_01.cc
+++ b/tests/trilinos_tpetra/sparse_matrix_set_01.cc
@@ -29,7 +29,8 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::SparseMatrix<double> &m)
+test(
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> &m)
 {
   // first set a few entries one-by-one
   for (unsigned int i = 0; i < m.m(); ++i)
@@ -39,9 +40,8 @@ test(LinearAlgebra::TpetraWrappers::SparseMatrix<double> &m)
 
   m.compress(VectorOperation::insert);
 
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double> m2(m.m(),
-                                                         m.n(),
-                                                         m.n() / 3 + 1);
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> m2(
+    m.m(), m.n(), m.n() / 3 + 1);
 
   // now set the same elements row-wise
   {
@@ -100,7 +100,9 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::SparseMatrix<double> m(5U, 5U, 3U);
+        LinearAlgebra::TpetraWrappers::SparseMatrix<double,
+                                                    MemorySpace::Default>
+          m(5U, 5U, 3U);
         test(m);
       }
     }
diff --git a/tests/trilinos_tpetra/sparse_matrix_set_02.cc b/tests/trilinos_tpetra/sparse_matrix_set_02.cc
index 0a294ebc7f..a69d6a7c54 100644
--- a/tests/trilinos_tpetra/sparse_matrix_set_02.cc
+++ b/tests/trilinos_tpetra/sparse_matrix_set_02.cc
@@ -31,9 +31,11 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::SparseMatrix<double> &m)
+test(
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> &m)
 {
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double> m2(m.m(), m.n(), 3U);
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> m2(
+    m.m(), m.n(), 3U);
 
   // first set a few entries one-by-one and
   // initialize the sparsity pattern for m2
@@ -98,7 +100,9 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::SparseMatrix<double> m(16U, 16U, 3U);
+        LinearAlgebra::TpetraWrappers::SparseMatrix<double,
+                                                    MemorySpace::Default>
+          m(16U, 16U, 3U);
         test(m);
       }
     }
diff --git a/tests/trilinos_tpetra/sparse_matrix_vector_01.cc b/tests/trilinos_tpetra/sparse_matrix_vector_01.cc
index 21bb646b88..7054b93320 100644
--- a/tests/trilinos_tpetra/sparse_matrix_vector_01.cc
+++ b/tests/trilinos_tpetra/sparse_matrix_vector_01.cc
@@ -28,12 +28,11 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v,
-     LinearAlgebra::TpetraWrappers::Vector<double> &w)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v,
+     LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &w)
 {
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double> m(v.size(),
-                                                        v.size(),
-                                                        v.size());
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> 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);
@@ -76,9 +75,9 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
-        LinearAlgebra::TpetraWrappers::Vector<double> w;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> w;
         w.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v, w);
       }
diff --git a/tests/trilinos_tpetra/sparse_matrix_vector_02.cc b/tests/trilinos_tpetra/sparse_matrix_vector_02.cc
index 20382f5c84..73bd4e1186 100644
--- a/tests/trilinos_tpetra/sparse_matrix_vector_02.cc
+++ b/tests/trilinos_tpetra/sparse_matrix_vector_02.cc
@@ -28,12 +28,11 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v,
-     LinearAlgebra::TpetraWrappers::Vector<double> &w)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v,
+     LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &w)
 {
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double> m(v.size(),
-                                                        v.size(),
-                                                        v.size());
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> 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);
@@ -76,9 +75,9 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
-        LinearAlgebra::TpetraWrappers::Vector<double> w;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> w;
         w.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v, w);
       }
diff --git a/tests/trilinos_tpetra/sparse_matrix_vector_03.cc b/tests/trilinos_tpetra/sparse_matrix_vector_03.cc
index 90e920f7f6..b19c9ae79e 100644
--- a/tests/trilinos_tpetra/sparse_matrix_vector_03.cc
+++ b/tests/trilinos_tpetra/sparse_matrix_vector_03.cc
@@ -28,12 +28,11 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v,
-     LinearAlgebra::TpetraWrappers::Vector<double> &w)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v,
+     LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &w)
 {
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double> m(v.size(),
-                                                        v.size(),
-                                                        v.size());
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> 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);
@@ -79,9 +78,9 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
-        LinearAlgebra::TpetraWrappers::Vector<double> w;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> w;
         w.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v, w);
       }
diff --git a/tests/trilinos_tpetra/sparse_matrix_vector_04.cc b/tests/trilinos_tpetra/sparse_matrix_vector_04.cc
index 5369704b5d..f99dfc872c 100644
--- a/tests/trilinos_tpetra/sparse_matrix_vector_04.cc
+++ b/tests/trilinos_tpetra/sparse_matrix_vector_04.cc
@@ -28,12 +28,11 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v,
-     LinearAlgebra::TpetraWrappers::Vector<double> &w)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v,
+     LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &w)
 {
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double> m(v.size(),
-                                                        v.size(),
-                                                        v.size());
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> 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);
@@ -78,9 +77,9 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
-        LinearAlgebra::TpetraWrappers::Vector<double> w;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> w;
         w.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v, w);
       }
diff --git a/tests/trilinos_tpetra/sparse_matrix_vector_05.cc b/tests/trilinos_tpetra/sparse_matrix_vector_05.cc
index 2e57681848..b18a4a6aec 100644
--- a/tests/trilinos_tpetra/sparse_matrix_vector_05.cc
+++ b/tests/trilinos_tpetra/sparse_matrix_vector_05.cc
@@ -28,12 +28,11 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v,
-     LinearAlgebra::TpetraWrappers::Vector<double> &w)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v,
+     LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &w)
 {
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double> m(v.size(),
-                                                        v.size(),
-                                                        v.size());
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> 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);
@@ -82,9 +81,9 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
-        LinearAlgebra::TpetraWrappers::Vector<double> w;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> w;
         w.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v, w);
       }
diff --git a/tests/trilinos_tpetra/sparse_matrix_vector_06.cc b/tests/trilinos_tpetra/sparse_matrix_vector_06.cc
index 62e3e07952..5f5a54980c 100644
--- a/tests/trilinos_tpetra/sparse_matrix_vector_06.cc
+++ b/tests/trilinos_tpetra/sparse_matrix_vector_06.cc
@@ -28,11 +28,10 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v)
 {
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double> m(v.size(),
-                                                        v.size(),
-                                                        v.size());
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> 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);
@@ -74,7 +73,7 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(30), MPI_COMM_WORLD);
         test(v);
       }
diff --git a/tests/trilinos_tpetra/sparse_matrix_vector_07.cc b/tests/trilinos_tpetra/sparse_matrix_vector_07.cc
index edaae70e3f..3037c1bef2 100644
--- a/tests/trilinos_tpetra/sparse_matrix_vector_07.cc
+++ b/tests/trilinos_tpetra/sparse_matrix_vector_07.cc
@@ -28,13 +28,12 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v,
-     LinearAlgebra::TpetraWrappers::Vector<double> &w,
-     LinearAlgebra::TpetraWrappers::Vector<double> &x)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v,
+     LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &w,
+     LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &x)
 {
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double> m(v.size(),
-                                                        v.size(),
-                                                        v.size());
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> 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);
@@ -84,11 +83,11 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
-        LinearAlgebra::TpetraWrappers::Vector<double> w;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> w;
         w.reinit(complete_index_set(100), MPI_COMM_WORLD);
-        LinearAlgebra::TpetraWrappers::Vector<double> x;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> x;
         x.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v, w, x);
       }
diff --git a/tests/trilinos_tpetra/sparse_matrix_vector_08.cc b/tests/trilinos_tpetra/sparse_matrix_vector_08.cc
index 47e28b6b4b..b0f1ab0301 100644
--- a/tests/trilinos_tpetra/sparse_matrix_vector_08.cc
+++ b/tests/trilinos_tpetra/sparse_matrix_vector_08.cc
@@ -30,9 +30,8 @@
 void
 test(Vector<double> &v, Vector<double> &w)
 {
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double> m(w.size(),
-                                                        v.size(),
-                                                        v.size());
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> m(
+    w.size(), v.size(), v.size());
   for (unsigned int i = 0; i < m.m(); ++i)
     for (unsigned int j = 0; j < m.n(); ++j)
       m.set(i, j, i + 2 * j);
diff --git a/tests/trilinos_tpetra/sparse_matrix_vector_14.cc b/tests/trilinos_tpetra/sparse_matrix_vector_14.cc
index f77ec86647..cce4fa7580 100644
--- a/tests/trilinos_tpetra/sparse_matrix_vector_14.cc
+++ b/tests/trilinos_tpetra/sparse_matrix_vector_14.cc
@@ -30,15 +30,14 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v,
-     LinearAlgebra::TpetraWrappers::Vector<double> &w)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v,
+     LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &w)
 {
   LinearAlgebra::ReadWriteVector<double> read_write_v(v.size());
   LinearAlgebra::ReadWriteVector<double> read_write_w(w.size());
 
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double> m(w.size(),
-                                                        v.size(),
-                                                        v.size());
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> m(
+    w.size(), v.size(), v.size());
   for (unsigned int i = 0; i < m.m(); ++i)
     for (unsigned int j = 0; j < m.n(); ++j)
       m.set(i, j, i + 2 * j);
@@ -91,10 +90,10 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v(complete_index_set(100),
-                                                        MPI_COMM_SELF);
-        LinearAlgebra::TpetraWrappers::Vector<double> w(complete_index_set(95),
-                                                        MPI_COMM_SELF);
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v(
+          complete_index_set(100), MPI_COMM_SELF);
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> w(
+          complete_index_set(95), MPI_COMM_SELF);
         test(v, w);
       }
     }
diff --git a/tests/trilinos_tpetra/sparse_matrix_vector_15.cc b/tests/trilinos_tpetra/sparse_matrix_vector_15.cc
index 4bb72978df..b4227d294b 100644
--- a/tests/trilinos_tpetra/sparse_matrix_vector_15.cc
+++ b/tests/trilinos_tpetra/sparse_matrix_vector_15.cc
@@ -30,15 +30,14 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v,
-     LinearAlgebra::TpetraWrappers::Vector<double> &w)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v,
+     LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &w)
 {
   LinearAlgebra::ReadWriteVector<double> read_write_v(v.size());
   LinearAlgebra::ReadWriteVector<double> read_write_w(w.size());
 
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double> m(v.size(),
-                                                        w.size(),
-                                                        w.size());
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> m(
+    v.size(), w.size(), w.size());
   for (unsigned int i = 0; i < m.m(); ++i)
     for (unsigned int j = 0; j < m.n(); ++j)
       m.set(i, j, i + 2 * j);
@@ -90,10 +89,10 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v(complete_index_set(95),
-                                                        MPI_COMM_SELF);
-        LinearAlgebra::TpetraWrappers::Vector<double> w(complete_index_set(100),
-                                                        MPI_COMM_SELF);
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v(
+          complete_index_set(95), MPI_COMM_SELF);
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> w(
+          complete_index_set(100), MPI_COMM_SELF);
         test(v, w);
       }
     }
diff --git a/tests/trilinos_tpetra/tpetra_vector_01.cc b/tests/trilinos_tpetra/tpetra_vector_01.cc
index 3504607d32..75d5895f8c 100644
--- a/tests/trilinos_tpetra/tpetra_vector_01.cc
+++ b/tests/trilinos_tpetra/tpetra_vector_01.cc
@@ -44,10 +44,10 @@ test()
     }
   parallel_partitioner_1.compress();
   parallel_partitioner_2.compress();
-  LinearAlgebra::TpetraWrappers::Vector<Number> a;
-  LinearAlgebra::TpetraWrappers::Vector<Number> b(parallel_partitioner_1,
-                                                  MPI_COMM_WORLD);
-  LinearAlgebra::TpetraWrappers::Vector<Number> c(b);
+  LinearAlgebra::TpetraWrappers::Vector<Number, MemorySpace::Default> a;
+  LinearAlgebra::TpetraWrappers::Vector<Number, MemorySpace::Default> b(
+    parallel_partitioner_1, MPI_COMM_WORLD);
+  LinearAlgebra::TpetraWrappers::Vector<Number, MemorySpace::Default> c(b);
 
   AssertThrow(a.size() == 0, ExcMessage("Vector has the wrong size."));
   AssertThrow(b.size() == 10, ExcMessage("Vector has the wrong size."));
diff --git a/tests/trilinos_tpetra/tpetra_vector_02.cc b/tests/trilinos_tpetra/tpetra_vector_02.cc
index 0072891a81..0a85342b8c 100644
--- a/tests/trilinos_tpetra/tpetra_vector_02.cc
+++ b/tests/trilinos_tpetra/tpetra_vector_02.cc
@@ -44,12 +44,12 @@ test()
     }
   parallel_partitioner_1.compress();
   parallel_partitioner_2.compress();
-  LinearAlgebra::TpetraWrappers::Vector<Number> a(parallel_partitioner_1,
-                                                  MPI_COMM_WORLD);
-  LinearAlgebra::TpetraWrappers::Vector<Number> b(parallel_partitioner_1,
-                                                  MPI_COMM_WORLD);
-  LinearAlgebra::TpetraWrappers::Vector<Number> c(parallel_partitioner_2,
-                                                  MPI_COMM_WORLD);
+  LinearAlgebra::TpetraWrappers::Vector<Number, MemorySpace::Default> a(
+    parallel_partitioner_1, MPI_COMM_WORLD);
+  LinearAlgebra::TpetraWrappers::Vector<Number, MemorySpace::Default> b(
+    parallel_partitioner_1, MPI_COMM_WORLD);
+  LinearAlgebra::TpetraWrappers::Vector<Number, MemorySpace::Default> c(
+    parallel_partitioner_2, MPI_COMM_WORLD);
 
   IndexSet read_write_index_set(10);
   if (rank == 0)
@@ -117,7 +117,7 @@ test()
     }
 
 
-  LinearAlgebra::TpetraWrappers::Vector<Number> d(a);
+  LinearAlgebra::TpetraWrappers::Vector<Number, MemorySpace::Default> d(a);
   a.add(2., b, 3., d);
   read_write_3.import_elements(a, VectorOperation::insert);
   if (rank == 0)
diff --git a/tests/trilinos_tpetra/tpetra_vector_03.cc b/tests/trilinos_tpetra/tpetra_vector_03.cc
index 29d705f544..43eb2a0b06 100644
--- a/tests/trilinos_tpetra/tpetra_vector_03.cc
+++ b/tests/trilinos_tpetra/tpetra_vector_03.cc
@@ -38,8 +38,8 @@ test()
   locally_owned.add_range(my_id * 2, my_id * 2 + 2);
   locally_owned.compress();
 
-  LinearAlgebra::TpetraWrappers::Vector<Number> v(locally_owned,
-                                                  MPI_COMM_WORLD);
+  LinearAlgebra::TpetraWrappers::Vector<Number, MemorySpace::Default> v(
+    locally_owned, MPI_COMM_WORLD);
 
   IndexSet     workaround_set(n_procs * 2);
   unsigned int my_first_index  = (my_id * 2 + 2) % (n_procs * 2);
diff --git a/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_01.cc b/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_01.cc
index 9c167c52a6..aeca862901 100644
--- a/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_01.cc
+++ b/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_01.cc
@@ -27,16 +27,19 @@
 void
 test()
 {
-  LinearAlgebra::TpetraWrappers::SparsityPattern sp(5, 5, 3);
+  LinearAlgebra::TpetraWrappers::SparsityPattern<MemorySpace::Default> sp(5,
+                                                                          5,
+                                                                          3);
   for (unsigned int i = 0; i < 5; ++i)
     for (unsigned int j = 0; j < 5; ++j)
       if ((i + 2 * j + 1) % 3 == 0)
         sp.add(i, j);
   sp.compress();
 
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double>                 m(sp);
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double>::const_iterator i =
-    m.begin();
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> m(
+    sp);
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default>::
+    const_iterator i = m.begin();
   deallog << i->value() << std::endl;
   ++i;
   deallog << i->value() << std::endl;
diff --git a/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_02.cc b/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_02.cc
index 0d4020bd1c..6e97f5d751 100644
--- a/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_02.cc
+++ b/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_02.cc
@@ -27,21 +27,24 @@
 void
 test()
 {
-  LinearAlgebra::TpetraWrappers::SparsityPattern sp(5, 5, 3);
+  LinearAlgebra::TpetraWrappers::SparsityPattern<MemorySpace::Default> sp(5,
+                                                                          5,
+                                                                          3);
   for (unsigned int i = 0; i < 5; ++i)
     for (unsigned int j = 0; j < 5; ++j)
       if (((i + 2 * j + 1) % 3 == 0) || (i == j))
         sp.add(i, j);
   sp.compress();
 
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double> m(sp);
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> m(
+    sp);
   for (unsigned int i = 0; i < 5; ++i)
     for (unsigned int j = 0; j < 5; ++j)
       if (((i + 2 * j + 1) % 3 == 0) || (i == j))
         m.set(i, j, i * j);
 
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double>::const_iterator i =
-    m.begin();
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default>::
+    const_iterator i = m.begin();
   for (; i != m.end(); ++i)
     {
       deallog << i->row() << ' ' << i->column() << ' ' << i->value()
diff --git a/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_03.cc b/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_03.cc
index 4ef6c64c24..bee0fe2616 100644
--- a/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_03.cc
+++ b/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_03.cc
@@ -27,20 +27,25 @@
 void
 test()
 {
-  LinearAlgebra::TpetraWrappers::SparsityPattern sp(5, 5, 3);
+  LinearAlgebra::TpetraWrappers::SparsityPattern<MemorySpace::Default> sp(5,
+                                                                          5,
+                                                                          3);
   for (unsigned int i = 0; i < 5; ++i)
     for (unsigned int j = 0; j < 5; ++j)
       if (((i + 2 * j + 1) % 3 == 0) || (i == j))
         sp.add(i, j);
   sp.compress();
 
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double> m(sp);
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> m(
+    sp);
   for (unsigned int i = 0; i < 5; ++i)
     for (unsigned int j = 0; j < 5; ++j)
       if (((i + 2 * j + 1) % 3 == 0) || (i == j))
         m.set(i, j, i * j);
 
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double>::iterator i = m.begin();
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double,
+                                              MemorySpace::Default>::iterator
+    i = m.begin();
   for (; i != m.end(); ++i)
     {
       deallog << i->row() << ' ' << i->column() << ' ' << i->value()
diff --git a/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_04.cc b/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_04.cc
index 60558dd201..8d4badda16 100644
--- a/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_04.cc
+++ b/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_04.cc
@@ -26,15 +26,20 @@
 void
 test()
 {
-  LinearAlgebra::TpetraWrappers::SparsityPattern sp(5, 5, 3);
+  LinearAlgebra::TpetraWrappers::SparsityPattern<MemorySpace::Default> sp(5,
+                                                                          5,
+                                                                          3);
   for (unsigned int i = 0; i < 5; ++i)
     for (unsigned int j = 0; j < 5; ++j)
       if (((i + 2 * j + 1) % 3 == 0) || (i == j))
         sp.add(i, j);
   sp.compress();
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double> m(sp);
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> m(
+    sp);
 
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double>::iterator i = m.begin();
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double,
+                                              MemorySpace::Default>::iterator
+    i = m.begin();
   for (; i != m.end(); ++i)
     i->value() = i->row() * i->column();
   m.compress(VectorOperation::insert);
diff --git a/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_05.cc b/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_05.cc
index 9470ad84f8..20332c4f4f 100644
--- a/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_05.cc
+++ b/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_05.cc
@@ -26,20 +26,25 @@
 void
 test()
 {
-  LinearAlgebra::TpetraWrappers::SparsityPattern sp(5, 5, 3);
+  LinearAlgebra::TpetraWrappers::SparsityPattern<MemorySpace::Default> sp(5,
+                                                                          5,
+                                                                          3);
   for (unsigned int i = 0; i < 5; ++i)
     for (unsigned int j = 0; j < 5; ++j)
       if (((i + 2 * j + 1) % 3 == 0) || (i == j))
         sp.add(i, j);
   sp.compress();
 
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double> m(sp);
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> m(
+    sp);
   for (unsigned int i = 0; i < 5; ++i)
     for (unsigned int j = 0; j < 5; ++j)
       if (((i + 2 * j + 1) % 3 == 0) || (i == j))
         m.set(i, j, 1.);
 
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double>::iterator i = m.begin();
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double,
+                                              MemorySpace::Default>::iterator
+    i = m.begin();
   for (; i != m.end(); ++i)
     i->value() += i->row() * i->column();
   m.compress(VectorOperation::insert);
diff --git a/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_06.cc b/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_06.cc
index e05fbd2e1f..3b974dadb0 100644
--- a/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_06.cc
+++ b/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_06.cc
@@ -26,20 +26,25 @@
 void
 test()
 {
-  LinearAlgebra::TpetraWrappers::SparsityPattern sp(5, 5, 3);
+  LinearAlgebra::TpetraWrappers::SparsityPattern<MemorySpace::Default> sp(5,
+                                                                          5,
+                                                                          3);
   for (unsigned int i = 0; i < 5; ++i)
     for (unsigned int j = 0; j < 5; ++j)
       if (((i + 2 * j + 1) % 3 == 0) || (i == j))
         sp.add(i, j);
   sp.compress();
 
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double> m(sp);
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> m(
+    sp);
   for (unsigned int i = 0; i < 5; ++i)
     for (unsigned int j = 0; j < 5; ++j)
       if (((i + 2 * j + 1) % 3 == 0) || (i == j))
         m.set(i, j, 1.);
 
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double>::iterator i = m.begin();
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double,
+                                              MemorySpace::Default>::iterator
+    i = m.begin();
   for (; i != m.end(); ++i)
     i->value() -= i->row() * i->column();
   m.compress(VectorOperation::insert);
diff --git a/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_07.cc b/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_07.cc
index b1ab103824..465369b044 100644
--- a/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_07.cc
+++ b/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_07.cc
@@ -26,20 +26,25 @@
 void
 test()
 {
-  LinearAlgebra::TpetraWrappers::SparsityPattern sp(5, 5, 3);
+  LinearAlgebra::TpetraWrappers::SparsityPattern<MemorySpace::Default> sp(5,
+                                                                          5,
+                                                                          3);
   for (unsigned int i = 0; i < 5; ++i)
     for (unsigned int j = 0; j < 5; ++j)
       if (((i + 2 * j + 1) % 3 == 0) || (i == j))
         sp.add(i, j);
   sp.compress();
 
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double> m(sp);
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> m(
+    sp);
   for (unsigned int i = 0; i < 5; ++i)
     for (unsigned int j = 0; j < 5; ++j)
       if (((i + 2 * j + 1) % 3 == 0) || (i == j))
         m.set(i, j, i * j);
 
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double>::iterator i = m.begin();
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double,
+                                              MemorySpace::Default>::iterator
+    i = m.begin();
   for (; i != m.end(); ++i)
     i->value() *= 2;
   m.compress(VectorOperation::insert);
diff --git a/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_08.cc b/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_08.cc
index af5b5a58ea..30a9ffb377 100644
--- a/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_08.cc
+++ b/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_08.cc
@@ -26,20 +26,25 @@
 void
 test()
 {
-  LinearAlgebra::TpetraWrappers::SparsityPattern sp(5, 5, 3);
+  LinearAlgebra::TpetraWrappers::SparsityPattern<MemorySpace::Default> sp(5,
+                                                                          5,
+                                                                          3);
   for (unsigned int i = 0; i < 5; ++i)
     for (unsigned int j = 0; j < 5; ++j)
       if (((i + 2 * j + 1) % 3 == 0) || (i == j))
         sp.add(i, j);
   sp.compress();
 
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double> m(sp);
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> m(
+    sp);
   for (unsigned int i = 0; i < 5; ++i)
     for (unsigned int j = 0; j < 5; ++j)
       if (((i + 2 * j + 1) % 3 == 0) || (i == j))
         m.set(i, j, i * j);
 
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double>::iterator i = m.begin();
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double,
+                                              MemorySpace::Default>::iterator
+    i = m.begin();
   for (; i != m.end(); ++i)
     i->value() /= 2;
   m.compress(VectorOperation::insert);
diff --git a/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_09.cc b/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_09.cc
index e5598b203a..761d04f536 100644
--- a/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_09.cc
+++ b/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_09.cc
@@ -31,17 +31,20 @@ test()
   // create a sparsity pattern with totally
   // empty lines (not even diagonals, since
   // not quadratic)
-  LinearAlgebra::TpetraWrappers::SparsityPattern sparsity(4, 5, 1);
+  LinearAlgebra::TpetraWrappers::SparsityPattern<MemorySpace::Default> sparsity(
+    4, 5, 1);
   sparsity.add(1, 1);
   sparsity.add(3, 1);
   sparsity.compress();
 
   // attach a sparse matrix to it
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double> A(sparsity);
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> A(
+    sparsity);
 
   // and loop over the elements of it
-  for (LinearAlgebra::TpetraWrappers::SparseMatrix<double>::const_iterator k =
-         A.begin();
+  for (LinearAlgebra::TpetraWrappers::
+         SparseMatrix<double, MemorySpace::Default>::const_iterator k =
+           A.begin();
        k != A.end();
        ++k)
     deallog << k->row() << ' ' << k->column() << ' ' << k->value() << std::endl;
diff --git a/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_10.cc b/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_10.cc
index 5ced6a843f..21f16b0dfe 100644
--- a/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_10.cc
+++ b/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_10.cc
@@ -33,17 +33,20 @@ test()
   // create a sparsity pattern with totally
   // empty lines (not even diagonals, since
   // not quadratic)
-  LinearAlgebra::TpetraWrappers::SparsityPattern sparsity(4, 5, 1);
+  LinearAlgebra::TpetraWrappers::SparsityPattern<MemorySpace::Default> sparsity(
+    4, 5, 1);
   sparsity.add(1, 1);
   sparsity.add(3, 1);
   sparsity.compress();
 
   // attach a sparse matrix to it
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double> A(sparsity);
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> A(
+    sparsity);
 
   // and loop over the elements of it
-  for (LinearAlgebra::TpetraWrappers::SparseMatrix<double>::const_iterator k =
-         A.begin();
+  for (LinearAlgebra::TpetraWrappers::
+         SparseMatrix<double, MemorySpace::Default>::const_iterator k =
+           A.begin();
        k != A.end();
        k++)
     deallog << k->row() << ' ' << k->column() << ' ' << k->value() << std::endl;
diff --git a/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_11.cc b/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_11.cc
index 2ae83bda57..ff77e782cd 100644
--- a/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_11.cc
+++ b/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_11.cc
@@ -30,17 +30,20 @@ test()
   // create a sparsity pattern with totally
   // empty lines (not even diagonals, since
   // not quadratic)
-  LinearAlgebra::TpetraWrappers::SparsityPattern sparsity(4, 5, 1);
+  LinearAlgebra::TpetraWrappers::SparsityPattern<MemorySpace::Default> sparsity(
+    4, 5, 1);
   sparsity.add(1, 1);
   sparsity.add(3, 1);
   sparsity.compress();
 
   // attach a sparse matrix to it
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double> A(sparsity);
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> A(
+    sparsity);
 
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double>::iterator k = A.begin(),
-                                                                j = std::next(
-                                                                  A.begin());
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double,
+                                              MemorySpace::Default>::iterator
+    k = A.begin(),
+    j = std::next(A.begin());
 
   AssertThrow(k < j, ExcInternalError());
   AssertThrow(j > k, ExcInternalError());
diff --git a/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_12.cc b/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_12.cc
index f9bcb8ca16..7832676093 100644
--- a/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_12.cc
+++ b/tests/trilinos_tpetra/trilinos_sparse_matrix_iterator_12.cc
@@ -30,17 +30,19 @@ test()
   // create a sparsity pattern with totally
   // empty lines (not even diagonals, since
   // not quadratic)
-  LinearAlgebra::TpetraWrappers::SparsityPattern sparsity(4, 5, 1);
+  LinearAlgebra::TpetraWrappers::SparsityPattern<MemorySpace::Default> sparsity(
+    4, 5, 1);
   sparsity.add(1, 1);
   sparsity.add(3, 1);
   sparsity.compress();
 
   // attach a sparse matrix to it
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double> A(sparsity);
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default> A(
+    sparsity);
 
-  LinearAlgebra::TpetraWrappers::SparseMatrix<double>::const_iterator
-    k = A.begin(),
-    j = std::next(A.begin());
+  LinearAlgebra::TpetraWrappers::SparseMatrix<double, MemorySpace::Default>::
+    const_iterator k = A.begin(),
+                   j = std::next(A.begin());
 
   AssertThrow(k < j, ExcInternalError());
   AssertThrow(j > k, ExcInternalError());
diff --git a/tests/trilinos_tpetra/vector_assign_01.cc b/tests/trilinos_tpetra/vector_assign_01.cc
index 43fa1cbaf8..f0e167e5e9 100644
--- a/tests/trilinos_tpetra/vector_assign_01.cc
+++ b/tests/trilinos_tpetra/vector_assign_01.cc
@@ -33,8 +33,8 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v,
-     LinearAlgebra::TpetraWrappers::Vector<double> &w)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v,
+     LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &w)
 {
   // set the first vector
   for (unsigned int i = 0; i < v.size(); ++i)
@@ -64,9 +64,9 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
-        LinearAlgebra::TpetraWrappers::Vector<double> w;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> w;
         w.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v, w);
       }
diff --git a/tests/trilinos_tpetra/vector_assign_02.cc b/tests/trilinos_tpetra/vector_assign_02.cc
index 9a25a63ccf..c351417b8a 100644
--- a/tests/trilinos_tpetra/vector_assign_02.cc
+++ b/tests/trilinos_tpetra/vector_assign_02.cc
@@ -32,8 +32,8 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v,
-     LinearAlgebra::TpetraWrappers::Vector<double> &w)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v,
+     LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &w)
 {
   // set the first vector
   for (unsigned int i = 0; i < v.size(); ++i)
@@ -63,9 +63,9 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
-        LinearAlgebra::TpetraWrappers::Vector<double> w;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> w;
         w.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v, w);
       }
diff --git a/tests/trilinos_tpetra/vector_equality_1.cc b/tests/trilinos_tpetra/vector_equality_1.cc
index 5b74a951b2..d2876eb603 100644
--- a/tests/trilinos_tpetra/vector_equality_1.cc
+++ b/tests/trilinos_tpetra/vector_equality_1.cc
@@ -28,8 +28,8 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v,
-     LinearAlgebra::TpetraWrappers::Vector<double> &w)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v,
+     LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &w)
 {
   // set only certain elements of each
   // vector
@@ -59,9 +59,9 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
-        LinearAlgebra::TpetraWrappers::Vector<double> w;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> w;
         w.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v, w);
       }
diff --git a/tests/trilinos_tpetra/vector_equality_2.cc b/tests/trilinos_tpetra/vector_equality_2.cc
index a7b43e2f9d..8b873db8f6 100644
--- a/tests/trilinos_tpetra/vector_equality_2.cc
+++ b/tests/trilinos_tpetra/vector_equality_2.cc
@@ -29,8 +29,8 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v,
-     LinearAlgebra::TpetraWrappers::Vector<double> &w)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v,
+     LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &w)
 {
   // set only certain elements of each
   // vector
@@ -62,9 +62,9 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
-        LinearAlgebra::TpetraWrappers::Vector<double> w;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> w;
         w.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v, w);
       }
diff --git a/tests/trilinos_tpetra/vector_equality_3.cc b/tests/trilinos_tpetra/vector_equality_3.cc
index 0e189172f7..9226e1c98b 100644
--- a/tests/trilinos_tpetra/vector_equality_3.cc
+++ b/tests/trilinos_tpetra/vector_equality_3.cc
@@ -28,8 +28,8 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v,
-     LinearAlgebra::TpetraWrappers::Vector<double> &w)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v,
+     LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &w)
 {
   // set only certain elements of each
   // vector
@@ -59,9 +59,9 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
-        LinearAlgebra::TpetraWrappers::Vector<double> w;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> w;
         w.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v, w);
       }
diff --git a/tests/trilinos_tpetra/vector_equality_4.cc b/tests/trilinos_tpetra/vector_equality_4.cc
index e86fc16318..1b910c76c7 100644
--- a/tests/trilinos_tpetra/vector_equality_4.cc
+++ b/tests/trilinos_tpetra/vector_equality_4.cc
@@ -28,8 +28,8 @@
 
 
 void
-test(LinearAlgebra::TpetraWrappers::Vector<double> &v,
-     LinearAlgebra::TpetraWrappers::Vector<double> &w)
+test(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v,
+     LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &w)
 {
   // set only certain elements of each
   // vector
@@ -61,9 +61,9 @@ main(int argc, char **argv)
   try
     {
       {
-        LinearAlgebra::TpetraWrappers::Vector<double> v;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
         v.reinit(complete_index_set(100), MPI_COMM_WORLD);
-        LinearAlgebra::TpetraWrappers::Vector<double> w;
+        LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> w;
         w.reinit(complete_index_set(100), MPI_COMM_WORLD);
         test(v, w);
       }
diff --git a/tests/trilinos_tpetra/vector_reference_01.cc b/tests/trilinos_tpetra/vector_reference_01.cc
index 27c67d6686..fbd7f8e46c 100644
--- a/tests/trilinos_tpetra/vector_reference_01.cc
+++ b/tests/trilinos_tpetra/vector_reference_01.cc
@@ -31,7 +31,7 @@
 void
 test()
 {
-  LinearAlgebra::TpetraWrappers::Vector<double> v;
+  LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
   v.reinit(complete_index_set(3), MPI_COMM_WORLD);
   v(0) = 0;
   v(1) = 1;
diff --git a/tests/trilinos_tpetra/vector_reinit.cc b/tests/trilinos_tpetra/vector_reinit.cc
index 24af4e3c73..7a3f4387cc 100644
--- a/tests/trilinos_tpetra/vector_reinit.cc
+++ b/tests/trilinos_tpetra/vector_reinit.cc
@@ -58,7 +58,7 @@ test()
 
   IndexSet locally_owned_dofs = dof_handler.locally_owned_dofs();
 
-  LinearAlgebra::TpetraWrappers::Vector<double> vector_Re;
+  LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> vector_Re;
   vector_Re.reinit(locally_owned_dofs, mpi_communicator);
 
   deallog << "OK" << std::endl;
diff --git a/tests/trilinos_tpetra/vector_swap_01.cc b/tests/trilinos_tpetra/vector_swap_01.cc
index 574b0deb7a..c7cad30bac 100644
--- a/tests/trilinos_tpetra/vector_swap_01.cc
+++ b/tests/trilinos_tpetra/vector_swap_01.cc
@@ -28,7 +28,7 @@
 #include "../tests.h"
 
 void
-print(LinearAlgebra::TpetraWrappers::Vector<double> &v)
+print(LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> &v)
 {
   deallog << "size= " << v.size() << " el(0)= " << v(0)
           << " l2norm()= " << v.l2_norm() << std::endl;
@@ -38,11 +38,11 @@ print(LinearAlgebra::TpetraWrappers::Vector<double> &v)
 void
 test()
 {
-  LinearAlgebra::TpetraWrappers::Vector<double> v;
+  LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> v;
   v.reinit(complete_index_set(5), MPI_COMM_WORLD);
   for (unsigned int i = 0; i < v.size(); ++i)
     v(i) = 1;
-  LinearAlgebra::TpetraWrappers::Vector<double> w;
+  LinearAlgebra::TpetraWrappers::Vector<double, MemorySpace::Default> w;
   w.reinit(complete_index_set(9), MPI_COMM_WORLD);
   for (unsigned int i = 0; i < w.size(); ++i)
     w(i) = 2;
-- 
2.39.5