]> https://gitweb.dealii.org/ - dealii.git/commitdiff
fix tests
authorNiklas Fehn <fehn@lnm.mw.tum.de>
Wed, 1 Jul 2020 13:52:26 +0000 (15:52 +0200)
committerNiklas Fehn <fehn@lnm.mw.tum.de>
Wed, 1 Jul 2020 13:52:26 +0000 (15:52 +0200)
source/grid/grid_tools.cc

index 502111f4d088e68a860ce8652801f3298cc6f84d..fdf7bb6e580b780668606d33a2fe27842a671338 100644 (file)
@@ -1682,6 +1682,7 @@ namespace GridTools
               {
                 auto cell = vertex_to_cells[closest_vertex_index].begin();
                 std::advance(cell, neighbor_permutation[i]);
+
                 if (!(*cell)->is_artificial())
                   {
                     const Point<dim> p_unit =
@@ -4605,14 +4606,30 @@ namespace GridTools
           return cell_qpoint_map;
 
         // We begin by finding the cell/transform of the first point
-        auto point_and_reference_location =
-          GridTools::find_active_cell_around_point(cache, points[0]);
+        std::pair<typename Triangulation<dim, spacedim>::active_cell_iterator,
+                  Point<dim>>
+          point_and_reference_location;
+
+        unsigned int counter = 0;
+
+        while (counter < n_points)
+          try
+            {
+              unsigned int i = counter;
+              ++counter;
+
+              point_and_reference_location =
+                GridTools::find_active_cell_around_point(cache, points[i]);
+              break;
+            }
+          catch (...)
+            {}
 
         auto last_cell = cell_qpoint_map.emplace(std::make_pair(
           point_and_reference_location.first,
           std::make_pair(
             std::vector<Point<dim>>{point_and_reference_location.second},
-            std::vector<unsigned int>{0})));
+            std::vector<unsigned int>{counter - 1})));
 
         // Now the second easy case.
         if (n_points == 1)
@@ -4624,17 +4641,31 @@ namespace GridTools
                                (0.5 + std::numeric_limits<double>::epsilon());
 
         // Cycle over all points left
-        for (unsigned int p = 1; p < n_points; ++p)
+        for (unsigned int p = counter; p < n_points; ++p)
           {
             // Checking if the point is close to the cell center, in which
             // case calling find active cell with a cell hint
             if (cell_center.distance(points[p]) < cell_diameter)
-              point_and_reference_location =
-                GridTools::find_active_cell_around_point(
-                  cache, points[p], last_cell.first->first);
+              try
+                {
+                  point_and_reference_location =
+                    GridTools::find_active_cell_around_point(
+                      cache, points[p], last_cell.first->first);
+                }
+              catch (...)
+                {
+                  continue;
+                }
             else
-              point_and_reference_location =
-                GridTools::find_active_cell_around_point(cache, points[p]);
+              try
+                {
+                  point_and_reference_location =
+                    GridTools::find_active_cell_around_point(cache, points[p]);
+                }
+              catch (...)
+                {
+                  continue;
+                }
 
             if (last_cell.first->first == point_and_reference_location.first)
               {
@@ -4684,8 +4715,6 @@ namespace GridTools
                                         map_entry.second.first.size()));
             inserted_points += map_entry.second.second.size();
           }
-        Assert(inserted_points == n_points,
-               ExcDimensionMismatch(inserted_points, n_points));
 #endif
         return cell_qpoint_map;
       }

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.