]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Convert vertex numbers in CellData to unsigned integers. Fix all the consequences.
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Sun, 10 Sep 2006 01:20:12 +0000 (01:20 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Sun, 10 Sep 2006 01:20:12 +0000 (01:20 +0000)
git-svn-id: https://svn.dealii.org/trunk@13877 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/grid/tria.h
deal.II/deal.II/source/grid/grid_in.cc
deal.II/deal.II/source/grid/grid_reordering.cc
deal.II/deal.II/source/grid/tria.cc
deal.II/deal.II/source/grid/tria_accessor.cc
deal.II/doc/news/changes.html

index 033a1430d39d1af07bc967e7001cb23a63f6bf75..f10f0335c21749414d8d44277471ee02a2617061 100644 (file)
@@ -58,7 +58,7 @@ struct CellData
                                      * Indices of the vertices of
                                      * this cell.
                                      */
-    int           vertices[GeometryInfo<dim>::vertices_per_cell];
+    unsigned int vertices[GeometryInfo<dim>::vertices_per_cell];
 
                                     /**
                                      * Material indicator of this
index 7d4b6ade0412fa105188f62380176f0632a3ca52..ee28d683c720f35c3dedc326c5e7271365abf570 100644 (file)
@@ -126,7 +126,7 @@ void GridIn<dim>::read_ucd (std::istream &in)
              {
                                                 // no such vertex index
                AssertThrow (false, ExcInvalidVertexIndex(cell, cells.back().vertices[i]));
-               cells.back().vertices[i] = -1;
+               cells.back().vertices[i] = deal_II_numbers::invalid_unsigned_int;
              };
        }
       else
@@ -152,7 +152,8 @@ void GridIn<dim>::read_ucd (std::istream &in)
                  AssertThrow (false,
                               ExcInvalidVertexIndex(cell,
                                                     subcelldata.boundary_lines.back().vertices[i]));
-                 subcelldata.boundary_lines.back().vertices[i] = -1;
+                 subcelldata.boundary_lines.back().vertices[i]
+                   = deal_II_numbers::invalid_unsigned_int;
                };
          }
        else
@@ -181,7 +182,8 @@ void GridIn<dim>::read_ucd (std::istream &in)
                    Assert (false,
                            ExcInvalidVertexIndex(cell,
                                                  subcelldata.boundary_quads.back().vertices[i]));
-                   subcelldata.boundary_quads.back().vertices[i] = -1;
+                   subcelldata.boundary_quads.back().vertices[i] =
+                     deal_II_numbers::invalid_unsigned_int;
                  };
              
            }
@@ -660,7 +662,7 @@ void GridIn<dim>::read_msh (std::istream &in)
              {
                                                 // no such vertex index
                AssertThrow (false, ExcInvalidVertexIndex(cell, cells.back().vertices[i]));
-               cells.back().vertices[i] = -1;
+               cells.back().vertices[i] = deal_II_numbers::invalid_unsigned_int;
              };
        }
       else
@@ -686,7 +688,8 @@ void GridIn<dim>::read_msh (std::istream &in)
                  AssertThrow (false,
                               ExcInvalidVertexIndex(cell,
                                                     subcelldata.boundary_lines.back().vertices[i]));
-                 subcelldata.boundary_lines.back().vertices[i] = -1;
+                 subcelldata.boundary_lines.back().vertices[i] =
+                   deal_II_numbers::invalid_unsigned_int;
                };
          }
        else
@@ -715,7 +718,8 @@ void GridIn<dim>::read_msh (std::istream &in)
                    Assert (false,
                            ExcInvalidVertexIndex(cell,
                                                  subcelldata.boundary_quads.back().vertices[i]));
-                   subcelldata.boundary_quads.back().vertices[i] = -1;
+                   subcelldata.boundary_quads.back().vertices[i] =
+                     deal_II_numbers::invalid_unsigned_int;
                  };
              
            }
index 1ffc1ca3e98a0a888331cfb8c155e1cead3ce651..0fcb4d0a9f33970287a2f7fc12d4212017fd9b87 100644 (file)
@@ -1451,7 +1451,7 @@ GridReordering<3>::invert_all_cells_of_negative_grid(
   const std::vector<Point<3> > &all_vertices,
   std::vector<CellData<3> > &cells)
 {
-  int vertices_lex[GeometryInfo<3>::vertices_per_cell];
+  unsigned int vertices_lex[GeometryInfo<3>::vertices_per_cell];
   unsigned int n_negative_cells=0;
   for (unsigned int cell_no=0; cell_no<cells.size(); ++cell_no)
     {
index 7c2c5625d95622b42ea1337675e63b78242407ab..e0b0d2c30290366230fb9d0cb6519f1f9ec16159 100644 (file)
@@ -480,8 +480,7 @@ void Triangulation<2>::create_triangulation (const std::vector<Point<2> >    &v,
   for (unsigned int cell=0; cell<cells.size(); ++cell)
     {
       for (unsigned int vertex=0; vertex<4; ++vertex)
-       if ( ! ((0<=cells[cell].vertices[vertex]) &&
-               (cells[cell].vertices[vertex]<static_cast<signed int>(vertices.size()))))
+       if ( ! (cells[cell].vertices[vertex]<vertices.size()))
          {
                                              // store the number of
                                              // vertices, as this
@@ -926,8 +925,7 @@ Triangulation<3>::create_triangulation (const std::vector<Point<3> >    &v,
                                       // check whether vertex indices
                                       // are valid ones
       for (unsigned int vertex=0; vertex<GeometryInfo<dim>::vertices_per_cell; ++vertex)
-       if (! ((0<=cells[cell].vertices[vertex]) &&
-              (cells[cell].vertices[vertex]<static_cast<signed int>(vertices.size()))))
+       if (! (cells[cell].vertices[vertex] < vertices.size()))
          {
                                              // store the number of
                                              // vertices, as this
index 430e97833cc120281b59c86c88357b9b63832e4d..7c447b85601eb4c3ada988499b95ec87129bdc8b 100644 (file)
@@ -1595,7 +1595,7 @@ Point<3> TriaObjectAccessor<3, 3>::barycenter () const
 template <>
 double TriaObjectAccessor<3, 3>::measure () const
 {
-  static int vertex_indices[GeometryInfo<3>::vertices_per_cell];
+  static unsigned int vertex_indices[GeometryInfo<3>::vertices_per_cell];
   for (unsigned int i=0; i<GeometryInfo<3>::vertices_per_cell; ++i)
     vertex_indices[i]=vertex_index(i);
 
index 841131cf8fc7b4343d902e7bc48ebabdffbdd304..2fed0e0fed9849a369c1689dbdff6505d70b32a3 100644 (file)
@@ -44,6 +44,14 @@ inconvenience this causes.
 
 
 <ol>
+  <li> <p>
+       Changed: Indices, such as vertex indices, are usually represented by
+       unsigned integers in deal.II. One place where we didn't do this was in
+       the <code>CellData</code> structure that can be used to describe cells
+       when building a triangulation. This has now been rectified. 
+       <br> 
+       (WB 2006/09/06)
+       </p>
 
   <li> <p>
        Changed: Lower dimensional objects have been removed from the

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.