]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Interpolate Stokes solution when the grid changes to reduce number of iterations...
authorkronbichler <kronbichler@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 3 Sep 2008 14:04:35 +0000 (14:04 +0000)
committerkronbichler <kronbichler@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 3 Sep 2008 14:04:35 +0000 (14:04 +0000)
git-svn-id: https://svn.dealii.org/trunk@16732 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/examples/step-31/step-31.cc
deal.II/lac/include/lac/trilinos_block_vector.h
deal.II/lac/include/lac/trilinos_precondition_amg.h
deal.II/lac/source/trilinos_precondition_amg.cc

index 29361942555e9d9ff8f686ffde30fa30e2fe0247..a7f8f66ec8370b0122313fb5b10bf6f3d18c3361 100644 (file)
@@ -1124,7 +1124,7 @@ BoussinesqFlowProblem<dim>::build_stokes_preconditioner ()
   DoFTools::extract_constant_modes (stokes_dof_handler, velocity_components, 
                                    null_space);
   Amg_preconditioner->initialize(stokes_preconditioner_matrix.block(0,0),
-                                true, true, null_space, false);
+                                true, true, 5e-2, null_space, false);
 
                                   // TODO: we could throw away the (0,0)
                                   // block here since things have been
@@ -1798,7 +1798,6 @@ void BoussinesqFlowProblem<dim>::solve ()
     SolverGMRES<TrilinosWrappers::BlockVector> gmres(solver_control,
       SolverGMRES<TrilinosWrappers::BlockVector >::AdditionalData(100));
 
-    //stokes_solution = 0;
     gmres.solve(stokes_matrix, stokes_solution, stokes_rhs, preconditioner);
 
     std::cout << "   "
@@ -1972,16 +1971,20 @@ void BoussinesqFlowProblem<dim>::refine_mesh (const unsigned int max_grid_level)
         cell != triangulation.end(); ++cell)
       cell->clear_refine_flag ();
   
-  std::vector<TrilinosWrappers::Vector> x_solution (2);
-  x_solution[0].reinit (temperature_solution);
-  x_solution[0] = temperature_solution;
-  x_solution[1].reinit (temperature_solution);
-  x_solution[1] = old_temperature_solution;
+  std::vector<TrilinosWrappers::Vector> x_temperature (2);
+  x_temperature[0].reinit (temperature_solution);
+  x_temperature[0] = temperature_solution;
+  x_temperature[1].reinit (temperature_solution);
+  x_temperature[1] = old_temperature_solution;
+  TrilinosWrappers::BlockVector x_stokes(2);
+  x_stokes = stokes_solution;
 
-  SolutionTransfer<dim,TrilinosWrappers::Vector> soltrans(temperature_dof_handler);
+  SolutionTransfer<dim,TrilinosWrappers::Vector> temperature_trans(temperature_dof_handler);
+  SolutionTransfer<dim,TrilinosWrappers::BlockVector> stokes_trans(stokes_dof_handler);
 
   triangulation.prepare_coarsening_and_refinement();
-  soltrans.prepare_for_coarsening_and_refinement(x_solution);
+  temperature_trans.prepare_for_coarsening_and_refinement(x_temperature);
+  stokes_trans.prepare_for_coarsening_and_refinement(x_stokes);
 
   triangulation.execute_coarsening_and_refinement ();
   setup_dofs ();
@@ -1989,10 +1992,12 @@ void BoussinesqFlowProblem<dim>::refine_mesh (const unsigned int max_grid_level)
   std::vector<TrilinosWrappers::Vector> tmp (2);
   tmp[0].reinit (temperature_solution);
   tmp[1].reinit (temperature_solution);
-  soltrans.interpolate(x_solution, tmp);
+  temperature_trans.interpolate(x_temperature, tmp);
 
   temperature_solution = tmp[0];
   old_temperature_solution = tmp[1];
+  
+  stokes_trans.interpolate (x_stokes, stokes_solution);
 
   rebuild_stokes_matrix         = true;
   rebuild_temperature_matrices  = true;
index b944b9800173acea15016fb99f8d6af8c9dc6e71..abca8be5e054086f50d1c8d7a5eb49bc3826a1d3 100644 (file)
@@ -291,7 +291,14 @@ namespace TrilinosWrappers
     BlockVector &
     BlockVector::operator = (const BlockVector &v)
     {
-      BaseClass::operator = (v);
+      if (this->n_blocks() != v.n_blocks())
+       reinit(v.n_blocks());
+
+      for (unsigned int i=0; i<this->n_blocks(); ++i)
+       this->block(i) = v.block(i);
+
+      collect_sizes();
+       
       return *this;
     }
 
@@ -338,6 +345,8 @@ namespace TrilinosWrappers
   
       for (unsigned int i=0;i<this->n_blocks();++i)
         block(i).reinit(v.block(i), fast);
+      
+      collect_sizes();
     }
 
 
index 93a5ce7be6b0e52fedd8d0c05f3c33607a6d5c1d..2c66315388608475628e479c9b05687b69f86a7f 100755 (executable)
@@ -100,6 +100,7 @@ namespace TrilinosWrappers
       void initialize (const SparseMatrix                    &matrix,
                       const bool                             elliptic = true,
                       const bool                             higher_order_elements = false,
+                      const double                           aggregation_threshold = 1e-4,
                       const std::vector<std::vector<bool> > &null_space = std::vector<std::vector<bool> > (),
                       const bool                            output_details = false);
 
@@ -117,6 +118,7 @@ namespace TrilinosWrappers
       void initialize (const ::dealii::SparseMatrix<double>  &deal_ii_sparse_matrix,
                       const bool                             elliptic = true,
                       const bool                             higher_order_elements = false,
+                      const double                           aggregation_threshold = 1e-4,
                       const std::vector<std::vector<bool> > &null_space = std::vector<std::vector<bool> > (),
                       const bool                            output_details = false);
 
index 2ac0fa1db29e548aca058d910cec313363b5e60c..95c98345fd9d8b6abd82adf42c457f0556ff308b 100755 (executable)
@@ -46,6 +46,7 @@ namespace TrilinosWrappers
   initialize (const SparseMatrix                    &matrix,
              const bool                             elliptic,
              const bool                             higher_order_elements,
+             const double                           aggregation_threshold,
              const std::vector<std::vector<bool> > &null_space,
              const bool                             output_details)
   {
@@ -68,7 +69,7 @@ namespace TrilinosWrappers
        parameter_list.set("aggregation: block scaling", true);
       }
   
-    parameter_list.set("aggregation: threshold", 1e-8);
+    parameter_list.set("aggregation: threshold", aggregation_threshold);
     
     if (output_details)
       parameter_list.set("ML output", 10);
@@ -114,6 +115,7 @@ namespace TrilinosWrappers
   initialize (const ::dealii::SparseMatrix<double>  &deal_ii_sparse_matrix,
              const bool                             elliptic,
              const bool                             higher_order_elements,
+             const double                           aggregation_threshold,
              const std::vector<std::vector<bool> > &null_space,
              const bool                             output_details)
   {
@@ -130,8 +132,8 @@ namespace TrilinosWrappers
     Matrix->reinit (*Map, deal_ii_sparse_matrix);
     Matrix->compress();
 
-    initialize (*Matrix, elliptic, higher_order_elements, null_space,
-               output_details);
+    initialize (*Matrix, elliptic, higher_order_elements, 
+               aggregation_threshold, null_space, output_details);
   }
   
   

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.