Mesh preparation fixes

examples/vector_fe/vector_fe_ex10/vector_fe_ex10.C

@@ -150,6 +150,9 @@ int main (int argc, char ** argv)
   const Real phi = infile("phi", 0.), theta = infile("theta", 0.), psi = infile("psi", 0.);
   GradDivExactSolution::RM(MeshTools::Modification::rotate(mesh, phi, theta, psi));
 
+  // Rotation can leave a mesh's caches unprepared
+  mesh.complete_preparation();
+
   // Print information about the mesh to the screen.
   mesh.print_info();
 

examples/vector_fe/vector_fe_ex3/vector_fe_ex3.C

@@ -113,6 +113,9 @@ int main (int argc, char ** argv)
   const Real phi = infile("phi", 0.), theta = infile("theta", 0.), psi = infile("psi", 0.);
   CurlCurlExactSolution::RM(MeshTools::Modification::rotate(mesh, phi, theta, psi));
 
+  // Rotation can leave a mesh's caches unprepared
+  mesh.complete_preparation();
+
   // Print information about the mesh to the screen.
   mesh.print_info();
 

src/mesh/boundary_info.C

@@ -2185,6 +2185,8 @@ void BoundaryInfo::renumber_id (boundary_id_type old_id,
     {
       _boundary_ids.erase(old_id);
       _boundary_ids.insert(new_id);
+      _global_boundary_ids.erase(old_id);
+      _global_boundary_ids.insert(new_id);
     }
 
   renumber_name(_ss_id_to_name, old_id, new_id);
@@ -2224,8 +2226,10 @@ void BoundaryInfo::renumber_side_id (boundary_id_type old_id,
           !_node_boundary_ids.count(old_id))
         {
           _boundary_ids.erase(old_id);
+          _global_boundary_ids.erase(old_id);
         }
       _boundary_ids.insert(new_id);
+      _global_boundary_ids.insert(new_id);
     }
 
   renumber_name(_ss_id_to_name, old_id, new_id);
@@ -2263,8 +2267,10 @@ void BoundaryInfo::renumber_edge_id (boundary_id_type old_id,
           !_node_boundary_ids.count(old_id))
         {
           _boundary_ids.erase(old_id);
+          _global_boundary_ids.erase(old_id);
         }
       _boundary_ids.insert(new_id);
+      _global_boundary_ids.insert(new_id);
     }
 
   renumber_name(_es_id_to_name, old_id, new_id);
@@ -2302,8 +2308,10 @@ void BoundaryInfo::renumber_shellface_id (boundary_id_type old_id,
           !_node_boundary_ids.count(old_id))
         {
           _boundary_ids.erase(old_id);
+          _global_boundary_ids.erase(old_id);
         }
       _boundary_ids.insert(new_id);
+      _global_boundary_ids.insert(new_id);
     }
 
   this->libmesh_assert_valid_multimaps();
@@ -2339,8 +2347,10 @@ void BoundaryInfo::renumber_node_id (boundary_id_type old_id,
           !_edge_boundary_ids.count(old_id))
         {
           _boundary_ids.erase(old_id);
+          _global_boundary_ids.erase(old_id);
         }
       _boundary_ids.insert(new_id);
+      _global_boundary_ids.insert(new_id);
     }
 
   renumber_name(_ns_id_to_name, old_id, new_id);

src/mesh/mesh_base.C

@@ -2072,19 +2072,31 @@ void MeshBase::detect_interior_parents()
   // This requires an inspection on every processor
   parallel_object_only();
 
-  // This requires up-to-date mesh dimensions in cache
-  libmesh_assert(_preparation.has_cached_elem_data);
+  // This requires up-to-date mesh dimensions, but if we don't have
+  // them cached then we can't update them without changing the mesh
+  // in unexpected ways that interfere with our tests of
+  // partially-prepared meshes in MeshTools::*valid_is_prepared
+  std::set<unsigned char> elem_dims_copy;
+  if (_preparation.has_cached_elem_data)
+    elem_dims_copy = this->elem_dimensions();
+  else
+    {
+      for (const auto & elem : this->active_element_ptr_range())
+        elem_dims_copy.insert(cast_int<unsigned char>(elem->dim()));
+      if (!this->is_serial())
+        this->comm().set_union(elem_dims_copy);
+    }
 
   // Early return if the mesh is empty or has elements of a single spatial dimension.
-  if (this->elem_dimensions().size() <= 1)
+  if (elem_dims_copy.size() <= 1)
     {
       _preparation.has_interior_parent_ptrs = true;
       return;
     }
 
   // Convenient elem_dimensions iterators
-  const auto dim_start = this->elem_dimensions().begin();
-  const auto dim_end = this->elem_dimensions().end();
+  const auto dim_start = elem_dims_copy.begin();
+  const auto dim_end = elem_dims_copy.end();
 
   // In this function we find only +1 dimensional interior parents,
   // (so, for a given element el, the interior parent p must satisfy p.dim() == el.dim() + 1).

src/mesh/mesh_modification.C

@@ -314,6 +314,10 @@ void MeshTools::Modification::redistribute (MeshBase & mesh,
 #endif
     }
 
+  // If we just moved a mesh in or out out of the X axis or XY plane
+  // then we might have changed its spatial_dimension()
+  mesh.unset_has_cached_elem_data();
+
   // We haven't changed any topology, but just changing geometry could
   // have invalidated a point locator.
   mesh.clear_point_locator();
@@ -331,6 +335,10 @@ void MeshTools::Modification::translate (MeshBase & mesh,
   for (auto & node : mesh.node_ptr_range())
     *node += p;
 
+  // If we just moved a mesh in or out out of the X axis or XY plane
+  // then we might have changed its spatial_dimension()
+  mesh.unset_has_cached_elem_data();
+
   // We haven't changed any topology, but just changing geometry could
   // have invalidated a point locator.
   mesh.clear_point_locator();
@@ -372,15 +380,16 @@ MeshTools::Modification::rotate (MeshBase & mesh,
 #if LIBMESH_DIM == 3
   const auto R = RealTensorValue::intrinsic_rotation_matrix(phi, theta, psi);
 
-  if (theta)
-    mesh.set_spatial_dimension(3);
-
   for (auto & node : mesh.node_ptr_range())
     {
       Point & pt = *node;
       pt = R * pt;
     }
 
+  // If we just moved a mesh in or out out of the X axis or XY plane
+  // then we might have changed its spatial_dimension()
+  mesh.unset_has_cached_elem_data();
+
   return R;
 
 #else
@@ -426,6 +435,10 @@ void MeshTools::Modification::scale (MeshBase & mesh,
   for (auto & node : mesh.node_ptr_range())
     (*node)(2) *= z_scale;
 
+  // If we just collapsed a manifold onto the X axis or XY plane
+  // then we might have changed its spatial_dimension()
+  mesh.unset_has_cached_elem_data();
+
   // We haven't changed any topology, but just changing geometry could
   // have invalidated a point locator.
   mesh.clear_point_locator();

tests/fe/fe_test.h

@@ -530,6 +530,7 @@ class FETestBase : public CppUnit::TestCase {
           }
       }
 
+    _mesh->complete_preparation();
     _es = std::make_unique<EquationSystems>(*_mesh);
     _sys = &(_es->add_system<System> ("SimpleSystem"));
     _sys->add_variable("u", order, family);