]> https://gitweb.dealii.org/ - dealii.git/commitdiff
pass a MGConstraints instead of the boundary_indices only in
authorBaerbel Jannsen <baerbel.janssen@gmail.com>
Wed, 25 Aug 2010 17:24:53 +0000 (17:24 +0000)
committerBaerbel Jannsen <baerbel.janssen@gmail.com>
Wed, 25 Aug 2010 17:24:53 +0000 (17:24 +0000)
build_matrices

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

deal.II/deal.II/include/multigrid/mg_transfer.h
deal.II/deal.II/source/multigrid/mg_transfer_prebuilt.cc
deal.II/examples/step-16/step-16.cc

index 7db0779cc55ce5e3ce02241fdafdf74f29d5138a..5ec3cf99ae95a4a77a0f5badf6556aa8aa74e938 100644 (file)
@@ -27,6 +27,7 @@
 #include <lac/vector_memory.h>
 
 #include <multigrid/mg_base.h>
+#include <multigrid/mg_constraints.h>
 #include <base/mg_level_object.h>
 
 
@@ -87,9 +88,7 @@ class MGTransferPrebuilt : public MGTransferBase<VECTOR>
                                      */
     template <int dim, int spacedim>
     void build_matrices (const MGDoFHandler<dim,spacedim> &mg_dof,
-    const std::vector<std::set<unsigned int> >&boundary_indices
-                        = std::vector<std::set<unsigned int> >()
-        );
+    const MGConstraints& mg_constraints);
 
     virtual void prolongate (const unsigned int    to_level,
                             VECTOR       &dst,
index abcf291e51c7db493033b09fdb0485e15914ac1a..b170691295e6e81eeb1d672d125b1ef36052a068 100644 (file)
@@ -155,9 +155,8 @@ MGTransferPrebuilt<VECTOR>::copy_to_mg (
 template <typename VECTOR>
 template <int dim, int spacedim>
 void MGTransferPrebuilt<VECTOR>::build_matrices (
-  const MGDoFHandler<dim,spacedim>           &mg_dof,
-  const std::vector<std::set<unsigned int> > &boundary_indices
-  )
+  const MGDoFHandler<dim,spacedim>  &mg_dof,
+  const MGConstraints               &mg_constraints)
 {
   const unsigned int n_levels      = mg_dof.get_tria().n_levels();
   const unsigned int dofs_per_cell = mg_dof.get_fe().dofs_per_cell;
@@ -287,12 +286,12 @@ void MGTransferPrebuilt<VECTOR>::build_matrices (
                                // impose boundary conditions
                                // but only in the column of
                                // the prolongation matrix
-  if (boundary_indices.size() != 0)
+  if (mg_constraints.get_boundary_indices().size() != 0)
     {
       std::vector<unsigned int> constrain_indices;
       for (int level=n_levels-2; level>=0; --level)
        {
-         if (boundary_indices[level].size() == 0)
+         if (mg_constraints.get_boundary_indices()[level].size() == 0)
            continue;
 
                                // need to delete all the columns in the
@@ -302,8 +301,9 @@ void MGTransferPrebuilt<VECTOR>::build_matrices (
                                // need to filter away.
          constrain_indices.resize (0);
          constrain_indices.resize (prolongation_matrices[level]->n(), 0);
-         std::set<unsigned int>::const_iterator dof = boundary_indices[level].begin(),
-           endd = boundary_indices[level].end();
+         std::set<unsigned int>::const_iterator dof 
+            = mg_constraints.get_boundary_indices()[level].begin(),
+           endd = mg_constraints.get_boundary_indices()[level].end();
          for (; dof != endd; ++dof)
            constrain_indices[*dof] = 1;
 
@@ -332,10 +332,6 @@ void MGTransferPrebuilt<VECTOR>::build_matrices (
                                // the std::vector<std::pair<uint,uint> > that
                                // only contains the active dofs on the
                                // levels.
-  interface_dofs.resize(mg_dof.get_tria().n_levels());
-  for(unsigned int l=0; l<mg_dof.get_tria().n_levels(); ++l)
-    interface_dofs[l].resize(mg_dof.n_dofs(l));
-  MGTools::extract_inner_interface_dofs (mg_dof, interface_dofs);
 
   copy_indices.resize(n_levels);
   std::vector<unsigned int> temp_copy_indices;
@@ -366,7 +362,7 @@ void MGTransferPrebuilt<VECTOR>::build_matrices (
 
          for (unsigned int i=0; i<dofs_per_cell; ++i)
           {
-           if(!interface_dofs[level][level_dof_indices[i]])
+           if(!mg_constraints.at_refinement_edge(level,level_dof_indices[i]))
              temp_copy_indices[level_dof_indices[i]] = global_dof_indices[i];
           }
        }
@@ -397,22 +393,22 @@ void MGTransferPrebuilt<VECTOR>::build_matrices (
 template
 void MGTransferPrebuilt<Vector<float> >::build_matrices<deal_II_dimension>
 (const MGDoFHandler<deal_II_dimension> &mg_dof,
- const std::vector<std::set<unsigned int> >&boundary_indices);
+ const MGConstraints& mg_constraints);
 
 template
 void MGTransferPrebuilt<Vector<double> >::build_matrices<deal_II_dimension>
 (const MGDoFHandler<deal_II_dimension> &mg_dof,
- const std::vector<std::set<unsigned int> >&boundary_indices);
+ const MGConstraints& mg_constraints);
 
 template
 void MGTransferPrebuilt<BlockVector<float> >::build_matrices<deal_II_dimension>
 (const MGDoFHandler<deal_II_dimension> &mg_dof,
- const std::vector<std::set<unsigned int> >&boundary_indices);
+ const MGConstraints& mg_constraints);
 
 template
 void MGTransferPrebuilt<BlockVector<double> >::build_matrices<deal_II_dimension>
 (const MGDoFHandler<deal_II_dimension> &mg_dof,
- const std::vector<std::set<unsigned int> >&boundary_indices);
+ const MGConstraints& mg_constraints);
 
 template void
 MGTransferPrebuilt<Vector<float> >::copy_to_mg (
index 78d5f2fdbe9d05e9fb696717e787cedbf443cda8..383d8a3feac7eeda6b99ffeed0cd1818f5d56c8b 100644 (file)
@@ -767,8 +767,11 @@ void LaplaceProblem<dim>::solve ()
   MGTransferPrebuilt<Vector<double> > 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());
+                                 // each level into account and needs to know about 
+                                 // the indices at the refinement egdes. The 
+                                 // <code>MGConstraints</code> knows about that so
+                                 // pass it as an argument.
+  mg_transfer.build_matrices(mg_dof_handler, mg_constraints);
 
   FullMatrix<double> coarse_matrix;
   coarse_matrix.copy_from (mg_matrices[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.