Possible undefined behavior when computing dual faces in the LFN function.

[Flo1314]

In the LFN function, the dual face is computed as follows (lines 1854-1884) in t8_forest_cxx.cxx:

/* We need to determine the dual face */
if (neigh_scheme->t8_element_level (ancestor) ==
    ts->t8_element_level (leaf)) {
  /* The ancestor is the same-level neighbor of leaf */
  if (!at_maxlevel) {
    /* its dual face is the face of the parent of the first neighbor leaf */
    *dual_faces[0] =
      neigh_scheme->t8_element_face_parent_face (neighbor_leafs[0],
                                                 *dual_faces[0]);
  }
}
else {
  /* The ancestor is the parent of the parent */
  T8_ASSERT (neigh_scheme->t8_element_level (ancestor) ==
             ts->t8_element_level (leaf) - 1);

  *dual_faces[0] =
    neigh_scheme->t8_element_face_parent_face (neighbor_leafs[0],
                                               *dual_faces[0]);
  if (!at_maxlevel) {
    /* We need to compute the dual face of the grandparent. */
    /* Construct the parent of the grand child */
    neigh_scheme->t8_element_parent (neighbor_leafs[0],
                                     neighbor_leafs[0]);
    /* Compute the face id of the parent's face */
    *dual_faces[0] =
      neigh_scheme->t8_element_face_parent_face (neighbor_leafs[0],
                                                 *dual_faces[0]);
  }
}

In the else case, if !maxlevel is true, then neigh_scheme->t8_element_face_parent_face will be executed twice in a row. Furthermore, this function may return -1 after its first execution.

Lets look inside the default_quad implementation of this function:

int
t8_default_scheme_quad_c::t8_element_face_parent_face (const t8_element_t
                                                       *elem, int face)
{
  const p4est_quadrant_t *q = (const p4est_quadrant_t *) elem;
  int                 child_id;

  T8_ASSERT (t8_element_is_valid (elem));
  if (q->level == 0) {
    return face;
  }
  /* Determine whether face is a subface of the parent.
   * This is the case if the child_id matches one of the faces corners */
  child_id = p4est_quadrant_child_id (q);
  if (child_id == p4est_face_corners[face][0]
      || child_id == p4est_face_corners[face][1]) {
    return face;
  }
  return -1;
}

Now, if face = -1 we run into problems when calling

child_id == p4est_face_corners[face][0] || child_id == p4est_face_corners[face][1]

I noticed this problem while using subelements and I am not sure if this is still a problem without them. However, this seems dangerous to me.

[holke]

This may be resolved with #1736, we need to check that.