]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Deprecate TrilinosWrappers::Vector and TrilinosWrappers::BlockVector and update the...
authorBruno Turcksin <bruno.turcksin@gmail.com>
Sat, 25 Apr 2015 00:17:56 +0000 (19:17 -0500)
committerBruno Turcksin <bruno.turcksin@gmail.com>
Tue, 5 May 2015 19:16:22 +0000 (14:16 -0500)
examples/step-31/step-31.cc
examples/step-32/step-32.cc
examples/step-41/step-41.cc
examples/step-43/step-43.cc
include/deal.II/lac/trilinos_block_vector.h
include/deal.II/lac/trilinos_vector.h

index 8d53da03f5b01928659dc954cc9bf1b48d04577f..d3ffb8ff99b377f9b640023c1f236b0e239352e4 100644 (file)
@@ -58,6 +58,7 @@
 // preconditioner classes that implement interfaces to the respective Trilinos
 // classes. In particular, we will need interfaces to the matrix and vector
 // classes based on Trilinos as well as Trilinos preconditioners:
+#include <deal.II/base/index_set.h>
 #include <deal.II/lac/trilinos_sparse_matrix.h>
 #include <deal.II/lac/trilinos_block_sparse_matrix.h>
 #include <deal.II/lac/trilinos_vector.h>
@@ -378,8 +379,8 @@ namespace Step31
         PreconditionerMp>         &Mpinv,
         const PreconditionerA                         &Apreconditioner);
 
-      void vmult (TrilinosWrappers::BlockVector       &dst,
-                  const TrilinosWrappers::BlockVector &src) const;
+      void vmult (TrilinosWrappers::MPI::BlockVector       &dst,
+                  const TrilinosWrappers::MPI::BlockVector &src) const;
 
     private:
       const SmartPointer<const TrilinosWrappers::BlockSparseMatrix> stokes_matrix;
@@ -387,7 +388,7 @@ namespace Step31
             PreconditionerMp > > m_inverse;
       const PreconditionerA &a_preconditioner;
 
-      mutable TrilinosWrappers::Vector tmp;
+      mutable TrilinosWrappers::MPI::Vector tmp;
     };
 
 
@@ -401,9 +402,19 @@ namespace Step31
       :
       stokes_matrix           (&S),
       m_inverse               (&Mpinv),
-      a_preconditioner        (Apreconditioner),
-      tmp                     (stokes_matrix->block(1,1).m())
-    {}
+      a_preconditioner        (Apreconditioner)
+    {
+      // When using a TrilinosWrappers::MPI::Vector or a
+      // TrilinosWrappers::MPI::BlockVector, the Vector is initialized using an
+      // IndexSet. IndexSet is used not only to resize the
+      // TrilinosWrappers::MPI::Vector but it also associates an index in the
+      // TrilinosWrappers::MPI::Vector with a degree of freedom (see step-40 for
+      // a more detailed explanation). This assocation is done by the add_range()
+      // function.
+      IndexSet tmp_index_set(stokes_matrix->block(1,1).m());
+      tmp_index_set.add_range(0,stokes_matrix->block(1,1).m());
+      tmp.reinit(tmp_index_set);
+    }
 
 
     // Next is the <code>vmult</code> function. We implement the action of
@@ -424,8 +435,8 @@ namespace Step31
     template <class PreconditionerA, class PreconditionerMp>
     void
     BlockSchurPreconditioner<PreconditionerA, PreconditionerMp>::
-    vmult (TrilinosWrappers::BlockVector       &dst,
-           const TrilinosWrappers::BlockVector &src) const
+    vmult (TrilinosWrappers::MPI::BlockVector       &dst,
+           const TrilinosWrappers::MPI::BlockVector &src) const
     {
       a_preconditioner.vmult (dst.block(0), src.block(0));
       stokes_matrix->block(1,0).residual(tmp, dst.block(0), src.block(1));
@@ -502,13 +513,13 @@ namespace Step31
     DoFHandler<dim>                     stokes_dof_handler;
     ConstraintMatrix                    stokes_constraints;
 
-    std::vector<types::global_dof_index> stokes_block_sizes;
+    std::vector<IndexSet>               stokes_block_sizes;
     TrilinosWrappers::BlockSparseMatrix stokes_matrix;
     TrilinosWrappers::BlockSparseMatrix stokes_preconditioner_matrix;
 
-    TrilinosWrappers::BlockVector       stokes_solution;
-    TrilinosWrappers::BlockVector       old_stokes_solution;
-    TrilinosWrappers::BlockVector       stokes_rhs;
+    TrilinosWrappers::MPI::BlockVector  stokes_solution;
+    TrilinosWrappers::MPI::BlockVector  old_stokes_solution;
+    TrilinosWrappers::MPI::BlockVector  stokes_rhs;
 
 
     const unsigned int                  temperature_degree;
@@ -520,10 +531,10 @@ namespace Step31
     TrilinosWrappers::SparseMatrix      temperature_stiffness_matrix;
     TrilinosWrappers::SparseMatrix      temperature_matrix;
 
-    TrilinosWrappers::Vector            temperature_solution;
-    TrilinosWrappers::Vector            old_temperature_solution;
-    TrilinosWrappers::Vector            old_old_temperature_solution;
-    TrilinosWrappers::Vector            temperature_rhs;
+    TrilinosWrappers::MPI::Vector       temperature_solution;
+    TrilinosWrappers::MPI::Vector       old_temperature_solution;
+    TrilinosWrappers::MPI::Vector       old_old_temperature_solution;
+    TrilinosWrappers::MPI::Vector       temperature_rhs;
 
 
     double                              time_step;
@@ -952,9 +963,12 @@ namespace Step31
     // Trilinos matrices store the sparsity pattern internally, there is no
     // need to keep the sparsity pattern around after the initialization of
     // the matrix.
+    stokes_block_sizes.clear();
     stokes_block_sizes.resize (2);
-    stokes_block_sizes[0] = n_u;
-    stokes_block_sizes[1] = n_p;
+    stokes_block_sizes[0].set_size(n_u);
+    stokes_block_sizes[1].set_size(n_p);
+    stokes_block_sizes[0].add_range(0,n_u);
+    stokes_block_sizes[1].add_range(0,n_p);
     {
       stokes_matrix.clear ();
 
@@ -1040,15 +1054,17 @@ namespace Step31
     // and $\mathbf u^{n-2}$, as well as for the temperatures $T^{n}$,
     // $T^{n-1}$ and $T^{n-2}$ (required for time stepping) and all the system
     // right hand sides to their correct sizes and block structure:
+    IndexSet temperature_partitioning (n_T);
+    temperature_partitioning.add_range(0,n_T);
     stokes_solution.reinit (stokes_block_sizes);
     old_stokes_solution.reinit (stokes_block_sizes);
     stokes_rhs.reinit (stokes_block_sizes);
 
-    temperature_solution.reinit (temperature_dof_handler.n_dofs());
-    old_temperature_solution.reinit (temperature_dof_handler.n_dofs());
-    old_old_temperature_solution.reinit (temperature_dof_handler.n_dofs());
+    temperature_solution.reinit (temperature_partitioning);
+    old_temperature_solution.reinit (temperature_partitioning);
+    old_old_temperature_solution.reinit (temperature_partitioning);
 
-    temperature_rhs.reinit (temperature_dof_handler.n_dofs());
+    temperature_rhs.reinit (temperature_partitioning);
   }
 
 
@@ -1793,9 +1809,9 @@ namespace Step31
       SolverControl solver_control (stokes_matrix.m(),
                                     1e-6*stokes_rhs.l2_norm());
 
-      SolverGMRES<TrilinosWrappers::BlockVector>
+      SolverGMRES<TrilinosWrappers::MPI::BlockVector>
       gmres (solver_control,
-             SolverGMRES<TrilinosWrappers::BlockVector >::AdditionalData(100));
+             SolverGMRES<TrilinosWrappers::MPI::BlockVector >::AdditionalData(100));
 
       for (unsigned int i=0; i<stokes_solution.size(); ++i)
         if (stokes_constraints.is_constrained(i))
@@ -1864,7 +1880,7 @@ namespace Step31
     // preconditioner (IC) as we also use for preconditioning the pressure
     // mass matrix solver. As a solver, we choose the conjugate gradient
     // method CG. As before, we tell the solver to use Trilinos vectors via
-    // the template argument <code>TrilinosWrappers::Vector</code>.  Finally,
+    // the template argument <code>TrilinosWrappers::MPI::Vector</code>.  Finally,
     // we solve, distribute the hanging node constraints and write out the
     // number of iterations.
     assemble_temperature_system (maximal_velocity);
@@ -1872,7 +1888,7 @@ namespace Step31
 
       SolverControl solver_control (temperature_matrix.m(),
                                     1e-8*temperature_rhs.l2_norm());
-      SolverCG<TrilinosWrappers::Vector> cg (solver_control);
+      SolverCG<TrilinosWrappers::MPI::Vector> cg (solver_control);
 
       TrilinosWrappers::PreconditionIC preconditioner;
       preconditioner.initialize (temperature_matrix);
@@ -2033,14 +2049,14 @@ namespace Step31
     // and temperature DoFHandler objects, by attaching them to the old dof
     // handlers. With this at place, we can prepare the triangulation and the
     // data vectors for refinement (in this order).
-    std::vector<TrilinosWrappers::Vector> x_temperature (2);
+    std::vector<TrilinosWrappers::MPI::Vector> x_temperature (2);
     x_temperature[0] = temperature_solution;
     x_temperature[1] = old_temperature_solution;
-    TrilinosWrappers::BlockVector x_stokes = stokes_solution;
+    TrilinosWrappers::MPI::BlockVector x_stokes = stokes_solution;
 
-    SolutionTransfer<dim,TrilinosWrappers::Vector>
+    SolutionTransfer<dim,TrilinosWrappers::MPI::Vector>
     temperature_trans(temperature_dof_handler);
-    SolutionTransfer<dim,TrilinosWrappers::BlockVector>
+    SolutionTransfer<dim,TrilinosWrappers::MPI::BlockVector>
     stokes_trans(stokes_dof_handler);
 
     triangulation.prepare_coarsening_and_refinement();
@@ -2062,7 +2078,7 @@ namespace Step31
     triangulation.execute_coarsening_and_refinement ();
     setup_dofs ();
 
-    std::vector<TrilinosWrappers::Vector> tmp (2);
+    std::vector<TrilinosWrappers::MPI::Vector> tmp (2);
     tmp[0].reinit (temperature_solution);
     tmp[1].reinit (temperature_solution);
     temperature_trans.interpolate(x_temperature, tmp);
index 2cc7273befa3c10e5074db966fc57c9353e392cb..dab38039d2c7c752f5f124e702be17407dced990 100644 (file)
@@ -865,18 +865,20 @@ namespace Step32
     // happen to use the same class names for %parallel and sequential data
     // structures, i.e., all matrices will actually be considered %parallel
     // below. On the other hand, for vectors, only those from namespace
-    // TrilinosWrappers::MPI are actually distributed. In particular, we will
-    // frequently have to query velocities and temperatures at arbitrary
-    // quadrature points; consequently, rather than importing ghost
-    // information of a vector whenever we need access to degrees of freedom
-    // that are relevant locally but owned by another processor, we solve
-    // linear systems in %parallel but then immediately initialize a vector
-    // including ghost entries of the solution for further processing. The
-    // various <code>*_solution</code> vectors are therefore filled
-    // immediately after solving their respective linear system in %parallel
-    // and will always contain values for all @ref GlossLocallyRelevantDof
-    // "locally relevant degrees of freedom"; the fully distributed vectors
-    // that we obtain from the solution process and that only ever contain the
+    // TrilinosWrappers::MPI are actually distributed (be aware that
+    // TrilinosWrappers::Vector and TrilinosWrappers::BlockVector are
+    // deprecated). In particular, we will frequently have to query velocities
+    // and temperatures at arbitrary quadrature points; consequently, rather
+    // than importing ghost information of a vector whenever we need access
+    // to degrees of freedom that are relevant locally but owned by another
+    // processor, we solve linear systems in %parallel but then immediately
+    // initialize a vector including ghost entries of the solution for further
+    // processing. The various <code>*_solution</code> vectors are therefore
+    // filled immediately after solving their respective linear system in
+    // %parallel and will always contain values for all
+    // @ref GlossLocallyRelevantDof "locally relevant degrees of freedom";
+    // the fully distributed vectors that we obtain from the solution process
+    // and that only ever contain the
     // @ref GlossLocallyOwnedDof "locally owned degrees of freedom" are
     // destroyed immediately after the solution process and after we have
     // copied the relevant values into the member variable vectors.
index 981df9e03e95af9eba9bc0439dea58293baf4942..8587e3cff89c73dd31fbf0c6c4694172e6d5d918 100644 (file)
@@ -94,15 +94,16 @@ namespace Step41
     DoFHandler<dim>      dof_handler;
     ConstraintMatrix     constraints;
     IndexSet             active_set;
+    IndexSet             solution_index_set;
 
     TrilinosWrappers::SparseMatrix system_matrix;
     TrilinosWrappers::SparseMatrix complete_system_matrix;
 
-    TrilinosWrappers::Vector       solution;
-    TrilinosWrappers::Vector       system_rhs;
-    TrilinosWrappers::Vector       complete_system_rhs;
-    TrilinosWrappers::Vector       diagonal_of_mass_matrix;
-    TrilinosWrappers::Vector       contact_force;
+    TrilinosWrappers::MPI::Vector  solution;
+    TrilinosWrappers::MPI::Vector  system_rhs;
+    TrilinosWrappers::MPI::Vector  complete_system_rhs;
+    TrilinosWrappers::MPI::Vector  diagonal_of_mass_matrix;
+    TrilinosWrappers::MPI::Vector  contact_force;
   };
 
 
@@ -255,10 +256,12 @@ namespace Step41
     system_matrix.reinit (dsp);
     complete_system_matrix.reinit (dsp);
 
-    solution.reinit (dof_handler.n_dofs());
-    system_rhs.reinit (dof_handler.n_dofs());
-    complete_system_rhs.reinit (dof_handler.n_dofs());
-    contact_force.reinit (dof_handler.n_dofs());
+    solution_index_set.set_size(dof_handler.n_dofs());
+    solution_index_set.add_range(0, dof_handler.n_dofs());
+    solution.reinit (solution_index_set);
+    system_rhs.reinit (solution_index_set);
+    complete_system_rhs.reinit (solution_index_set);
+    contact_force.reinit (solution_index_set);
 
     // The only other thing to do here is to compute the factors in the $B$
     // matrix which is used to scale the residual. As discussed in the
@@ -268,7 +271,7 @@ namespace Step41
     TrilinosWrappers::SparseMatrix mass_matrix;
     mass_matrix.reinit (dsp);
     assemble_mass_matrix_diagonal (mass_matrix);
-    diagonal_of_mass_matrix.reinit (dof_handler.n_dofs());
+    diagonal_of_mass_matrix.reinit (solution_index_set);
     for (unsigned int j=0; j<solution.size (); j++)
       diagonal_of_mass_matrix (j) = mass_matrix.diag_element (j);
   }
@@ -300,7 +303,7 @@ namespace Step41
     const unsigned int        n_q_points    = quadrature_formula.size();
 
     FullMatrix<double>        cell_matrix (dofs_per_cell, dofs_per_cell);
-    TrilinosWrappers::Vector  cell_rhs (dofs_per_cell);
+    Vector<double>            cell_rhs (dofs_per_cell);
 
     std::vector<types::global_dof_index> local_dof_indices (dofs_per_cell);
 
@@ -432,7 +435,7 @@ namespace Step41
 
     const double penalty_parameter = 100.0;
 
-    TrilinosWrappers::Vector lambda (dof_handler.n_dofs());
+    TrilinosWrappers::MPI::Vector lambda (solution_index_set);
     complete_system_matrix.residual (lambda,
                                      solution, complete_system_rhs);
     contact_force.ratio (lambda, diagonal_of_mass_matrix);
@@ -556,7 +559,7 @@ namespace Step41
     std::cout << "   Solving system..." << std::endl;
 
     ReductionControl                    reduction_control (100, 1e-12, 1e-3);
-    SolverCG<TrilinosWrappers::Vector>  solver (reduction_control);
+    SolverCG<TrilinosWrappers::MPI::Vector>  solver (reduction_control);
     TrilinosWrappers::PreconditionAMG   precondition;
     precondition.initialize (system_matrix);
 
index 88c273c1df1743e4b98969313c39fd13ecfeac99..2f429b1b6e65de74daf35d5de16d2decd185f3f2 100644 (file)
@@ -33,6 +33,7 @@
 #include <deal.II/base/function.h>
 #include <deal.II/base/tensor_function.h>
 #include <deal.II/base/std_cxx11/shared_ptr.h>
+#include <deal.II/base/index_set.h>
 
 #include <deal.II/lac/full_matrix.h>
 #include <deal.II/lac/solver_gmres.h>
@@ -421,8 +422,8 @@ namespace Step43
         PreconditionerMp>         &Mpinv,
         const PreconditionerA                         &Apreconditioner);
 
-      void vmult (TrilinosWrappers::BlockVector       &dst,
-                  const TrilinosWrappers::BlockVector &src) const;
+      void vmult (TrilinosWrappers::MPI::BlockVector       &dst,
+                  const TrilinosWrappers::MPI::BlockVector &src) const;
 
     private:
       const SmartPointer<const TrilinosWrappers::BlockSparseMatrix> darcy_matrix;
@@ -430,7 +431,7 @@ namespace Step43
             PreconditionerMp > > m_inverse;
       const PreconditionerA &a_preconditioner;
 
-      mutable TrilinosWrappers::Vector tmp;
+      mutable TrilinosWrappers::MPI::Vector tmp;
     };
 
 
@@ -444,15 +445,18 @@ namespace Step43
       :
       darcy_matrix            (&S),
       m_inverse               (&Mpinv),
-      a_preconditioner        (Apreconditioner),
-      tmp                     (darcy_matrix->block(1,1).m())
-    {}
+      a_preconditioner        (Apreconditioner)
+    {
+      IndexSet tmp_index_set(darcy_matrix->block(1,1).m());
+      tmp_index_set.add_range(0,darcy_matrix->block(1,1).m()); 
+      tmp.reinit(tmp_index_set);
+    }
 
 
     template <class PreconditionerA, class PreconditionerMp>
     void BlockSchurPreconditioner<PreconditionerA, PreconditionerMp>::vmult (
-      TrilinosWrappers::BlockVector       &dst,
-      const TrilinosWrappers::BlockVector &src) const
+      TrilinosWrappers::MPI::BlockVector       &dst,
+      const TrilinosWrappers::MPI::BlockVector &src) const
     {
       a_preconditioner.vmult (dst.block(0), src.block(0));
       darcy_matrix->block(1,0).residual(tmp, dst.block(0), src.block(1));
@@ -546,33 +550,35 @@ namespace Step43
     FESystem<dim>                        darcy_fe;
     DoFHandler<dim>                      darcy_dof_handler;
     ConstraintMatrix                     darcy_constraints;
+    std::vector<IndexSet>                darcy_index_set;
 
     ConstraintMatrix                     darcy_preconditioner_constraints;
 
     TrilinosWrappers::BlockSparseMatrix  darcy_matrix;
     TrilinosWrappers::BlockSparseMatrix  darcy_preconditioner_matrix;
 
-    TrilinosWrappers::BlockVector        darcy_solution;
-    TrilinosWrappers::BlockVector        darcy_rhs;
+    TrilinosWrappers::MPI::BlockVector   darcy_solution;
+    TrilinosWrappers::MPI::BlockVector   darcy_rhs;
 
-    TrilinosWrappers::BlockVector        last_computed_darcy_solution;
-    TrilinosWrappers::BlockVector        second_last_computed_darcy_solution;
+    TrilinosWrappers::MPI::BlockVector   last_computed_darcy_solution;
+    TrilinosWrappers::MPI::BlockVector   second_last_computed_darcy_solution;
 
 
     const unsigned int                   saturation_degree;
     FE_Q<dim>                            saturation_fe;
     DoFHandler<dim>                      saturation_dof_handler;
     ConstraintMatrix                     saturation_constraints;
+    IndexSet                             saturation_index_set;
 
     TrilinosWrappers::SparseMatrix       saturation_matrix;
 
 
-    TrilinosWrappers::Vector             saturation_solution;
-    TrilinosWrappers::Vector             old_saturation_solution;
-    TrilinosWrappers::Vector             old_old_saturation_solution;
-    TrilinosWrappers::Vector             saturation_rhs;
+    TrilinosWrappers::MPI::Vector        saturation_solution;
+    TrilinosWrappers::MPI::Vector        old_saturation_solution;
+    TrilinosWrappers::MPI::Vector        old_old_saturation_solution;
+    TrilinosWrappers::MPI::Vector        saturation_rhs;
 
-    TrilinosWrappers::Vector             saturation_matching_last_computed_darcy_solution;
+    TrilinosWrappers::MPI::Vector        saturation_matching_last_computed_darcy_solution;
 
     const double                         saturation_refinement_threshold;
 
@@ -802,33 +808,34 @@ namespace Step43
       saturation_matrix.reinit (dsp);
     }
 
-    darcy_solution.reinit (2);
-    darcy_solution.block(0).reinit (n_u);
-    darcy_solution.block(1).reinit (n_p);
+    darcy_index_set.clear();
+    darcy_index_set.resize(2);
+    darcy_index_set[0].set_size(n_u);
+    darcy_index_set[1].set_size(n_p);
+    darcy_index_set[0].add_range(0,n_u);
+    darcy_index_set[1].add_range(0,n_p);
+    darcy_solution.reinit (darcy_index_set);
     darcy_solution.collect_sizes ();
 
-    last_computed_darcy_solution.reinit (2);
-    last_computed_darcy_solution.block(0).reinit (n_u);
-    last_computed_darcy_solution.block(1).reinit (n_p);
+    last_computed_darcy_solution.reinit (darcy_index_set);
     last_computed_darcy_solution.collect_sizes ();
 
-    second_last_computed_darcy_solution.reinit (2);
-    second_last_computed_darcy_solution.block(0).reinit (n_u);
-    second_last_computed_darcy_solution.block(1).reinit (n_p);
+    second_last_computed_darcy_solution.reinit (darcy_index_set);
     second_last_computed_darcy_solution.collect_sizes ();
 
-    darcy_rhs.reinit (2);
-    darcy_rhs.block(0).reinit (n_u);
-    darcy_rhs.block(1).reinit (n_p);
+    darcy_rhs.reinit (darcy_index_set);
     darcy_rhs.collect_sizes ();
 
-    saturation_solution.reinit (n_s);
-    old_saturation_solution.reinit (n_s);
-    old_old_saturation_solution.reinit (n_s);
+    saturation_index_set.clear();
+    saturation_index_set.set_size(n_s);
+    saturation_index_set.add_range(0,n_s);
+    saturation_solution.reinit (saturation_index_set);
+    old_saturation_solution.reinit (saturation_index_set);
+    old_old_saturation_solution.reinit (saturation_index_set);
 
-    saturation_matching_last_computed_darcy_solution.reinit (n_s);
+    saturation_matching_last_computed_darcy_solution.reinit (saturation_index_set);
 
-    saturation_rhs.reinit (n_s);
+    saturation_rhs.reinit (saturation_index_set);
   }
 
 
@@ -1536,9 +1543,9 @@ namespace Step43
           SolverControl solver_control (darcy_matrix.m(),
                                         1e-16*darcy_rhs.l2_norm());
 
-          SolverGMRES<TrilinosWrappers::BlockVector>
+          SolverGMRES<TrilinosWrappers::MPI::BlockVector>
           gmres (solver_control,
-                 SolverGMRES<TrilinosWrappers::BlockVector >::AdditionalData(100));
+                 SolverGMRES<TrilinosWrappers::MPI::BlockVector >::AdditionalData(100));
 
           for (unsigned int i=0; i<darcy_solution.size(); ++i)
             if (darcy_constraints.is_constrained(i))
@@ -1631,7 +1638,7 @@ namespace Step43
 
       SolverControl solver_control (saturation_matrix.m(),
                                     1e-16*saturation_rhs.l2_norm());
-      SolverCG<TrilinosWrappers::Vector> cg (solver_control);
+      SolverCG<TrilinosWrappers::MPI::Vector> cg (solver_control);
 
       TrilinosWrappers::PreconditionIC preconditioner;
       preconditioner.initialize (saturation_matrix);
@@ -1673,7 +1680,7 @@ namespace Step43
       FEValues<dim> fe_values (saturation_fe, quadrature_formula, update_gradients);
       std::vector<Tensor<1,dim> > grad_saturation (1);
 
-      TrilinosWrappers::Vector extrapolated_saturation_solution (saturation_solution);
+      TrilinosWrappers::MPI::Vector extrapolated_saturation_solution (saturation_solution);
       if (timestep_number != 0)
         extrapolated_saturation_solution.sadd ((1. + time_step/old_time_step),
                                                time_step/old_time_step, old_saturation_solution);
@@ -1713,18 +1720,18 @@ namespace Step43
     triangulation.prepare_coarsening_and_refinement ();
 
     {
-      std::vector<TrilinosWrappers::Vector> x_saturation (3);
+      std::vector<TrilinosWrappers::MPI::Vector> x_saturation (3);
       x_saturation[0] = saturation_solution;
       x_saturation[1] = old_saturation_solution;
       x_saturation[2] = saturation_matching_last_computed_darcy_solution;
 
-      std::vector<TrilinosWrappers::BlockVector> x_darcy (2);
+      std::vector<TrilinosWrappers::MPI::BlockVector> x_darcy (2);
       x_darcy[0] = last_computed_darcy_solution;
       x_darcy[1] = second_last_computed_darcy_solution;
 
-      SolutionTransfer<dim,TrilinosWrappers::Vector> saturation_soltrans(saturation_dof_handler);
+      SolutionTransfer<dim,TrilinosWrappers::MPI::Vector> saturation_soltrans(saturation_dof_handler);
 
-      SolutionTransfer<dim,TrilinosWrappers::BlockVector> darcy_soltrans(darcy_dof_handler);
+      SolutionTransfer<dim,TrilinosWrappers::MPI::BlockVector> darcy_soltrans(darcy_dof_handler);
 
 
       triangulation.prepare_coarsening_and_refinement();
@@ -1735,7 +1742,7 @@ namespace Step43
       triangulation.execute_coarsening_and_refinement ();
       setup_dofs ();
 
-      std::vector<TrilinosWrappers::Vector> tmp_saturation (3);
+      std::vector<TrilinosWrappers::MPI::Vector> tmp_saturation (3);
       tmp_saturation[0].reinit (saturation_solution);
       tmp_saturation[1].reinit (saturation_solution);
       tmp_saturation[2].reinit (saturation_solution);
@@ -1745,7 +1752,7 @@ namespace Step43
       old_saturation_solution = tmp_saturation[1];
       saturation_matching_last_computed_darcy_solution = tmp_saturation[2];
 
-      std::vector<TrilinosWrappers::BlockVector> tmp_darcy (2);
+      std::vector<TrilinosWrappers::MPI::BlockVector> tmp_darcy (2);
       tmp_darcy[0].reinit (darcy_solution);
       tmp_darcy[1].reinit (darcy_solution);
       darcy_soltrans.interpolate(x_darcy, tmp_darcy);
index 90c2d703ff6c231a556b218a727d31c34f456f40..f9eb32422f855ab553d31467e11854e7831c014d 100644 (file)
@@ -59,6 +59,8 @@ namespace TrilinosWrappers
    * block vector class do only work in case the program is run on only one
    * processor, since the Trilinos matrices are inherently parallel.
    *
+   * This class is deprecated use TrilinosWrappers::MPI::BlockVector instead.
+   *
    * @ingroup Vectors
    * @ingroup TrilinosWrappers @see
    * @ref GlossBlockLA "Block (linear algebra)"
@@ -92,14 +94,14 @@ namespace TrilinosWrappers
     /**
      * Default constructor. Generate an empty vector without any blocks.
      */
-    BlockVector ();
+    BlockVector () DEAL_II_DEPRECATED;
 
     /**
      * Constructor. Generate a block vector with as many blocks as there are
      * entries in Input_Maps.  For this non-distributed vector, the %parallel
      * partitioning is not used, just the global size of the partitioner.
      */
-    explicit BlockVector (const std::vector<Epetra_Map> &partitioner);
+    explicit BlockVector (const std::vector<Epetra_Map> &partitioner) DEAL_II_DEPRECATED;
 
     /**
      * Constructor. Generate a block vector with as many blocks as there are
@@ -107,26 +109,26 @@ namespace TrilinosWrappers
      * partitioning is not used, just the global size of the partitioner.
      */
     explicit BlockVector (const std::vector<IndexSet> &partitioner,
-                          const MPI_Comm              &communicator = MPI_COMM_WORLD);
+                          const MPI_Comm              &communicator = MPI_COMM_WORLD) DEAL_II_DEPRECATED;
 
     /**
      * Copy-Constructor. Set all the properties of the non-%parallel vector to
      * those of the given %parallel vector and import the elements.
      */
-    BlockVector (const MPI::BlockVector &V);
+    BlockVector (const MPI::BlockVector &V) DEAL_II_DEPRECATED;
 
     /**
      * Copy-Constructor. Set all the properties of the vector to those of the
      * given input vector and copy the elements.
      */
-    BlockVector (const BlockVector  &V);
+    BlockVector (const BlockVector  &V) DEAL_II_DEPRECATED;
 
     /**
      * Creates a block vector consisting of <tt>num_blocks</tt> components,
      * but there is no content in the individual components and the user has
      * to fill appropriate data using a reinit of the blocks.
      */
-    explicit BlockVector (const size_type num_blocks);
+    explicit BlockVector (const size_type num_blocks) DEAL_II_DEPRECATED;
 
     /**
      * Constructor. Set the number of blocks to <tt>n.size()</tt> and
@@ -134,7 +136,7 @@ namespace TrilinosWrappers
      *
      * References BlockVector.reinit().
      */
-    explicit BlockVector (const std::vector<size_type> &N);
+    explicit BlockVector (const std::vector<size_type> &N) DEAL_II_DEPRECATED;
 
     /**
      * Constructor. Set the number of blocks to <tt>n.size()</tt>. Initialize
@@ -147,7 +149,7 @@ namespace TrilinosWrappers
     template <typename InputIterator>
     BlockVector (const std::vector<size_type> &n,
                  const InputIterator           first,
-                 const InputIterator           end);
+                 const InputIterator           end) DEAL_II_DEPRECATED;
 
     /**
      * Destructor. Clears memory
index 542473c1891ba32ec06d1c2477595fa65adb5e00..241a334d84e045fa4a7eb8adf8518e3af61c651d 100644 (file)
@@ -710,6 +710,8 @@ namespace TrilinosWrappers
    * in order to be able to access all elements in the vector or to apply
    * certain deal.II functions.
    *
+   * This class is deprecated use TrilinosWrappers::MPI::Vector instead.
+   *
    * @ingroup TrilinosWrappers
    * @ingroup Vectors
    * @author Martin Kronbichler, 2008
@@ -739,12 +741,12 @@ namespace TrilinosWrappers
      * function <tt>reinit()</tt> will have to give the vector the correct
      * size.
      */
-    Vector ();
+    Vector () DEAL_II_DEPRECATED;
 
     /**
      * This constructor takes as input the number of elements in the vector.
      */
-    explicit Vector (const size_type n);
+    explicit Vector (const size_type n) DEAL_II_DEPRECATED;
 
     /**
      * This constructor takes as input the number of elements in the vector.
@@ -755,7 +757,7 @@ namespace TrilinosWrappers
      * ignored, the only thing that matters is the size of the index space
      * described by this argument.
      */
-    explicit Vector (const Epetra_Map &partitioning);
+    explicit Vector (const Epetra_Map &partitioning) DEAL_II_DEPRECATED;
 
     /**
      * This constructor takes as input the number of elements in the vector.
@@ -767,20 +769,20 @@ namespace TrilinosWrappers
      * size of the index space described by this argument.
      */
     explicit Vector (const IndexSet &partitioning,
-                     const MPI_Comm &communicator = MPI_COMM_WORLD);
+                     const MPI_Comm &communicator = MPI_COMM_WORLD) DEAL_II_DEPRECATED;
 
     /**
      * This constructor takes a (possibly parallel) Trilinos Vector and
      * generates a localized version of the whole content on each processor.
      */
-    explicit Vector (const VectorBase &V);
+    explicit Vector (const VectorBase &V) DEAL_II_DEPRECATED;
 
     /**
      * Copy-constructor from deal.II vectors. Sets the dimension to that of
      * the given vector, and copies all elements.
      */
     template <typename Number>
-    explicit Vector (const dealii::Vector<Number> &v);
+    explicit Vector (const dealii::Vector<Number> &v) DEAL_II_DEPRECATED;
 
     /**
      * Reinit function that resizes the vector to the size specified by

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

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.