1 /*============================================================================
2  * Management of mesh groups.
3  *===========================================================================*/
4 
5 /*
6   This file is part of Code_Saturne, a general-purpose CFD tool.
7 
8   Copyright (C) 1998-2021 EDF S.A.
9 
10   This program is free software; you can redistribute it and/or modify it under
11   the terms of the GNU General Public License as published by the Free Software
12   Foundation; either version 2 of the License, or (at your option) any later
13   version.
14 
15   This program is distributed in the hope that it will be useful, but WITHOUT
16   ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
17   FOR A PARTICULAR PURPOSE.  See the GNU General Public License for more
18   details.
19 
20   You should have received a copy of the GNU General Public License along with
21   this program; if not, write to the Free Software Foundation, Inc., 51 Franklin
22   Street, Fifth Floor, Boston, MA 02110-1301, USA.
23 */
24 
25 /*----------------------------------------------------------------------------*/
26 
27 #include "cs_defs.h"
28 
29 /*----------------------------------------------------------------------------
30  * Standard C library headers
31  *---------------------------------------------------------------------------*/
32 
33 #include <stdio.h>
34 #include <string.h>
35 #include <assert.h>
36 #include <math.h>
37 #include <float.h>
38 
39 /*----------------------------------------------------------------------------
40  *  Local headers
41  *---------------------------------------------------------------------------*/
42 
43 #include "bft_mem.h"
44 #include "bft_error.h"
45 #include "bft_printf.h"
46 
47 #include "fvm_io_num.h"
48 #include "fvm_periodicity.h"
49 
50 #include "cs_order.h"
51 #include "cs_search.h"
52 #include "cs_sort.h"
53 #include "cs_join_perio.h"
54 #include "cs_join_post.h"
55 #include "cs_join_util.h"
56 
57 /*----------------------------------------------------------------------------
58  *  Header for the current file
59  *---------------------------------------------------------------------------*/
60 
61 #include "cs_mesh_group.h"
62 
63 /*---------------------------------------------------------------------------*/
64 
65 BEGIN_C_DECLS
66 
67 /*! \cond DOXYGEN_SHOULD_SKIP_THIS */
68 
69 /*============================================================================
70  * Static global variables
71  *===========================================================================*/
72 
73 /*=============================================================================
74  * Local Structure Definitions
75  *===========================================================================*/
76 
77 /*============================================================================
78  * Private function definitions
79  *===========================================================================*/
80 
81 #if defined(HAVE_MPI)
82 
83 /*----------------------------------------------------------------------------*/
84 /*!
85  * \brief Compare 2 family definitions
86  *
87  * \param[in]  n0    number of values in family 0
88  * \param[in]  n1    number of values in family 1
89  * \param[in]  val0  family 0 values
90  * \param[in]  val1  family 1 values
91  *
92  * \return  -1 if family 0 is smaller than family 1, 0 if families are equal,
93  *   1 otherwise
94  */
95 /*----------------------------------------------------------------------------*/
96 
97 inline static int
_sync_compare_families(int n0,int n1,const int val0[],const int val1[])98 _sync_compare_families(int        n0,
99                        int        n1,
100                        const int  val0[],
101                        const int  val1[])
102 {
103   int i;
104   int retval = 0;
105 
106   for (i = 0; i < n0 && i < n1; i++) {
107     if (val0[i] < val1[i]) {
108       retval = -1;
109       break;
110     }
111     else if (val0[i] > val1[i]) {
112       retval = 1;
113       break;
114     }
115   }
116 
117   if (retval == 0) {
118     if (n0 < n1)
119       retval = -1;
120     else if (n0 > n1)
121       retval = 1;
122   }
123 
124   return retval;
125 }
126 
127 /*----------------------------------------------------------------------------*/
128 /*!
129  * \brief Synchronize family combinations globally.
130  *
131  * \param[in, out]  n_fam       number of families
132  * \param[in, out]  family_idx  index of element families
133  * \param[in, out]  family      element family ancestors
134  *
135  * \return  renumbering from previous to synchronized family combinations
136  */
137 /*----------------------------------------------------------------------------*/
138 
139 static int *
_sync_family_combinations(int * n_fam,int ** family_idx,int ** family)140 _sync_family_combinations(int   *n_fam,
141                           int  **family_idx,
142                           int  **family)
143 {
144   int i;
145   int sizes[2];
146   int rank_mult, rank_id, n_ranks;
147 
148   int _n_fam = *n_fam;
149   int *_family_idx = NULL;
150   int *_family = NULL;
151 
152   int *renum = NULL;
153 
154   BFT_MALLOC(_family_idx, _n_fam+1, int);
155   memcpy(_family_idx, *family_idx, (_n_fam+1)*sizeof(int));
156 
157   BFT_MALLOC(_family, _family_idx[_n_fam], int);
158   if (_family != NULL)
159     memcpy(_family, *family, _family_idx[_n_fam]*sizeof(int));
160 
161   /* Merge all info by stages up to rank 0 */
162 
163   for (rank_id = cs_glob_rank_id, n_ranks = cs_glob_n_ranks, rank_mult = 1;
164        n_ranks > 1;
165        rank_id /= 2, n_ranks = (n_ranks+1)/2, rank_mult *= 2) {
166 
167     /* Even ranks receive and merge */
168 
169     if (rank_id %2 == 0 && rank_id + 1 < n_ranks) {
170 
171       MPI_Status status;
172       int recv_rank = (rank_id + 1)*rank_mult;
173 
174       MPI_Recv(sizes, 2, MPI_INT, recv_rank, 1, cs_glob_mpi_comm, &status);
175 
176       if (sizes[0] > 0) {
177 
178         int j, k;
179         int n_fam_new = 0;
180         int *idx_recv = NULL, *fam_recv = NULL;
181         int *idx_new = NULL, *fam_new = NULL;
182 
183         BFT_MALLOC(idx_recv, sizes[0] + 1, int);
184         MPI_Recv(idx_recv, sizes[0] + 1, MPI_INT, recv_rank, 2,
185                  cs_glob_mpi_comm, &status);
186         BFT_MALLOC(fam_recv, sizes[1], int);
187         MPI_Recv(fam_recv, sizes[1], MPI_INT, recv_rank, 3,
188                  cs_glob_mpi_comm, &status);
189 
190         /* Merge data */
191 
192         BFT_MALLOC(idx_new, _n_fam + sizes[0] + 1, int);
193         BFT_MALLOC(fam_new, _family_idx[_n_fam] + sizes[1], int);
194 
195         idx_new[0] = 0;
196         i = 0; j = 0;
197 
198         while (i < _n_fam && j < sizes[0]) {
199           int cmp;
200           int n0 = _family_idx[i+1] - _family_idx[i];
201           int n1 = idx_recv[j+1] - idx_recv[j];
202           cmp = _sync_compare_families(n0,
203                                        n1,
204                                        _family + _family_idx[i],
205                                        fam_recv + idx_recv[j]);
206           if (cmp <= 0) {
207             for (k = 0; k < n0; k++)
208               fam_new[idx_new[n_fam_new] + k] = _family[_family_idx[i] + k];
209             idx_new[n_fam_new + 1] = idx_new[n_fam_new] + n0;
210             n_fam_new += 1;
211             i += 1;
212             if (cmp == 0)
213               j += 1;
214           }
215           else if (cmp > 0) {
216             for (k = 0; k < n1; k++)
217               fam_new[idx_new[n_fam_new] + k] = fam_recv[idx_recv[j] + k];
218             idx_new[n_fam_new + 1] = idx_new[n_fam_new] + n1;
219             n_fam_new += 1;
220             j += 1;
221           }
222         }
223 
224         while (i < _n_fam) {
225           int n0 = _family_idx[i+1] - _family_idx[i];
226           for (k = 0; k < n0; k++)
227             fam_new[idx_new[n_fam_new] + k] = _family[_family_idx[i] + k];
228           idx_new[n_fam_new + 1] = idx_new[n_fam_new] + n0;
229           n_fam_new += 1;
230           i += 1;
231         }
232         while (j < sizes[0]) {
233           int n1 = idx_recv[j+1] - idx_recv[j];
234           for (k = 0; k < n1; k++)
235             fam_new[idx_new[n_fam_new] + k] = fam_recv[idx_recv[j] + k];
236           idx_new[n_fam_new + 1] = idx_new[n_fam_new] + n1;
237           n_fam_new += 1;
238           j += 1;
239         }
240 
241         BFT_FREE(fam_recv);
242         BFT_FREE(idx_recv);
243 
244         BFT_REALLOC(idx_new, n_fam_new + 1, int);
245         BFT_REALLOC(fam_new, idx_new[n_fam_new], int);
246 
247         BFT_FREE(_family_idx);
248         BFT_FREE(_family);
249 
250         _family_idx = idx_new;
251         _family = fam_new;
252         _n_fam = n_fam_new;
253 
254       } /* if (sizes[0] > 0) */
255 
256     }
257 
258     /* Odd ranks send once, then are finished for the merge step */
259 
260     else if (rank_id % 2 == 1) {
261 
262       int send_rank = (rank_id-1)*rank_mult;
263 
264       sizes[0] = _n_fam;
265       sizes[1] = _family_idx[_n_fam];
266       MPI_Send(sizes, 2, MPI_INT, send_rank, 1, cs_glob_mpi_comm);
267 
268       if (sizes[0] > 0) {
269         MPI_Send(_family_idx, sizes[0] + 1, MPI_INT, send_rank, 2,
270                  cs_glob_mpi_comm);
271         MPI_Send(_family, sizes[1], MPI_INT, send_rank, 3, cs_glob_mpi_comm);
272       }
273 
274       break;
275 
276     }
277 
278   }
279 
280   /* Now rank 0 broadcasts */
281 
282   sizes[0] = _n_fam;
283   sizes[1] = _family_idx[_n_fam];
284   MPI_Bcast(sizes, 2, MPI_INT, 0, cs_glob_mpi_comm);
285 
286   _n_fam = sizes[0];
287 
288   if (cs_glob_rank_id != 0) {
289     BFT_REALLOC(_family_idx, sizes[0] + 1, int);
290     BFT_REALLOC(_family, sizes[1], int);
291   }
292 
293   MPI_Bcast(_family_idx, sizes[0] + 1, MPI_INT, 0, cs_glob_mpi_comm);
294   MPI_Bcast(_family, sizes[1], MPI_INT, 0, cs_glob_mpi_comm);
295 
296   /* Finally generate renumbering array */
297 
298   BFT_MALLOC(renum, *n_fam, int);
299 
300   for (i = 0; i < *n_fam; i++) {
301 
302     int start_id, end_id, mid_id;
303     int cmp_ret = 1;
304     int n1 = (*family_idx)[i+1] - (*family_idx)[i];
305 
306     /* Use binary search to find entry */
307 
308     start_id = 0;
309     end_id = _n_fam - 1;
310     mid_id = start_id + ((end_id -start_id) / 2);
311 
312     while (start_id <= end_id) {
313       int n0 = _family_idx[mid_id+1] - _family_idx[mid_id];
314       cmp_ret = _sync_compare_families(n0,
315                                        n1,
316                                        _family + _family_idx[mid_id],
317                                        (*family) + (*family_idx)[i]);
318       if (cmp_ret < 0)
319         start_id = mid_id + 1;
320       else if (cmp_ret > 0)
321         end_id = mid_id - 1;
322       else
323         break;
324       mid_id = start_id + ((end_id -start_id) / 2);
325     }
326 
327     assert(cmp_ret == 0);
328 
329     renum[i] = mid_id;
330   }
331 
332   BFT_FREE(*family_idx);
333   BFT_FREE(*family);
334 
335   *n_fam = _n_fam;
336   *family_idx = _family_idx;
337   *family = _family;
338 
339   return renum;
340 }
341 
342 #endif /* defined(HAVE_MPI) */
343 
344 /*----------------------------------------------------------------------------*/
345 /*!
346  * \brief Descend binary tree for the ordering of a mesh's groups.
347  *
348  * \param[in, out]  mesh   pointer to mesh structure
349  * \param[in]       level  level of the binary tree to descend
350  * \param[in]       n      number of groups in the binary tree to descend
351  * \param[in, out]  order  pre-allocated ordering table
352  */
353 /*----------------------------------------------------------------------------*/
354 
355 inline static void
_groups_descend_tree(const cs_mesh_t * mesh,size_t level,const size_t n,int order[])356 _groups_descend_tree(const cs_mesh_t  *mesh,
357                      size_t            level,
358                      const size_t      n,
359                      int               order[])
360 {
361   size_t i_save, i1, i2, lv_cur;
362 
363   i_save = (size_t)(order[level]);
364 
365   while (level <= (n/2)) {
366 
367     lv_cur = (2*level) + 1;
368 
369     if (lv_cur < n - 1) {
370 
371       i1 = (size_t)(order[lv_cur+1]);
372       i2 = (size_t)(order[lv_cur]);
373 
374       if (strcmp(mesh->group + mesh->group_idx[i1],
375                  mesh->group + mesh->group_idx[i2]) > 0)
376         lv_cur++;
377     }
378 
379     if (lv_cur >= n) break;
380 
381     i1 = i_save;
382     i2 = (size_t)(order[lv_cur]);
383 
384     if (strcmp(mesh->group + mesh->group_idx[i1],
385                mesh->group + mesh->group_idx[i2]) >= 0)
386       break;
387 
388     order[level] = order[lv_cur];
389     level = lv_cur;
390   }
391 
392   order[level] = i_save;
393 }
394 
395 /*----------------------------------------------------------------------------*/
396 /*!
397  * \brief Order mesh groups.
398  *
399  * \param[in, out]  mesh   pointer to mesh structure
400  * \param[out]      order  pre-allocated ordering table
401  */
402 /*----------------------------------------------------------------------------*/
403 
404 static void
_order_groups(const cs_mesh_t * mesh,int order[])405 _order_groups(const cs_mesh_t  *mesh,
406               int               order[])
407 {
408   int    o_save;
409   size_t i;
410   size_t n = mesh->n_groups;
411 
412   /* Initialize ordering array */
413 
414   for (i = 0; i < n; i++)
415     order[i] = i;
416 
417   if (n < 2)
418     return;
419 
420   /* Create binary tree */
421 
422   i = (n / 2);
423   do {
424     i--;
425     _groups_descend_tree(mesh, i, n, order);
426   } while (i > 0);
427 
428   /* Sort binary tree */
429 
430   for (i = n - 1; i > 0; i--) {
431     o_save   = order[0];
432     order[0] = order[i];
433     order[i] = o_save;
434     _groups_descend_tree(mesh, 0, i, order);
435   }
436 }
437 
438 /*----------------------------------------------------------------------------*/
439 /*!
440  * \brief Add mesh group definition if needed.
441  *
442  * \param[in, out]  mesh  pointer to mesh structure to modify
443  * \param[in]       name  group name to add
444  *
445  * \return  id of group
446  */
447 /*----------------------------------------------------------------------------*/
448 
449 static int
_add_group(cs_mesh_t * mesh,const char * name)450 _add_group(cs_mesh_t   *mesh,
451            const char  *name)
452 {
453   /* Check if already present */
454 
455   for (int i = 0; i < mesh->n_groups; i++) {
456     const char *g_cur = mesh->group + mesh->group_idx[i];
457     int d = strcmp(g_cur, name);
458       if (d == 0)
459         return i;
460   }
461 
462   /* Add in case of empty groups */
463 
464   if (mesh->n_groups == 0) {
465     mesh->n_groups = 1;
466     BFT_MALLOC(mesh->group_idx, 2, int);
467     mesh->group_idx[0] = 0;
468     mesh->group_idx[1] = strlen(name) + 1;
469     BFT_MALLOC(mesh->group, mesh->group_idx[1], char);
470     strcpy(mesh->group, name);
471     return 0;
472   }
473 
474   /* Add to end then renumber */
475 
476   else {
477 
478     size_t l = strlen(name) + 1;
479     int n_groups_o = mesh->n_groups;
480     mesh->n_groups += 1;
481     BFT_REALLOC(mesh->group_idx, mesh->n_groups + 1, int);
482     BFT_REALLOC(mesh->group,
483                 mesh->group_idx[n_groups_o] + l,
484                 char);
485     strcpy(mesh->group + mesh->group_idx[n_groups_o], name);
486     mesh->group_idx[mesh->n_groups] = mesh->group_idx[n_groups_o] + l;
487     cs_mesh_group_clean(mesh);
488 
489     for (int i = 0; i < mesh->n_groups; i++) {
490       const char *g_cur = mesh->group + mesh->group_idx[i];
491       int d = strcmp(g_cur, name);
492       if (d == 0)
493         return i;
494     }
495 
496   }
497 
498   /* We should not arrive here */
499 
500   assert(0);
501   return -1;
502 }
503 
504 /*----------------------------------------------------------------------------*/
505 /*!
506  * \brief Define a group class for a group with a given name.
507  *
508  * If a group class with that group an no other is already present,
509  * its id is returned.
510  *
511  * \param[in, out]  mesh   pointer to mesh structure to modify
512  * \param[in]       name   group name to add
513  *
514  * \return          id of matching or added group class
515  */
516 /*----------------------------------------------------------------------------*/
517 
518 static int
_add_gc(cs_mesh_t * mesh,const char * name)519 _add_gc(cs_mesh_t   *mesh,
520         const char  *name)
521 {
522   int g_id = _add_group(mesh, name);
523   int g_id_cmp = -g_id - 1;
524 
525   /* Check if a group class definition with only the selected group
526      exists */
527 
528   int gc_id = -1;
529   int gc_id_shift = 0;
530 
531   if (mesh->n_max_family_items > 1) {
532     for (int i = 0; i < mesh->n_families; i++) {
533       if (   mesh->family_item[i] == g_id_cmp
534           && mesh->family_item[mesh->n_families + i] == 0) {
535         gc_id = i;
536         break;
537       }
538     }
539   }
540 
541   else if (mesh->n_max_family_items == 1) {
542     for (int i = 0; i < mesh->n_families; i++) {
543       if (mesh->family_item[i] == g_id_cmp) {
544         gc_id = i;
545         break;
546       }
547     }
548   }
549 
550   /* Add a family if needed */
551 
552   if (gc_id < 0) {
553 
554     int n_f_prv = mesh->n_families;
555     int *f_prv = NULL;
556 
557     if (n_f_prv*mesh->n_max_family_items > 0) {
558       BFT_MALLOC(f_prv, n_f_prv * mesh->n_max_family_items, int);
559       memcpy(f_prv,
560              mesh->family_item,
561              (n_f_prv * mesh->n_max_family_items)*sizeof(int));
562     }
563 
564     gc_id = mesh->n_families + gc_id_shift;
565 
566     mesh->n_families += 1;
567 
568     if (mesh->n_max_family_items == 0) {
569       mesh->n_max_family_items = 1;
570       BFT_REALLOC(mesh->family_item,
571                   mesh->n_families*mesh->n_max_family_items,
572                   int);
573       for (int i = 0; i < mesh->n_families; i++)
574         mesh->family_item[i] = 0;
575       mesh->family_item[gc_id] = g_id_cmp;
576     }
577     else {
578       BFT_REALLOC(mesh->family_item,
579                   mesh->n_families*mesh->n_max_family_items,
580                   int);
581       for (int j = 0; j < mesh->n_max_family_items; j++) {
582         for (int i = 0; i < n_f_prv; i++)
583           mesh->family_item[mesh->n_families*j+i] = f_prv[n_f_prv*j + i];
584       }
585       mesh->family_item[gc_id - gc_id_shift] = g_id_cmp;
586       for (int j = 1; j < mesh->n_max_family_items; j++)
587         mesh->family_item[mesh->n_families*j + gc_id - gc_id_shift] = 0;
588     }
589 
590     BFT_FREE(f_prv);
591   }
592 
593   return gc_id;
594 }
595 
596 /*----------------------------------------------------------------------------*/
597 /*!
598  * \brief Assign a given group to mesh entities, removing those entities
599  *        from previous groups if present.
600  *
601  * The group is created if necessary.
602  *
603  * \param[in, out]  mesh             pointer to mesh structure to modify
604  * \param[in]       name             group name to assign to selected elements
605  * \param[in]       n_selected_elts  number of selected elements
606  * \param[in]       selected_elt_id  selected element ids (size: n_selected_elts)
607  * \param[in, out]  gc_id            element group class ids (size: n_elts)
608  */
609 /*----------------------------------------------------------------------------*/
610 
611 static void
_mesh_group_set(cs_mesh_t * mesh,const char * name,cs_lnum_t n_selected_elts,const cs_lnum_t selected_elt_id[],int gc_id[])612 _mesh_group_set(cs_mesh_t        *mesh,
613                 const char       *name,
614                 cs_lnum_t         n_selected_elts,
615                 const cs_lnum_t   selected_elt_id[],
616                 int               gc_id[])
617 {
618   int _gc_id = _add_gc(mesh, name);
619 
620   for (cs_lnum_t i = 0; i < n_selected_elts; i++)
621     gc_id[selected_elt_id[i]] = _gc_id + 1;
622 
623   if (mesh->class_defs != NULL)
624     cs_mesh_update_selectors(mesh);
625 }
626 
627 /*----------------------------------------------------------------------------*/
628 /*!
629  * \brief Add selected mesh entities to a given group.
630  *
631  * The group is created if necessary.
632  *
633  * \param[in, out]  mesh             pointer to mesh structure to modify
634  * \param[in]       name             group name to assign to selected elements
635  * \param[in]       n_elts           number of elements in mesh entity
636  * \param[in]       n_selected_elts  number of selected elements
637  * \param[in]       selected_elt_id  selected element ids (size: n_selected_elts)
638  * \param[in, out]  gc_id            element group class ids (size: n_elts)
639  */
640 /*----------------------------------------------------------------------------*/
641 
642 static void
_mesh_group_add(cs_mesh_t * mesh,const char * name,cs_lnum_t n_elts,cs_lnum_t n_selected_elts,const cs_lnum_t selected_elt_id[],int gc_id[])643 _mesh_group_add(cs_mesh_t        *mesh,
644                 const char       *name,
645                 cs_lnum_t         n_elts,
646                 cs_lnum_t         n_selected_elts,
647                 const cs_lnum_t   selected_elt_id[],
648                 int               gc_id[])
649 {
650   int _gc_id = _add_gc(mesh, name);
651 
652   int null_family = 0;
653   if (mesh->n_families > 0) {
654     if (mesh->family_item[0] == 0)
655       null_family = 1;
656   }
657 
658   /* Build index on entities (previous group class for elements
659      not selected, previous + new for those selected) */
660 
661   cs_lnum_t *gc_tmp_idx = NULL;
662   int *gc_tmp = NULL;
663 
664   BFT_MALLOC(gc_tmp_idx, n_elts + 1, cs_lnum_t);
665   gc_tmp_idx[0] = 0;
666 
667   for (cs_lnum_t i = 0; i < n_elts; i++)
668     gc_tmp_idx[i+1] = 1;
669   for (cs_lnum_t i = 0; i < n_selected_elts; i++) {
670     cs_lnum_t j= selected_elt_id[i];
671     if (gc_id[j] != null_family)
672       gc_tmp_idx[j+1] += 1;
673   }
674 
675   /* Transform count to index */
676 
677   for (cs_lnum_t i = 0; i < n_elts; i++)
678     gc_tmp_idx[i+1] += gc_tmp_idx[i];
679 
680   /* Now assign multiple group classes */
681 
682   BFT_MALLOC(gc_tmp, gc_tmp_idx[n_elts], int);
683   for (cs_lnum_t i = 0; i < n_elts; i++)
684     gc_tmp[gc_tmp_idx[i]] = gc_id[i];
685   for (cs_lnum_t i = 0; i < n_selected_elts; i++) {
686     cs_lnum_t j = selected_elt_id[i];
687     if (gc_id[j] != null_family)
688       gc_tmp[gc_tmp_idx[j] + 1] = _gc_id + 1;
689     else
690       gc_tmp[gc_tmp_idx[j]] = _gc_id + 1;
691   }
692 
693   /* Merge definitions */
694 
695   cs_mesh_group_combine_classes(mesh, n_elts, gc_tmp_idx, gc_tmp, gc_id);
696 
697   /* Cleanup */
698 
699   BFT_FREE(gc_tmp_idx);
700   BFT_FREE(gc_tmp);
701 
702   if (mesh->class_defs != NULL)
703     cs_mesh_update_selectors(mesh);
704 }
705 
706 /*! (DOXYGEN_SHOULD_SKIP_THIS) \endcond */
707 
708 /*============================================================================
709  * Public function definitions
710  *===========================================================================*/
711 
712 /*----------------------------------------------------------------------------*/
713 /*!
714  * \brief Clean mesh group definitions.
715  *
716  * \param[in]  mesh  pointer to mesh structure to modify
717  */
718 /*----------------------------------------------------------------------------*/
719 
720 void
cs_mesh_group_clean(cs_mesh_t * mesh)721 cs_mesh_group_clean(cs_mesh_t  *mesh)
722 {
723   int i;
724   size_t j;
725   int n_groups = 0;
726   size_t size_tot = 0;
727   char *g_prev = NULL, *g_cur = NULL, *g_lst = NULL;
728   int *order = NULL, *renum = NULL;
729 
730   if (mesh->n_groups < 1)
731     return;
732 
733   /* Order group names */
734 
735   BFT_MALLOC(renum, mesh->n_groups, int);
736   BFT_MALLOC(order, mesh->n_groups, int);
737 
738   _order_groups(mesh, order);
739 
740   /* Build compact group copy */
741 
742   BFT_MALLOC(g_lst, mesh->group_idx[mesh->n_groups], char);
743 
744   g_cur = mesh->group + mesh->group_idx[order[0]];
745   g_prev = g_cur;
746   n_groups += 1;
747   strcpy(g_lst, g_cur);
748   size_tot += strlen(g_cur) + 1;
749   g_lst[size_tot - 1] = '\0';
750   renum[order[0]] = 0;
751 
752   for (i = 1; i < mesh->n_groups; i++) {
753     g_cur = mesh->group + mesh->group_idx[order[i]];
754     if (strcmp(g_cur, g_prev) != 0) {
755       g_prev = g_cur;
756       strcpy(g_lst + size_tot, g_cur);
757       n_groups += 1;
758       size_tot += strlen(g_cur) + 1;
759       g_lst[size_tot - 1] = '\0';
760     }
761     renum[order[i]] = n_groups - 1;
762   }
763 
764   BFT_FREE(order);
765 
766   BFT_REALLOC(mesh->group_idx, n_groups + 1, int);
767   BFT_REALLOC(mesh->group, size_tot, char);
768 
769   mesh->n_groups = n_groups;
770   memcpy(mesh->group, g_lst, size_tot);
771 
772   mesh->group_idx[0] = 0;
773   for (i = 0; i < mesh->n_groups; i++) {
774     j = strlen(mesh->group + mesh->group_idx[i]) + 1;
775     mesh->group_idx[i + 1] = mesh->group_idx[i] + j;
776   }
777 
778   BFT_FREE(g_lst);
779 
780   /* Now renumber groups in group class description */
781 
782   size_tot = mesh->n_families * mesh->n_max_family_items;
783 
784   for (j = 0; j < size_tot; j++) {
785     int gc_id = mesh->family_item[j];
786     if (gc_id < 0)
787       mesh->family_item[j] = - renum[-gc_id - 1] - 1;
788   }
789 
790   BFT_FREE(renum);
791 
792   /* Remove empty group if present (it should appear first now) */
793 
794   if (mesh->n_groups > 1) {
795 
796     if ((mesh->group_idx[1] - mesh->group_idx[0]) == 1) {
797 
798       size_t new_lst_size = (  mesh->group_idx[mesh->n_groups]
799                              - mesh->group_idx[1]);
800       for (i = 0; i < mesh->n_groups; i++)
801         mesh->group_idx[i] = mesh->group_idx[i+1] - 1;
802       mesh->n_groups -= 1;
803       memmove(mesh->group, mesh->group+1, new_lst_size);
804 
805       BFT_REALLOC(mesh->group_idx, mesh->n_groups + 1, int);
806       BFT_REALLOC(mesh->group, new_lst_size, char);
807 
808       for (j = 0; j < size_tot; j++) {
809         int gc_id = mesh->family_item[j];
810         if (gc_id < 0)
811           mesh->family_item[j] += 1;
812       }
813 
814     }
815   }
816 
817 }
818 
819 /*----------------------------------------------------------------------------*/
820 /*!
821  * \brief Combine mesh group classes.
822  *
823  * \param[in, out]  mesh          pointer to mesh structure to modify
824  * \param[in]       n_elts        number of local elements
825  * \param[in]       gc_id_idx     element group class index (size: n_elts +1)
826  * \param[in]       gc_id         input element group classes
827  *                                (size: gc_id_idx[n_elts])
828  * \param[in]       gc_id_merged  output element group classes (size: n_elts)
829  *
830  * \return  array of new element group class ids
831  */
832 /*----------------------------------------------------------------------------*/
833 
834 void
cs_mesh_group_combine_classes(cs_mesh_t * mesh,cs_lnum_t n_elts,cs_lnum_t gc_id_idx[],int gc_id[],int gc_id_merged[])835 cs_mesh_group_combine_classes(cs_mesh_t   *mesh,
836                               cs_lnum_t    n_elts,
837                               cs_lnum_t    gc_id_idx[],
838                               int          gc_id[],
839                               int          gc_id_merged[])
840 {
841   int   n_fam = 0;
842   int  *_gc_id_idx = NULL;
843   int  *_gc_id = NULL;
844 
845   /* First, build families locally */
846 
847   if (n_elts > 0) {
848 
849     cs_lnum_t i, j, j_prev, n_prev;
850     cs_lnum_t *order = NULL;
851     cs_gnum_t *tmp_gc_id = NULL;
852     const cs_lnum_t n_gc_id_values = gc_id_idx[n_elts];
853 
854     /* Build ordering of elements by associated families */
855 
856     BFT_MALLOC(tmp_gc_id, n_gc_id_values, cs_gnum_t);
857 
858     for (i = 0; i < n_gc_id_values; i++)
859       tmp_gc_id[i] = gc_id[i] + 1;
860 
861     order = cs_order_gnum_i(NULL, tmp_gc_id, gc_id_idx, n_elts);
862 
863     BFT_FREE(tmp_gc_id);
864 
865     /* Build new family array, merging identical definitions */
866 
867     BFT_MALLOC(_gc_id_idx, n_elts + 1, int);
868     BFT_MALLOC(_gc_id, gc_id_idx[n_elts], int);
869 
870     _gc_id_idx[0] = 0;
871 
872     j_prev = -1;
873     n_prev = -1;
874 
875     for (i = 0; i < n_elts; i++) {
876       cs_lnum_t k, l, n;
877       bool is_same = true;
878       j = order[i];
879       n = gc_id_idx[j+1] - gc_id_idx[j];
880       if (n != n_prev)
881         is_same = false;
882       else {
883         for (k = gc_id_idx[j], l = gc_id_idx[j_prev];
884              k < gc_id_idx[j + 1];
885              k++, l++) {
886           if (gc_id[k] != gc_id[l])
887             is_same = false;
888         }
889       }
890       if (is_same)
891         gc_id_merged[j] = gc_id_merged[j_prev];
892       else if (n == 0)
893         gc_id_merged[j] = 0;
894       else if (n == 1)
895         gc_id_merged[j] = gc_id[gc_id_idx[j]];
896       else {
897         gc_id_merged[j] = mesh->n_families + 1 + n_fam;
898         for (k = gc_id_idx[j], l = _gc_id_idx[n_fam];
899              k < gc_id_idx[j + 1];
900              k++, l++)
901           _gc_id[l] = gc_id[k];
902         _gc_id_idx[n_fam+1] = _gc_id_idx[n_fam] + n;
903         n_fam += 1;
904       }
905       j_prev = j;
906       n_prev = n;
907     }
908 
909     BFT_FREE(order);
910 
911     BFT_REALLOC(_gc_id_idx, n_fam + 1, int);
912     BFT_REALLOC(_gc_id, _gc_id_idx[n_fam], int);
913 
914   }
915   else {
916     BFT_MALLOC(_gc_id_idx, 1, int);
917     _gc_id_idx[0] = 0;
918   }
919 
920 #if defined(HAVE_MPI)
921 
922   if (cs_glob_n_ranks > 1) {
923     cs_lnum_t i;
924     int *renum = _sync_family_combinations(&n_fam,
925                                            &_gc_id_idx,
926                                            &_gc_id);
927 
928     for (i = 0; i < n_elts; i++) {
929       int j = gc_id_merged[i] - mesh->n_families - 1;
930       if (j >= 0)
931         gc_id_merged[i] = mesh->n_families + renum[j] + 1;
932     }
933 
934     BFT_FREE(renum);
935   }
936 
937 #endif
938 
939   /* Update mesh definition */
940 
941   {
942     int i, j, k, l;
943     int n_max_family_items = 0;
944 
945     for (i = 0; i < n_fam; i++) {
946       int n_family_items = 0;
947       for (j = _gc_id_idx[i]; j < _gc_id_idx[i+1]; j++) {
948         int f_id = _gc_id[j] - 1;
949         for (k = 0; k < mesh->n_max_family_items; k++) {
950           if (mesh->family_item[mesh->n_families*k + f_id] != 0)
951             n_family_items++;
952         }
953       }
954       if (n_family_items > n_max_family_items)
955         n_max_family_items = n_family_items;
956     }
957 
958     /* Increase maximum number of definitions and pad it necessary */
959 
960     if (n_max_family_items > mesh->n_max_family_items) {
961       BFT_REALLOC(mesh->family_item,
962                   mesh->n_families*n_max_family_items,
963                   int);
964       for (i = mesh->n_max_family_items;
965            i < n_max_family_items;
966            i++) {
967         for (j = 0; j < mesh->n_families; j++)
968           mesh->family_item[mesh->n_families*i + j] = 0;
969       }
970       mesh->n_max_family_items = n_max_family_items;
971     }
972 
973     /* Increase number of families */
974 
975     mesh->n_families += n_fam;
976 
977     BFT_REALLOC(mesh->family_item,
978                 mesh->n_families * mesh->n_max_family_items,
979                 int);
980     for (j = mesh->n_max_family_items - 1; j > 0; j--) {
981       for (i = mesh->n_families - n_fam - 1; i > -1; i--)
982         mesh->family_item[mesh->n_families*j + i]
983           = mesh->family_item[(mesh->n_families - n_fam)*j + i];
984     }
985     for (i = mesh->n_families - n_fam, j = 0; i < mesh->n_families; i++, j++) {
986       int n_family_items = 0;
987       for (k = _gc_id_idx[j]; k < _gc_id_idx[j+1]; k++) {
988         int f_id = _gc_id[k] - 1;
989         for (l = 0; l < mesh->n_max_family_items; l++) {
990           if (mesh->family_item[mesh->n_families*l + f_id] != 0) {
991             mesh->family_item[mesh->n_families*n_family_items + i]
992               = mesh->family_item[mesh->n_families*l + f_id];
993             n_family_items++;
994           }
995         }
996       }
997       for (k = n_family_items; k < mesh->n_max_family_items; k++)
998         mesh->family_item[mesh->n_families*k + i] = 0;
999     }
1000 
1001   }
1002 
1003   BFT_FREE(_gc_id_idx);
1004   BFT_FREE(_gc_id);
1005 }
1006 
1007 /*----------------------------------------------------------------------------*/
1008 /*!
1009  * \brief Assign a given group to cells, removing those entities
1010  *        from previous groups if present.
1011  *
1012  * The group is created if necessary.
1013  *
1014  * \param[in, out]  mesh              pointer to mesh structure to modify
1015  * \param[in]       name              group name to assign to selected cells
1016  * \param[in]       n_selected_cells  number of selected cells
1017  * \param[in]       selected_cell_id  selected cell ids (size: n_selected_cells)
1018  */
1019 /*----------------------------------------------------------------------------*/
1020 
1021 void
cs_mesh_group_cells_set(cs_mesh_t * mesh,const char * name,cs_lnum_t n_selected_cells,const cs_lnum_t selected_cell_id[])1022 cs_mesh_group_cells_set(cs_mesh_t        *mesh,
1023                         const char       *name,
1024                         cs_lnum_t         n_selected_cells,
1025                         const cs_lnum_t   selected_cell_id[])
1026 {
1027   _mesh_group_set(mesh,
1028                   name,
1029                   n_selected_cells,
1030                   selected_cell_id,
1031                   mesh->cell_family);
1032 }
1033 
1034 /*----------------------------------------------------------------------------*/
1035 /*!
1036  * \brief Assign a given group to interior faces, removing those entities
1037  *        from previous groups if present.
1038  *
1039  * The group is created if necessary.
1040  *
1041  * \param[in, out]  mesh              pointer to mesh structure to modify
1042  * \param[in]       name              group name to assign to selected faces
1043  * \param[in]       n_selected_faces  number of selected faces
1044  * \param[in]       selected_face_id  selected face ids (size: n_selected_faces)
1045  */
1046 /*----------------------------------------------------------------------------*/
1047 
1048 void
cs_mesh_group_i_faces_set(cs_mesh_t * mesh,const char * name,cs_lnum_t n_selected_faces,const cs_lnum_t selected_face_id[])1049 cs_mesh_group_i_faces_set(cs_mesh_t        *mesh,
1050                           const char       *name,
1051                           cs_lnum_t         n_selected_faces,
1052                           const cs_lnum_t   selected_face_id[])
1053 {
1054   _mesh_group_set(mesh,
1055                   name,
1056                   n_selected_faces,
1057                   selected_face_id,
1058                   mesh->i_face_family);
1059 }
1060 
1061 /*----------------------------------------------------------------------------*/
1062 /*!
1063  * \brief Assign a given group to boundary faces, removing those entities
1064  *        from previous groups if present.
1065  *
1066  * The group is created if necessary.
1067  *
1068  * \param[in, out]  mesh              pointer to mesh structure to modify
1069  * \param[in]       name              group name to assign to selected faces
1070  * \param[in]       n_selected_faces  number of selected faces
1071  * \param[in]       selected_face_id  selected face ids (size: n_selected_faces)
1072  */
1073 /*----------------------------------------------------------------------------*/
1074 
1075 void
cs_mesh_group_b_faces_set(cs_mesh_t * mesh,const char * name,cs_lnum_t n_selected_faces,const cs_lnum_t selected_face_id[])1076 cs_mesh_group_b_faces_set(cs_mesh_t        *mesh,
1077                           const char       *name,
1078                           cs_lnum_t         n_selected_faces,
1079                           const cs_lnum_t   selected_face_id[])
1080 {
1081   _mesh_group_set(mesh,
1082                   name,
1083                   n_selected_faces,
1084                   selected_face_id,
1085                   mesh->b_face_family);
1086 }
1087 
1088 /*----------------------------------------------------------------------------*/
1089 /*!
1090  * \brief Add selected cells to a given group.
1091  *
1092  * The group is created if necessary.
1093  *
1094  * \param[in, out]  mesh              pointer to mesh structure to modify
1095  * \param[in]       name              group name to assign to selected cells
1096  * \param[in]       n_selected_cells  number of selected cells
1097  * \param[in]       selected_cell_id  selected cell ids (size: n_selected_cells)
1098  */
1099 /*----------------------------------------------------------------------------*/
1100 
1101 void
cs_mesh_group_cells_add(cs_mesh_t * mesh,const char * name,cs_lnum_t n_selected_cells,const cs_lnum_t selected_cell_id[])1102 cs_mesh_group_cells_add(cs_mesh_t        *mesh,
1103                         const char       *name,
1104                         cs_lnum_t         n_selected_cells,
1105                         const cs_lnum_t   selected_cell_id[])
1106 {
1107   _mesh_group_add(mesh,
1108                   name,
1109                   mesh->n_cells,
1110                   n_selected_cells,
1111                   selected_cell_id,
1112                   mesh->cell_family);
1113 }
1114 
1115 /*----------------------------------------------------------------------------*/
1116 /*!
1117  * \brief Add selected interior faces to a given group.
1118  *
1119  * The group is created if necessary.
1120  *
1121  * \param[in, out]  mesh              pointer to mesh structure to modify
1122  * \param[in]       name              group name to assign to selected faces
1123  * \param[in]       n_selected_faces  number of selected faces
1124  * \param[in]       selected_face_id  selected face ids (size: n_selected_faces)
1125  */
1126 /*----------------------------------------------------------------------------*/
1127 
1128 void
cs_mesh_group_i_faces_add(cs_mesh_t * mesh,const char * name,cs_lnum_t n_selected_faces,const cs_lnum_t selected_face_id[])1129 cs_mesh_group_i_faces_add(cs_mesh_t        *mesh,
1130                           const char       *name,
1131                           cs_lnum_t         n_selected_faces,
1132                           const cs_lnum_t   selected_face_id[])
1133 {
1134   _mesh_group_add(mesh,
1135                   name,
1136                   mesh->n_i_faces,
1137                   n_selected_faces,
1138                   selected_face_id,
1139                   mesh->i_face_family);
1140 }
1141 
1142 /*----------------------------------------------------------------------------*/
1143 /*!
1144  * \brief Add selected boundary faces to a given group.
1145  *
1146  * The group is created if necessary.
1147  *
1148  * \param[in, out]  mesh              pointer to mesh structure to modify
1149  * \param[in]       name              group name to assign to selected faces
1150  * \param[in]       n_selected_faces  number of selected faces
1151  * \param[in]       selected_face_id  selected face ids (size: n_selected_faces)
1152  */
1153 /*----------------------------------------------------------------------------*/
1154 
1155 void
cs_mesh_group_b_faces_add(cs_mesh_t * mesh,const char * name,cs_lnum_t n_selected_faces,const cs_lnum_t selected_face_id[])1156 cs_mesh_group_b_faces_add(cs_mesh_t        *mesh,
1157                           const char       *name,
1158                           cs_lnum_t         n_selected_faces,
1159                           const cs_lnum_t   selected_face_id[])
1160 {
1161   _mesh_group_add(mesh,
1162                   name,
1163                   mesh->n_b_faces,
1164                   n_selected_faces,
1165                   selected_face_id,
1166                   mesh->b_face_family);
1167 }
1168 
1169 /*---------------------------------------------------------------------------*/
1170 
1171 END_C_DECLS
1172