]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Provide subdomain ids to triangulation cells.
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 29 Oct 2001 08:43:35 +0000 (08:43 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 29 Oct 2001 08:43:35 +0000 (08:43 +0000)
git-svn-id: https://svn.dealii.org/trunk@5166 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/grid/tria_accessor.h
deal.II/deal.II/include/grid/tria_accessor.templates.h
deal.II/deal.II/include/grid/tria_levels.h
deal.II/deal.II/source/grid/tria_accessor.cc
tests/deal.II/Makefile
tests/deal.II/subdomain_ids.cc [new file with mode: 0644]
tests/deal.II/subdomain_ids.checked [new file with mode: 0644]

index 5385ceff5e76133f10eeadcf1f9c5022271e062c..20271723070a5d6d32d625790c9206e234c349a3 100644 (file)
@@ -14,6 +14,7 @@
 #define __deal2__tria_accessor_h
 
 
+#include <base/config.h>
 #include <base/exceptions.h>
 #include <grid/tria_iterator_base.h>
 
@@ -1966,14 +1967,28 @@ class CellAccessor :  public TriaObjectAccessor<dim,dim>
     face (const unsigned int i) const;
     
                                     /**
-                                     * Return the material id of this cell.
+                                     * Return the material id of this
+                                     * cell.
                                      */
     unsigned char material_id () const;
 
                                     /**
-                                     * Set the material id of this cell.
+                                     * Set the material id of this
+                                     * cell.
                                      */
-    void set_material_id (const unsigned char) const;
+    void set_material_id (const unsigned char new_material_id) const;
+
+                                    /**
+                                     * Return the subdomain id of
+                                     * this cell.
+                                     */
+    unsigned int subdomain_id () const;
+
+                                    /**
+                                     * Set the subdomain id of this
+                                     * cell.
+                                     */
+    void set_subdomain_id (const unsigned int new_subdomain_id) const;
 
                                     /**
                                      * Test whether the cell has children
@@ -2056,6 +2071,7 @@ template <> bool CellAccessor<3>::at_boundary () const;
 template <> unsigned char CellAccessor<3>::material_id () const;
 template <> void CellAccessor<3>::set_material_id (const unsigned char mat_id) const;
 template <> bool CellAccessor<3>::point_inside (const Point<3> &) const;
+
 template <> bool CellAccessor<1>::has_boundary_lines () const;
 
 template <> double TriaObjectAccessor<2, 2>::measure () const;
index 367fdb0599a3111730f1dfb2f8c0efddb4886c50..1fd585ed08032ef57ac4a85ca419f4495dbcec1e 100644 (file)
@@ -10,6 +10,7 @@
 #define __deal2__tria_accessor_templates_h
 
 
+#include <base/config.h>
 #include <grid/tria.h>
 #include <grid/tria_levels.h>
 #include <grid/tria_iterator.h>
index 4c3c454c477bf77bd95aa32b39ac62a79b2f4838..8e91ee985bd2064ac2a64bb6c5715b9644eb3080 100644 (file)
@@ -14,6 +14,7 @@
 #define __deal2__tria_levels_h
 
 
+#include <base/config.h>
 #include <vector>
 #include <grid/tria_line.h>
 #include <grid/tria_quad.h>
@@ -127,6 +128,18 @@ class TriangulationLevel<0>
                                      *  to the mother cell of this cell).
                                      */
     std::vector<std::pair<int,int> > neighbors;
+
+                                    /**
+                                     * One integer per cell to store
+                                     * which subdomain it belongs
+                                     * to. This field is most often
+                                     * used in parallel computations,
+                                     * where it denotes which
+                                     * processor shall work on the
+                                     * cells with a given subdomain
+                                     * number.
+                                     */
+    std::vector<unsigned int> subdomain_ids;
     
                                     /**
                                      *  Reserve enough space to accomodate
index c1c548281d1e9985287ab19ed6cf166c99b54e14..9701c1dd08e8ddfe14a111fdb5218681306504e1 100644 (file)
@@ -1562,7 +1562,7 @@ void CellAccessor<1>::set_material_id (const unsigned char mat_id) const
 {
   Assert (used(), TriaAccessor<1>::ExcCellNotUsed());
   tria->levels[present_level]->lines.material_id[present_index]
-    = mat_id;                                           
+    = mat_id;
 };
 
 
@@ -1705,6 +1705,26 @@ bool CellAccessor<3>::point_inside (const Point<3> &) const
 /*------------------------ Functions: CellAccessor<dim> -----------------------*/
 
 
+template <int dim>
+unsigned int CellAccessor<dim>::subdomain_id () const
+{
+  Assert (used(), typename TriaAccessor<dim>::ExcCellNotUsed());
+  return tria->levels[present_level]->subdomain_ids[present_index];
+};
+
+
+
+template <int dim>
+void
+CellAccessor<dim>::set_subdomain_id (const unsigned int new_subdomain_id) const
+{
+  Assert (used(), typename TriaAccessor<dim>::ExcCellNotUsed());
+  tria->levels[present_level]->subdomain_ids[present_index]
+    = new_subdomain_id;
+};
+
+
+
 template <int dim>
 void CellAccessor<dim>::set_neighbor (const unsigned int i,
                                      const TriaIterator<dim,CellAccessor<dim> > &pointer) const
index f421960a69493c55f2d27c9c228f18df699ea301..27474438a4903df55c91f487200d7a210fbb2dd0 100644 (file)
@@ -45,12 +45,14 @@ filtered_matrix.exe          : filtered_matrix.go          $(lib-1d) $(lib-2d) $
 boundaries.exe               : boundaries.go               $(lib-1d) $(lib-2d) $(lib-3d) $(libraries)
 sparsity_pattern.exe         : sparsity_pattern.go         $(lib-1d) $(lib-2d) $(lib-3d) $(libraries)
 grid_tools.exe               : grid_tools.go               $(lib-1d) $(lib-2d) $(lib-3d) $(libraries)
+subdomain_ids.exe            : subdomain_ids.go            $(lib-1d) $(lib-2d) $(lib-3d) $(libraries)
 
 tests = grid_test grid_transform dof_test data_out derivatives gradients \
        constraints block_matrices second_derivatives derivative_approximation \
         matrices error_estimator intergrid_constraints intergrid_map \
         wave-test-3 dof_renumbering support_point_map filtered_matrix \
-       boundaries sparsity_pattern grid_tools
+       boundaries sparsity_pattern grid_tools subdomain_ids
+
 ############################################################
 
 
diff --git a/tests/deal.II/subdomain_ids.cc b/tests/deal.II/subdomain_ids.cc
new file mode 100644 (file)
index 0000000..be7d3ec
--- /dev/null
@@ -0,0 +1,185 @@
+//----------------------------  subdomain_ids.cc  ---------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2001 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.
+//
+//----------------------------  subdomain_ids.cc  ---------------------------
+
+// check some things about subdomain_ids
+
+
+#include <base/logstream.h>
+#include <grid/tria.h>
+#include <grid/grid_generator.h>
+#include <grid/tria_iterator.h>
+#include <dofs/dof_handler.h>
+#include <dofs/dof_accessor.h>
+#include <fe/fe_dgq.h>
+#include <fe/fe_q.h>
+#include <dofs/dof_tools.h>
+
+#include <fstream>
+#include <algorithm>
+#include <cmath>
+
+
+std::ofstream logfile("subdomain_ids.output");
+
+
+DeclException2 (ExcNumberMismatch,
+               int, int,
+               << "The numbers " << arg1 << " and " << arg2
+               << " should be equation, but are not.");
+
+
+template <int dim>
+void test ()
+{
+  Triangulation<dim> tria;
+  GridGenerator::hyper_cube(tria, -1, 1);
+  tria.refine_global (3);
+
+                                  // we now have a number of cells,
+                                  // flag them with some subdomain
+                                  // ids based on their position, in
+                                  // particular we take the quadrant
+                                  // (octant)
+  typename Triangulation<dim>::active_cell_iterator
+    cell = tria.begin_active (),
+    endc = tria.end ();
+  for (; cell!=endc; ++cell)
+    {
+      unsigned int subdomain = 0;
+      for (unsigned int d=0; d<dim; ++d)
+       if (cell->center()(d) > 0)
+         subdomain |= (1<<d);
+      Assert (subdomain < (1<<dim), ExcInternalError());
+
+      cell->set_subdomain_id (subdomain);
+    };
+
+                                  // refine twice (for 3d only once)
+                                  // and count again the numbers of
+                                  // cells in each subdomain
+  tria.refine_global ((dim != 3) ? 2 : 1);
+  if (true)
+    {
+      cell = tria.begin_active();
+      endc = tria.end();
+      std::vector<unsigned int> subdomain_cells(1<<dim, 0);
+      for (; cell!=endc; ++cell)
+       {
+         Assert (cell->subdomain_id() < (1<<dim), ExcInternalError());
+         ++subdomain_cells[cell->subdomain_id()];
+       };
+      for (unsigned int i=0; i<(1<<dim); ++i)
+       Assert (subdomain_cells[i] == tria.n_active_cells()/(1<<dim),
+               ExcNumberMismatch(subdomain_cells[i],
+                                 tria.n_active_cells()/(1<<dim)));
+      deallog << "Check 1 (dim=" << dim << ") ok" << std::endl;
+    };
+  
+                                  // coarsen once and check again
+  if (true)
+    {
+      cell = tria.begin_active();
+      endc = tria.end();
+      for (; cell!=endc; ++cell)
+       cell->set_coarsen_flag ();
+      tria.execute_coarsening_and_refinement ();
+      
+      cell = tria.begin_active();
+      endc = tria.end();
+      std::vector<unsigned int> subdomain_cells(1<<dim, 0);
+      for (; cell!=endc; ++cell)
+       {
+         Assert (cell->subdomain_id() < (1<<dim), ExcInternalError());
+         ++subdomain_cells[cell->subdomain_id()];
+       };
+      for (unsigned int i=0; i<(1<<dim); ++i)
+       Assert (subdomain_cells[i] == tria.n_active_cells()/(1<<dim),
+               ExcNumberMismatch(subdomain_cells[i],
+                                 tria.n_active_cells()/(1<<dim)));
+      deallog << "Check 2 (dim=" << dim << ") ok" << std::endl;
+    };
+
+                                  // check 3: assign DQ2 elements to
+                                  // cells and count degrees of
+                                  // freedom associated with each
+                                  // subdomain
+  if (true)
+    {
+      FE_DGQ<dim> fe(2);
+      DoFHandler<dim> dof_handler (tria);
+      dof_handler.distribute_dofs (fe);
+      std::vector<bool> selected_dofs (dof_handler.n_dofs());
+      for (unsigned int subdomain=0; subdomain<(1<<dim); ++subdomain)
+       {
+         DoFTools::extract_subdomain_dofs (dof_handler, subdomain,
+                                           selected_dofs);
+         Assert (static_cast<unsigned int>(std::count (selected_dofs.begin(),
+                                                       selected_dofs.end(),
+                                                       true))
+                 ==
+                 dof_handler.n_dofs() / (1<<dim),
+                 ExcNumberMismatch(std::count (selected_dofs.begin(),
+                                               selected_dofs.end(),
+                                               true),
+                                   dof_handler.n_dofs() / (1<<dim)));
+       }
+      deallog << "Check 3 (dim=" << dim << ") ok" << std::endl;      
+    };
+
+
+                                  // check 4: check again for
+                                  // continuous elements. note that
+                                  // the number of dofs here is
+                                  // different, since dofs can be on
+                                  // several subdomain at once
+  if (true)
+    {
+      FE_Q<dim> fe(1);
+      DoFHandler<dim> dof_handler (tria);
+      dof_handler.distribute_dofs (fe);
+
+      const unsigned int cells_per_direction
+       = static_cast<unsigned int>(rint(pow(tria.n_active_cells(), 1./dim)));
+      
+      std::vector<bool> selected_dofs (dof_handler.n_dofs());
+      for (unsigned int subdomain=0; subdomain<(1<<dim); ++subdomain)
+       {
+         DoFTools::extract_subdomain_dofs (dof_handler, subdomain,
+                                           selected_dofs);
+         Assert (static_cast<unsigned int>(std::count (selected_dofs.begin(),
+                                                       selected_dofs.end(),
+                                                       true))
+                 ==
+                 pow(cells_per_direction/2+1,dim),
+                 ExcNumberMismatch(std::count (selected_dofs.begin(),
+                                               selected_dofs.end(),
+                                               true),
+                                   pow(cells_per_direction/2+1,dim)));
+       }
+      deallog << "Check 4 (dim=" << dim << ") ok" << std::endl;      
+    };
+};
+
+
+int main ()
+{
+  logfile.precision(4);
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+
+  test<1> ();
+  test<2> ();
+  test<3> ();
+  
+  return 0;
+};
diff --git a/tests/deal.II/subdomain_ids.checked b/tests/deal.II/subdomain_ids.checked
new file mode 100644 (file)
index 0000000..b8cc67c
--- /dev/null
@@ -0,0 +1,13 @@
+
+DEAL::Check 1 (dim=1) ok
+DEAL::Check 2 (dim=1) ok
+DEAL::Check 3 (dim=1) ok
+DEAL::Check 4 (dim=1) ok
+DEAL::Check 1 (dim=2) ok
+DEAL::Check 2 (dim=2) ok
+DEAL::Check 3 (dim=2) ok
+DEAL::Check 4 (dim=2) ok
+DEAL::Check 1 (dim=3) ok
+DEAL::Check 2 (dim=3) ok
+DEAL::Check 3 (dim=3) ok
+DEAL::Check 4 (dim=3) ok

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.