Make minor improvements and fix bugs

- module -M: fixed -part for periodic meshes, fixed compatibility with Gmsh
  4.13.1.

VERSIONS

@@ -1,5 +1,6 @@
-New in 4.10.1-2 (09 Oct 2024):
-- module -M: Fixed compatibility with Gmsh 4.13.1.
+New in 4.10.1 (09 Oct 2024):
+- module -M: fixed -part for periodic meshes, fixed compatibility with Gmsh
+  4.13.1.
 
 New in 4.10.0 (01 Oct 2024):
 - general: switched "active" and "passive" for orientation convention (see

doc/conf.py

@@ -11,8 +11,8 @@
 import sphinx_rtd_theme
 
 project = u'Neper'
-version = u'4.10.1-2'
-release = u'4.10.1-2'
+version = u'4.10.1'
+release = u'4.10.1'
 author = u'Romain Quey'
 copyright = u'Romain Quey'
 language = 'en'

src/CMakeLists.txt

@@ -7,7 +7,7 @@ if(POLICY CMP0077)
   cmake_policy(SET CMP0077 NEW)
 endif()
 
-set(NEPER_VERSION \"4.10.1-2\")
+set(NEPER_VERSION \"4.10.1\")
 project(neper)
 
 if(CMAKE_CXX_COMPILER_VERSION VERSION_LESS 4.8.1)

src/neper_m/nem_order/nem_order2.c

@@ -82,6 +82,8 @@ nem_order_periodic (struct TESS Tess, struct NODES *pNodes, struct MESH *Mesh)
   if (Tess.Dim == 3)
     nem_order_periodic_faces (Tess, pNodes, Mesh);
 
+  neut_nodes_init_nodeslave (pNodes);
+
   return;
 }
 

src/neper_m/nem_order/nem_order3.c

@@ -40,7 +40,7 @@ nem_order_periodic_edges (struct TESS Tess, struct NODES *pNodes,
   int i, j, k, l, nodeqty_tot, nodeqty_o1, nodes[2], nodemasters[2], node,
     nodemaster;
   int **fir = NULL;
-  int *node_slave = ut_alloc_1d_int ((*pNodes).NodeQty);
+  int *node_slave = ut_alloc_1d_int ((*pNodes).NodeQty + 1);
   int *shift = ut_alloc_1d_int (3);
 
   neut_elt_orderarrays (Mesh[1].EltType, Mesh[1].Dimension, &fir,

src/neper_m/nem_partition/nem_partition2.c

@@ -102,7 +102,7 @@ nem_partition_prep (struct NODES *pNodes, struct MESH *Mesh,
     neut_mesh_init_nodeelts (Mesh + dim, (*pNodes).NodeQty);
 
   SCOTCH_meshInit (pSCMesh);
-  neut_mesh_scotchmesh (Mesh[dim], (*pNodes).NodeQty, pSCMesh);
+  neut_mesh_scotchmesh (*pNodes, Mesh[dim], pSCMesh);
 
   return;
 }

src/neut/neut_mesh/neut_mesh_fprintf/neut_mesh_fprintf_scotch/neut_mesh_fprintf_scotch1.c

@@ -17,7 +17,7 @@ neut_mesh_fprintf_scotch (char *body, struct NODES Nodes, struct MESH Mesh)
 
   SCOTCH_meshInit (&SCMesh);
 
-  neut_mesh_scotchmesh (Mesh, Nodes.NodeQty, &SCMesh);
+  neut_mesh_scotchmesh (Nodes, Mesh, &SCMesh);
 
   file = ut_file_open (scm, "w");
   SCOTCH_meshSave (&SCMesh, file);

src/neut/neut_mesh/neut_mesh_op/neut_mesh_op.c

@@ -294,47 +294,43 @@ neut_mesh_memcpy (struct MESH Old, struct MESH *pNew)
 
 #ifdef HAVE_LIBSCOTCH
 void
-neut_mesh_scotchmesh (struct MESH Mesh, int vnodnbr, SCOTCH_Mesh * pSCMesh)
+neut_mesh_scotchmesh (struct NODES Nodes, struct MESH Mesh, SCOTCH_Mesh * pSCMesh)
 {
-  int i, id, status;
-  int velmbas, vnodbas, velmnbr, vertnbr, edgenbr;
-  int *verttab, *edgetab;
-  int eltnodeqty =
-    neut_elt_nodeqty (Mesh.EltType, Mesh.Dimension, Mesh.EltOrder);
+  int i, j, status;
+  int verqty, edgeqty, nodeqty, eltqty;
+  int *vers = NULL, *edges = NULL, *nodes = NULL, *elts = NULL;
+  int eltnodeqty = neut_elt_nodeqty (Mesh.EltType, Mesh.Dimension, Mesh.EltOrder);
 
-  velmbas = 1;                  /* element base number */
-  vnodbas = 1 + Mesh.EltQty;    /* node base number */
-  velmnbr = Mesh.EltQty;
+  verqty = Mesh.EltQty + Nodes.NodeQty;
+  edgeqty = 2 * Mesh.EltQty * eltnodeqty;
 
-  vertnbr = velmnbr + vnodnbr;
-  edgenbr = 2 * velmnbr * eltnodeqty;
-
-  verttab = ut_alloc_1d_int (vertnbr + 1);
-  edgetab = ut_alloc_1d_int (edgenbr + 1);
+  vers = ut_alloc_1d_int (verqty + 1);
 
   /* recording elements */
-  id = 0;
-  for (i = 1; i <= velmnbr; i++)
+  edgeqty = 0;
+  for (i = 1; i <= Mesh.EltQty; i++)
   {
-    verttab[i - 1] = id + 1;
-    ut_array_1d_int_memcpy (Mesh.EltNodes[i], eltnodeqty, edgetab + id);
-    ut_array_1d_int_addval (edgetab + id, eltnodeqty, vnodbas - 1,
-                            edgetab + id);
-    id += eltnodeqty;
+    vers[i - 1] = edgeqty + 1;
+
+    neut_mesh_elt_nodes_per (Nodes, Mesh, i, &nodes, &nodeqty);
+    for (j = 0; j < nodeqty; j++)
+      ut_array_1d_int_list_addval_nocheck (&edges, &edgeqty, nodes[j]);
   }
+  ut_array_1d_int_addval (edges, edgeqty, Mesh.EltQty, edges);
 
   /* recording nodes */
-  for (i = 1; i <= vnodnbr; i++)
+  for (i = 1; i <= Nodes.NodeQty; i++)
   {
-    verttab[i + Mesh.EltQty - 1] = id + 1;
-    ut_array_1d_int_memcpy (Mesh.NodeElts[i] + 1, Mesh.NodeElts[i][0],
-                            edgetab + id);
-    id += Mesh.NodeElts[i][0];
+    vers[i + Mesh.EltQty - 1] = edgeqty + 1;
+
+    neut_mesh_node_elts_per (Nodes, Mesh, i, &elts, &eltqty);
+    for (j = 0; j < eltqty; j++)
+      ut_array_1d_int_list_addval_nocheck (&edges, &edgeqty, elts[j]);
   }
-  verttab[vertnbr] = edgenbr + 1;
+  vers[verqty] = edgeqty + 1;
 
-  SCOTCH_meshBuild (pSCMesh, velmbas, vnodbas, velmnbr, vnodnbr, verttab,
-                    NULL, NULL, NULL, NULL, edgenbr, edgetab);
+  SCOTCH_meshBuild (pSCMesh, 1, Mesh.EltQty + 1, Mesh.EltQty, Nodes.NodeQty, vers,
+                    NULL, NULL, NULL, NULL, edgeqty, edges);
   /*
      file = ut_file_open ("dd", "w");
      SCOTCH_meshSave (pSCMesh, file);

src/neut/neut_mesh/neut_mesh_op/neut_mesh_op.h

@@ -71,7 +71,7 @@ extern "C"
 ///
 ///
 #ifdef HAVE_LIBSCOTCH
-  extern void neut_mesh_scotchmesh (struct MESH, int, SCOTCH_Mesh *);
+  extern void neut_mesh_scotchmesh (struct NODES Nodes, struct MESH Mesh, SCOTCH_Mesh *);
 #endif /* HAVE_LIBSCOTCH */
 
 /// \brief Initialize NodeElts of a MESH structure.

src/neut/neut_mesh/neut_mesh_topo/neut_mesh_topo.c

@@ -2504,3 +2504,63 @@ neut_mesh_order1nodes_order2node (struct MESH Mesh, int node1, int node2,
 
   return;
 }
+
+void
+neut_mesh_elt_nodes_per (struct NODES Nodes, struct MESH Mesh, int elt, int **pnodes, int *pnodeqty)
+{
+  int i, j, node, eltnodeqty = neut_elt_nodeqty (Mesh.EltType, Mesh.Dimension, Mesh.EltOrder);
+
+  (*pnodeqty) = eltnodeqty;
+  (*pnodes) = ut_alloc_1d_int (*pnodeqty);
+
+  ut_array_1d_int_memcpy (Mesh.EltNodes[elt], eltnodeqty, *pnodes);
+
+  // periodic
+  if (Nodes.PerNodeQty > 0)
+  {
+    for (i = 0; i < eltnodeqty; i++)
+    {
+      node = Mesh.EltNodes[elt][i];
+
+      // is slave
+      if (Nodes.PerNodeMaster[node] > 0)
+        ut_array_1d_int_list_addval_nocheck (pnodes, pnodeqty, Nodes.PerNodeMaster[node]);
+
+      // is master
+      else if (Nodes.PerNodeSlaveQty[node] > 0)
+        for (j = 1; j <= Nodes.PerNodeSlaveQty[node]; j++)
+          ut_array_1d_int_list_addval_nocheck (pnodes, pnodeqty, Nodes.PerNodeSlaveNb[node][j]);
+    }
+  }
+
+  return;
+}
+
+void
+neut_mesh_node_elts_per (struct NODES Nodes, struct MESH Mesh, int node, int **pelts, int *peltqty)
+{
+  int nodeeltqty, *nodeelts = NULL;
+
+  (*peltqty) = Mesh.NodeElts[node][0];
+  (*pelts) = ut_alloc_1d_int (*peltqty);
+
+  ut_array_1d_int_memcpy (Mesh.NodeElts[node] + 1, Mesh.NodeElts[node][0], *pelts);
+
+  // periodic
+  if (Nodes.PerNodeQty > 0)
+  {
+    // is slave
+    if (Nodes.PerNodeMaster[node] > 0)
+      ut_array_1d_int_listpair_merge (*pelts, *peltqty, Mesh.NodeElts[Nodes.PerNodeMaster[node]] + 1,
+                                      Mesh.NodeElts[Nodes.PerNodeMaster[node]][0], pelts, peltqty);
+
+    // is master
+    else if (Nodes.PerNodeSlaveQty[node] > 0)
+    {
+      neut_mesh_nodes_allelts (Mesh, Nodes.PerNodeSlaveNb[node] + 1, Nodes.PerNodeSlaveQty[node], &nodeelts, &nodeeltqty);
+      ut_array_1d_int_listpair_merge (*pelts, *peltqty, nodeelts, nodeeltqty, pelts, peltqty);
+    }
+  }
+
+  return;
+}

src/neut/neut_mesh/neut_mesh_topo/neut_mesh_topo.h

@@ -395,6 +395,9 @@ extern void neut_mesh_elt_domface (struct TESS Tess, struct MESH Mesh2D, int elt
 
 extern void neut_mesh_order1nodes_order2node (struct MESH Mesh, int node1, int node2, int *pnode);
 
+extern void neut_mesh_elt_nodes_per (struct NODES Nodes, struct MESH Mesh, int elt, int **pnodes, int *pnodeqty);
+extern void neut_mesh_node_elts_per (struct NODES Nodes, struct MESH Mesh, int node, int **pelts, int *peltqty);
+
 #endif /* NEUT_MESH_TOPO_H */
 
 #ifdef __cplusplus

src/neut/neut_nodes/neut_nodes.c

@@ -316,16 +316,15 @@ neut_nodes_memcpy (struct NODES Nodes, struct NODES *pNodes2)
     if (Nodes.PerNodeSlaveNb)
     {
       (*pNodes2).PerNodeSlaveNb = ut_alloc_1d_pint ((*pNodes2).NodeQty + 1);
-      for (i = 1; i <= (*pNodes2).PerNodeQty; i++)
-      {
-        id = (*pNodes2).PerNodeNb[i];
-        (*pNodes2).PerNodeSlaveNb[id] =
-          ut_alloc_1d_int ((*pNodes2).PerNodeSlaveQty[id] + 1);
-
-        ut_array_1d_int_memcpy (Nodes.PerNodeSlaveNb[id] + 1,
-                                (*pNodes2).PerNodeSlaveQty[id],
-                                (*pNodes2).PerNodeSlaveNb[id] + 1);
-      }
+      for (i = 1; i <= (*pNodes2).NodeQty; i++)
+        if ((*pNodes2).PerNodeSlaveQty[i])
+        {
+          (*pNodes2).PerNodeSlaveNb[i] = ut_alloc_1d_int ((*pNodes2).PerNodeSlaveQty[i] + 1);
+
+          ut_array_1d_int_memcpy (Nodes.PerNodeSlaveNb[i] + 1,
+                                  (*pNodes2).PerNodeSlaveQty[i],
+                                  (*pNodes2).PerNodeSlaveNb[i] + 1);
+        }
     }
   }
 
@@ -629,21 +628,69 @@ neut_nodes_wbary (struct NODES Nodes, int *nodes, double *nodeweights,
 }
 
 /* node_nbs[...] = new pos */
-/* RAM could be improved by more elegant exchange of lines (not through
- * a big copy coo array */
 void
 neut_nodes_renumber_switch (struct NODES *pNodes, int *node_nbs)
 {
-  int i;
-  double **coo = ut_alloc_2d ((*pNodes).NodeQty + 1, 3);
+  int i, j;
+  struct NODES Nodes2;
 
-  for (i = 1; i <= (*pNodes).NodeQty; i++)
-    ut_array_1d_memcpy ((*pNodes).NodeCoo[i], 3, coo[i]);
+  neut_nodes_set_zero (&Nodes2);
+  neut_nodes_memcpy (*pNodes, &Nodes2);
 
   for (i = 1; i <= (*pNodes).NodeQty; i++)
-    ut_array_1d_memcpy (coo[i], 3, (*pNodes).NodeCoo[node_nbs[i]]);
+  {
+    ut_array_1d_memcpy (Nodes2.NodeCoo[i], 3, (*pNodes).NodeCoo[node_nbs[i]]);
+    if ((*pNodes).NodeCl)
+      (*pNodes).NodeCl[node_nbs[i]] = Nodes2.NodeCl[i];
+  }
 
-  ut_free_2d (&coo, (*pNodes).NodeQty);
+  if ((*pNodes).PerNodeQty > 0)
+  {
+    for (i = 1; i <= (*pNodes).PerNodeQty; i++)
+      (*pNodes).PerNodeNb[i] = node_nbs[Nodes2.PerNodeNb[i]];
+
+    ut_array_1d_int_zero ((*pNodes).PerNodeMaster + 1, (*pNodes).NodeQty);
+    ut_array_2d_int_zero ((*pNodes).PerNodeShift + 1, (*pNodes).NodeQty, 3);
+    ut_array_1d_int_zero ((*pNodes).PerNodeSlaveQty + 1, (*pNodes).NodeQty);
+    for (i = 1; i <= (*pNodes).NodeQty; i++)
+      ut_free_1d_int ((*pNodes).PerNodeSlaveNb + i);
+
+    for (i = 1; i <= (*pNodes).NodeQty; i++)
+    {
+      (*pNodes).PerNodeMaster[node_nbs[i]] = node_nbs[Nodes2.PerNodeMaster[i]];
+      ut_array_1d_int_memcpy (Nodes2.PerNodeShift[i], 3, (*pNodes).PerNodeShift[node_nbs[i]]);
+      (*pNodes).PerNodeSlaveQty[node_nbs[i]] = Nodes2.PerNodeSlaveQty[i];
+      if ((*pNodes).PerNodeSlaveQty[node_nbs[i]])
+      {
+        (*pNodes).PerNodeSlaveNb[node_nbs[i]] = ut_realloc_1d_int ((*pNodes).PerNodeSlaveNb[node_nbs[i]],
+                                                                   (*pNodes).PerNodeSlaveQty[node_nbs[i]] + 1);
+        for (j = 1; j <= Nodes2.PerNodeSlaveQty[i]; j++)
+          (*pNodes).PerNodeSlaveNb[node_nbs[i]][j] = node_nbs[Nodes2.PerNodeSlaveNb[i][j]];
+      }
+    }
+  }
+
+  if ((*pNodes).DupNodeQty > 0)
+  {
+    for (i = 1; i <= (*pNodes).DupNodeQty; i++)
+      (*pNodes).DupNodeNb[i] = node_nbs[i];
+
+    for (i = 1; i <= (*pNodes).NodeQty; i++)
+    {
+      (*pNodes).DupNodeMaster[node_nbs[i]] = node_nbs[Nodes2.DupNodeMaster[i]];
+      (*pNodes).DupNodeSeed[node_nbs[i]] = Nodes2.DupNodeSeed[i];
+      (*pNodes).DupNodeSlaveQty[node_nbs[i]] = Nodes2.DupNodeSlaveQty[i];
+      if ((*pNodes).DupNodeSlaveQty[node_nbs[i]])
+      {
+        (*pNodes).DupNodeSlaveNb[node_nbs[i]] = ut_realloc_1d_int ((*pNodes).DupNodeSlaveNb[node_nbs[i]],
+                                                                    (*pNodes).DupNodeSlaveQty[node_nbs[i]] + 1);
+        for (j = 1; j <= Nodes2.DupNodeSlaveQty[i]; j++)
+          (*pNodes).DupNodeSlaveNb[node_nbs[i]][j] = node_nbs[Nodes2.DupNodeSlaveNb[i][j]];
+      }
+    }
+  }
+
+  neut_nodes_free (&Nodes2);
 
   return;
 }
@@ -716,6 +763,11 @@ neut_nodes_init_nodeslave (struct NODES *pNodes)
 {
   int i, master, slave;
 
+  if ((*pNodes).PerNodeSlaveQty)
+    ut_free_1d_int (&(*pNodes).PerNodeSlaveQty);
+  if ((*pNodes).PerNodeSlaveNb)
+    ut_free_2d_int (&(*pNodes).PerNodeSlaveNb, (*pNodes).PerNodeQty);
+
   (*pNodes).PerNodeSlaveQty = ut_alloc_1d_int ((*pNodes).NodeQty + 1);
   (*pNodes).PerNodeSlaveNb = ut_alloc_1d_pint ((*pNodes).NodeQty + 1);