]> https://gitweb.dealii.org/ - dealii.git/commitdiff
changed the call of the constructor in MGTransferPrebuilt due to changes
authorBaerbel Jannsen <baerbel.janssen@gmail.com>
Thu, 26 Aug 2010 15:50:29 +0000 (15:50 +0000)
committerBaerbel Jannsen <baerbel.janssen@gmail.com>
Thu, 26 Aug 2010 15:50:29 +0000 (15:50 +0000)
in the library

git-svn-id: https://svn.dealii.org/trunk@21731 0785d39b-7218-0410-832d-ea1e28bc413d

tests/multigrid/mg_output.cc
tests/multigrid/mg_renumbered_01.cc
tests/multigrid/mg_renumbered_03.cc
tests/multigrid/transfer_system_adaptive_02.cc

index 72aaf40a216b194fdad14d931e722b6cbcfe4fe9..408b43fecb3bdd06c21ba17cda09e3aa4971a0ac 100644 (file)
@@ -231,10 +231,17 @@ void check_simple(const FiniteElement<dim>& fe)
   for (unsigned int level=0;level<tr.n_levels();++level)
     DoFRenumbering::component_wise (mgdof_renumbered, level, block_component);
 
-  MGTransferPrebuilt<Vector<double> > transfer(hnc);
+
+  MGConstrainedDoFs mg_constrained_dofs;
+  mg_constrained_dofs.initialize(mgdof);
+
+  MGConstrainedDoFs mg_constrained_dofs_renumbered;
+  mg_constrained_dofs_renumbered.initialize(mgdof_renumbered);
+
+  MGTransferPrebuilt<Vector<double> > transfer(hnc, mg_constrained_dofs);
   transfer.build_matrices(mgdof);
 
-  MGTransferPrebuilt<Vector<double> > transfer_renumbered(hnc_renumbered);
+  MGTransferPrebuilt<Vector<double> > transfer_renumbered(hnc_renumbered, mg_constrained_dofs_renumbered);
   transfer_renumbered.build_matrices(mgdof_renumbered);
 
   Vector<double> u(mgdof.n_dofs());
index f6c60157907f62d4db5b40df6f41a92653d1b408..bc7115297a7b0f1e497f027d821e56999a4233a7 100644 (file)
@@ -322,10 +322,7 @@ void LaplaceProblem<dim>::test ()
   typename FunctionMap<dim>::type      dirichlet_boundary;
   ZeroFunction<dim>                    dirichlet_bc(fe.n_components());
   dirichlet_boundary[0] =             &dirichlet_bc;
-  MGTools::make_boundary_list (mg_dof_handler, dirichlet_boundary,
-                              boundary_indices);
-  MGTools::make_boundary_list (mg_dof_handler_renumbered, dirichlet_boundary,
-                              boundary_indices_renumbered);
+
   const unsigned int min_l = mg_matrices.get_minlevel();
   const unsigned int max_l = mg_matrices.get_maxlevel();
   for(unsigned int l=min_l; l<max_l; ++l)
@@ -337,10 +334,19 @@ void LaplaceProblem<dim>::test ()
   GrowingVectorMemory<>   vector_memory;
   GrowingVectorMemory<>   vector_memory_renumbered;
 
-  MGTransferPrebuilt<Vector<double> > mg_transfer;
-  mg_transfer.build_matrices(mg_dof_handler, boundary_indices);
-  MGTransferPrebuilt<Vector<double> > mg_transfer_renumbered;
-  mg_transfer_renumbered.build_matrices(mg_dof_handler_renumbered, boundary_indices_renumbered);
+  MGConstrainedDoFs mg_constrained_dofs;
+  mg_constrained_dofs.initialize(mg_dof_handler, dirichlet_boundary);
+  MGConstrainedDoFs mg_constrained_dofs_renumbered;
+  mg_constrained_dofs_renumbered.initialize(mg_dof_handler_renumbered, dirichlet_boundary);
+
+  ConstraintMatrix constraints;
+  constraints.close();
+
+  MGTransferPrebuilt<Vector<double> > mg_transfer (constraints, mg_constrained_dofs);
+  mg_transfer.build_matrices(mg_dof_handler);
+  MGTransferPrebuilt<Vector<double> > mg_transfer_renumbered (constraints, mg_constrained_dofs_renumbered);
+  mg_transfer_renumbered.build_matrices(mg_dof_handler_renumbered);
+
   FullMatrix<double> coarse_matrix;
   coarse_matrix.copy_from (mg_matrices[0]);
   MGCoarseGridHouseholder<double, Vector<double> > mg_coarse;
index d3b13db33538a44b1393fc7bfa3c868229a88ed3..709599ae678a09bf9f581267adbb84405bbd6d68 100644 (file)
@@ -313,14 +313,19 @@ void LaplaceProblem<dim>::test ()
   typename FunctionMap<dim>::type      dirichlet_boundary;
   ZeroFunction<dim>                    dirichlet_bc(fe.n_components());
   dirichlet_boundary[0] =             &dirichlet_bc;
-  MGTools::make_boundary_list (mg_dof_handler, dirichlet_boundary,
-                              boundary_indices);
-  MGTools::make_boundary_list (mg_dof_handler_renumbered, dirichlet_boundary,
-                              boundary_indices_renumbered);
-  MGTransferPrebuilt<Vector<double> > mg_transfer;
-  mg_transfer.build_matrices(mg_dof_handler, boundary_indices);
-  MGTransferPrebuilt<Vector<double> > mg_transfer_renumbered;
-  mg_transfer_renumbered.build_matrices(mg_dof_handler_renumbered, boundary_indices_renumbered);
+
+  MGConstrainedDoFs mg_constrained_dofs;
+  mg_constrained_dofs.initialize(mg_dof_handler, dirichlet_boundary);
+  MGConstrainedDoFs mg_constrained_dofs_renumbered;
+  mg_constrained_dofs_renumbered.initialize(mg_dof_handler_renumbered, dirichlet_boundary);
+
+  ConstraintMatrix constraints;
+  constraints.close();
+
+  MGTransferPrebuilt<Vector<double> > mg_transfer (constraints, mg_constrained_dofs);
+  mg_transfer.build_matrices(mg_dof_handler);
+  MGTransferPrebuilt<Vector<double> > mg_transfer_renumbered (constraints, mg_constrained_dofs_renumbered);
+  mg_transfer_renumbered.build_matrices(mg_dof_handler_renumbered);
 
   Vector<double> test;
   test.reinit(mg_dof_handler.n_dofs());
index f749a8532b1bba604b8279dd8f0cd40cfc1f60c7..aa7b0e46f0679a83dab5e12be30978309f11e8f8 100644 (file)
@@ -148,11 +148,12 @@ void check (const FiniteElement<dim>& fe)
   ZeroFunction<dim>                    dirichlet_bc(fe.n_components());
   dirichlet_boundary[3] =             &dirichlet_bc;
 
-  MGTools::make_boundary_list (mg_dof_handler, dirichlet_boundary,
-                              boundary_indices);
+  MGConstrainedDoFs mg_constrained_dofs;
+  mg_constrained_dofs.initialize(mg_dof_handler, dirichlet_boundary);
 
-  MGTransferPrebuilt<Vector<double> > transfer;
-  transfer.build_matrices(mg_dof_handler, boundary_indices);
+  ConstraintMatrix constraints;
+  MGTransferPrebuilt<Vector<double> > transfer(constraints, mg_constrained_dofs);
+  transfer.build_matrices(mg_dof_handler);
 
   FullMatrix<double> prolong_0_1 (mg_dof_handler.n_dofs(1),
                                  mg_dof_handler.n_dofs(0));

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.