From: Baerbel Jannsen Date: Wed, 25 Aug 2010 16:54:29 +0000 (+0000) Subject: changed the interface of MGConstraints. Made the member variables X-Git-Tag: v8.0.0~5686 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=61c76669246b25f627c004b320a4d37f7e4cd522;p=dealii.git changed the interface of MGConstraints. Made the member variables private and wrote functions to acess them. Changes made in step-16 accordingly. git-svn-id: https://svn.dealii.org/trunk@21712 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/multigrid/mg_constraints.h b/deal.II/deal.II/include/multigrid/mg_constraints.h index 4eaef26ea2..4796aeeea0 100644 --- a/deal.II/deal.II/include/multigrid/mg_constraints.h +++ b/deal.II/deal.II/include/multigrid/mg_constraints.h @@ -106,7 +106,7 @@ class MGConstraints : public Subscriptor * this level and cells on the level * below). */ - const std::vector > & + const std::vector > & get_refinement_edge_indices () const; /** @@ -116,7 +116,7 @@ class MGConstraints : public Subscriptor * get_boundary_indices() and * get_refinement_edge_indices(). */ - const std::vector > & + const std::vector > & get_refinement_edge_boundary_indices () const; private: @@ -234,6 +234,26 @@ MGConstraints::at_refinement_edge_boundary (const unsigned int level, return refinement_edge_boundary_indices[level][index]; } +inline +const std::vector > & +MGConstraints::get_boundary_indices () const +{ + return boundary_indices; +} + +inline +const std::vector > & +MGConstraints::get_refinement_edge_indices () const +{ + return refinement_edge_indices; +} + +inline +const std::vector > & +MGConstraints::get_refinement_edge_boundary_indices () const +{ + return refinement_edge_boundary_indices; +} DEAL_II_NAMESPACE_CLOSE diff --git a/deal.II/examples/step-16/step-16.cc b/deal.II/examples/step-16/step-16.cc index d7f0072846..78d5f2fdbe 100644 --- a/deal.II/examples/step-16/step-16.cc +++ b/deal.II/examples/step-16/step-16.cc @@ -75,6 +75,7 @@ // operator (solver or preconditioner). #include #include +#include #include #include #include @@ -122,6 +123,15 @@ class LaplaceProblem SparsityPattern sparsity_pattern; SparseMatrix system_matrix; + // We need an additional object for the + // hanging nodes constraints. They are + // handed to the transfer object in the + // multigrid. Since we call a compress + // inside the multigrid these constraints + // are not allowed to be inhomogeneous so + // we store them in different ConstraintMatrix + // objects. + ConstraintMatrix hanging_node_constraints; ConstraintMatrix constraints; Vector solution; @@ -129,13 +139,19 @@ class LaplaceProblem const unsigned int degree; - // The following three objects are the + // The following four objects are the // only additional member variables, - // compared to step-6. They represent the + // compared to step-6. They first three + // represent the // operators that act on individual // levels of the multilevel hierarchy, // rather than on the finest mesh as do - // the objects above. + // the objects above while the last object + // stores information about the boundary + // indices on each level and information + // about indices lying on a refinement + // edge between two different refinement + // levels. // // To facilitate having objects on each // level of a multilevel hierarchy, @@ -163,6 +179,7 @@ class LaplaceProblem MGLevelObject mg_sparsity_patterns; MGLevelObject > mg_matrices; MGLevelObject > mg_interface_matrices; + MGConstraints mg_constraints; }; @@ -284,9 +301,7 @@ void LaplaceProblem::setup_system () sparsity_pattern.reinit (mg_dof_handler.n_dofs(), mg_dof_handler.n_dofs(), mg_dof_handler.max_couplings_between_dofs()); - DoFTools::make_sparsity_pattern ( - static_cast&>(mg_dof_handler), - sparsity_pattern); + DoFTools::make_sparsity_pattern (mg_dof_handler, sparsity_pattern); solution.reinit (mg_dof_handler.n_dofs()); system_rhs.reinit (mg_dof_handler.n_dofs()); @@ -305,25 +320,39 @@ void LaplaceProblem::setup_system () // matrix. To this end, we here do not just // compute the constraints do to hanging // nodes, but also due to zero boundary - // conditions. Both kinds of constraints - // can be put into the same object - // (constraints), and we will + // conditions. We will // use this set of constraints later on to // help us copy local contributions // correctly into the global linear system // right away, without the need for a later // clean-up stage: constraints.clear (); + hanging_node_constraints.clear (); + DoFTools::make_hanging_node_constraints (mg_dof_handler, hanging_node_constraints); DoFTools::make_hanging_node_constraints (mg_dof_handler, constraints); - VectorTools::interpolate_boundary_values (mg_dof_handler, - 0, - ZeroFunction(), + + typename FunctionMap::type dirichlet_boundary; + ZeroFunction homogeneous_dirichlet_bc (1); + dirichlet_boundary[0] = &homogeneous_dirichlet_bc; + MappingQ1 mapping; + VectorTools::interpolate_boundary_values (mapping, mg_dof_handler, + dirichlet_boundary, constraints); constraints.close (); + hanging_node_constraints.close (); constraints.condense (sparsity_pattern); sparsity_pattern.compress(); system_matrix.reinit (sparsity_pattern); + // The multigrid constraints have to be + // initialized. They need to know about + // the boundary values as well, so we + // pass the dirichlet_boundary + // here as well. + mg_constraints.clear(); + mg_constraints.initialize(mg_dof_handler, dirichlet_boundary); + + // Now for the things that concern the // multigrid data structures. First, we // resize the multi-level objects to hold @@ -522,43 +551,23 @@ void LaplaceProblem::assemble_multigrid () // i.e. vectors of booleans each element of // which indicates whether the // corresponding degree of freedom index is - // an interface DoF or not: - std::vector > interface_dofs; - std::vector > boundary_interface_dofs; - for (unsigned int level = 0; level - (mg_dof_handler.n_dofs(level))); - boundary_interface_dofs.push_back (std::vector - (mg_dof_handler.n_dofs(level))); - } - MGTools::extract_inner_interface_dofs (mg_dof_handler, - interface_dofs, - boundary_interface_dofs); + // an interface DoF or not. The MGConstraints + // already computed the information for us + // when we called initialize in setup_system(). + std::vector > interface_dofs + = mg_constraints.get_refinement_edge_indices (); + std::vector > boundary_interface_dofs + = mg_constraints.get_refinement_edge_boundary_indices (); // The indices just identified will later - // be used to impose zero boundary - // conditions for the operator that we will - // apply on each level. On the other hand, + // be used to decide where the assembled value + // has to be added into on each level. + // On the other hand, // we also have to impose zero boundary // conditions on the external boundary of - // each level. So let's identify these - // nodes as well (this time as a set of - // degrees of freedom, rather than a - // boolean mask; the reason for this being - // that we will not need fast tests whether - // a certain degree of freedom is in the - // boundary list, though we will need such - // access for the interface degrees of - // freedom further down below): - typename FunctionMap::type dirichlet_boundary; - ZeroFunction homogeneous_dirichlet_bc (1); - dirichlet_boundary[0] = &homogeneous_dirichlet_bc; - - std::vector boundary_indices (triangulation.n_levels()); - MGTools::make_boundary_list (mg_dof_handler, dirichlet_boundary, - boundary_indices); - + // each level. But this the MGConstraints + // knows it. So we simply ask for them by calling + // get_boundary_indices (). // The third step is to construct // constraints on all those degrees of // freedom: their value should be zero @@ -580,7 +589,7 @@ void LaplaceProblem::assemble_multigrid () for (unsigned int level=0; level::assemble_multigrid () // the loop in // assemble_system, with two // exceptions: (i) we don't need a right - // han side, and more significantly (ii) we + // hand side, and more significantly (ii) we // don't just loop over all active cells, // but in fact all cells, active or // not. Consequently, the correct iterator @@ -751,8 +760,15 @@ void LaplaceProblem::assemble_multigrid () template void LaplaceProblem::solve () { - MGTransferPrebuilt > mg_transfer(constraints); - mg_transfer.build_matrices(mg_dof_handler); + + // Create the object that deals with the transfer + // between different refinement levels. We need to + // pass it the hanging node constraints. + MGTransferPrebuilt > mg_transfer(hanging_node_constraints); + // Now the prolongation matrix has to be built. + // This matrix needs to take the boundary values on + // each level into account. + mg_transfer.build_matrices(mg_dof_handler, mg_constraints.get_boundary_indices()); FullMatrix coarse_matrix; coarse_matrix.copy_from (mg_matrices[0]);