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