Merge pull request #1173 from DLR-AMR/feature-multiple_attributes_gho…

…st_commit_from_stash

Feature multiple attributes ghost commit from stash

src/CMakeLists.txt

@@ -94,7 +94,7 @@ target_sources( T8 PRIVATE
     t8_cmesh/t8_cmesh_vtk_reader.cxx 
     t8_cmesh/t8_cmesh_save.cxx 
     t8_cmesh/t8_cmesh_netcdf.c 
-    t8_cmesh/t8_cmesh_trees.c 
+    t8_cmesh/t8_cmesh_trees.cxx 
     t8_cmesh/t8_cmesh_commit.cxx 
     t8_cmesh/t8_cmesh_partition.cxx 
     t8_cmesh/t8_cmesh_copy.c 

src/Makefile.am

@@ -114,7 +114,7 @@ libt8_compiled_sources = \
   src/t8_cmesh/t8_cmesh_vtk_reader.cxx \
   src/t8_cmesh/t8_cmesh_save.cxx \
   src/t8_cmesh/t8_cmesh_netcdf.c \
-  src/t8_cmesh/t8_cmesh_trees.c src/t8_cmesh/t8_cmesh_commit.cxx \
+  src/t8_cmesh/t8_cmesh_trees.cxx src/t8_cmesh/t8_cmesh_commit.cxx \
   src/t8_cmesh/t8_cmesh_partition.cxx\
   src/t8_cmesh/t8_cmesh_copy.c src/t8_data/t8_shmem.c \
   src/t8_cmesh/t8_cmesh_geometry.cxx \

src/t8_cmesh/t8_cmesh_commit.cxx

@@ -84,7 +84,7 @@ t8_cmesh_set_shmem_type (sc_MPI_Comm comm)
 }
 
 static void
-t8_cmesh_add_attributes (const t8_cmesh_t cmesh, sc_hash_t *ghost_ids, size_t *attribute_data_offset)
+t8_cmesh_add_attributes (const t8_cmesh_t cmesh, sc_hash_t *ghost_ids)
 {
   t8_stash_attribute_struct_t *attribute;
   const t8_stash_t stash = cmesh->stash;
@@ -94,6 +94,7 @@ t8_cmesh_add_attributes (const t8_cmesh_t cmesh, sc_hash_t *ghost_ids, size_t *a
 
   temp_facejoin = T8_ALLOC_ZERO (t8_ghost_facejoin_t, 1);
 
+  t8_locidx_t ghosts_inserted = 0;
   ltree = -1;
   for (si = 0, sj = 0; si < stash->attributes.elem_count; si++, sj++) {
     attribute = (t8_stash_attribute_struct_t *) sc_array_index (&stash->attributes, si);
@@ -112,9 +113,13 @@ t8_cmesh_add_attributes (const t8_cmesh_t cmesh, sc_hash_t *ghost_ids, size_t *a
       T8_ASSERT (ghost_ids != NULL);
       temp_facejoin->ghost_id = attribute->id;
       if (sc_hash_lookup (ghost_ids, temp_facejoin, (void ***) &facejoin_pp)) {
+        T8_ASSERT ((t8_locidx_t) sj == (t8_locidx_t) (*facejoin_pp)->attr_id);
+        if (sj == 0) {
+          ghosts_inserted++;
+        }
         /* attribute is on a ghost tree */
-        t8_cmesh_trees_add_ghost_attribute (cmesh->trees, 0, attribute, (*facejoin_pp)->local_id,
-                                            (*facejoin_pp)->attr_id, attribute_data_offset);
+        t8_cmesh_trees_add_ghost_attribute (cmesh->trees, attribute, (*facejoin_pp)->local_id, ghosts_inserted,
+                                            (*facejoin_pp)->attr_id);
         (*facejoin_pp)->attr_id++;
       }
     }
@@ -167,7 +172,7 @@ t8_cmesh_commit_replicated_new (t8_cmesh_t cmesh)
     t8_stash_attribute_sort (cmesh->stash);
     cmesh->num_trees = cmesh->num_local_trees = num_trees;
     cmesh->first_tree = 0;
-    t8_cmesh_add_attributes (cmesh, NULL, NULL);
+    t8_cmesh_add_attributes (cmesh, NULL);
 
     /* Set all face connections */
     t8_cmesh_trees_set_all_boundary (cmesh, cmesh->trees);
@@ -198,7 +203,6 @@ t8_cmesh_commit_partitioned_new (t8_cmesh_t cmesh, sc_MPI_Comm comm)
   t8_cghost_t ghost1;
   int F;
   size_t si;
-  size_t attribute_data_offset;
 
 #if T8_ENABLE_DEBUG
   sc_flopinfo_t fi, snapshot;
@@ -375,7 +379,7 @@ t8_cmesh_commit_partitioned_new (t8_cmesh_t cmesh, sc_MPI_Comm comm)
         ghost1->att_offset += attribute->attr_size;
       }
     }
-    attribute_data_offset = t8_cmesh_trees_finish_part (cmesh->trees, 0);
+    t8_cmesh_trees_finish_part (cmesh->trees, 0);
     t8_cmesh_trees_set_all_boundary (cmesh, cmesh->trees);
 
     /* Go through all face_neighbour entries and parse every
@@ -473,7 +477,7 @@ t8_cmesh_commit_partitioned_new (t8_cmesh_t cmesh, sc_MPI_Comm comm)
      *       counting the attributes per tree. */
 
     t8_stash_attribute_sort (cmesh->stash);
-    t8_cmesh_add_attributes (cmesh, ghost_ids, &attribute_data_offset);
+    t8_cmesh_add_attributes (cmesh, ghost_ids);
 
     /* compute global number of trees. id1 serves as buffer since
      * global number and local number have different datatypes */

src/t8_cmesh/t8_cmesh_examples.cxx

@@ -371,7 +371,7 @@ t8_cmesh_new_pyramid (sc_MPI_Comm comm)
   /* Use linear geometry */
   t8_cmesh_register_geometry<t8_geometry_linear> (cmesh, 3);
   t8_cmesh_set_tree_class (cmesh, 0, T8_ECLASS_PYRAMID);
-  t8_cmesh_set_tree_vertices (cmesh, 0, vertices, 15);
+  t8_cmesh_set_tree_vertices (cmesh, 0, vertices, 5);
   t8_cmesh_commit (cmesh, comm);
   return cmesh;
 }

src/t8_cmesh/t8_cmesh_trees.c → src/t8_cmesh/t8_cmesh_trees.cxx

@@ -25,6 +25,9 @@
  * TODO: document this file
  */
 
+#include <vector>
+#include <algorithm>
+#include <iostream>
 #include "t8_cmesh_stash.h"
 #include "t8_cmesh_trees.h"
 
@@ -210,6 +213,46 @@ t8_cmesh_trees_start_part (const t8_cmesh_trees_t trees, const int proc, const t
   part->first_ghost_id = lfirst_ghost;
 }
 
+/* Helper struct for sorting the number of ghost attributes by global id.
+ * In order to sort them, we need the part ghost id to access the global id.
+ * Thus, we store both the part id and the number of attributes. */
+typedef struct
+{
+  t8_locidx_t part_ghost_id;
+  t8_gloidx_t global_id;
+  int num_attributes;
+  int attribute_offset;
+} t8_part_ghost_id_and_num_atts;
+
+/* Compare function for t8_part_ghost_id_and_num_atts to compare by global id.
+*
+* Return True if global_id of if_A < global_id of id_B
+* Return False otherwise
+* */
+bool
+t8_compare_id_and_att_by_global_id (t8_part_ghost_id_and_num_atts &id_A, t8_part_ghost_id_and_num_atts &id_B)
+{
+  return id_A.global_id < id_B.global_id;
+}
+
+/* Compare function for t8_part_ghost_id_and_num_atts to compare by local id.
+*
+* Return True if local id of if_A < local id of id_B
+* Return False otherwise
+* */
+bool
+t8_compare_id_and_att_by_part_id (t8_part_ghost_id_and_num_atts &id_A, t8_part_ghost_id_and_num_atts &id_B)
+{
+  return id_A.part_ghost_id < id_B.part_ghost_id;
+}
+
+static std::vector<t8_part_ghost_id_and_num_atts>
+t8_cmesh_allocate_ghost_num_atts_array (t8_locidx_t part_num_ghosts)
+{
+  std::vector<t8_part_ghost_id_and_num_atts> num_attributes_of_ghosts (part_num_ghosts);
+  return num_attributes_of_ghosts;
+}
+
 /* After all classes of trees and ghosts have been set and after the
  * number of tree attributes  was set and their total size (per tree)
  * stored temporarily in the att_offset variable
@@ -218,21 +261,17 @@ t8_cmesh_trees_start_part (const t8_cmesh_trees_t trees, const int proc, const t
 /* The workflow can be: call start_part, set tree and ghost classes maually, call
  * init_attributes, call finish_part, successively call add_attributes
  * and also set all face neighbors (TODO: write function)*/
-size_t
+void
 t8_cmesh_trees_finish_part (const t8_cmesh_trees_t trees, const int proc)
 {
   t8_part_tree_t part;
   t8_ctree_t tree;
   t8_cghost_t ghost;
-  size_t tree_attr_bytes;
-  size_t ghost_attr_bytes;
-  size_t face_neigh_bytes;     /* count the total number of bytes needed for face_neighbor information */
-  size_t first_face;           /* offset of the first face neighbor information */
-  size_t first_tree;           /* offset of the first tree */
-  size_t first_ghost;          /* offset of the first ghost */
-  size_t temp_offset;          /* offset of the currently looked at tree/ghost */
-  size_t num_tree_attributes;  /* total number of tree attributes */
-  size_t num_ghost_attributes; /* total number of ghost attributes */
+  size_t tree_attr_data_bytes;
+  size_t ghost_attr_data_bytes;
+  size_t face_neigh_bytes;    /* count the total number of bytes needed for face_neighbor information */
+  size_t temp_offset;         /* offset of the currently looked at tree/ghost */
+  size_t num_tree_attributes; /* total number of tree attributes */
   t8_attribute_info_struct_t *attr;
   t8_locidx_t it;
 #ifndef SC_ENABLE_REALLOC
@@ -243,16 +282,22 @@ t8_cmesh_trees_finish_part (const t8_cmesh_trees_t trees, const int proc)
   part = t8_cmesh_trees_get_part (trees, proc);
   T8_ASSERT (part != NULL);
 
-  num_tree_attributes = num_ghost_attributes = 0;
-  tree_attr_bytes = ghost_attr_bytes = face_neigh_bytes = 0;
+  num_tree_attributes = 0;
+  tree_attr_data_bytes = ghost_attr_data_bytes = face_neigh_bytes = 0;
   /* The offset of the first tree */
-  first_tree = 0;
+  const size_t first_tree = 0;
   /* The offset of the first ghost */
-  first_ghost = first_tree + part->num_trees * sizeof (t8_ctree_struct_t);
+  const size_t first_ghost = first_tree + part->num_trees * sizeof (t8_ctree_struct_t);
   /* The offset of the first ghost face */
-  first_face = first_ghost + part->num_ghosts * sizeof (t8_cghost_struct_t);
+  const size_t first_face = first_ghost + part->num_ghosts * sizeof (t8_cghost_struct_t);
+
+  /* First pass through ghosts to set the face neighbor offsets. */
+  /* Additionally, we need to sort the number of attributes according to the global id in order
+   * to properly compute the attribute info offsets later. */
+
+  std::vector<t8_part_ghost_id_and_num_atts> num_attributes_of_ghosts
+    = t8_cmesh_allocate_ghost_num_atts_array (part->num_ghosts);
 
-  /* First pass through ghosts to set the face neighbor offsets */
   temp_offset = first_ghost;
   for (it = 0; it < part->num_ghosts; it++) {
     ghost = t8_part_tree_get_ghost (part, it + part->first_ghost_id);
@@ -264,7 +309,27 @@ t8_cmesh_trees_finish_part (const t8_cmesh_trees_t trees, const int proc)
     face_neigh_bytes += T8_ADD_PADDING (face_neigh_bytes);
     T8_ASSERT (face_neigh_bytes % T8_PADDING_SIZE == 0);
     temp_offset += sizeof (t8_cghost_struct_t);
+
+    /* Add the number of attributes of this ghost to the array. */
+    t8_part_ghost_id_and_num_atts &lid_and_num_atts_entry = num_attributes_of_ghosts[it];
+    lid_and_num_atts_entry.part_ghost_id = it;
+    lid_and_num_atts_entry.num_attributes = ghost->num_attributes;
+    lid_and_num_atts_entry.global_id = ghost->treeid;
+    lid_and_num_atts_entry.attribute_offset = 0;
+  }
+  /* We now sort the array of num attributes by global id */
+  std::sort (num_attributes_of_ghosts.begin (), num_attributes_of_ghosts.end (), t8_compare_id_and_att_by_global_id);
+
+  if (part->num_ghosts > 0) {
+    num_attributes_of_ghosts[0].attribute_offset = 0;
   }
+  for (t8_locidx_t ighost = 1; ighost < part->num_ghosts; ++ighost) {
+    // Build the sum: 0 a a+b a+b+c a+b+c+d
+    num_attributes_of_ghosts[ighost].attribute_offset
+      = num_attributes_of_ghosts[ighost - 1].attribute_offset + num_attributes_of_ghosts[ighost - 1].num_attributes;
+  }
+  /* We now need to sort the offset sum by local id again. */
+  std::sort (num_attributes_of_ghosts.begin (), num_attributes_of_ghosts.end (), t8_compare_id_and_att_by_part_id);
 
   /* First pass through trees to set the face neighbor offsets */
   temp_offset = 0;
@@ -280,41 +345,51 @@ t8_cmesh_trees_finish_part (const t8_cmesh_trees_t trees, const int proc)
 
   /* Second pass through trees to set attribute offsets */
   temp_offset = 0;
+  size_t next_tree_offset = 0;  // Compute the offset of the next tree
   for (it = 0; it < part->num_trees; it++) {
     tree = t8_part_tree_get_tree (part, it + part->first_tree_id);
-    tree_attr_bytes += tree->att_offset; /* att_offset temporarily stored the total size of the attributes */
+    tree_attr_data_bytes += tree->att_offset; /* att_offset temporarily stored the total size of the attributes */
     /* The att_offset of the tree is the first_face plus the number of attribute
      * bytes used by previous trees minus the temp_offset */
     tree->att_offset
       = first_face - temp_offset + face_neigh_bytes + num_tree_attributes * sizeof (t8_attribute_info_struct_t);
+    T8_ASSERT (it == 0 || tree->att_offset == next_tree_offset);
     num_tree_attributes += tree->num_attributes;
     temp_offset += sizeof (t8_ctree_struct_t);
+    next_tree_offset
+      = tree->att_offset + tree->num_attributes * sizeof (t8_attribute_info_struct_t) - sizeof (t8_ctree_struct_t);
   }
-  tree_attr_bytes += num_tree_attributes * sizeof (t8_attribute_info_struct_t);
+  const size_t tree_attr_total_bytes = tree_attr_data_bytes + num_tree_attributes * sizeof (t8_attribute_info_struct_t);
 
   /* Second pass through ghosts to set attribute offsets */
-  temp_offset = first_ghost;
+  temp_offset = 0;
+  size_t num_ghost_attributes = 0; /* total number of ghost attributes */
+  /* To get the offset of the first ghost attribute info, we have to add the bytes
+   * of all tree attributes. */
+  size_t first_ghost_offset = next_tree_offset + tree_attr_data_bytes;
   for (it = 0; it < part->num_ghosts; it++) {
     ghost = t8_part_tree_get_ghost (part, it + part->first_ghost_id);
-    ghost_attr_bytes += ghost->att_offset; /* att_offset temporarily stored the total size of the attributes */
-    /* The att_offset of the tree is the first_face plus the number of attribute
-     * bytes used by previous trees minus the temp_offset */
-    ghost->att_offset = first_face - temp_offset + face_neigh_bytes + tree_attr_bytes
-                        + num_ghost_attributes * sizeof (t8_attribute_info_struct_t);
-    num_ghost_attributes += ghost->num_attributes;
+    ghost_attr_data_bytes += ghost->att_offset; /* att_offset temporarily stored the total size of the attributes */
+    /* The att_offset of the ghost is the offset of the first ghost + the attribute
+     * offset of this ghost minus the size of all previous ghosts. */
+    ghost->att_offset = first_ghost_offset
+                        + num_attributes_of_ghosts[it].attribute_offset * sizeof (t8_attribute_info_struct_t)
+                        - temp_offset;
     temp_offset += sizeof (t8_cghost_struct_t);
+    num_ghost_attributes += ghost->num_attributes;
   }
-  ghost_attr_bytes += num_ghost_attributes * sizeof (t8_attribute_info_struct_t);
-  size_t attr_bytes = tree_attr_bytes + ghost_attr_bytes;
+  const size_t ghost_attr_total_bytes
+    = ghost_attr_data_bytes + num_ghost_attributes * sizeof (t8_attribute_info_struct_t);
+  const size_t attr_total_bytes = tree_attr_total_bytes + ghost_attr_total_bytes;
 
   /* Done setting all tree and ghost offsets */
-  /* Allocate memory, first_face + attr_bytes + face_neigh_bytes gives the new total byte count */
+  /* Allocate memory, first_face + attr_total_bytes + face_neigh_bytes gives the new total byte count */
 #ifdef SC_ENABLE_REALLOC
   /* Since we use realloc and padding, memcmp will not work if we don't set everything to zero, solved with memset */
-  SC_REALLOC (part->first_tree, char, first_face + attr_bytes + +face_neigh_bytes);
-  memset (part->first_tree + first_face, 0, attr_bytes + face_neigh_bytes)
+  SC_REALLOC (part->first_tree, char, first_face + attr_total_bytes + face_neigh_bytes);
+  memset (part->first_tree + first_face, 0, attr_total_bytes + face_neigh_bytes)
 #else
-  temp = T8_ALLOC_ZERO (char, first_face + attr_bytes + face_neigh_bytes);
+  temp = T8_ALLOC_ZERO (char, first_face + attr_total_bytes + face_neigh_bytes);
   memcpy (temp, part->first_tree, first_face);
   T8_FREE (part->first_tree);
   part->first_tree = temp;
@@ -324,7 +399,11 @@ t8_cmesh_trees_finish_part (const t8_cmesh_trees_t trees, const int proc)
     attr = (t8_attribute_info_struct_t *) (part->first_tree + first_face + face_neigh_bytes);
     attr->attribute_offset = num_tree_attributes * sizeof (t8_attribute_info_struct_t);
   }
-  return num_ghost_attributes * sizeof (t8_attribute_info_struct_t);
+  if (num_ghost_attributes > 0) {
+    ghost = t8_part_tree_get_ghost (part, 0);
+    attr = (t8_attribute_info_struct_t *) (part->first_tree + first_face + face_neigh_bytes + tree_attr_total_bytes);
+    attr->attribute_offset = num_ghost_attributes * sizeof (t8_attribute_info_struct_t);
+  }
 }
 
 void
@@ -608,7 +687,7 @@ t8_cmesh_trees_copy_toproc (t8_cmesh_trees_t trees_dest, const t8_cmesh_trees_t
 
 void
 t8_cmesh_trees_init_attributes (const t8_cmesh_trees_t trees, const t8_locidx_t ltree_id, const size_t num_attributes,
-                                const size_t attr_bytes)
+                                const size_t attr_info_bytes)
 {
   int proc;
   t8_ctree_t tree;
@@ -619,7 +698,7 @@ t8_cmesh_trees_init_attributes (const t8_cmesh_trees_t trees, const t8_locidx_t
   T8_ASSERT (proc >= 0 && proc < t8_cmesh_trees_get_num_procs (trees));
   tree = t8_part_tree_get_tree (t8_cmesh_trees_get_part (trees, proc), ltree_id);
 
-  tree->att_offset = attr_bytes; /* This is only temporary until t8_cmesh_trees_finish_part
+  tree->att_offset = attr_info_bytes; /* This is only temporary until t8_cmesh_trees_finish_part
                                            is called */
   tree->num_attributes = num_attributes;
 }
@@ -676,9 +755,9 @@ t8_cmesh_trees_add_attribute (const t8_cmesh_trees_t trees, const int proc, cons
 }
 
 void
-t8_cmesh_trees_add_ghost_attribute (const t8_cmesh_trees_t trees, const int proc,
-                                    const t8_stash_attribute_struct_t *attr, const t8_locidx_t local_ghost_id,
-                                    const size_t index, size_t *attribute_data_offset)
+t8_cmesh_trees_add_ghost_attribute (const t8_cmesh_trees_t trees, const t8_stash_attribute_struct_t *attr,
+                                    const t8_locidx_t local_ghost_id, const t8_locidx_t ghosts_inserted,
+                                    const size_t index)
 {
   t8_part_tree_t part;
   t8_cghost_t ghost;
@@ -690,11 +769,10 @@ t8_cmesh_trees_add_ghost_attribute (const t8_cmesh_trees_t trees, const int proc
   T8_ASSERT (attr->attr_data != NULL || attr->attr_size == 0);
   T8_ASSERT (attr->id >= 0);
 
-  part = t8_cmesh_trees_get_part (trees, proc);
+  part = t8_cmesh_trees_get_part (trees, 0);
   ghost = t8_part_tree_get_ghost (part, local_ghost_id);
 
   attr_info = T8_GHOST_ATTR_INFO (ghost, index);
-  attr_info->attribute_offset = *attribute_data_offset;
   new_attr_data = T8_GHOST_ATTR (ghost, attr_info);
 
   memcpy (new_attr_data, attr->attr_data, attr->attr_size);
@@ -704,9 +782,19 @@ t8_cmesh_trees_add_ghost_attribute (const t8_cmesh_trees_t trees, const int proc
   attr_info->package_id = attr->package_id;
   attr_info->attribute_size = attr->attr_size;
 
-  *attribute_data_offset += attr->attr_size;
-  if (index == (size_t) ghost->num_attributes - 1) {
-    *attribute_data_offset -= ghost->num_attributes * sizeof (t8_attribute_info_struct_t);
+  /* If we are not yet at the last attribute of the part,
+   * get next attribute and set its offset*/
+  if (!(index == (size_t) ghost->num_attributes - 1 && part->num_ghosts == ghosts_inserted)) {
+    /* Store offset of current attribute */
+    const size_t offset = attr_info->attribute_offset;
+    attr_info = attr_info + 1;
+    attr_info->attribute_offset = offset + attr->attr_size;
+    /* if the current attribute was the last attribute of the tree
+     * the next attribute offset must be corrected by the size of
+     * the attribute infos of the current tree */
+    if (index == (size_t) ghost->num_attributes - 1) {
+      attr_info->attribute_offset -= ghost->num_attributes * sizeof (t8_attribute_info_struct_t);
+    }
   }
 }