From aa3af9899efb9e1b2868c965ad70423122150c9d Mon Sep 17 00:00:00 2001
From: bangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Date: Sun, 11 Sep 2011 05:10:03 +0000
Subject: [PATCH] Add a function cell->is_locally_owned.

git-svn-id: https://svn.dealii.org/trunk@24303 0785d39b-7218-0410-832d-ea1e28bc413d
---
 deal.II/doc/news/changes.h                    |  11 ++
 deal.II/include/deal.II/grid/tria_accessor.h  |  25 ++++
 .../deal.II/grid/tria_accessor.templates.h    |  23 ++++
 tests/mpi/is_locally_owned.cc                 | 114 ++++++++++++++++++
 .../mpi/is_locally_owned/ncpu_10/cmp/generic  |  11 ++
 tests/mpi/is_locally_owned/ncpu_4/cmp/generic |  11 ++
 6 files changed, 195 insertions(+)
 create mode 100644 tests/mpi/is_locally_owned.cc
 create mode 100644 tests/mpi/is_locally_owned/ncpu_10/cmp/generic
 create mode 100644 tests/mpi/is_locally_owned/ncpu_4/cmp/generic

diff --git a/deal.II/doc/news/changes.h b/deal.II/doc/news/changes.h
index 2055fe1710..a5248326e6 100644
--- a/deal.II/doc/news/changes.h
+++ b/deal.II/doc/news/changes.h
@@ -283,6 +283,17 @@ and DoF handlers embedded in higher dimensional space have been added.
 <h3>Specific improvements</h3>
 
 <ol>
+<li> New: The CellAccessor::is_locally_owned function is a shortcut
+for the <code>!cell-@>is_ghost() && !cell-@>is_artificial()</code> pattern
+found in many places when dealing with distributed meshes.
+<br>
+(Wolfgang Bangerth, 2011/09/10)
+
+<li> New: The GridGenerator::torus function and TorusBoundary class can
+create and describe the surface of a torus.
+<br>
+(Daniel Castanon Quiroz, 2011/09/08)
+
 <li> New: The GridGenerator::torus function and TorusBoundary class can
 create and describe the surface of a torus.
 <br>
diff --git a/deal.II/include/deal.II/grid/tria_accessor.h b/deal.II/include/deal.II/grid/tria_accessor.h
index e5ba1977bc..5d182ff8e8 100644
--- a/deal.II/include/deal.II/grid/tria_accessor.h
+++ b/deal.II/include/deal.II/grid/tria_accessor.h
@@ -2719,6 +2719,23 @@ class CellAccessor :  public TriaAccessor<dim,dim,spacedim>
 				      */
     bool active () const;
 
+    /**
+     * Return whether this cell is owned by the current processor
+     * or is owned by another processor. The function always returns
+     * true if applied to an object of type dealii::Triangulation,
+     * but may yield false if the triangulation is of type
+     * parallel::distributed::Triangulation.
+     *
+     * See the @ref GlossGhostCell
+     * "glossary" and the @ref
+     * distributed module for more
+     * information.
+     * 
+     * @post The returned value is equal to <code>!is_ghost() &&
+     * !is_artificial()</code>.
+     **/
+    bool is_locally_owned () const;
+    
 				     /**
 				      * Return whether this cell
 				      * exists in the global mesh but
@@ -2741,6 +2758,10 @@ class CellAccessor :  public TriaAccessor<dim,dim,spacedim>
 				      * "glossary" and the @ref
 				      * distributed module for more
 				      * information.
+				      *
+				      * @post The returned value is equal to
+				      * <code>!is_locally_owned() &&
+				      * !is_artificial()</code>.
 				      */
     bool is_ghost () const;
 
@@ -2775,6 +2796,10 @@ class CellAccessor :  public TriaAccessor<dim,dim,spacedim>
 				      * GlossArtificialCell "glossary"
 				      * and the @ref distributed
 				      * module for more information.
+				      *
+				      * @post The returned value is equal to
+				      * <code>!is_ghost() &&
+				      * !is_artificial()</code>.
 				      */
     bool is_artificial () const;
 
diff --git a/deal.II/include/deal.II/grid/tria_accessor.templates.h b/deal.II/include/deal.II/grid/tria_accessor.templates.h
index eb0879bbb8..967c3c5639 100644
--- a/deal.II/include/deal.II/grid/tria_accessor.templates.h
+++ b/deal.II/include/deal.II/grid/tria_accessor.templates.h
@@ -2806,6 +2806,29 @@ CellAccessor<dim,spacedim>::active () const
 
 
 
+template <int dim, int spacedim>
+inline
+bool
+CellAccessor<dim,spacedim>::is_locally_owned () const
+{
+#ifndef DEAL_II_USE_P4EST
+  return true;
+#else
+  const types::subdomain_id_t subdomain = this->subdomain_id();
+  if (subdomain == types::artificial_subdomain_id)
+    return false;
+
+  const parallel::distributed::Triangulation<dim,spacedim> *pdt
+    = dynamic_cast<const parallel::distributed::Triangulation<dim,spacedim> *>(this->tria);
+
+  if (pdt == 0)
+    return true;
+  else
+    return (subdomain == pdt->locally_owned_subdomain());
+#endif
+}
+
+
 template <int dim, int spacedim>
 inline
 bool
diff --git a/tests/mpi/is_locally_owned.cc b/tests/mpi/is_locally_owned.cc
new file mode 100644
index 0000000000..d48d7bc973
--- /dev/null
+++ b/tests/mpi/is_locally_owned.cc
@@ -0,0 +1,114 @@
+//---------------------------------------------------------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2009, 2010, 2011 by the deal.II authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//---------------------------------------------------------------------------
+
+
+// check CellAccessor::is_locally_owned
+
+#include "../tests.h"
+#include "coarse_grid_common.h"
+#include <deal.II/base/logstream.h>
+#include <deal.II/base/tensor.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/distributed/tria.h>
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_out.h>
+#include <deal.II/grid/grid_tools.h>
+#include <deal.II/base/utilities.h>
+
+
+#include <fstream>
+
+
+template<int dim>
+void test()
+{
+  unsigned int myid = Utilities::System::get_this_mpi_process (MPI_COMM_WORLD);
+  unsigned int numprocs = Utilities::System::get_n_mpi_processes (MPI_COMM_WORLD);
+
+  if (Utilities::System::get_this_mpi_process (MPI_COMM_WORLD) == 0)
+    deallog << "hyper_cube" << std::endl;
+
+  parallel::distributed::Triangulation<dim> tr(MPI_COMM_WORLD);
+
+  GridGenerator::subdivided_hyper_cube(tr, 3);
+
+  typename Triangulation<dim,dim>::active_cell_iterator cell;
+
+  for (cell = tr.begin_active();
+       cell != tr.end();
+       ++cell)
+    {
+      if (cell->is_locally_owned())
+	{
+	  if (myid==0)
+	    deallog << cell << ": locally owned"
+		    << std::endl;
+	  Assert (!cell->is_ghost() && !cell->is_artificial(),
+		  ExcInternalError());
+	}
+      else if (cell->is_ghost())
+	{
+	  if (myid==0)
+	    deallog << cell << ": ghost"
+		    << std::endl;
+	  Assert (!cell->is_locally_owned() && !cell->is_artificial(),
+		  ExcInternalError());
+	}
+      else if (cell->is_artificial())
+	{
+	  if (myid==0)
+	    deallog << cell << ": artificial"
+		    << std::endl;
+	  Assert (!cell->is_locally_owned() && !cell->is_ghost(),
+		  ExcInternalError());
+	}
+      else
+	Assert (false, ExcInternalError());
+    }
+}
+
+
+int main(int argc, char *argv[])
+{
+#ifdef DEAL_II_COMPILER_SUPPORTS_MPI
+  MPI_Init (&argc,&argv);
+#else
+  (void)argc;
+  (void)argv;
+#endif
+
+  unsigned int myid = Utilities::System::get_this_mpi_process (MPI_COMM_WORLD);
+
+
+  deallog.push(Utilities::int_to_string(myid));
+
+  if (myid == 0)
+    {
+      std::ofstream logfile(output_file_for_mpi("is_locally_owned").c_str());
+      deallog.attach(logfile);
+      deallog.depth_console(0);
+      deallog.threshold_double(1.e-10);
+
+      deallog.push("2d");
+      test<2>();
+      deallog.pop();
+    }
+  else
+    test<2>();
+
+
+#ifdef DEAL_II_COMPILER_SUPPORTS_MPI
+  MPI_Finalize();
+#endif
+}
diff --git a/tests/mpi/is_locally_owned/ncpu_10/cmp/generic b/tests/mpi/is_locally_owned/ncpu_10/cmp/generic
new file mode 100644
index 0000000000..57dd7499f6
--- /dev/null
+++ b/tests/mpi/is_locally_owned/ncpu_10/cmp/generic
@@ -0,0 +1,11 @@
+
+DEAL:0:2d::hyper_cube
+DEAL:0:2d::0.0: artificial
+DEAL:0:2d::0.1: artificial
+DEAL:0:2d::0.2: artificial
+DEAL:0:2d::0.3: artificial
+DEAL:0:2d::0.4: artificial
+DEAL:0:2d::0.5: artificial
+DEAL:0:2d::0.6: artificial
+DEAL:0:2d::0.7: artificial
+DEAL:0:2d::0.8: artificial
diff --git a/tests/mpi/is_locally_owned/ncpu_4/cmp/generic b/tests/mpi/is_locally_owned/ncpu_4/cmp/generic
new file mode 100644
index 0000000000..fdbd4c193f
--- /dev/null
+++ b/tests/mpi/is_locally_owned/ncpu_4/cmp/generic
@@ -0,0 +1,11 @@
+
+DEAL:0:2d::hyper_cube
+DEAL:0:2d::0.0: locally owned
+DEAL:0:2d::0.1: locally owned
+DEAL:0:2d::0.2: ghost
+DEAL:0:2d::0.3: ghost
+DEAL:0:2d::0.4: ghost
+DEAL:0:2d::0.5: ghost
+DEAL:0:2d::0.6: artificial
+DEAL:0:2d::0.7: artificial
+DEAL:0:2d::0.8: artificial
-- 
2.39.5