]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Fix MGTransfer::copy_from_mg_add
authorkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 1 May 2007 23:02:40 +0000 (23:02 +0000)
committerkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 1 May 2007 23:02:40 +0000 (23:02 +0000)
git-svn-id: https://svn.dealii.org/trunk@14656 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/multigrid/mg_transfer.h
deal.II/deal.II/include/multigrid/mg_transfer.templates.h
deal.II/deal.II/source/multigrid/mg_transfer_prebuilt.cc

index 1500a037dd95f236406105f6de5ec0bf4c81caa0..9bb9fc1f8e5e130135f7b00ba9d3a31eea492bce 100644 (file)
@@ -167,6 +167,15 @@ class MGTransferPrebuilt : public MGTransferBase<VECTOR>
                                      * child cell, i.e. the fine level.
                                      */
     std::vector<boost::shared_ptr<SparseMatrix<double> > > prolongation_matrices;    
+    
+                                    /**
+                                     * Mapping for the
+                                     * <tt>copy_to/from_mg</tt>-functions.
+                                     * The data is first the global
+                                     * index, then the level index.
+                                    */
+    std::vector<std::map<unsigned int, unsigned int> >
+    copy_indices;
 };
 
 
index bef63bb1283318d48ef4e4889c20e2a6040bdf11..a09d83816645e17a28efa145a9c4bb83596364ea 100644 (file)
@@ -29,76 +29,29 @@ DEAL_II_NAMESPACE_OPEN
 
 /* --------------------- MGTransferPrebuilt -------------- */
 
+// Simplify some things below
+typedef std::map<unsigned int, unsigned int>::const_iterator IT;
 
 
 template <class VECTOR>
 template <int dim, class InVector>
 void
 MGTransferPrebuilt<VECTOR>::copy_to_mg (
-  const MGDoFHandler<dim>        &mg_dof_handler,
-  MGLevelObject<VECTOR> &dst,
-  const InVector                 &src) const
+  const MGDoFHandler<dim>mg_dof_handler,
+  MGLevelObject<VECTOR>dst,
+  const InVectorsrc) const
 {
-  const FiniteElement<dim>& fe = mg_dof_handler.get_fe();
-  
-  const unsigned int dofs_per_cell = fe.dofs_per_cell;
-  
-                                  // set the elements of the vectors
-                                  // on all levels to zero
-  unsigned int minlevel = dst.get_minlevel();
-  unsigned int maxlevel = dst.get_maxlevel();
-
   MGTools::reinit_vector(mg_dof_handler, dst);
-
-  Assert(sizes.size()==mg_dof_handler.get_tria().n_levels(),
-        ExcMatricesNotBuilt());
-  
-  std::vector<unsigned int> global_dof_indices (dofs_per_cell);
-  std::vector<unsigned int> level_dof_indices  (dofs_per_cell);
-  
-                                  // traverse the grid top-down
-                                  // (i.e. starting with the most
-                                  // refined grid). this way, we can
-                                  // always get that part of one
-                                  // level of the output vector which
-                                  // corresponds to a region which is
-                                  // more refined, by restriction of
-                                  // the respective vector on the
-                                  // next finer level, which we then
-                                  // already have built.
-  for (int level=maxlevel; level>=static_cast<signed int>(minlevel); --level)
+  bool first = true;
+  for (unsigned int level=mg_dof_handler.get_tria().n_levels();level != 0;)
     {
-      typename MGDoFHandler<dim>::active_cell_iterator
-       level_cell = mg_dof_handler.begin_active(level);
-      const typename MGDoFHandler<dim>::active_cell_iterator
-       level_end  = mg_dof_handler.end_active(level);
-
-                                      // Compute coarse level right hand side
-                                      // by restricting from fine level.
-      for (; level_cell!=level_end; ++level_cell)
-       {
-                                          // get the dof numbers of
-                                          // this cell for the global
-                                          // and the level-wise
-                                          // numbering
-         level_cell->get_dof_indices(global_dof_indices);
-         level_cell->get_mg_dof_indices (level_dof_indices);
-
-                                          // transfer the global
-                                          // defect in the vector
-                                          // into the level-wise one
-         for (unsigned int i=0; i<dofs_per_cell; ++i)
-           dst[level](level_dof_indices[i]) = src(global_dof_indices[i]);
-       }
-                                      // for that part of the level
-                                      // which is further refined:
-                                      // get the defect by
-                                      // restriction of the defect on
-                                      // one level higher
-      if (static_cast<unsigned int>(level) < maxlevel)
-       {
-         restrict_and_add (level+1, dst[level], dst[level+1]);
-       }
+      --level;
+      for (IT i= copy_indices[level].begin();
+          i != copy_indices[level].end();++i)
+       dst[level](i->second) = src(i->first);
+      if (!first)
+       restrict_and_add (level+1, dst[level], dst[level+1]);
+      first = false;
     }
 }
 
@@ -112,39 +65,10 @@ MGTransferPrebuilt<VECTOR>::copy_from_mg(
   OutVector&                     dst,
   const MGLevelObject<VECTOR>& src) const
 {
-  const FiniteElement<dim>& fe = mg_dof_handler.get_fe();
-  
-  const unsigned int dofs_per_cell = fe.dofs_per_cell;
-  std::vector<unsigned int> global_dof_indices (dofs_per_cell);
-  std::vector<unsigned int> level_dof_indices (dofs_per_cell);
-
-  typename MGDoFHandler<dim>::active_cell_iterator
-    level_cell = mg_dof_handler.begin_active();
-  const typename MGDoFHandler<dim>::active_cell_iterator
-    endc = mg_dof_handler.end();
-
-                                  // traverse all cells and copy the
-                                  // data appropriately to the output
-                                  // vector
-
-                                  // Note that the level is
-                                  // monotonuosly increasing
-  for (; level_cell != endc; ++level_cell)
-    {
-      const unsigned int level = level_cell->level();
-      
-                                      // get the dof numbers of
-                                      // this cell for the global
-                                      // and the level-wise
-                                      // numbering
-      level_cell->get_dof_indices (global_dof_indices);
-      level_cell->get_mg_dof_indices(level_dof_indices);
-
-                                      // copy level-wise data to
-                                      // global vector
-      for (unsigned int i=0; i<dofs_per_cell; ++i)
-       dst(global_dof_indices[i]) = src[level](level_dof_indices[i]);
-    }
+  for (unsigned int level=0;level<mg_dof_handler.get_tria().n_levels();++level)
+    for (IT i= copy_indices[level].begin();
+        i != copy_indices[level].end();++i)
+      dst(i->first) = src[level](i->second);
 }
 
 
@@ -157,38 +81,10 @@ MGTransferPrebuilt<VECTOR>::copy_from_mg_add (
   OutVector                            &dst,
   const MGLevelObject<VECTOR> &src) const
 {
-  const FiniteElement<dim>& fe = mg_dof_handler.get_fe();
-  
-  const unsigned int dofs_per_cell = fe.dofs_per_cell;
-  std::vector<unsigned int> global_dof_indices (dofs_per_cell);
-  std::vector<unsigned int> level_dof_indices (dofs_per_cell);
-
-  typename MGDoFHandler<dim>::active_cell_iterator
-    level_cell = mg_dof_handler.begin_active();
-  const typename MGDoFHandler<dim>::active_cell_iterator
-    endc = mg_dof_handler.end();
-
-                                  // traverse all cells and copy the
-                                  // data appropriately to the output
-                                  // vector
-
-                                  // Note that the level is
-                                  // monotonuosly increasing
-  for (; level_cell != endc; ++level_cell)
-    {
-      const unsigned int level = level_cell->level();
-      
-                                      // get the dof numbers of
-                                      // this cell for the global
-                                      // and the level-wise
-                                      // numbering
-      level_cell->get_dof_indices (global_dof_indices);
-      level_cell->get_mg_dof_indices(level_dof_indices);
-                                      // copy level-wise data to
-                                      // global vector
-      for (unsigned int i=0; i<dofs_per_cell; ++i)
-       dst(global_dof_indices[i]) += src[level](level_dof_indices[i]);
-    }
+  for (unsigned int level=0;level<mg_dof_handler.get_tria().n_levels();++level)
+    for (IT i= copy_indices[level].begin();
+        i != copy_indices[level].end();++i)
+      dst(i->first) += src[level](i->second);
 }
 
 
index 291d56830198f066c6794a58a53127a5b58d0ec4..9ab4e8f59997f767b92a8ba4fd04854359743da2 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 2003, 2004, 2005, 2006 by the deal.II authors
+//    Copyright (C) 2003, 2004, 2005, 2006, 2007 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -149,9 +149,39 @@ void MGTransferPrebuilt<number>::build_matrices (
                      prolongation_matrices[level]->set (dof_indices_child[i],
                                                         dof_indices_parent[j],
                                                         prolongation(i,j));
-             };
-         };
-    };
+             }
+         }
+    }
+
+  copy_indices.resize(mg_dof.get_tria().n_levels());
+  std::vector<unsigned int> global_dof_indices (dofs_per_cell);
+  std::vector<unsigned int> level_dof_indices  (dofs_per_cell);
+  for (int level=mg_dof.get_tria().n_levels()-1; level>=0; --level)
+    {
+      typename MGDoFHandler<dim>::active_cell_iterator
+       level_cell = mg_dof.begin_active(level);
+      const typename MGDoFHandler<dim>::active_cell_iterator
+       level_end  = mg_dof.end_active(level);
+      
+                                      // Compute coarse level right hand side
+                                      // by restricting from fine level.
+      for (; level_cell!=level_end; ++level_cell)
+       {
+         DoFObjectAccessor<dim, DoFHandler<dim> >& global_cell = *level_cell;
+                                          // get the dof numbers of
+                                          // this cell for the global
+                                          // and the level-wise
+                                          // numbering
+         global_cell.get_dof_indices(global_dof_indices);
+         level_cell->get_mg_dof_indices (level_dof_indices);
+         
+         for (unsigned int i=0; i<dofs_per_cell; ++i)
+           {
+             copy_indices[level].insert(
+               std::make_pair(global_dof_indices[i], level_dof_indices[i]));
+           }
+       }
+    }
 }
 
 

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.