From 3a7ab47ca58a63fb03452be9b924d31e108d6874 Mon Sep 17 00:00:00 2001 From: Baerbel Jannsen Date: Fri, 27 Aug 2010 20:21:04 +0000 Subject: [PATCH] changed to MGConstrainedDoFs git-svn-id: https://svn.dealii.org/trunk@21758 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/examples/step-42/step-42.cc | 50 +++++++++++++---------------- 1 file changed, 23 insertions(+), 27 deletions(-) diff --git a/deal.II/examples/step-42/step-42.cc b/deal.II/examples/step-42/step-42.cc index d5bcacb38b..abf30b6f53 100644 --- a/deal.II/examples/step-42/step-42.cc +++ b/deal.II/examples/step-42/step-42.cc @@ -60,6 +60,7 @@ #include #include #include +#include #include #include #include @@ -151,6 +152,7 @@ class StokesProblem MGLevelObject > mg_matrices; MGLevelObject > mg_interface_matrices; + MGConstrainedDoFs mg_constrained_dofs; std::vector > mg_dofs_per_component; std::vector::type> > mg_A_preconditioner; @@ -439,15 +441,25 @@ void StokesProblem::setup_dofs () { constraints.clear (); + typename FunctionMap::type dirichlet_boundary; + ZeroFunction homogeneous_dirichlet_bc (dim+1); //TODO: go back to BoundaryValues + + dirichlet_boundary[0] = &homogeneous_dirichlet_bc; + MappingQ1 mapping; + std::vector component_mask (dim+1, true); component_mask[dim] = false; - VectorTools::interpolate_boundary_values (dof_handler, - 0, - ZeroFunction(dim+1), //TODO: go back to BoundaryValues - constraints, - component_mask); + VectorTools::interpolate_boundary_values (mapping, + dof_handler, + dirichlet_boundary, + constraints, + component_mask); + DoFTools::make_hanging_node_constraints (dof_handler, constraints); + + mg_constrained_dofs.clear(); + mg_constrained_dofs.initialize(dof_handler, dirichlet_boundary); } constraints.close (); @@ -647,33 +659,17 @@ void StokesProblem::assemble_multigrid () std::vector div_phi_u (dofs_per_cell); std::vector phi_p (dofs_per_cell); - std::vector > interface_dofs; - std::vector > boundary_interface_dofs; - for (unsigned int level = 0; level tmp (dof_handler.n_dofs(level)); - interface_dofs.push_back (tmp); - boundary_interface_dofs.push_back (tmp); - } - MGTools::extract_inner_interface_dofs (dof_handler, - interface_dofs, boundary_interface_dofs); - - typename FunctionMap::type dirichlet_boundary; - BoundaryValues dirichlet_bc; - dirichlet_boundary[0] = &dirichlet_bc; - - std::vector > boundary_indices(triangulation.n_levels()); - std::vector component_mask (dim+1, true); - component_mask[dim] = false; - MGTools::make_boundary_list (dof_handler, dirichlet_boundary, - boundary_indices, component_mask); + std::vector > interface_dofs + = mg_constrained_dofs.get_refinement_edge_indices (); + std::vector > boundary_interface_dofs + = mg_constrained_dofs.get_refinement_edge_boundary_indices (); std::vector boundary_constraints (triangulation.n_levels()); std::vector boundary_interface_constraints (triangulation.n_levels()); for (unsigned int level=0; level::solve () GrowingVectorMemory > mg_vector_memory; - MGTransferPrebuilt > mg_transfer(constraints); + MGTransferPrebuilt > mg_transfer(constraints, mg_constrained_dofs); std::vector block_component (dim+1,0); block_component[dim] = 1; mg_transfer.set_component_to_block_map (block_component); -- 2.39.5