Merge pull request #1205 from DLR-AMR/enhancement-geometry_checking

Add geometry checking to cmesh commit

example/geometry/t8_example_geometries.cxx

@@ -139,6 +139,22 @@ struct t8_geometry_sincos: public t8_geometry
     return 0;
   }
 
+  /**
+   * Check for compatibility of the currently loaded tree with the geometry.
+   * Only quad elements are supported by this geometry.
+   */
+  bool
+  t8_geom_check_tree_compatibility () const
+  {
+    if (active_tree_class != T8_ECLASS_QUAD) {
+      t8_productionf (
+        "t8_geometry_sincos is not compatible with tree type %s\n It is only compatible with quad elements.\n",
+        t8_eclass_to_string[active_tree_class]);
+      return false;
+    }
+    return true;
+  }
+
   /**
    * Get the type of this geometry.
    * \return The type.
@@ -203,6 +219,22 @@ struct t8_geometry_moebius: public t8_geometry_with_vertices
     SC_ABORT_NOT_REACHED ();
   }
 
+  /**
+   * Check for compatibility of the currently loaded tree with the geometry.
+   * Only quad elements are supported by this geometry.
+   */
+  bool
+  t8_geom_check_tree_compatibility () const
+  {
+    if (active_tree_class != T8_ECLASS_QUAD) {
+      t8_productionf (
+        "t8_geometry_moebius is not compatible with tree type %s\n It is only compatible with quad elements.\n",
+        t8_eclass_to_string[active_tree_class]);
+      return false;
+    }
+    return true;
+  }
+
   /**
    * Get the type of this geometry.
    * \return The type.
@@ -273,6 +305,22 @@ struct t8_geometry_cylinder: public t8_geometry
     return 0;
   }
 
+  /**
+   * Check for compatibility of the currently loaded tree with the geometry.
+   * Only quad elements are supported by this geometry.
+   */
+  bool
+  t8_geom_check_tree_compatibility () const
+  {
+    if (active_tree_class != T8_ECLASS_QUAD) {
+      t8_productionf (
+        "t8_geometry_cylinder is not compatible with tree type %s\n It is only compatible with quad elements.\n",
+        t8_eclass_to_string[active_tree_class]);
+      return false;
+    }
+    return true;
+  }
+
   /**
    * Get the type of this geometry.
    * \return The type.
@@ -341,6 +389,22 @@ struct t8_geometry_circle: public t8_geometry_with_vertices
     SC_ABORT_NOT_REACHED ();
   }
 
+  /**
+   * Check for compatibility of the currently loaded tree with the geometry.
+   * Only quad elements are supported by this geometry.
+   */
+  bool
+  t8_geom_check_tree_compatibility () const
+  {
+    if (active_tree_class != T8_ECLASS_QUAD) {
+      t8_productionf (
+        "t8_geometry_circle is not compatible with tree type %s\n It is only compatible with quad elements.\n",
+        t8_eclass_to_string[active_tree_class]);
+      return false;
+    }
+    return true;
+  }
+
   /**
    * Get the type of this geometry.
    * \return The type.
@@ -434,6 +498,22 @@ struct t8_geometry_moving: public t8_geometry
     return 0;
   }
 
+  /**
+   * Check for compatibility of the currently loaded tree with the geometry.
+   * Only quad elements are supported by this geometry.
+   */
+  bool
+  t8_geom_check_tree_compatibility () const
+  {
+    if (active_tree_class != T8_ECLASS_QUAD) {
+      t8_productionf (
+        "t8_geometry_moving is not compatible with tree type %s\n It is only compatible with quad elements.\n",
+        t8_eclass_to_string[active_tree_class]);
+      return false;
+    }
+    return true;
+  }
+
   /**
    * Get the type of this geometry.
    * \return The type.
@@ -504,6 +584,22 @@ struct t8_geometry_cube_zdistorted: public t8_geometry
     return 0;
   }
 
+  /**
+   * Check for compatibility of the currently loaded tree with the geometry.
+   * Only hex elements are supported by this geometry.
+   */
+  bool
+  t8_geom_check_tree_compatibility () const
+  {
+    if (active_tree_class != T8_ECLASS_HEX) {
+      t8_productionf (
+        "t8_geometry_cube_zdistorted is not compatible with tree type %s\n It is only compatible with hex elements.\n",
+        t8_eclass_to_string[active_tree_class]);
+      return false;
+    }
+    return true;
+  }
+
   /**
    * Get the type of this geometry.
    * \return The type.
@@ -653,7 +749,7 @@ t8_analytic_geom (int level, t8_example_geom_type geom_type)
   case T8_GEOM_ANALYTIC_QUAD_TO_SPHERE:
     t8_global_productionf ("Wrapping a quad around a sphere.\n");
     t8_cmesh_register_geometry<t8_geometry_analytic> (cmesh, 3, "geom_quad_to_sphere", quad_to_sphere_callback, nullptr,
-                                                      nullptr, nullptr, nullptr);
+                                                      nullptr, nullptr, nullptr, nullptr);
     t8_cmesh_set_tree_class (cmesh, 0, T8_ECLASS_QUAD);
     t8_cmesh_set_join (cmesh, 0, 0, 1, 0, 0);
 

src/t8_cmesh.h

@@ -104,17 +104,6 @@ t8_cmesh_is_committed (const t8_cmesh_t cmesh);
  */
 int
 t8_cmesh_validate_geometry (const t8_cmesh_t cmesh);
-
-/** After a cmesh is committed, check whether all trees in a cmesh do have positive volume.
- * Returns true if all trees have positive volume.
- * \param [in]  cmesh           This cmesh is examined. May be NULL.
- * \return                      True if \a cmesh is not NULL and all trees for
- *                              which \ref  t8_cmesh_set_tree_vertices
- *                              was called, do have positive geometric volume.
- *                              False otherwise.
- */
-int
-t8_cmesh_no_negative_volume (t8_cmesh_t cmesh);
 #endif
 
 /** Given a set of vertex coordinates for a tree of a given eclass.

src/t8_cmesh/t8_cmesh.cxx

@@ -124,15 +124,45 @@ t8_cmesh_is_committed (const t8_cmesh_t cmesh)
   return 1;
 }
 
-#ifdef T8_ENABLE_DEBUG
+#if T8_ENABLE_DEBUG
 int
 t8_cmesh_validate_geometry (const t8_cmesh_t cmesh)
 {
+  /* After a cmesh is committed, check whether all trees in a cmesh are compatible
+ * with their geometry and if they have positive volume.
+ * Returns true if all trees are valid. Returns also true if no geometries are
+ * registered yet, since the validity computation depends on the used geometry.
+ */
+
   /* Geometry handler is not constructed yet */
   if (cmesh->geometry_handler == NULL) {
-    return 1;
+    return true;
   }
-  return t8_cmesh_no_negative_volume (cmesh);
+  if (cmesh == NULL) {
+    return true;
+  }
+  if (cmesh->geometry_handler->get_num_geometries () > 0) {
+    /* Iterate over all trees, get their vertices and check the volume */
+    for (t8_locidx_t itree = 0; itree < cmesh->num_local_trees; itree++) {
+      /* Check if tree and geometry are compatible. */
+      const int geometry_compatible
+        = cmesh->geometry_handler->tree_compatible_with_geom (cmesh, t8_cmesh_get_global_id (cmesh, itree));
+      if (!geometry_compatible) {
+        t8_debugf ("Detected incompatible geometry for tree %li\n", (long) itree);
+        return false;
+      }
+      if (geometry_compatible) {
+        /* Check for negative volume. This only makes sense if the geometry is valid for the tree. */
+        const int negative_volume
+          = cmesh->geometry_handler->tree_negative_volume (cmesh, t8_cmesh_get_global_id (cmesh, itree));
+        if (negative_volume) {
+          t8_debugf ("Detected negative volume in tree %li\n", (long) itree);
+          return false;
+        }
+      }
+    }
+  }
+  return true;
 }
 #endif /* T8_ENABLE_DEBUG */
 
@@ -400,7 +430,6 @@ t8_cmesh_get_tree_vertices (const t8_cmesh_t cmesh, const t8_locidx_t ltreeid)
 {
   T8_ASSERT (t8_cmesh_is_committed (cmesh));
   T8_ASSERT (t8_cmesh_treeid_is_local_tree (cmesh, ltreeid) || t8_cmesh_treeid_is_ghost (cmesh, ltreeid));
-
   return (double *) t8_cmesh_get_attribute (cmesh, t8_get_package_id (), T8_CMESH_VERTICES_ATTRIBUTE_KEY, ltreeid);
 }
 
@@ -582,39 +611,6 @@ t8_cmesh_tree_vertices_negative_volume (const t8_eclass_t eclass, const double *
   return eclass == T8_ECLASS_TET ? sc_prod > 0 : sc_prod < 0;
 }
 
-#ifdef T8_ENABLE_DEBUG
-/* After a cmesh is committed, check whether all trees in a cmesh do have positive volume.
- * Returns true if all trees have positive volume. Returns also true if no geometries are
- * registered yet, since the volume computation depends on the used geometry.
- */
-int
-t8_cmesh_no_negative_volume (const t8_cmesh_t cmesh)
-{
-  bool res = false;
-
-  if (cmesh == NULL) {
-    return 0;
-  }
-  if (cmesh->geometry_handler == NULL) {
-    return 0;
-  }
-  if (cmesh->geometry_handler->get_num_geometries () > 0) {
-    /* Iterate over all trees, get their vertices and check the volume */
-    for (t8_locidx_t itree = 0; itree < cmesh->num_local_trees; itree++) {
-      const int ret = cmesh->geometry_handler->tree_negative_volume (cmesh, t8_cmesh_get_global_id (cmesh, itree));
-      if (ret) {
-        t8_debugf ("Detected negative volume in tree %li\n", (long) itree);
-      }
-      res |= ret; /* res is true if one ret value is true */
-    }
-    return !res;
-  }
-  else {
-    return true;
-  }
-}
-#endif
-
 void
 t8_cmesh_set_tree_vertices (t8_cmesh_t cmesh, const t8_gloidx_t gtree_id, const double *vertices,
                             const int num_vertices)

src/t8_cmesh/t8_cmesh_types.h

@@ -59,9 +59,9 @@ typedef struct t8_cprofile t8_cprofile_t; /* Defined below */
   +T8_ECLASS_MAX_EDGES /* Used to store which face is linked to which surface */
 #define T8_CMESH_CAD_FACE_PARAMETERS_ATTRIBUTE_KEY \
   T8_CMESH_CAD_FACE_ATTRIBUTE_KEY + 1 /* Used to store face parameters */
-#define T8_CMESH_LAGRANGE_POLY_DEGREE T8_CMESH_CAD_FACE_PARAMETERS_ATTRIBUTE_KEY + T8_ECLASS_MAX_FACES
+#define T8_CMESH_LAGRANGE_POLY_DEGREE_KEY T8_CMESH_CAD_FACE_PARAMETERS_ATTRIBUTE_KEY + T8_ECLASS_MAX_FACES
 #define T8_CMESH_NEXT_POSSIBLE_KEY \
-  T8_CMESH_LAGRANGE_POLY_DEGREE + 1 /* The next free value for a t8code attribute key */
+  T8_CMESH_LAGRANGE_POLY_DEGREE_KEY + 1 /* The next free value for a t8code attribute key */
 
 /** This structure holds the connectivity data of the coarse mesh.
  *  It can either be replicated, then each process stores a copy of the whole

src/t8_geometry/t8_geometry_base.cxx

@@ -62,3 +62,25 @@ t8_geom_get_type (const t8_geometry_c *geom)
 
   return geom->t8_geom_get_type ();
 }
+
+/* Load the id and class of the newly active tree to the active_tree and active_tree_class variable. */
+void
+t8_geometry::t8_geom_load_tree_data (const t8_cmesh_t cmesh, const t8_gloidx_t gtreeid)
+{
+  /* Set active id and eclass */
+  const t8_locidx_t ltreeid = t8_cmesh_get_local_id (cmesh, gtreeid);
+  active_tree = gtreeid;
+  const t8_locidx_t num_local_trees = t8_cmesh_get_num_local_trees (cmesh);
+  if (0 <= ltreeid && ltreeid < num_local_trees) {
+    active_tree_class = t8_cmesh_get_tree_class (cmesh, ltreeid);
+  }
+  else {
+    active_tree_class = t8_cmesh_get_ghost_class (cmesh, ltreeid - num_local_trees);
+  }
+
+  /* Check whether we support this class */
+  T8_ASSERT (active_tree_class == T8_ECLASS_VERTEX || active_tree_class == T8_ECLASS_TRIANGLE
+             || active_tree_class == T8_ECLASS_TET || active_tree_class == T8_ECLASS_QUAD
+             || active_tree_class == T8_ECLASS_HEX || active_tree_class == T8_ECLASS_LINE
+             || active_tree_class == T8_ECLASS_PRISM || active_tree_class == T8_ECLASS_PYRAMID);
+}

src/t8_geometry/t8_geometry_base.hxx

@@ -98,14 +98,14 @@ struct t8_geometry
 
   /** Update a possible internal data buffer for per tree data.
    * This function is called before the first coordinates in a new tree are
-   * evaluated. You can use it for example to load the vertex coordinates of the 
-   * tree into an internal buffer (as is done in the linear geometry).
+   * evaluated.
+   * In this base implementation we use it to load the treeid and class
+   * to the internal member variables \a active_tree and \a active_tree_class.
    * \param [in]  cmesh      The cmesh.
    * \param [in]  gtreeid    The global tree.
    */
   virtual void
-  t8_geom_load_tree_data (t8_cmesh_t cmesh, t8_gloidx_t gtreeid)
-    = 0;
+  t8_geom_load_tree_data (const t8_cmesh_t cmesh, const t8_gloidx_t gtreeid);
 
   /** Query whether a batch of points lies inside an element. 
  * \param [in]      forest      The forest.
@@ -131,8 +131,8 @@ struct t8_geometry
   };
 
   /**
-   * Check if  the currently active tree has a negative volume
-   * \return                True (non-zero) if the currently loaded tree has a negative volume. 0 otherwise.  
+   * Check if the currently active tree has a negative volume.
+   * \return                True if the currently loaded tree has a negative volume.
    */
   virtual bool
   t8_geom_tree_negative_volume () const
@@ -142,6 +142,17 @@ struct t8_geometry
     return 0;
   };
 
+  /**
+   * Check for compatibility of the currently loaded tree with the geometry.
+   * If the geometry has limitations these can be checked here.
+   * This includes for example if only specific tree types or dimensions are supported.
+   * If all trees are supported, this function should return true.
+   * \return                True if the geometry is compatible with the tree.
+   */
+  virtual bool
+  t8_geom_check_tree_compatibility () const
+    = 0;
+
   /**
    * Get the dimension of this geometry.
    * \return The dimension.
@@ -177,14 +188,11 @@ struct t8_geometry
     = 0;
 
  protected:
-  int dimension;
-  /**< The dimension of reference space for which this is a geometry. */
-
-  std::string name;
-  /**< The name of this geometry. */
-
-  size_t hash;
-  /**< The hash of the name of this geometry. */
+  int dimension;                 /**< The dimension of reference space for which this is a geometry. */
+  std::string name;              /**< The name of this geometry. */
+  size_t hash;                   /**< The hash of the name of this geometry. */
+  t8_gloidx_t active_tree;       /**< The tree of which currently vertices are loaded. */
+  t8_eclass_t active_tree_class; /**< The class of the currently active tree. */
 };
 
 T8_EXTERN_C_END ();