From 33c37e368ca147292644aae667cf16fefd41e44b Mon Sep 17 00:00:00 2001 From: kronbichler Date: Fri, 23 Jan 2009 17:12:30 +0000 Subject: [PATCH] A few small changes that make the program a bit faster. git-svn-id: https://svn.dealii.org/trunk@18281 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/examples/step-31/doc/results.dox | 6 +-- deal.II/examples/step-31/step-31.cc | 62 +++++++++++++----------- 2 files changed, 36 insertions(+), 32 deletions(-) diff --git a/deal.II/examples/step-31/doc/results.dox b/deal.II/examples/step-31/doc/results.dox index b54c469206..4a70bc0f8c 100644 --- a/deal.II/examples/step-31/doc/results.dox +++ b/deal.II/examples/step-31/doc/results.dox @@ -69,7 +69,7 @@ Timestep 0: t=0 Timestep 1: t=0.0610352 Assembling... Solving... - 30 GMRES iterations for Stokes subsystem. + 33 GMRES iterations for Stokes subsystem. Time step: 0.0610352 16 CG iterations for temperature. Temperature range: -0.0285735 0.153308 @@ -226,10 +226,10 @@ Timestep 0: t=0 Timestep 1: t=0.325521 Assembling... Solving... - 60 GMRES iterations for Stokes subsystem. + 55 GMRES iterations for Stokes subsystem. Time step: 0.325521 26 CG iterations for temperature. - Temperature range: -0.488042 1.02058 + Temperature range: -0.488053 1.02071 ... diff --git a/deal.II/examples/step-31/step-31.cc b/deal.II/examples/step-31/step-31.cc index 0f89ce6740..77d795f16f 100644 --- a/deal.II/examples/step-31/step-31.cc +++ b/deal.II/examples/step-31/step-31.cc @@ -347,7 +347,7 @@ namespace LinearSolvers vmult (VectorType &dst, const VectorType &src) const { - SolverControl solver_control (src.size(), 5e-8*src.l2_norm()); + SolverControl solver_control (src.size(), 1e-7*src.l2_norm()); SolverCG cg (solver_control); dst = 0; @@ -1112,7 +1112,7 @@ void BoussinesqFlowProblem::setup_dofs () temperature_constraints); temperature_constraints.close (); } - + std::vector stokes_dofs_per_block (2); DoFTools::count_dofs_per_block (stokes_dof_handler, stokes_dofs_per_block, stokes_sub_blocks); @@ -1481,33 +1481,36 @@ BoussinesqFlowProblem::build_stokes_preconditioner () // need to tell the AMG setup that we use // quadratic basis functions for the // velocity matrix (this implies more - // nonzero elements in the matrix, so that - // a more rubust algorithm needs to be - // chosen internally). Moreover, we want to - // be able to control how the coarsening - // structure is build up. The way AMG does - // this is to look which matrix entries are - // of similar size as the diagonal entry in - // order to algebraically build a - // coarse-grid structure. By setting the - // parameter + // nonzero elements in the matrix, so + // that a more rubust algorithm needs to + // be chosen internally). Moreover, we + // want to be able to control how the + // coarsening structure is build up. The + // way the Trilinos smoothed aggregation + // AMG does this is to look which matrix + // entries are of similar size as the + // diagonal entry in order to + // algebraically build a coarse-grid + // structure. By setting the parameter // aggregation_threshold to - // 0.05, we specify that all entries that - // are more than five precent of size of - // some diagonal pivots in that row should - // form one coarse grid point. This - // parameter is rather ad-hoc, and some - // fine-tuning of it can influence the - // performance of the preconditioner. As a - // rule of thumb, larger values of - // aggregation_threshold will - // decrease the number of iterations, but - // increase the costs per iteration. A look - // at the Trilinos documentation will - // provide more information on these - // parameters. With this data set, we then - // initialize the preconditioner with the - // matrix we want it to apply to. + // 0.02, we specify that all entries that + // are more than two precent of size of + // some diagonal pivots in that row + // should form one coarse grid + // point. This parameter is rather + // ad-hoc, and some fine-tuning of it can + // influence the performance of the + // preconditioner. As a rule of thumb, + // larger values of + // aggregation_threshold + // will decrease the number of + // iterations, but increase the costs per + // iteration. A look at the Trilinos + // documentation will provide more + // information on these parameters. With + // this data set, we then initialize the + // preconditioner with the matrix we want + // it to apply to. // // Finally, we also initialize the // preconditioner for the inversion of @@ -1529,7 +1532,8 @@ BoussinesqFlowProblem::build_stokes_preconditioner () // object. amg_data.elliptic = true; amg_data.higher_order_elements = true; - amg_data.aggregation_threshold = 5e-2; + amg_data.smoother_sweeps = 2; + amg_data.aggregation_threshold = 0.02; Amg_preconditioner->initialize(stokes_preconditioner_matrix.block(0,0), amg_data); -- 2.39.5