Skip to content

Commit

Permalink
Merge pull request #863 from DLR-AMR/fix_pyramid-face-corner
Browse files Browse the repository at this point in the history
Fix pyramid face corner
  • Loading branch information
lukasdreyer authored Feb 7, 2024
2 parents 02aa85e + e05cba3 commit fb4728a
Show file tree
Hide file tree
Showing 8 changed files with 200 additions and 9 deletions.
23 changes: 23 additions & 0 deletions src/t8_forest/t8_forest.c
Original file line number Diff line number Diff line change
Expand Up @@ -1061,6 +1061,29 @@ t8_forest_get_local_id (const t8_forest_t forest, const t8_gloidx_t gtreeid)
return -1;
}
}
t8_locidx_t
t8_forest_get_local_or_ghost_id (const t8_forest_t forest, const t8_gloidx_t gtreeid)
{
t8_gloidx_t ltreeid;
T8_ASSERT (t8_forest_is_committed (forest));
T8_ASSERT (0 <= gtreeid && gtreeid < t8_forest_get_num_global_trees (forest));

/* If the tree is local then its local id is the global id minus the
* first global tree id on this forest. If this number is not in the
* range of local tree ids then the tree is not local. */
/* we use a gloidx for ltreeid to prevent overflow and false positives */
ltreeid = gtreeid - t8_forest_get_first_local_tree_id (forest);
/* Check if this tree is a local tree and if so return its local id */
if (0 <= ltreeid && ltreeid < t8_forest_get_num_local_trees (forest)) {
return ltreeid;
}
else {
t8_locidx_t ghost_id = t8_forest_ghost_get_ghost_treeid (forest, gtreeid);
if (ghost_id >= 0)
return t8_forest_get_num_local_trees (forest) + ghost_id;
return -1;
}
}

t8_locidx_t
t8_forest_ltreeid_to_cmesh_ltreeid (t8_forest_t forest, t8_locidx_t ltreeid)
Expand Down
21 changes: 17 additions & 4 deletions src/t8_forest/t8_forest_cxx.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -2001,10 +2001,10 @@ t8_forest_element_half_face_neighbors (t8_forest_t forest, t8_locidx_t ltreeid,
}

void
t8_forest_leaf_face_neighbors (t8_forest_t forest, t8_locidx_t ltreeid, const t8_element_t *leaf,
t8_element_t **pneighbor_leaves[], int face, int *dual_faces[], int *num_neighbors,
t8_locidx_t **pelement_indices, t8_eclass_scheme_c **pneigh_scheme,
int forest_is_balanced)
t8_forest_leaf_face_neighbors_ext (t8_forest_t forest, t8_locidx_t ltreeid, const t8_element_t *leaf,
t8_element_t **pneighbor_leaves[], int face, int *dual_faces[], int *num_neighbors,
t8_locidx_t **pelement_indices, t8_eclass_scheme_c **pneigh_scheme,
int forest_is_balanced, t8_gloidx_t *gneigh_tree)
{
t8_eclass_t neigh_class, eclass;
t8_gloidx_t gneigh_treeid;
Expand Down Expand Up @@ -2056,6 +2056,9 @@ t8_forest_leaf_face_neighbors (t8_forest_t forest, t8_locidx_t ltreeid, const t8
gneigh_treeid = t8_forest_element_half_face_neighbors (forest, ltreeid, leaf, neighbor_leaves, neigh_scheme, face,
num_children_at_face, *dual_faces);
}
if (gneigh_tree) {
*gneigh_tree = gneigh_treeid;
}
if (gneigh_treeid < 0) {
/* There exists no face neighbor across this face, we return with this info */
neigh_scheme->t8_element_destroy (num_children_at_face, neighbor_leaves);
Expand Down Expand Up @@ -2232,6 +2235,16 @@ t8_forest_leaf_face_neighbors (t8_forest_t forest, t8_locidx_t ltreeid, const t8
}
}

void
t8_forest_leaf_face_neighbors (t8_forest_t forest, t8_locidx_t ltreeid, const t8_element_t *leaf,
t8_element_t **pneighbor_leaves[], int face, int *dual_faces[], int *num_neighbors,
t8_locidx_t **pelement_indices, t8_eclass_scheme_c **pneigh_scheme,
int forest_is_balanced)
{
t8_forest_leaf_face_neighbors_ext (forest, ltreeid, leaf, pneighbor_leaves, face, dual_faces, num_neighbors,
pelement_indices, pneigh_scheme, forest_is_balanced, NULL);
}

void
t8_forest_print_all_leaf_neighbors (t8_forest_t forest)
{
Expand Down
21 changes: 21 additions & 0 deletions src/t8_forest/t8_forest_general.h
Original file line number Diff line number Diff line change
Expand Up @@ -448,6 +448,21 @@ t8_forest_get_eclass (const t8_forest_t forest, const t8_locidx_t ltreeid);
t8_locidx_t
t8_forest_get_local_id (const t8_forest_t forest, const t8_gloidx_t gtreeid);

/** Given a global tree id compute the forest local id of this tree.
* If the tree is a local tree, then the local id is between 0 and the number
* of local trees. If the tree is a ghost, then the local id is between num_local_trees and
* num_local_trees + num_ghost_trees.
* If the tree is neither a local tree nor a ghost tree, a negative number is returned.
* \param [in] forest The forest.
* \param [in] gtreeid The global id of a tree.
* \return The tree's local id in \a forest, if it is a local tree.
* num_local_trees + the ghosts id, if it is a ghost tree.
* A negative number if not.
* \see https://github.com/DLR-AMR/t8code/wiki/Tree-indexing for more details about tree indexing.
*/
t8_locidx_t
t8_forest_get_local_or_ghost_id (const t8_forest_t forest, const t8_gloidx_t gtreeid);

/** Given the local id of a tree in a forest, compute the tree's local id in the associated cmesh.
* \param [in] forest The forest.
* \param [in] ltreeid The local id of a tree or ghost in the forest.
Expand Down Expand Up @@ -506,6 +521,12 @@ t8_forest_leaf_face_neighbors (t8_forest_t forest, t8_locidx_t ltreeid, const t8
t8_locidx_t **pelement_indices, t8_eclass_scheme_c **pneigh_scheme,
int forest_is_balanced);

void
t8_forest_leaf_face_neighbors_ext (t8_forest_t forest, t8_locidx_t ltreeid, const t8_element_t *leaf,
t8_element_t **pneighbor_leaves[], int face, int *dual_faces[], int *num_neighbors,
t8_locidx_t **pelement_indices, t8_eclass_scheme_c **pneigh_scheme,
int forest_is_balanced, t8_gloidx_t *gneigh_tree);

/** Exchange ghost information of user defined element data.
* \param[in] forest The forest. Must be committed.
* \param[in] element_data An array of length num_local_elements + num_ghosts
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -1401,7 +1401,7 @@ t8_dpyramid_get_face_corner (const t8_dpyramid_t *pyra, int face, int corner)
return t8_dtet_face_corner[face][corner];
}
else {
int corner_number = t8_dpyramid_face_corner[face][corner];
const int corner_number = t8_dpyramid_face_corner[pyra->pyramid.type - T8_DPYRAMID_FIRST_TYPE][face][corner];
T8_ASSERT (0 <= corner_number && corner_number < T8_DPYRAMID_FACES);
return corner_number;
}
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -128,10 +128,19 @@ const int t8_dpyramid_tritype_rootface_to_face[2][4] = {
{ 2, 0, 2, 0 },
{ 1, 0, 1, 0 } };

const int t8_dpyramid_face_corner[5][4] = {
const int t8_dpyramid_face_corner[2][5][4] = {
{
{ 0, 2, 4, -1 },
{ 1, 3, 4, -1 },
{ 0, 1, 4, -1 },
{ 2, 3, 4, -1 },
{ 0, 1, 2, 3 } };
{ 0, 1, 2, 3 }
},{
{ 2, 3, 4, -1 },
{ 0, 1, 4, -1 },
{ 1, 3, 4, -1 },
{ 0, 2, 4, -1 },
{ 0, 1, 2, 3 }
}
};
/* clang-format on */
Original file line number Diff line number Diff line change
Expand Up @@ -111,8 +111,8 @@ extern const int t8_dpyramid_tritype_rootface_to_face[2][4];
extern const int t8_dtet_type_cid_to_pyramid_parenttype[6][8];

/** The corner numbers of the face of a pyramid
* corner_number = A(face, corner_number_at_face)
* corner_number = A(pyramid_type, face, corner_number_at_face)
*/
extern const int t8_dpyramid_face_corner[5][4];
extern const int t8_dpyramid_face_corner[2][5][4];

#endif // T8_DPYRAMID_CONNECTIVITY_H
10 changes: 10 additions & 0 deletions test/Makefile.am
Original file line number Diff line number Diff line change
Expand Up @@ -52,6 +52,7 @@ t8code_googletest_programs = \
test/t8_data/t8_gtest_shmem \
test/t8_forest/t8_gtest_half_neighbors \
test/t8_forest/t8_gtest_find_owner \
test/t8_forest/t8_gtest_forest_face_normal \
test/t8_schemes/t8_gtest_face_descendant \
test/t8_geometry/t8_gtest_point_inside \
test/t8_forest/t8_gtest_user_data \
Expand Down Expand Up @@ -228,6 +229,10 @@ test_t8_forest_t8_gtest_find_owner_SOURCES = \
test/t8_gtest_main.cxx \
test/t8_forest/t8_gtest_find_owner.cxx

test_t8_forest_t8_gtest_forest_face_normal_SOURCES = \
test/t8_gtest_main.cxx \
test/t8_forest/t8_gtest_forest_face_normal.cxx

test_t8_schemes_t8_gtest_face_descendant_SOURCES = \
test/t8_gtest_main.cxx \
test/t8_schemes/t8_gtest_face_descendant.cxx
Expand Down Expand Up @@ -448,6 +453,10 @@ test_t8_forest_t8_gtest_find_owner_LDADD = $(t8_gtest_target_ld_add)
test_t8_forest_t8_gtest_find_owner_LDFLAGS = $(t8_gtest_target_ld_flags)
test_t8_forest_t8_gtest_find_owner_CPPFLAGS = $(t8_gtest_target_cpp_flags)

test_t8_forest_t8_gtest_forest_face_normal_LDADD = $(t8_gtest_target_ld_add)
test_t8_forest_t8_gtest_forest_face_normal_LDFLAGS = $(t8_gtest_target_ld_flags)
test_t8_forest_t8_gtest_forest_face_normal_CPPFLAGS = $(t8_gtest_target_cpp_flags)

test_t8_schemes_t8_gtest_face_descendant_LDADD = $(t8_gtest_target_ld_add)
test_t8_schemes_t8_gtest_face_descendant_LDFLAGS = $(t8_gtest_target_ld_flags)
test_t8_schemes_t8_gtest_face_descendant_CPPFLAGS = $(t8_gtest_target_cpp_flags)
Expand Down Expand Up @@ -550,6 +559,7 @@ test_t8_gtest_vtk_linkage_CPPFLAGS += $(t8_gtest_target_mpi_cpp_flags)
test_t8_data_t8_gtest_shmem_CPPFLAGS += $(t8_gtest_target_mpi_cpp_flags)
test_t8_forest_t8_gtest_half_neighbors_CPPFLAGS += $(t8_gtest_target_mpi_cpp_flags)
test_t8_forest_t8_gtest_find_owner_CPPFLAGS += $(t8_gtest_target_mpi_cpp_flags)
test_t8_forest_t8_gtest_forest_face_normal_CPPFLAGS += $(t8_gtest_target_mpi_cpp_flags)
test_t8_schemes_t8_gtest_face_descendant_CPPFLAGS += $(t8_gtest_target_mpi_cpp_flags)
test_t8_geometry_t8_gtest_point_inside_CPPFLAGS += $(t8_gtest_target_mpi_cpp_flags)
test_t8_forest_t8_gtest_user_data_CPPFLAGS += $(t8_gtest_target_mpi_cpp_flags)
Expand Down
115 changes: 115 additions & 0 deletions test/t8_forest/t8_gtest_forest_face_normal.cxx
Original file line number Diff line number Diff line change
@@ -0,0 +1,115 @@
/*
This file is part of t8code.
t8code is a C library to manage a collection (a forest) of multiple
connected adaptive space-trees of general element classes in parallel.
Copyright (C) 2023 the developers
t8code is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or
(at your option) any later version.
t8code is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with t8code; if not, write to the Free Software Foundation, Inc.,
51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
*/

#include <sc/src/sc_functions.h>
#include <gtest/gtest.h>
#include <t8_eclass.h>
#include <t8_schemes/t8_default/t8_default_cxx.hxx>
#include <t8_schemes/t8_default/t8_default_pyramid/t8_dpyramid_bits.h>
#include <t8_cmesh/t8_cmesh_examples.h>
#include <t8_forest/t8_forest_general.h>
#include <t8_forest/t8_forest_geometrical.h>
#include <test/t8_gtest_macros.hxx>

/**
* This file tests the face normal computation of elements.
*/

class class_forest_face_normal: public testing::TestWithParam<std::tuple<t8_eclass_t, int>> {
protected:
void
SetUp () override
{
eclass = std::get<0> (GetParam ());
level = std::get<1> (GetParam ());
scheme = t8_scheme_new_default_cxx ();
t8_cmesh_t cmesh = t8_cmesh_new_hypercube (eclass, sc_MPI_COMM_WORLD, 0, 0, 0);
const int do_face_ghost = 1;
forest = t8_forest_new_uniform (cmesh, scheme, level, do_face_ghost, sc_MPI_COMM_WORLD);
}
void
TearDown () override
{
t8_forest_unref (&forest);
}
t8_forest_t forest;
t8_scheme_cxx *scheme;
t8_eclass_t eclass;
int level;
};

TEST_P (class_forest_face_normal, back_and_forth)
{
/** Iterate over all elements of a uniformly refined forest. For all faceneighbors elements, if they are local,
* check if their facenormal is the negative of the corresponding facenormal of the neighbor elements.
*/
const t8_locidx_t local_num_trees = t8_forest_get_num_local_trees (forest);
/* Iterate over all elements. */
for (t8_locidx_t itree = 0; itree < local_num_trees; itree++) {
const t8_locidx_t tree_elements = t8_forest_get_tree_num_elements (forest, itree);
const t8_eclass_t tree_eclass = t8_forest_get_tree_class (forest, itree);
ASSERT_EQ (eclass, tree_eclass);
const t8_eclass_scheme_c *escheme = t8_forest_get_eclass_scheme (forest, tree_eclass);
for (t8_locidx_t ielement = 0; ielement < tree_elements; ielement++) {
const t8_element_t *element = t8_forest_get_element_in_tree (forest, itree, ielement);
const int num_faces = escheme->t8_element_num_faces (element);
for (int iface = 0; iface < num_faces; iface++) {
/* Compute facenormal */
double face_normal[3] = { 0, 0, 0 };
t8_forest_element_face_normal (forest, itree, element, iface, face_normal);

/* Get all face neighbors */

t8_element_t **neighbors;
int num_neighbors;
const int forest_is_balanced = 1;
t8_eclass_scheme_c *neigh_scheme;
int *dual_faces;
t8_locidx_t *neigh_ids;

t8_gloidx_t gneightree;
t8_forest_leaf_face_neighbors_ext (forest, itree, element, &neighbors, iface, &dual_faces, &num_neighbors,
&neigh_ids, &neigh_scheme, forest_is_balanced, &gneightree);

/* Iterate and compute their facenormal */
for (int ineigh = 0; ineigh < num_neighbors; ineigh++) {
t8_locidx_t lneightree = t8_forest_get_local_or_ghost_id (forest, gneightree);
double neigh_face_normal[3];
t8_forest_element_face_normal (forest, lneightree, neighbors[ineigh], dual_faces[ineigh], neigh_face_normal);
EXPECT_NEAR (face_normal[0], -neigh_face_normal[0], 10 * T8_PRECISION_EPS);
EXPECT_NEAR (face_normal[1], -neigh_face_normal[1], 10 * T8_PRECISION_EPS);
EXPECT_NEAR (face_normal[2], -neigh_face_normal[2], 10 * T8_PRECISION_EPS);
}

if (num_neighbors > 0) {
neigh_scheme->t8_element_destroy (num_neighbors, neighbors);
T8_FREE (neigh_ids);
T8_FREE (neighbors);
T8_FREE (dual_faces);
}
}
}
}
}

INSTANTIATE_TEST_SUITE_P (t8_gtest_forest_face_normal, class_forest_face_normal,
testing::Combine (AllEclasses, testing::Range (0, 2)));

0 comments on commit fb4728a

Please sign in to comment.