From: kronbichler Date: Mon, 1 Sep 2008 12:46:00 +0000 (+0000) Subject: Added the function DoFTools::extract_constant_modes(). This simplifies the preconditi... X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=8cae5ec8c2687805e6bd38ba191108003541cb7e;p=dealii-svn.git Added the function DoFTools::extract_constant_modes(). This simplifies the preconditioner initialization in step-31. git-svn-id: https://svn.dealii.org/trunk@16701 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/dofs/dof_tools.h b/deal.II/deal.II/include/dofs/dof_tools.h index cd67ee3180..1a1982eccd 100644 --- a/deal.II/deal.II/include/dofs/dof_tools.h +++ b/deal.II/deal.II/include/dofs/dof_tools.h @@ -1008,6 +1008,70 @@ class DoFTools const unsigned int subdomain_id, std::vector &selected_dofs); + /** + * Extract a vector that + * represents the constant + * modes of the DoFHandler + * for the components + * chosen by component_select. + * The constant modes on a + * discretization are the + * null space of a Laplace + * operator on the selected + * components with Neumann + * boundary conditions + * applied. The null space is + * a necessary ingredient for + * obtaining a good AMG + * preconditioner when using + * the class + * TrilinosWrappers::PreconditionAMG. + * Since the ML AMG package + * only works on algebraic + * properties of the respective + * matrix, it has no chance + * to detect whether the matrix + * comes from a scalar or a + * vector valued problem. However, + * a near null space supplies + * exactly the needed information + * about these components. + * The null space will consist + * of as many vectors as there + * are true arguments in + * component_select, + * say d, + * each of which will be one + * in one component and + * zero in all others. We store + * the null space contiguously + * in the output vector + * constant_modes + * to enable the internal + * Trilinos structures the + * access to it. This means that + * the length of the vector + * will be d times + * the number of degrees of + * freedom in the selected + * components. Note that any + * matrix associated with this + * null space must have been + * constructed from the + * same argument as in + * component_select. + * + * The main reason for this + * program is the use of the + * null space with the + * AMG preconditioner. + */ + template + static void + extract_constant_modes (const DH &dof_handler, + const std::vector &component_select, + std::vector &constant_modes); + /** * For each active cell of a DoFHandler * or hp::DoFHandler, extract the active diff --git a/deal.II/deal.II/source/dofs/dof_tools.cc b/deal.II/deal.II/source/dofs/dof_tools.cc index ffb981acd2..19372180ac 100644 --- a/deal.II/deal.II/source/dofs/dof_tools.cc +++ b/deal.II/deal.II/source/dofs/dof_tools.cc @@ -3613,7 +3613,60 @@ DoFTools::extract_subdomain_dofs (const DH &dof_handler, - +template +void +DoFTools::extract_constant_modes (const DH &dof_handler, + const std::vector &component_select, + std::vector &constant_modes) +{ + const unsigned int n_components = dof_handler.get_fe().n_components(); + Assert (n_components == component_select.size(), + ExcDimensionMismatch(n_components, + component_select.size())); + + // First count the number of dofs + // in the current component. + unsigned int n_components_selected = 0; + std::vector component_list (n_components, 0); + for (unsigned int d=0; d dofs_per_block(2); + count_dofs_per_block(dof_handler, dofs_per_block, component_list); + + const unsigned int n_u = dofs_per_block[1]; + std::vector selection_dof_list (dof_handler.n_dofs(), false); + std::vector temporary_dof_list (dof_handler.n_dofs(), false); + extract_dofs (dof_handler, component_select, selection_dof_list); + + constant_modes.resize (n_components_selected * n_u, 0.); + + for (unsigned int component=0, component_used=0; + component < n_components; ++component, ++component_used) + if (component_select[component]) + { + std::vector selection_mask (n_components, false); + selection_mask[component] = true; + extract_dofs (dof_handler, selection_mask, temporary_dof_list); + + unsigned int counter = 0; + for (unsigned int i=0; i void DoFTools::get_active_fe_indices (const DH &dof_handler, @@ -5550,6 +5603,13 @@ DoFTools::extract_subdomain_dofs > const unsigned int subdomain_id, std::vector &selected_dofs); +template +void +DoFTools::extract_constant_modes > +(const DoFHandler &dof_handler, + const std::vector &selected_components, + std::vector &constant_modes); + template void DoFTools::get_active_fe_indices > diff --git a/deal.II/examples/step-31/step-31.cc b/deal.II/examples/step-31/step-31.cc index 49c44bc325..c2ac126e02 100644 --- a/deal.II/examples/step-31/step-31.cc +++ b/deal.II/examples/step-31/step-31.cc @@ -1188,42 +1188,12 @@ BoussinesqFlowProblem::build_stokes_preconditioner () Amg_preconditioner = boost::shared_ptr (new TrilinosWrappers::PreconditionAMG()); - - const unsigned int n_u = stokes_preconditioner_matrix.block(0,0).m(); - std::vector null_space (dim * n_u, 0.); - - std::vector precondition_dof_list (stokes_dof_handler.n_dofs(), false); - - for (unsigned int component=0; component < dim; ++component) - { - std::vector precondition_mask (dim + 1, false); - precondition_mask[component] = true; - DoFTools::extract_dofs (stokes_dof_handler, precondition_mask, - precondition_dof_list); - - // TODO: The current implementation - // assumes that we are working on - // the first components of a system when - // writing into the null vector. - // Change this to the general case, - // probably use something similar as - // for block vectors. - unsigned int counter = 0; - for (unsigned int i=0; i null_space; + std::vector velocity_components (dim+1,true); + velocity[dim] = false; + DoFTools::extract_constant_modes (stokes_dof_handler, velocity_components, + null_space); Amg_preconditioner->initialize(stokes_preconditioner_matrix.block(0,0), true, true, null_space, dim, false);