]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Make step-32 compliant with new TrilinosWrappers structure.
authorkronbichler <kronbichler@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 16 Sep 2008 14:06:41 +0000 (14:06 +0000)
committerkronbichler <kronbichler@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 16 Sep 2008 14:06:41 +0000 (14:06 +0000)
git-svn-id: https://svn.dealii.org/trunk@16835 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/examples/step-32/step-32.cc
deal.II/lac/include/lac/trilinos_block_sparse_matrix.h
deal.II/lac/include/lac/trilinos_block_vector.h
deal.II/lac/include/lac/trilinos_sparse_matrix.h
deal.II/lac/include/lac/trilinos_vector_base.h
deal.II/lac/source/trilinos_block_vector.cc
deal.II/lac/source/trilinos_sparse_matrix.cc
deal.II/lac/source/trilinos_vector_base.cc

index c0c13e2f4a5ce8d11037bdc578f3f703e08a0249..b72f043873e6ebcc00e44b8446979fe0b26a6bf3 100644 (file)
@@ -189,8 +189,8 @@ namespace LinearSolvers
                     const Preconditioner &preconditioner);
 
 
-      void vmult (TrilinosWrappers::Vector       &dst,
-                 const TrilinosWrappers::Vector &src) const;
+      void vmult (TrilinosWrappers::MPI::Vector       &dst,
+                 const TrilinosWrappers::MPI::Vector &src) const;
 
     private:
       const SmartPointer<const Matrix> matrix;
@@ -210,11 +210,11 @@ namespace LinearSolvers
 
   template <class Matrix, class Preconditioner>
   void InverseMatrix<Matrix,Preconditioner>::vmult (
-                               TrilinosWrappers::Vector       &dst,
-                               const TrilinosWrappers::Vector &src) const
+                               TrilinosWrappers::MPI::Vector       &dst,
+                               const TrilinosWrappers::MPI::Vector &src) const
   {
     SolverControl solver_control (src.size(), 1e-6*src.l2_norm());
-    SolverCG<TrilinosWrappers::Vector> cg (solver_control);
+    SolverCG<TrilinosWrappers::MPI::Vector> cg (solver_control);
 
     dst = 0;
 
@@ -228,6 +228,8 @@ namespace LinearSolvers
       }
   }
 
+
+
   template <class PreconditionerA, class PreconditionerMp>
   class BlockSchurPreconditioner : public Subscriptor
   {
@@ -237,14 +239,16 @@ namespace LinearSolvers
        const InverseMatrix<TrilinosWrappers::SparseMatrix,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;
       const SmartPointer<const InverseMatrix<TrilinosWrappers::SparseMatrix,
                                             PreconditionerMp > > m_inverse;
       const PreconditionerA &a_preconditioner;
+      mutable TrilinosWrappers::MPI::Vector tmp;
+
 };
 
 
@@ -257,16 +261,18 @@ namespace LinearSolvers
                  :
                  stokes_matrix           (&S),
                  m_inverse               (&Mpinv),
-                 a_preconditioner        (Apreconditioner)
+                 a_preconditioner        (Apreconditioner),
+                 tmp                     (stokes_matrix->block(1,1).matrix->RowMap())
   {}
 
+
+
   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));
-    TrilinosWrappers::Vector tmp (src.block(1), true);
     stokes_matrix->block(1,0).residual(tmp, dst.block(0), src.block(1));
     tmp *= -1;
     m_inverse->vmult (dst.block(1), tmp);
@@ -332,8 +338,8 @@ class BoussinesqFlowProblem
     TrilinosWrappers::BlockSparseMatrix stokes_matrix;
     TrilinosWrappers::BlockSparseMatrix stokes_preconditioner_matrix;
 
-    TrilinosWrappers::BlockVector       stokes_solution;
-    TrilinosWrappers::BlockVector       stokes_rhs;
+    TrilinosWrappers::MPI::BlockVector  stokes_solution;
+    TrilinosWrappers::MPI::BlockVector  stokes_rhs;
 
 
     const unsigned int                  temperature_degree;    
@@ -346,10 +352,10 @@ class BoussinesqFlowProblem
     TrilinosWrappers::SparseMatrix      temperature_stiffness_matrix;
     TrilinosWrappers::SparseMatrix      temperature_matrix;
 
-    TrilinosWrappers::Vector            temperature_solution;
+    TrilinosWrappers::MPI::Vector       temperature_solution;
     TrilinosWrappers::Vector            old_temperature_solution;
     TrilinosWrappers::Vector            old_old_temperature_solution;
-    TrilinosWrappers::Vector            temperature_rhs;
+    TrilinosWrappers::MPI::Vector       temperature_rhs;
 
 
     double time_step;
@@ -450,16 +456,13 @@ BoussinesqFlowProblem<dim>::get_extrapolated_temperature_range () const
                            update_values);
   std::vector<double> old_temperature_values(n_q_points);
   std::vector<double> old_old_temperature_values(n_q_points);
-  
-  Vector<double> localized_old_temperature_solution (old_temperature_solution);
-  Vector<double> old_localized_old_temperature_solution (old_old_temperature_solution);
 
   double min_temperature = (1. + time_step/old_time_step) *
                           old_temperature_solution.linfty_norm()
                           +
                           time_step/old_time_step *
                           old_old_temperature_solution.linfty_norm(),
-        max_temperature = -min_temperature;
+         max_temperature = -min_temperature;
 
   typename DoFHandler<dim>::active_cell_iterator
     cell = temperature_dof_handler.begin_active(),
@@ -468,8 +471,8 @@ BoussinesqFlowProblem<dim>::get_extrapolated_temperature_range () const
     if (cell->subdomain_id() == (unsigned int)trilinos_communicator.MyPID())
       {
        fe_values.reinit (cell);
-       fe_values.get_function_values (localized_old_temperature_solution, old_temperature_values);
-       fe_values.get_function_values (old_localized_old_temperature_solution, old_old_temperature_values);
+       fe_values.get_function_values (old_temperature_solution, old_temperature_values);
+       fe_values.get_function_values (old_old_temperature_solution, old_old_temperature_values);
 
        for (unsigned int q=0; q<n_q_points; ++q)
          {
@@ -620,31 +623,28 @@ void BoussinesqFlowProblem<dim>::setup_dofs ()
 
     BlockSparsityPattern stokes_sparsity_pattern (2,2);
 
-    //if (trilinos_communicator.MyPID() == 0)
-      {
-       BlockCompressedSetSparsityPattern csp (2,2);
+    BlockCompressedSetSparsityPattern csp (2,2);
 
-       csp.block(0,0).reinit (n_u, n_u);
-       csp.block(0,1).reinit (n_u, n_p);
-       csp.block(1,0).reinit (n_p, n_u);
-       csp.block(1,1).reinit (n_p, n_p);
+    csp.block(0,0).reinit (n_u, n_u);
+    csp.block(0,1).reinit (n_u, n_p);
+    csp.block(1,0).reinit (n_p, n_u);
+    csp.block(1,1).reinit (n_p, n_p);
 
-       csp.collect_sizes ();
+    csp.collect_sizes ();
 
-       Table<2,DoFTools::Coupling> coupling (dim+1, dim+1);
+    Table<2,DoFTools::Coupling> coupling (dim+1, dim+1);
 
-       for (unsigned int c=0; c<dim+1; ++c)
-         for (unsigned int d=0; d<dim+1; ++d)
-           if (! ((c==dim) && (d==dim)))
-             coupling[c][d] = DoFTools::always;
-           else
-             coupling[c][d] = DoFTools::none;
+    for (unsigned int c=0; c<dim+1; ++c)
+      for (unsigned int d=0; d<dim+1; ++d)
+       if (! ((c==dim) && (d==dim)))
+         coupling[c][d] = DoFTools::always;
+       else
+         coupling[c][d] = DoFTools::none;
 
-       DoFTools::make_sparsity_pattern (stokes_dof_handler, coupling, csp);
-       stokes_constraints.condense (csp);
+    DoFTools::make_sparsity_pattern (stokes_dof_handler, coupling, csp);
+    stokes_constraints.condense (csp);
 
-       stokes_sparsity_pattern.copy_from (csp);
-      }
+    stokes_sparsity_pattern.copy_from (csp);
 
     stokes_matrix.reinit (stokes_partitioner, stokes_sparsity_pattern);
   }
@@ -656,33 +656,30 @@ void BoussinesqFlowProblem<dim>::setup_dofs ()
 
     BlockSparsityPattern stokes_preconditioner_sparsity_pattern (2,2);
     
-    //if (trilinos_communicator.MyPID() == 0)
-      {
-       BlockCompressedSetSparsityPattern csp (2,2);
-
-       csp.block(0,0).reinit (n_u, n_u);
-       csp.block(0,1).reinit (n_u, n_p);
-       csp.block(1,0).reinit (n_p, n_u);
-       csp.block(1,1).reinit (n_p, n_p);
-         
-       csp.collect_sizes ();
-
-       Table<2,DoFTools::Coupling> coupling (dim+1, dim+1);
-       for (unsigned int c=0; c<dim+1; ++c)
-         for (unsigned int d=0; d<dim+1; ++d)
-           if (c == d)
-             coupling[c][d] = DoFTools::always;
-           else
-             coupling[c][d] = DoFTools::none;
+    BlockCompressedSetSparsityPattern csp (2,2);
 
-       DoFTools::make_sparsity_pattern (stokes_dof_handler, coupling, csp);
-       stokes_constraints.condense (csp);
-       
-       stokes_preconditioner_sparsity_pattern.copy_from (csp);
-      }
+    csp.block(0,0).reinit (n_u, n_u);
+    csp.block(0,1).reinit (n_u, n_p);
+    csp.block(1,0).reinit (n_p, n_u);
+    csp.block(1,1).reinit (n_p, n_p);
+  
+    csp.collect_sizes ();
+
+    Table<2,DoFTools::Coupling> coupling (dim+1, dim+1);
+    for (unsigned int c=0; c<dim+1; ++c)
+      for (unsigned int d=0; d<dim+1; ++d)
+       if (c == d)
+         coupling[c][d] = DoFTools::always;
+       else
+         coupling[c][d] = DoFTools::none;
+
+    DoFTools::make_sparsity_pattern (stokes_dof_handler, coupling, csp);
+    stokes_constraints.condense (csp);
+
+    stokes_preconditioner_sparsity_pattern.copy_from (csp);
 
     stokes_preconditioner_matrix.reinit (stokes_partitioner,
-                                        stokes_preconditioner_sparsity_pattern);
+                                        stokes_preconditioner_sparsity_pattern);
 
   }
 
@@ -694,21 +691,18 @@ void BoussinesqFlowProblem<dim>::setup_dofs ()
 
     SparsityPattern temperature_sparsity_pattern;
 
-    //if (trilinos_communicator.MyPID() == 0)
-      {
-       CompressedSetSparsityPattern csp (n_T, n_T);
-       DoFTools::make_sparsity_pattern (temperature_dof_handler, csp);
-       temperature_constraints.condense (csp);
+    CompressedSetSparsityPattern csp (n_T, n_T);
+    DoFTools::make_sparsity_pattern (temperature_dof_handler, csp);
+    temperature_constraints.condense (csp);
 
-       temperature_sparsity_pattern.copy_from (csp);
-      }
+    temperature_sparsity_pattern.copy_from (csp);
 
     temperature_matrix.reinit (temperature_partitioner,
-                               temperature_sparsity_pattern);
+                              temperature_sparsity_pattern);
     temperature_mass_matrix.reinit (temperature_partitioner,
                                    temperature_sparsity_pattern);
     temperature_stiffness_matrix.reinit (temperature_partitioner,
-                                         temperature_sparsity_pattern);
+                                        temperature_sparsity_pattern);
   }
 
   stokes_solution.reinit (stokes_partitioner);
@@ -717,7 +711,6 @@ void BoussinesqFlowProblem<dim>::setup_dofs ()
   temperature_solution.reinit (temperature_partitioner);
   old_temperature_solution.reinit (temperature_partitioner);
   old_old_temperature_solution.reinit (temperature_partitioner);
-
   temperature_rhs.reinit (temperature_partitioner);
 }
 
@@ -875,8 +868,6 @@ void BoussinesqFlowProblem<dim>::assemble_stokes_system ()
   const FEValuesExtractors::Vector velocities (0);
   const FEValuesExtractors::Scalar pressure (dim);
 
-  Vector<double> localized_old_temperature_solution (old_temperature_solution);
-
   typename DoFHandler<dim>::active_cell_iterator
     cell = stokes_dof_handler.begin_active(),
     endc = stokes_dof_handler.end();
@@ -892,7 +883,7 @@ void BoussinesqFlowProblem<dim>::assemble_stokes_system ()
        local_matrix = 0;
        local_rhs = 0;
   
-       temperature_fe_values.get_function_values (localized_old_temperature_solution, 
+       temperature_fe_values.get_function_values (old_temperature_solution, 
                                                   old_temperature_values);
   
        for (unsigned int q=0; q<n_q_points; ++q)
@@ -1116,9 +1107,7 @@ void BoussinesqFlowProblem<dim>::assemble_temperature_system ()
     global_T_range = get_extrapolated_temperature_range();
   const double global_Omega_diameter = GridTools::diameter (triangulation);
 
-  const Vector<double> localized_old_temperature_solution (old_temperature_solution);
-  const Vector<double> localized_old_old_temperature_solution (old_old_temperature_solution);
-  const BlockVector<double> localized_stokes_solution (stokes_solution);
+  const TrilinosWrappers::BlockVector localized_stokes_solution (stokes_solution);
 
   typename DoFHandler<dim>::active_cell_iterator
     cell = temperature_dof_handler.begin_active(),
@@ -1134,19 +1123,19 @@ void BoussinesqFlowProblem<dim>::assemble_temperature_system ()
        temperature_fe_values.reinit (cell);
        stokes_fe_values.reinit (stokes_cell);
   
-       temperature_fe_values.get_function_values (localized_old_temperature_solution,
+       temperature_fe_values.get_function_values (old_temperature_solution,
                                                   old_temperature_values);
-       temperature_fe_values.get_function_values (localized_old_old_temperature_solution,
+       temperature_fe_values.get_function_values (old_old_temperature_solution,
                                                   old_old_temperature_values);
   
-       temperature_fe_values.get_function_gradients (localized_old_temperature_solution,
+       temperature_fe_values.get_function_gradients (old_temperature_solution,
                                                      old_temperature_grads);
-       temperature_fe_values.get_function_gradients (localized_old_old_temperature_solution,
+       temperature_fe_values.get_function_gradients (old_old_temperature_solution,
                                                      old_old_temperature_grads);
        
-       temperature_fe_values.get_function_hessians (localized_old_temperature_solution,
+       temperature_fe_values.get_function_hessians (old_temperature_solution,
                                                     old_temperature_hessians);
-       temperature_fe_values.get_function_hessians (localized_old_old_temperature_solution,
+       temperature_fe_values.get_function_hessians (old_old_temperature_solution,
                                                     old_old_temperature_hessians);
        
        temperature_right_hand_side.value_list (temperature_fe_values.get_quadrature_points(),
@@ -1266,8 +1255,8 @@ void BoussinesqFlowProblem<dim>::solve ()
     SolverControl solver_control (stokes_matrix.m(),
                                  1e-6*stokes_rhs.l2_norm());
 
-    SolverGMRES<TrilinosWrappers::BlockVector> gmres(solver_control,
-      SolverGMRES<TrilinosWrappers::BlockVector >::AdditionalData(100));
+    SolverGMRES<TrilinosWrappers::MPI::BlockVector> gmres(solver_control,
+      SolverGMRES<TrilinosWrappers::MPI::BlockVector >::AdditionalData(100));
 
     gmres.solve(stokes_matrix, stokes_solution, stokes_rhs, preconditioner);
 
@@ -1276,7 +1265,7 @@ void BoussinesqFlowProblem<dim>::solve ()
          << " GMRES iterations for Stokes subsystem."
          << std::endl;
 
-    BlockVector<double> localized_stokes_solution (stokes_solution);
+    TrilinosWrappers::BlockVector localized_stokes_solution (stokes_solution);
     stokes_constraints.distribute (localized_stokes_solution);
     stokes_solution = localized_stokes_solution;
   }
@@ -1295,7 +1284,7 @@ void BoussinesqFlowProblem<dim>::solve ()
 
     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::PreconditionSSOR preconditioner (temperature_matrix,
                                                       1.2);
@@ -1303,7 +1292,7 @@ void BoussinesqFlowProblem<dim>::solve ()
              temperature_rhs,
              preconditioner);
 
-    Vector<double> localized_temperature_solution (temperature_solution);
+    TrilinosWrappers::Vector localized_temperature_solution (temperature_solution);
     temperature_constraints.distribute (localized_temperature_solution);
     temperature_solution = localized_temperature_solution;
 
@@ -1346,8 +1335,8 @@ void BoussinesqFlowProblem<dim>::output_results ()  const
          ExcInternalError());
   
   Vector<double> joint_solution (joint_dof_handler.n_dofs());
-  BlockVector<double> localized_stokes_solution (stokes_solution);
-  Vector<double> localized_temperature_solution (temperature_solution);
+  TrilinosWrappers::BlockVector localized_stokes_solution (stokes_solution);
+  TrilinosWrappers::Vector localized_temperature_solution (temperature_solution);
 
   {
     std::vector<unsigned int> local_joint_dof_indices (joint_fe.dofs_per_cell);
@@ -1396,25 +1385,28 @@ void BoussinesqFlowProblem<dim>::output_results ()  const
 
   DataOut<dim> data_out;
 
-  data_out.attach_dof_handler (joint_dof_handler);
+  if (trilinos_communicator.MyPID() == 0)
+    {
+      data_out.attach_dof_handler (joint_dof_handler);
 
-  std::vector<DataComponentInterpretation::DataComponentInterpretation>
-    data_component_interpretation
-    (dim+2, DataComponentInterpretation::component_is_scalar);
-  for (unsigned int i=0; i<dim; ++i)
-    data_component_interpretation[i]
-      = DataComponentInterpretation::component_is_part_of_vector;
+      std::vector<DataComponentInterpretation::DataComponentInterpretation>
+       data_component_interpretation
+       (dim+2, DataComponentInterpretation::component_is_scalar);
+      for (unsigned int i=0; i<dim; ++i)
+       data_component_interpretation[i]
+         = DataComponentInterpretation::component_is_part_of_vector;
 
-  data_out.add_data_vector (joint_solution, joint_solution_names,
-                           DataOut<dim>::type_dof_data,
-                           data_component_interpretation);
-  data_out.build_patches (std::min(stokes_degree, temperature_degree));
+      data_out.add_data_vector (joint_solution, joint_solution_names,
+                               DataOut<dim>::type_dof_data,
+                               data_component_interpretation);
+      data_out.build_patches (std::min(stokes_degree, temperature_degree));
 
-  std::ostringstream filename;
-  filename << "solution-" << Utilities::int_to_string(timestep_number, 4) << ".vtk";
+      std::ostringstream filename;
+      filename << "solution-" << Utilities::int_to_string(timestep_number, 4) << ".vtk";
 
-  std::ofstream output (filename.str().c_str());
-  data_out.write_vtk (output);
+      std::ofstream output (filename.str().c_str());
+      data_out.write_vtk (output);
+    }
 }
 
 
@@ -1425,7 +1417,7 @@ void BoussinesqFlowProblem<dim>::refine_mesh (const unsigned int max_grid_level)
 {
   Vector<float> estimated_error_per_cell (triangulation.n_active_cells());
 
-  Vector<double> localized_temperature_solution (temperature_solution);
+  TrilinosWrappers::Vector localized_temperature_solution (temperature_solution);
 
   KellyErrorEstimator<dim>::estimate (temperature_dof_handler,
                                      QGauss<dim-1>(temperature_degree+1),
@@ -1442,11 +1434,11 @@ void BoussinesqFlowProblem<dim>::refine_mesh (const unsigned int max_grid_level)
         cell != triangulation.end(); ++cell)
       cell->clear_refine_flag ();
   
-  std::vector<Vector<double> > x_solution (2);
+  std::vector<TrilinosWrappers::Vector > x_solution (2);
   x_solution[0] = temperature_solution;
   x_solution[1] = old_temperature_solution;
 
-  SolutionTransfer<dim,Vector<double> > soltrans(temperature_dof_handler);
+  SolutionTransfer<dim,TrilinosWrappers::Vector > soltrans(temperature_dof_handler);
 
   triangulation.prepare_coarsening_and_refinement();
   soltrans.prepare_for_coarsening_and_refinement(x_solution);
@@ -1454,7 +1446,7 @@ void BoussinesqFlowProblem<dim>::refine_mesh (const unsigned int max_grid_level)
   triangulation.execute_coarsening_and_refinement ();
   setup_dofs ();
 
-  std::vector<Vector<double> > tmp (2);
+  std::vector<TrilinosWrappers::Vector > tmp (2);
   tmp[0] = temperature_solution;
   tmp[1] = temperature_solution;
   soltrans.interpolate(x_solution, tmp);
index e50998fc2feab3dcbd0e0c4ae7c76349c38773b7..546378454d1a48d4a950b08d8161d75e4c8801ab 100644 (file)
@@ -211,6 +211,21 @@ namespace TrilinosWrappers
                                        */
       void compress ();
 
+                                      /**
+                                       * Returns the state of the
+                                       * matrix, i.e., whether
+                                       * compress() needs to be called
+                                       * after an operation requiring
+                                       * data exchange. Does only
+                                       * return non-true values when
+                                       * used in <tt>debug</tt> mode,
+                                       * since it is quite expensive to
+                                       * keep track of all operations
+                                       * that lead to the need for
+                                       * compress().
+                                       */
+      bool is_compressed () const;
+
                                        /**
                                         * This function collects the
                                         * sizes of the sub-objects and
@@ -221,10 +236,16 @@ namespace TrilinosWrappers
                                         * subobjects. You *must* call
                                         * this function each time after
                                         * you have changed the size of
-                                        * the sub-objects. Note that 
-                                       * this is a collective 
-                                       * operation, i.e., it needs to
-                                       * be called on all MPI processes.
+                                        * the sub-objects. Note that
+                                        * this is a collective
+                                        * operation, i.e., it needs to
+                                        * be called on all MPI
+                                        * processes. This command
+                                        * internally calls the method
+                                        * <tt>compress()</tt>, so you
+                                        * don't need to call that
+                                        * function in case you use
+                                        * <tt>collect_sizes()</tt>.
                                         */
       void collect_sizes ();
 
@@ -452,6 +473,24 @@ namespace TrilinosWrappers
 
 // ------------- inline and template functions -----------------
 
+  inline
+  bool
+  BlockSparseMatrix::is_compressed () const
+  {
+    bool compressed = true;
+    for (unsigned int row=0; row<n_block_rows(); ++row)
+      for (unsigned int col=0; col<n_block_cols(); ++col)
+       if (block(row, col).is_compressed() == false)
+         {
+           compressed = false;
+           break;
+         }
+
+    return compressed;
+  }
+
+
+
   inline
   void
   BlockSparseMatrix::vmult (MPI::BlockVector       &dst,
index e4f1d5a946fc45b5e2420d230430d41a0370b235..f861f868070ff8aa3cee6d56b118822e33fb895b 100644 (file)
@@ -147,6 +147,14 @@ namespace TrilinosWrappers
         BlockVector &
          operator = (const BlockVector &V);
 
+                                       /**
+                                       * Copy operator for arguments of
+                                       * the localized Trilinos vector
+                                       * type.
+                                       */
+        BlockVector &
+         operator = (const ::dealii::TrilinosWrappers::BlockVector &V);
+
                                        /**
                                        * Another copy function. This
                                        * one takes a deal.II block
@@ -235,6 +243,22 @@ namespace TrilinosWrappers
                                          */
        void compress ();
 
+                                        /**
+                                         * Returns the state of the
+                                         * matrix, i.e., whether
+                                         * compress() needs to be
+                                         * called after an operation
+                                         * requiring data
+                                         * exchange. Does only return
+                                         * non-true values when used in
+                                         * <tt>debug</tt> mode, since
+                                         * it is quite expensive to
+                                         * keep track of all operations
+                                         * that lead to the need for
+                                         * compress().
+                                         */
+       bool is_compressed () const;
+
                                          /**
                                           * Swap the contents of this
                                           * vector and the other vector
@@ -319,6 +343,22 @@ namespace TrilinosWrappers
 
 
 
+    inline
+    bool
+    BlockVector::is_compressed () const
+    {
+      bool compressed = true;
+      for (unsigned int row=0; row<n_blocks(); ++row)
+       if (block(row).is_compressed() == false)
+         {
+           compressed = false;
+           break;
+         }
+
+      return compressed;
+    }
+
+
 
 /**
  * Global function which overloads the default implementation
index 4df28a109d0f39ce7d4371593b67d0b791dbc8d7..cfe80656069a5353013b1a9f6228ad7ab23958db 100755 (executable)
@@ -638,6 +638,21 @@ namespace TrilinosWrappers
                                         */
       void compress ();
 
+                                      /**
+                                       * Returns the state of the
+                                       * matrix, i.e., whether
+                                       * compress() needs to be called
+                                       * after an operation requiring
+                                       * data exchange. Does only
+                                       * return non-true values when
+                                       * used in <tt>debug</tt> mode,
+                                       * since it is quite expensive to
+                                       * keep track of all operations
+                                       * that lead to the need for
+                                       * compress().
+                                       */
+      bool is_compressed () const;
+
                                        /**
                                         * This operator assigns a scalar
                                         * to a matrix. Since this does
@@ -1384,6 +1399,13 @@ namespace TrilinosWrappers
                                         */
       Epetra_CombineMode last_action;
 
+                                      /**
+                                       * A boolean variable to hold
+                                       * information on whether the
+                                       * vector is compressed or not.
+                                       */
+      bool compressed;
+
     public:
                                        /**
                                         * A sparse matrix object in
@@ -1633,6 +1655,15 @@ namespace TrilinosWrappers
             (index < static_cast<unsigned int>(end)));
   }
 
+
+
+  inline
+  bool
+  SparseMatrix::is_compressed () const
+  {
+    return compressed;
+  }
+
 #endif // DOXYGEN      
 }
 
index fa8a0fdb0c13538873c89cca09ada454fc893ef1..a24525734d7bbcbb25fbee795da211a4c769d066 100644 (file)
@@ -1,5 +1,5 @@
 //---------------------------------------------------------------------------
-//    $Id: trilinos_vector.h 16818 2008-09-12 14:30:24Z kronbichler $
+//    $Id: trilinos_vector_base.h 16818 2008-09-12 14:30:24Z kronbichler $
 //    Version: $Name$
 //
 //    Copyright (C) 2008 by the deal.II authors
@@ -313,9 +313,12 @@ namespace TrilinosWrappers
                                        * compress() has already been
                                        * called after an operation
                                        * requiring data exchange. Does
-                                       * only return non-false values
+                                       * only return non-true values
                                        * when used in <tt>debug</tt>
-                                       * mode.
+                                       * mode, since it is quite
+                                       * expensive to keep track of all
+                                       * operations that lead to a
+                                       * need for compress().
                                        */
       bool is_compressed () const;
 
index 864f0b40dee1b724d141687c38ae1cba4531497f..a401e49e3d67cf5c283414cbc514b2225f293c75 100644 (file)
@@ -62,6 +62,20 @@ namespace TrilinosWrappers
 
 
 
+    BlockVector &
+    BlockVector::operator = (const ::dealii::TrilinosWrappers::BlockVector &v)
+    {
+      Assert (n_blocks() == v.n_blocks(),
+             ExcDimensionMismatch(n_blocks(),v.n_blocks()));
+
+      for (unsigned int i=0; i<this->n_blocks(); ++i)
+       this->components[i] = v.block(i);
+
+      return *this;
+    }
+
+
+
     BlockVector::~BlockVector ()
     {}
 
index 5881615a484b8caf36ad6148a59d7c531c02dab3..38988ac32b3bc01d6d6f83f98956e0f1b5fa2295 100755 (executable)
@@ -11,7 +11,6 @@
 //
 //---------------------------------------------------------------------------
 
-#include <lac/trilinos_vector_base.h>
 #include <lac/trilinos_sparse_matrix.h>
 
 #ifdef DEAL_II_USE_TRILINOS
@@ -86,6 +85,7 @@ namespace TrilinosWrappers
 #endif
                  col_map (row_map),
                  last_action (Zero),
+                 compressed (true),
                  matrix (std::auto_ptr<Epetra_FECrsMatrix>
                                (new Epetra_FECrsMatrix(Copy, row_map, 0)))
   {}
@@ -96,6 +96,7 @@ namespace TrilinosWrappers
                   row_map (InputMap),
                  col_map (row_map),
                  last_action (Zero),
+                 compressed (true),
                  matrix (std::auto_ptr<Epetra_FECrsMatrix>
                                (new Epetra_FECrsMatrix(Copy, row_map, 
                                        int(n_max_entries_per_row), false)))
@@ -107,6 +108,7 @@ namespace TrilinosWrappers
                   row_map (InputMap),
                  col_map (row_map),
                  last_action (Zero),
+                 compressed (true),
                  matrix (std::auto_ptr<Epetra_FECrsMatrix>
                    (new Epetra_FECrsMatrix(Copy, row_map, 
                      (int*)const_cast<unsigned int*>(&(n_entries_per_row[0])),
@@ -120,6 +122,7 @@ namespace TrilinosWrappers
                   row_map (InputRowMap),
                   col_map (InputColMap),
                  last_action (Zero),
+                 compressed (true),
                  matrix (std::auto_ptr<Epetra_FECrsMatrix>
                                (new Epetra_FECrsMatrix(Copy, row_map, col_map, 
                                        int(n_max_entries_per_row), false)))
@@ -132,6 +135,7 @@ namespace TrilinosWrappers
                   row_map (InputRowMap),
                   col_map (InputColMap),
                  last_action (Zero),
+                 compressed (true),
                  matrix (std::auto_ptr<Epetra_FECrsMatrix>
                    (new Epetra_FECrsMatrix(Copy, row_map, col_map, 
                      (int*)const_cast<unsigned int*>(&(n_entries_per_row[0])),
@@ -270,10 +274,12 @@ namespace TrilinosWrappers
   SparseMatrix::reinit (const SparseMatrix &sparse_matrix)
   {
     matrix.reset();
+
     row_map = sparse_matrix.row_map;
     col_map = sparse_matrix.col_map;
     matrix = std::auto_ptr<Epetra_FECrsMatrix>(new Epetra_FECrsMatrix(
                                                   *sparse_matrix.matrix));
+    compressed = true;
   }
 
 
@@ -354,6 +360,7 @@ namespace TrilinosWrappers
                                   &values[0], &row_indices[0]);
       }
 
+    compress();
   }
 
 
@@ -376,6 +383,7 @@ namespace TrilinosWrappers
 
     matrix = std::auto_ptr<Epetra_FECrsMatrix> 
              (new Epetra_FECrsMatrix(Copy, row_map, 0));
+    compressed = true;
   }
 
 
@@ -394,6 +402,8 @@ namespace TrilinosWrappers
 
     ierr = matrix->OptimizeStorage ();
     AssertThrow (ierr == 0, ExcTrilinosError(ierr));
+
+    compressed = true;
   }
 
 
@@ -437,6 +447,11 @@ namespace TrilinosWrappers
     if (last_action != Insert)
       last_action = Insert;
 
+#ifdef DEBUG
+    if (in_local_range (i) == false)
+      compressed = false;
+#endif
+      
     int trilinos_i = i;
     int trilinos_j = j;
 
@@ -484,6 +499,11 @@ namespace TrilinosWrappers
     if (value == 0)
       return;
 
+#ifdef DEBUG
+    if (in_local_range (i) == false)
+      compressed = false;
+#endif
+      
     int trilinos_i = i;
     int trilinos_j = j;
 
index 4d841242dd58cc9e84b08755afe0c83a8ecff33c..dda97c4bca6e287d1f73cf590c80d276a8113d58 100644 (file)
@@ -1,5 +1,5 @@
 //---------------------------------------------------------------------------
-//    $Id: trilinos_vector.cc 16818 2008-09-12 14:30:24Z kronbichler $
+//    $Id: trilinos_vector_base.cc 16818 2008-09-12 14:30:24Z kronbichler $
 //    Version: $Name$
 //
 //    Copyright (C) 2008 by the deal.II authors

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.