From: kronbichler Date: Wed, 3 Sep 2008 14:04:35 +0000 (+0000) Subject: Interpolate Stokes solution when the grid changes to reduce number of iterations... X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=41b84443a0c57ecc735926e1dd6764546ce790e2;p=dealii-svn.git Interpolate Stokes solution when the grid changes to reduce number of iterations even during the first step after remeshing (is it necessary?). Added a parameter aggregation_threshold to the Trilinos AMG preconditioner to allow control over the aggregation in the matrix. git-svn-id: https://svn.dealii.org/trunk@16732 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/examples/step-31/step-31.cc b/deal.II/examples/step-31/step-31.cc index 2936194255..a7f8f66ec8 100644 --- a/deal.II/examples/step-31/step-31.cc +++ b/deal.II/examples/step-31/step-31.cc @@ -1124,7 +1124,7 @@ BoussinesqFlowProblem::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::solve () SolverGMRES gmres(solver_control, SolverGMRES::AdditionalData(100)); - //stokes_solution = 0; gmres.solve(stokes_matrix, stokes_solution, stokes_rhs, preconditioner); std::cout << " " @@ -1972,16 +1971,20 @@ void BoussinesqFlowProblem::refine_mesh (const unsigned int max_grid_level) cell != triangulation.end(); ++cell) cell->clear_refine_flag (); - std::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 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 soltrans(temperature_dof_handler); + SolutionTransfer temperature_trans(temperature_dof_handler); + SolutionTransfer 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::refine_mesh (const unsigned int max_grid_level) std::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; diff --git a/deal.II/lac/include/lac/trilinos_block_vector.h b/deal.II/lac/include/lac/trilinos_block_vector.h index b944b98001..abca8be5e0 100644 --- a/deal.II/lac/include/lac/trilinos_block_vector.h +++ b/deal.II/lac/include/lac/trilinos_block_vector.h @@ -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; in_blocks(); ++i) + this->block(i) = v.block(i); + + collect_sizes(); + return *this; } @@ -338,6 +345,8 @@ namespace TrilinosWrappers for (unsigned int i=0;in_blocks();++i) block(i).reinit(v.block(i), fast); + + collect_sizes(); } diff --git a/deal.II/lac/include/lac/trilinos_precondition_amg.h b/deal.II/lac/include/lac/trilinos_precondition_amg.h index 93a5ce7be6..2c66315388 100755 --- a/deal.II/lac/include/lac/trilinos_precondition_amg.h +++ b/deal.II/lac/include/lac/trilinos_precondition_amg.h @@ -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 > &null_space = std::vector > (), const bool output_details = false); @@ -117,6 +118,7 @@ namespace TrilinosWrappers void initialize (const ::dealii::SparseMatrix &deal_ii_sparse_matrix, const bool elliptic = true, const bool higher_order_elements = false, + const double aggregation_threshold = 1e-4, const std::vector > &null_space = std::vector > (), const bool output_details = false); diff --git a/deal.II/lac/source/trilinos_precondition_amg.cc b/deal.II/lac/source/trilinos_precondition_amg.cc index 2ac0fa1db2..95c98345fd 100755 --- a/deal.II/lac/source/trilinos_precondition_amg.cc +++ b/deal.II/lac/source/trilinos_precondition_amg.cc @@ -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 > &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 &deal_ii_sparse_matrix, const bool elliptic, const bool higher_order_elements, + const double aggregation_threshold, const std::vector > &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); }