]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Use geometry information instead of hand-crafted values.
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Wed, 6 May 1998 14:10:07 +0000 (14:10 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Wed, 6 May 1998 14:10:07 +0000 (14:10 +0000)
git-svn-id: https://svn.dealii.org/trunk@255 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/source/dofs/dof_handler.cc
deal.II/deal.II/source/grid/tria.cc
deal.II/deal.II/source/grid/tria_accessor.cc
deal.II/deal.II/source/numerics/data_io.cc

index 3fa5fda2a755ddd113a2dafb0c2540f93f4ace04..1f325e92709d2dacef0e0c71a79010cce72448f6 100644 (file)
@@ -6,6 +6,7 @@
 #include <grid/tria_accessor.h>
 #include <grid/tria_iterator.h>
 #include <grid/tria.h>
+#include <grid/geometry_info.h>
 #include <fe/fe.h>
 #include <lac/dsmatrix.h>
 #include <map>
@@ -665,7 +666,7 @@ int DoFHandler<1>::distribute_dofs_on_cell (active_cell_iterator   &cell,
                                            unsigned int            next_free_dof) {
 
                                   // distribute dofs of vertices
-  for (unsigned int v=0; v<2; ++v)
+  for (unsigned int v=0; v<GeometryInfo<1>::vertices_per_cell; ++v)
     {
       cell_iterator neighbor = cell->neighbor(v);
 
@@ -716,7 +717,7 @@ int DoFHandler<2>::distribute_dofs_on_cell (active_cell_iterator   &cell,
                                            unsigned int            next_free_dof) {
   if (selected_fe->dofs_per_vertex > 0)
                                     // number dofs on vertices
-    for (unsigned int vertex=0; vertex<4; ++vertex)
+    for (unsigned int vertex=0; vertex<GeometryInfo<2>::vertices_per_cell; ++vertex)
                                       // check whether dofs for this
                                       // vertex have been distributed
                                       // (only check the first dof)
@@ -725,7 +726,7 @@ int DoFHandler<2>::distribute_dofs_on_cell (active_cell_iterator   &cell,
          cell->set_vertex_dof_index (vertex, d, next_free_dof++);
     
                                   // for the four sides
-  for (unsigned int side=0; side<4; ++side)
+  for (unsigned int side=0; side<GeometryInfo<2>::faces_per_cell; ++side)
     {
       line_iterator line = cell->line(side);
 
@@ -1105,11 +1106,11 @@ void DoFHandler<2>::make_sparsity_pattern (dSMatrixStruct &sparsity) const {
       int dof_number=0;
 
                                       // fill dof indices on vertices
-      for (unsigned int vertex=0; vertex<4; ++vertex)
+      for (unsigned int vertex=0; vertex<GeometryInfo<2>::vertices_per_cell; ++vertex)
        for (unsigned int d=0; d<selected_fe->dofs_per_vertex; ++d)
          dofs_on_this_cell[dof_number++] = cell->vertex_dof_index (vertex,d);
 
-      for (unsigned int line=0; line<4; ++line)
+      for (unsigned int line=0; line<GeometryInfo<2>::faces_per_cell; ++line)
        for (unsigned int d=0; d<selected_fe->dofs_per_line; ++d)
          dofs_on_this_cell[dof_number++] = cell->line(line)->dof_index (d);
       
@@ -1139,7 +1140,7 @@ void DoFHandler<dim>::make_transfer_matrix (const DoFHandler<dim> &transfer_from
                                   // assert for once at the beginning the
                                   // the matrices have the correct sizes
 #ifdef DEBUG
-  for (unsigned int c=0; c<(1<<dim); ++c)
+  for (unsigned int c=0; c<GeometryInfo<dim>::children_per_cell; ++c)
     {
       Assert (selected_fe->prolongate(c).m() == selected_fe->total_dofs,
              ExcMatrixHasWrongSize(selected_fe->prolongate(c).m()));
@@ -1185,7 +1186,7 @@ void DoFHandler<dim>::transfer_cell (const typename DoFHandler<dim>::cell_iterat
                                     dSMatrixStruct      &transfer_pattern) const {
   if (!new_cell->active() && !old_cell->active())
                                     // both cells have children; go deeper
-    for (unsigned int child=0; child<(1<<dim); ++child)
+    for (unsigned int child=0; child<GeometryInfo<dim>::children_per_cell; ++child)
       transfer_cell (old_cell->child(child), new_cell->child(child),
                     transfer_pattern);
   else
@@ -1208,8 +1209,8 @@ void DoFHandler<dim>::transfer_cell (const typename DoFHandler<dim>::cell_iterat
       if (!new_cell->active() && old_cell->active())
                                         // new cell has children, old one has not
        {
-         cell_iterator child[1<<dim];
-         for (unsigned int c=0; c<(1<<dim); ++c) 
+         cell_iterator child[GeometryInfo<dim>::children_per_cell];
+         for (unsigned int c=0; c<GeometryInfo<dim>::children_per_cell; ++c) 
            {
              child[c] = new_cell->child(c);
              Assert (child[c]->active(),
@@ -1223,7 +1224,7 @@ void DoFHandler<dim>::transfer_cell (const typename DoFHandler<dim>::cell_iterat
          Assert (old_dof_indices.size() == selected_fe->total_dofs,
                  ExcInternalError ());
 
-         for (unsigned int c=0; c<(1<<dim); ++c)
+         for (unsigned int c=0; c<GeometryInfo<dim>::children_per_cell; ++c)
            {
                                               // numbers of child dofs
              vector<int> child_dof_indices;
@@ -1240,8 +1241,8 @@ void DoFHandler<dim>::transfer_cell (const typename DoFHandler<dim>::cell_iterat
            };
        } else {
                                           // old cell has children, new one has not
-         cell_iterator child[1<<dim];
-         for (unsigned int c=0; c<(1<<dim); ++c) 
+         cell_iterator child[GeometryInfo<dim>::children_per_cell];
+         for (unsigned int c=0; c<GeometryInfo<dim>::children_per_cell; ++c) 
            {
              child[c] = old_cell->child(c);
              Assert (child[c]->active(),
@@ -1255,7 +1256,7 @@ void DoFHandler<dim>::transfer_cell (const typename DoFHandler<dim>::cell_iterat
          Assert (new_dof_indices.size() == selected_fe->total_dofs,
                  ExcInternalError ());
          
-         for (unsigned int c=0; c<(1<<dim); ++c)
+         for (unsigned int c=0; c<GeometryInfo<dim>::children_per_cell; ++c)
            {
                                               // numbers of child dofs
              vector<int> child_dof_indices;
@@ -1281,7 +1282,7 @@ void DoFHandler<dim>::transfer_cell (const typename DoFHandler<dim>::cell_iterat
                                     dSMatrix            &transfer_matrix) const {
   if (!new_cell->active() && !old_cell->active())
                                     // both cells have children; go deeper
-    for (unsigned int child=0; child<(1<<dim); ++child)
+    for (unsigned int child=0; child<GeometryInfo<dim>::children_per_cell; ++child)
       transfer_cell (old_cell->child(child), new_cell->child(child),
                     transfer_matrix);
   else
@@ -1304,8 +1305,8 @@ void DoFHandler<dim>::transfer_cell (const typename DoFHandler<dim>::cell_iterat
       if (!new_cell->active() && old_cell->active())
                                         // new cell has children, old one has not
        {
-         cell_iterator child[1<<dim];
-         for (unsigned int c=0; c<(1<<dim); ++c) 
+         cell_iterator child[GeometryInfo<dim>::children_per_cell];
+         for (unsigned int c=0; c<GeometryInfo<dim>::children_per_cell; ++c) 
            {
              child[c] = new_cell->child(c);
              Assert (child[c]->active(),
@@ -1319,7 +1320,7 @@ void DoFHandler<dim>::transfer_cell (const typename DoFHandler<dim>::cell_iterat
          Assert (old_dof_indices.size() == selected_fe->total_dofs,
                  ExcInternalError ());
 
-         for (unsigned int c=0; c<(1<<dim); ++c)
+         for (unsigned int c=0; c<GeometryInfo<dim>::children_per_cell; ++c)
            {
                                               // numbers of child dofs
              vector<int> child_dof_indices;
@@ -1337,8 +1338,8 @@ void DoFHandler<dim>::transfer_cell (const typename DoFHandler<dim>::cell_iterat
            };
        } else {
                                           // old cell has children, new one has not
-         cell_iterator child[1<<dim];
-         for (unsigned int c=0; c<(1<<dim); ++c) 
+         cell_iterator child[GeometryInfo<dim>::children_per_cell];
+         for (unsigned int c=0; c<GeometryInfo<dim>::children_per_cell; ++c) 
            {
              child[c] = old_cell->child(c);
              Assert (child[c]->active(),
@@ -1352,7 +1353,7 @@ void DoFHandler<dim>::transfer_cell (const typename DoFHandler<dim>::cell_iterat
          Assert (new_dof_indices.size() == selected_fe->total_dofs,
                  ExcInternalError ());
          
-         for (unsigned int c=0; c<(1<<dim); ++c)
+         for (unsigned int c=0; c<GeometryInfo<dim>::children_per_cell; ++c)
            {
                                               // numbers of child dofs
              vector<int> child_dof_indices;
index acc8ae2a49ac2036ff453245e31a94b601ee26ce..0a2daf6ca7154bef72eeb3078567d0ee57f7a749 100644 (file)
@@ -3,6 +3,7 @@
 #include <grid/tria_boundary.h>
 #include <grid/tria_iterator.h>
 #include <grid/tria.h>
+#include <grid/geometry_info.h>
 #include <basic/magic_numbers.h>
 #include <basic/data_io.h>
 #include <lac/dvector.h>
@@ -109,7 +110,7 @@ void Triangulation<1>::create_triangulation (const vector<Point<1> >    &v,
                                   // for all lines
   for (; line!=end(); ++line)
                                     // for each of the two vertices
-    for (unsigned int vertex=0; vertex<(1<<dim); ++vertex)
+    for (unsigned int vertex=0; vertex<GeometryInfo<dim>::vertices_per_cell; ++vertex)
                                       // if first cell adjacent to this
                                       // vertex is the present one, then
                                       // the neighbor is the second adjacent
@@ -210,10 +211,14 @@ void Triangulation<2>::create_triangulation (const vector<Point<2> >    &v,
            swap (cells[cell].vertices[0], cells[cell].vertices[2]);
            swap (cells[cell].vertices[1], cells[cell].vertices[3]);
                                             // remake lines
-           line_vertices[0] = make_pair (cells[cell].vertices[0], cells[cell].vertices[1]);
-           line_vertices[1] = make_pair (cells[cell].vertices[1], cells[cell].vertices[2]);
-           line_vertices[2] = make_pair (cells[cell].vertices[0], cells[cell].vertices[3]);
-           line_vertices[3] = make_pair (cells[cell].vertices[3], cells[cell].vertices[2]);
+           line_vertices[0]
+             = make_pair (cells[cell].vertices[0], cells[cell].vertices[1]);
+           line_vertices[1]
+             = make_pair (cells[cell].vertices[1], cells[cell].vertices[2]);
+           line_vertices[2]
+             = make_pair (cells[cell].vertices[0], cells[cell].vertices[3]);
+           line_vertices[3]
+             = make_pair (cells[cell].vertices[3], cells[cell].vertices[2]);
                                             // allow for only one such
                                             // rotation
            break;
@@ -1654,7 +1659,7 @@ unsigned int Triangulation<dim>::max_adjacent_cells () const {
                                   // used on level 0
   unsigned int max_vertex_index = 0;
   for (; cell!=endc; ++cell)
-    for (unsigned vertex=0; vertex<(1<<dim); ++vertex)
+    for (unsigned vertex=0; vertex<GeometryInfo<dim>::vertices_per_cell; ++vertex)
       if (cell->vertex_index(vertex) > (signed int)max_vertex_index)
        max_vertex_index = cell->vertex_index(vertex);
 
@@ -1665,10 +1670,10 @@ unsigned int Triangulation<dim>::max_adjacent_cells () const {
                                   // touch a vertex's usage count everytime
                                   // we find an adjacent element
   for (cell=begin(); cell!=endc; ++cell)
-    for (unsigned vertex=0; vertex<(1<<dim); ++vertex)
+    for (unsigned vertex=0; vertex<GeometryInfo<dim>::vertices_per_cell; ++vertex)
       ++usage_count[cell->vertex_index(vertex)];
 
-  return max (static_cast<unsigned int>(1<<dim),
+  return max (GeometryInfo<dim>::vertices_per_cell,
              static_cast<unsigned int>(*(max_element (usage_count.begin(),
                                                       usage_count.end()))));
 };
@@ -1701,7 +1706,7 @@ void Triangulation<dim>::print_gnuplot (ostream &out, const unsigned int level)
 
 
 void Triangulation<1>::print_gnuplot (ostream &out,
-                                     const TriaDimensionInfo<1>::active_cell_iterator & cell) const {
+                                     const active_cell_iterator & cell) const {
   out << cell->vertex(0) << " " << cell->level() << endl
       << cell->vertex(1) << " " << cell->level() << endl
       << endl;
@@ -1710,7 +1715,7 @@ void Triangulation<1>::print_gnuplot (ostream &out,
 
 
 void Triangulation<2>::print_gnuplot (ostream &out,
-                                     const TriaDimensionInfo<2>::active_cell_iterator & cell) const {
+                                     const active_cell_iterator & cell) const {
   out << cell->vertex(0) << " " << cell->level() << endl
       << cell->vertex(1) << " " << cell->level() << endl
       << cell->vertex(2) << " " << cell->level() << endl
@@ -1792,12 +1797,14 @@ void Triangulation<1>::execute_refinement () {
                                       // 2*flagged_cells that will be created
                                       // on that level
       levels[level+1]->
-       TriangulationLevel<0>::reserve_space (used_cells+2*flagged_cells, 1);
+       TriangulationLevel<0>::reserve_space (used_cells+
+                                             GeometryInfo<1>::children_per_cell *
+                                             flagged_cells, 1);
                                       // reserve space for 2*flagged_cells
                                       // new lines on the next higher
                                       // level
       levels[level+1]->
-       TriangulationLevel<1>::reserve_space (2*flagged_cells);
+       TriangulationLevel<1>::reserve_space (GeometryInfo<1>::children_per_cell*flagged_cells);
       
       needed_vertices += flagged_cells;
     };
@@ -2517,7 +2524,7 @@ void Triangulation<dim>::prepare_refinement () {
       for (; cell != endc; --cell)
        if (cell->refine_flag_set() == true)
                                           // loop over neighbors of cell
-         for (unsigned int i=0; i<(2*dim); ++i)
+         for (unsigned int i=0; i<GeometryInfo<dim>::faces_per_cell; ++i)
            if (cell->neighbor(i).state() == valid)
              if ((cell->neighbor_level(i) == cell->level()-1)
                  &&
index cd609f7d7c19333227fe3f6b156198177e05a3ca..0992b5c01771a99cd2fb4efaedcd22c69ea43758 100644 (file)
@@ -4,6 +4,7 @@
 #include <grid/tria_iterator.h>
 #include <grid/tria_accessor.h>
 #include <grid/tria_iterator.templates.h>
+#include <grid/geometry_info.h>
 
 /* Note: explicit instantiations at the end of the different sections!       */
 
@@ -595,18 +596,20 @@ CellAccessor<2>::face (const unsigned int i) const {
 
 template <int dim>
 int CellAccessor<dim>::neighbor_index (const unsigned int i) const {
-  Assert (i<2*dim,
+  Assert (i<GeometryInfo<dim>::faces_per_cell,
          typename TriaSubstructAccessor<dim>::ExcInvalidNeighbor(i));
-  return tria->levels[present_level]->neighbors[present_index*2*dim+i].second;
+  return tria->levels[present_level]->
+    neighbors[present_index*GeometryInfo<dim>::faces_per_cell+i].second;
 };
 
 
 
 template <int dim>
 int CellAccessor<dim>::neighbor_level (const unsigned int i) const {
-  Assert (i<2*dim,
+  Assert (i<GeometryInfo<dim>::faces_per_cell,
          typename TriaSubstructAccessor<dim>::ExcInvalidNeighbor(i));
-  return tria->levels[present_level]->neighbors[present_index*2*dim+i].first;
+  return tria->levels[present_level]->
+    neighbors[present_index*GeometryInfo<dim>::faces_per_cell+i].first;
 };
 
 
@@ -614,10 +617,13 @@ int CellAccessor<dim>::neighbor_level (const unsigned int i) const {
 template <int dim>
 void CellAccessor<dim>::set_neighbor (const unsigned int i,
                                      const TriaIterator<dim,CellAccessor<dim> > &pointer) const {
-  Assert (i<2*dim, typename TriaSubstructAccessor<dim>::ExcInvalidNeighbor(i));
-  tria->levels[present_level]->neighbors[present_index*2*dim+i].first
+  Assert (i<GeometryInfo<dim>::faces_per_cell,
+         typename TriaSubstructAccessor<dim>::ExcInvalidNeighbor(i));
+  tria->levels[present_level]->
+    neighbors[present_index*GeometryInfo<dim>::faces_per_cell+i].first
     = pointer.accessor.present_level;
-  tria->levels[present_level]->neighbors[present_index*2*dim+i].second
+  tria->levels[present_level]->
+    neighbors[present_index*GeometryInfo<dim>::faces_per_cell+i].second
     = pointer.accessor.present_index;
 };
 
@@ -626,8 +632,9 @@ void CellAccessor<dim>::set_neighbor (const unsigned int i,
 template <int dim>
 bool CellAccessor<dim>::at_boundary (const unsigned int i) const {
   Assert (used(), typename TriaSubstructAccessor<dim>::ExcCellNotUsed());
-  Assert (i<2*dim,
-         typename TriaSubstructAccessor<dim>::ExcInvalidIndex (i,0,2*dim-1));
+  Assert (i<GeometryInfo<dim>::faces_per_cell,
+         typename TriaSubstructAccessor<dim>::
+         ExcInvalidIndex (i,0,GeometryInfo<dim>::faces_per_cell-1));
   
   return (neighbor_index(i) == -1);
 };
index c5194bb06d8970403c07e49899934b24ed952117..a6531dd25344b8aecb502e2a4f0bbee683455b1b 100644 (file)
@@ -5,6 +5,7 @@
 #include <grid/dof_accessor.h>
 #include <grid/tria_iterator.h>
 #include <grid/tria.h>
+#include <grid/geometry_info.h>
 #include <fe/fe.h>
 
 #include <map>
@@ -119,13 +120,13 @@ void DataIn<dim>::read_ucd (istream &in) {
        {
                                           // allocate and read indices
          cells.push_back (CellData<dim>());
-         for (unsigned int i=0; i<(1<<dim); ++i)
+         for (unsigned int i=0; i<GeometryInfo<dim>::vertices_per_cell; ++i)
            in >> cells.back().vertices[i];
          cells.back().material_id = material_id;
 
                                           // transform from ucd to
                                           // consecutive numbering
-         for (unsigned int i=0; i<(1<<dim); ++i)
+         for (unsigned int i=0; i<GeometryInfo<dim>::vertices_per_cell; ++i)
            if (vertex_indices.find (cells.back().vertices[i]) != vertex_indices.end())
                                               // vertex with this index exists
              cells.back().vertices[i] = vertex_indices[cells.back().vertices[i]];
@@ -246,7 +247,8 @@ void DataOut<dim>::write_ucd (ostream &out) const {
                                       // use
       vector<bool> is_vertex_dof (dofs->n_dofs(), false);
       for (cell=dofs->begin_active(); cell!=endc; ++cell)
-       for (unsigned int vertex=0; vertex<(1<<dim); ++vertex) 
+       for (unsigned int vertex=0; vertex<GeometryInfo<dim>::vertices_per_cell;
+            ++vertex) 
          is_vertex_dof[cell->vertex_dof_index(vertex,0)] = true;
 
       n_vertex_dofs = 0;
@@ -295,7 +297,8 @@ void DataOut<dim>::write_ucd (ostream &out) const {
       vector<bool> already_written (dofs->n_dofs(), false);
                                       // write vertices
       for (cell=dofs->begin_active(); cell!=endc; ++cell)
-       for (unsigned int vertex=0; vertex<(1<<dim); ++vertex) 
+       for (unsigned int vertex=0; vertex<GeometryInfo<dim>::vertices_per_cell;
+            ++vertex) 
          if (already_written[cell->vertex_dof_index(vertex,0)]==false)
            {
              out << cell->vertex_dof_index(vertex,0) // vertex index
@@ -329,7 +332,8 @@ void DataOut<dim>::write_ucd (ostream &out) const {
              default:
                    Assert (false, ExcNotImplemented());
            };
-         for (unsigned int vertex=0; vertex<(1<<dim); ++vertex)
+         for (unsigned int vertex=0; vertex<GeometryInfo<dim>::vertices_per_cell;
+              ++vertex)
            out << cell->vertex_dof_index(vertex,0) << ' ';
          out << endl;
        };
@@ -360,7 +364,8 @@ void DataOut<dim>::write_ucd (ostream &out) const {
       vector<bool> already_written (dofs->n_dofs(), false);
                                       // write vertices
       for (cell=dofs->begin_active(); cell!=endc; ++cell)
-       for (unsigned int vertex=0; vertex<(1<<dim); ++vertex) 
+       for (unsigned int vertex=0; vertex<GeometryInfo<dim>::vertices_per_cell;
+            ++vertex) 
          if (already_written[cell->vertex_dof_index(vertex,0)]==false)
            {
              out << cell->vertex_dof_index(vertex,0) // vertex index
@@ -429,7 +434,7 @@ void DataOut<dim>::write_ucd_faces (ostream &out,
            default:
                  Assert (false, ExcNotImplemented());
          };
-       for (unsigned int vertex=0; vertex<(1<<(dim-1)); ++vertex)
+       for (unsigned int vertex=0; vertex<GeometryInfo<dim>::vertices_per_face; ++vertex)
          out << face->vertex_dof_index(vertex,0) << ' ';
        out << endl;
 

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.