From 1722a701a0c2ad1faa3a36394f781b306db9298d Mon Sep 17 00:00:00 2001
From: "denis.davydov" <denis.davydov@0785d39b-7218-0410-832d-ea1e28bc413d>
Date: Fri, 21 Mar 2014 13:44:45 +0000
Subject: [PATCH] removed marking of artifical cells; Introduced bones of the
 new policy class; put back assert in get_subdomain_association.

git-svn-id: https://svn.dealii.org/branches/branch_sharedtria@32680 0785d39b-7218-0410-832d-ea1e28bc413d
---
 deal.II/include/deal.II/distributed/tria.h    |  5 +--
 .../include/deal.II/dofs/dof_handler_policy.h | 19 ++++++++++
 deal.II/source/distributed/tria.cc            | 36 -------------------
 deal.II/source/dofs/dof_handler_policy.cc     | 34 ++++++++++++++++++
 deal.II/source/dofs/dof_tools.cc              |  7 ++--
 5 files changed, 59 insertions(+), 42 deletions(-)

diff --git a/deal.II/include/deal.II/distributed/tria.h b/deal.II/include/deal.II/distributed/tria.h
index 2e3bd67ed6..a54ae63267 100644
--- a/deal.II/include/deal.II/distributed/tria.h
+++ b/deal.II/include/deal.II/distributed/tria.h
@@ -939,6 +939,9 @@ namespace parallel
     class Triangulation : public dealii::Triangulation<dim,spacedim>
     {
     public:
+      typedef typename dealii::Triangulation<dim,spacedim>::active_cell_iterator active_cell_iterator;
+      typedef typename dealii::Triangulation<dim,spacedim>::cell_iterator        cell_iterator;
+
       /**
        * Constructor.
        */
@@ -967,8 +970,6 @@ namespace parallel
       MPI_Comm mpi_communicator;
       types::subdomain_id my_subdomain;
       types::subdomain_id num_subdomains;
-      
-      void mark_artificial();
     };
   }
 }
diff --git a/deal.II/include/deal.II/dofs/dof_handler_policy.h b/deal.II/include/deal.II/dofs/dof_handler_policy.h
index 42eecf8afa..16ef98281b 100644
--- a/deal.II/include/deal.II/dofs/dof_handler_policy.h
+++ b/deal.II/include/deal.II/dofs/dof_handler_policy.h
@@ -127,6 +127,25 @@ namespace internal
                        dealii::DoFHandler<dim,spacedim> &dof_handler) const;
       };
 
+      template <int dim, int spacedim>
+      class ParallelShared : public Sequential<dim,spacedim>
+      {
+      public:
+    	  virtual
+    	  NumberCache
+    	  distribute_dofs (dealii::DoFHandler<dim,spacedim> &dof_handler) const;
+
+    	  virtual
+    	  void
+    	  distribute_mg_dofs (dealii::DoFHandler<dim,spacedim> &dof_handler,
+    	                      std::vector<NumberCache> &number_caches) const;
+
+    	  virtual
+    	  NumberCache
+    	  renumber_dofs (const std::vector<types::global_dof_index>  &new_numbers,
+    	                 dealii::DoFHandler<dim,spacedim> &dof_handler) const;
+      };
+
 
       /**
        * This class implements the
diff --git a/deal.II/source/distributed/tria.cc b/deal.II/source/distributed/tria.cc
index ef29d175da..75aeab6d02 100644
--- a/deal.II/source/distributed/tria.cc
+++ b/deal.II/source/distributed/tria.cc
@@ -4151,7 +4151,6 @@ namespace parallel
     Triangulation<dim,spacedim>::execute_coarsening_and_refinement () {
     	  dealii::Triangulation<dim,spacedim>::execute_coarsening_and_refinement ();
     	  dealii::GridTools::partition_triangulation (num_subdomains, *this);
-    	  mark_artificial();
     }
     
     template <int dim, int spacedim>
@@ -4171,7 +4170,6 @@ namespace parallel
           AssertThrow (false, ExcInternalError());
         }
       dealii::GridTools::partition_triangulation (num_subdomains, *this);
-      mark_artificial();
     }
 
 
@@ -4181,40 +4179,6 @@ namespace parallel
     {
       return mpi_communicator;
     }
-    
-    template <int dim, int spacedim>
-    void
-    Triangulation<dim,spacedim>::mark_artificial() 
-    {
-    
-    	std::vector<bool > marked_vertices(this->n_vertices(),false);
-        for (typename dealii::Triangulation<dim,spacedim>::active_cell_iterator
-                   cell = this->begin_active();
-                   cell != this->end(); ++cell)
-            if (cell->subdomain_id() == this->locally_owned_subdomain())
-                for (unsigned int v=0; v<GeometryInfo<dim>::vertices_per_cell; ++v)
-                    marked_vertices[cell->vertex_index(v)] = true;
-                    
-                    
-                    
-        for (typename dealii::Triangulation<dim,spacedim>::active_cell_iterator
-                   cell = this->begin_active();
-                   cell != this->end(); ++cell)
-            if (cell->subdomain_id() != this->locally_owned_subdomain()) {
-            	bool is_ghost = false;
-            	for (unsigned int v=0; v<GeometryInfo<dim>::vertices_per_cell; ++v)
-            		if (marked_vertices[cell->vertex_index(v)] == true)
-            		{
-            			is_ghost = true;
-            			break;
-        			}
-        		if (is_ghost)
-        			continue;
-        			
-        		cell->set_subdomain_id(numbers::artificial_subdomain_id);	
-        	}		    
-    }
-    
 
   }
 }
diff --git a/deal.II/source/dofs/dof_handler_policy.cc b/deal.II/source/dofs/dof_handler_policy.cc
index 3af161ec67..30926ebeba 100644
--- a/deal.II/source/dofs/dof_handler_policy.cc
+++ b/deal.II/source/dofs/dof_handler_policy.cc
@@ -963,7 +963,41 @@ namespace internal
         return number_cache;
       }
 
+      /* --------------------- class ParallelSequential ---------------- */
+      template <int dim, int spacedim>
+      NumberCache
+      ParallelShared<dim,spacedim>::
+      distribute_dofs (DoFHandler<dim,spacedim> &dof_handler) const
+      {
+    	  NumberCache number_cache = Sequential<dim,spacedim>::distribute_dofs (dof_handler);
+    	  //correct number_cache:
+
 
+    	  return number_cache;
+      }
+
+      template <int dim, int spacedim>
+      void
+      ParallelShared<dim,spacedim>::
+      distribute_mg_dofs (DoFHandler<dim,spacedim> &dof_handler,
+                          std::vector<NumberCache> &number_caches) const
+      {
+    	  NumberCache number_cache = Sequential<dim,spacedim>:: distribute_mg_dofs (dof_handler, number_caches);
+    	  Assert(false,ExcMessage("Not implemented"));
+    	  return number_cache;
+      }
+
+      template <int dim, int spacedim>
+      NumberCache
+      ParallelShared<dim,spacedim>::
+      renumber_dofs (const std::vector<types::global_dof_index> &new_numbers,
+                     dealii::DoFHandler<dim,spacedim> &dof_handler) const
+      {
+    	  NumberCache number_cache = Sequential<dim,spacedim>::renumber_dofs (new_numbers,dof_handler);
+    	  //correct number_cache:
+
+    	  return number_cache;
+      }
 
       /* --------------------- class ParallelDistributed ---------------- */
 
diff --git a/deal.II/source/dofs/dof_tools.cc b/deal.II/source/dofs/dof_tools.cc
index 1402055720..3ea5499334 100644
--- a/deal.II/source/dofs/dof_tools.cc
+++ b/deal.II/source/dofs/dof_tools.cc
@@ -1096,10 +1096,9 @@ namespace DoFTools
     endc = dof_handler.end();
     for (; cell!=endc; ++cell)
       {
-        //can comment it since the case of distributed triangulations is handled above.
-        //Assert (cell->is_artificial() == false,
-        //        ExcMessage ("You can't call this function for meshes that "
-        //                    "have artificial cells."));
+        Assert (cell->is_artificial() == false,
+                ExcMessage ("You can't call this function for meshes that "
+                            "have artificial cells."));
 
         const types::subdomain_id subdomain_id = cell->subdomain_id();
         const unsigned int dofs_per_cell = cell->get_fe().dofs_per_cell;
-- 
2.39.5