1 /*============================================================================
2 * \file Define cs_mesh_builder_t fields from cs_mesh_t fields.
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 <assert.h>
34 #include <stdio.h>
35 #include <stdlib.h>
36 #include <string.h>
37 #include <limits.h>
38 #include <math.h>
39 #include <assert.h>
40
41 /*----------------------------------------------------------------------------
42 * Local headers
43 *----------------------------------------------------------------------------*/
44
45 #include "bft_mem.h"
46 #include "bft_printf.h"
47
48 #include "fvm_io_num.h"
49
50 #include "cs_base.h"
51 #include "cs_interface.h"
52 #include "cs_io.h"
53 #include "cs_mesh.h"
54 #include "cs_mesh_builder.h"
55 #include "cs_order.h"
56
57 /*----------------------------------------------------------------------------
58 * Header for the current file
59 *----------------------------------------------------------------------------*/
60
61 #include "cs_mesh_to_builder.h"
62
63 /*----------------------------------------------------------------------------*/
64
65 BEGIN_C_DECLS
66
67 /*! \cond DOXYGEN_SHOULD_SKIP_THIS */
68
69 /*=============================================================================
70 * Local Macro Definitions
71 *============================================================================*/
72
73 /*============================================================================
74 * Local structure definitions
75 *============================================================================*/
76
77 /*=============================================================================
78 * Local Macro definitions
79 *============================================================================*/
80
81 /*============================================================================
82 * Static global variables
83 *============================================================================*/
84
85 /*============================================================================
86 * Private function definitions
87 *============================================================================*/
88
89 #if defined(HAVE_MPI)
90
91 /*----------------------------------------------------------------------------
92 * Define mesh builder's face vertices arrays in parallel mode.
93 *
94 * parameters:
95 * mesh <-- pointer to mesh structure
96 * mb <-> pointer to mesh builder
97 * d <-> face part to block distribution helper
98 *----------------------------------------------------------------------------*/
99
100 static void
_mesh_to_builder_face_vertices_g(const cs_mesh_t * mesh,cs_mesh_builder_t * mb,cs_part_to_block_t * d)101 _mesh_to_builder_face_vertices_g(const cs_mesh_t *mesh,
102 cs_mesh_builder_t *mb,
103 cs_part_to_block_t *d)
104 {
105 cs_lnum_t i, j, k;
106 cs_gnum_t block_size, n_block_faces, n_face_vertices;
107
108 cs_lnum_t *face_vtx_idx = NULL;
109 cs_gnum_t *face_vtx_g = NULL;
110
111 const cs_lnum_t n_i_faces = mesh->n_i_faces;
112 const cs_lnum_t n_b_faces = mesh->n_b_faces;
113 const cs_lnum_t n_faces = mesh->n_i_faces + mesh->n_b_faces;
114
115 const cs_datatype_t gnum_type
116 = (sizeof(cs_gnum_t) == 8) ? CS_UINT64 : CS_UINT32;
117
118 /* Face -> vertex connectivity */
119 /*-----------------------------*/
120
121 BFT_FREE(mb->face_vertices_idx);
122 BFT_FREE(mb->face_vertices);
123
124 BFT_MALLOC(mb->face_vertices_idx,
125 (mb->face_bi.gnum_range[1] - mb->face_bi.gnum_range[0]) + 1,
126 cs_lnum_t);
127
128 BFT_MALLOC(face_vtx_idx, n_faces + 1, cs_lnum_t);
129
130 face_vtx_idx[0] = 0;
131 for (i = 0; i < n_i_faces; i++) {
132 n_face_vertices = mesh->i_face_vtx_idx[i+1] - mesh->i_face_vtx_idx[i];
133 face_vtx_idx[i+1] = face_vtx_idx[i] + n_face_vertices;
134 }
135 for (i = 0, j = n_i_faces; i < n_b_faces; i++, j++) {
136 n_face_vertices = mesh->b_face_vtx_idx[i+1] - mesh->b_face_vtx_idx[i];
137 face_vtx_idx[j+1] = face_vtx_idx[j] + n_face_vertices;
138 }
139
140 cs_part_to_block_copy_index(d, face_vtx_idx, mb->face_vertices_idx);
141
142 /* Build connectivity */
143
144 BFT_MALLOC(face_vtx_g,
145 mesh->i_face_vtx_connect_size + mesh->b_face_vtx_connect_size,
146 cs_gnum_t);
147
148 k = 0;
149 for (i = 0; i < n_i_faces; i++) {
150 for (j = mesh->i_face_vtx_idx[i]; j < mesh->i_face_vtx_idx[i+1]; j++)
151 face_vtx_g[k++]
152 = mesh->global_vtx_num[mesh->i_face_vtx_lst[j]];
153 }
154 for (i = 0; i < n_b_faces; i++) {
155 for (j = mesh->b_face_vtx_idx[i]; j < mesh->b_face_vtx_idx[i+1]; j++)
156 face_vtx_g[k++]
157 = mesh->global_vtx_num[mesh->b_face_vtx_lst[j]];
158 }
159
160 n_block_faces = (mb->face_bi.gnum_range[1] - mb->face_bi.gnum_range[0]);
161 block_size = mb->face_vertices_idx[n_block_faces];
162
163 BFT_MALLOC(mb->face_vertices, block_size, cs_gnum_t);
164
165 cs_part_to_block_copy_indexed(d,
166 gnum_type,
167 face_vtx_idx,
168 face_vtx_g,
169 mb->face_vertices_idx,
170 mb->face_vertices);
171
172 BFT_FREE(face_vtx_g);
173 BFT_FREE(face_vtx_idx);
174 }
175
176 /*----------------------------------------------------------------------------
177 * Write face->vertices connectivity in parallel.
178 *
179 * parameters:
180 * mesh <-- pointer to mesh structure
181 * mb <-> mesh builder
182 * transfer <-- if true, mesh transferred to builder; if false, builder
183 * is a temporary copy
184 * pp_out <-- pointer to output file
185 *----------------------------------------------------------------------------*/
186
187 static void
_write_face_vertices_g(const cs_mesh_t * mesh,cs_mesh_builder_t * mb,bool transfer,cs_io_t * pp_out)188 _write_face_vertices_g(const cs_mesh_t *mesh,
189 cs_mesh_builder_t *mb,
190 bool transfer,
191 cs_io_t *pp_out)
192 {
193 cs_lnum_t i;
194 cs_gnum_t block_size, g_vtx_connect_size, n_block_faces, n_face_vertices;
195 cs_gnum_t idx_range[4];
196
197 cs_gnum_t *face_vtx_idx_g = NULL;
198
199 const cs_gnum_t n_g_faces = mesh->n_g_i_faces + mesh->n_g_b_faces;
200 const cs_datatype_t gnum_type
201 = (sizeof(cs_gnum_t) == 8) ? CS_UINT64 : CS_UINT32;
202
203 /* Face -> vertex connectivity */
204 /*-----------------------------*/
205
206 /* Copy index from block to global values */
207
208 n_block_faces = (mb->face_bi.gnum_range[1] - mb->face_bi.gnum_range[0]);
209 BFT_MALLOC(face_vtx_idx_g, n_block_faces+ 1, cs_gnum_t);
210
211 block_size = mb->face_vertices_idx[n_block_faces];
212 MPI_Scan(&block_size, face_vtx_idx_g, 1, CS_MPI_GNUM, MPI_SUM,
213 cs_glob_mpi_comm);
214 face_vtx_idx_g[0] -= block_size;
215 face_vtx_idx_g[0] += 1;
216
217 for (i = 0; i < (cs_lnum_t)n_block_faces; i++) {
218 n_face_vertices = mb->face_vertices_idx[i+1] - mb->face_vertices_idx[i];
219 face_vtx_idx_g[i+1] = face_vtx_idx_g[i] + n_face_vertices;
220 }
221
222 idx_range[0] = mb->face_bi.gnum_range[0];
223 idx_range[1] = mb->face_bi.gnum_range[1];
224 if (mb->face_bi.gnum_range[0] >= n_g_faces) {
225 idx_range[0] += 1;
226 idx_range[1] += 1;
227 }
228 else if (mb->face_bi.gnum_range[1] >= n_g_faces + 1)
229 idx_range[1] += 1;
230
231 /* Set idx_range and compute global size for next write
232 with indexed values */
233
234 idx_range[2] = face_vtx_idx_g[0];
235 idx_range[3] = face_vtx_idx_g[0] + block_size;
236
237 MPI_Allreduce(face_vtx_idx_g + n_block_faces, &g_vtx_connect_size, 1,
238 CS_MPI_GNUM, MPI_MAX, cs_glob_mpi_comm);
239 g_vtx_connect_size -= 1;
240
241 /* Now write buffer */
242
243 cs_io_write_block_buffer("face_vertices_index",
244 n_g_faces + 1,
245 idx_range[0],
246 idx_range[1],
247 2, /* location_id, */
248 1, /* index id */
249 1, /* n_location_vals */
250 gnum_type,
251 face_vtx_idx_g,
252 pp_out);
253
254 BFT_FREE(face_vtx_idx_g);
255
256 /* Build connectivity */
257
258 if (transfer == true)
259 cs_io_write_block("face_vertices",
260 g_vtx_connect_size,
261 idx_range[2],
262 idx_range[3],
263 2, /* location_id, */
264 1, /* index id */
265 1, /* n_location_vals */
266 gnum_type,
267 mb->face_vertices,
268 pp_out);
269 else
270 cs_io_write_block_buffer("face_vertices",
271 g_vtx_connect_size,
272 idx_range[2],
273 idx_range[3],
274 2, /* location_id, */
275 1, /* index id */
276 1, /* n_location_vals */
277 gnum_type,
278 mb->face_vertices,
279 pp_out);
280 }
281
282 /*----------------------------------------------------------------------------
283 * Transfer or save mesh data to builder in parallel.
284 *
285 * parameters:
286 * mesh <-> pointer to mesh structure
287 * mb <-> mesh builder
288 * transfer <-- if true, data is transferred from mesh to builder;
289 * if false, builder fields are only used as a temporary
290 * arrays.
291 * pp_out <-> optional output file, or NULL
292 *----------------------------------------------------------------------------*/
293
294 static void
_mesh_to_builder_g(cs_mesh_t * mesh,cs_mesh_builder_t * mb,bool transfer,cs_io_t * pp_out)295 _mesh_to_builder_g(cs_mesh_t *mesh,
296 cs_mesh_builder_t *mb,
297 bool transfer,
298 cs_io_t *pp_out)
299 {
300 cs_lnum_t i, j;
301
302 cs_part_to_block_t *d = NULL;
303
304 int *face_gc_id = NULL;
305 cs_gnum_t *cell_gnum = NULL, *face_gnum = NULL;
306 cs_gnum_t *face_cell_g = NULL;
307
308 const cs_lnum_t n_i_faces = mesh->n_i_faces;
309 const cs_lnum_t n_b_faces = mesh->n_b_faces;
310 const cs_lnum_t n_faces = mesh->n_i_faces + mesh->n_b_faces;
311
312 const cs_datatype_t real_type
313 = (sizeof(cs_real_t) == 8) ? CS_DOUBLE : CS_FLOAT;
314 const cs_datatype_t gnum_type
315 = (sizeof(cs_gnum_t) == 8) ? CS_UINT64 : CS_UINT32;
316 const cs_datatype_t int_type
317 = (sizeof(int) == 8) ? CS_INT64 : CS_INT32;
318
319 /* Distribute cell group class info to blocks */
320 /*---------------------------------------------*/
321
322 BFT_MALLOC(mb->cell_gc_id,
323 (mb->cell_bi.gnum_range[1] - mb->cell_bi.gnum_range[0]),
324 int);
325
326 d = cs_part_to_block_create_by_gnum(cs_glob_mpi_comm,
327 mb->cell_bi,
328 mesh->n_cells,
329 mesh->global_cell_num);
330
331 cs_part_to_block_copy_array(d,
332 int_type,
333 1,
334 mesh->cell_family,
335 mb->cell_gc_id);
336
337 cs_part_to_block_destroy(&d);
338
339 if (transfer == true)
340 BFT_FREE(mesh->cell_family);
341
342 /* Build global face part to block distribution structures */
343 /*---------------------------------------------------------*/
344
345 BFT_MALLOC(face_gnum, n_faces, cs_gnum_t);
346
347 for (i = 0; i < n_i_faces; i++)
348 face_gnum[i] = mesh->global_i_face_num[i];
349 for (i = 0, j = n_i_faces; i < n_b_faces; i++, j++)
350 face_gnum[j] = mesh->global_b_face_num[i] + mesh->n_g_i_faces;
351
352 if (transfer == true) {
353 BFT_FREE(mesh->global_i_face_num);
354 BFT_FREE(mesh->global_b_face_num);
355 }
356
357 d = cs_part_to_block_create_by_gnum(cs_glob_mpi_comm,
358 mb->face_bi,
359 n_faces,
360 face_gnum);
361 cs_part_to_block_transfer_gnum(d, face_gnum);
362 face_gnum = NULL;
363
364 /* Face -> cell connectivity */
365 /*---------------------------*/
366
367 /* Build global cell numbering including parallel halos,
368 except for periodic values */
369
370 cell_gnum = cs_mesh_get_cell_gnum(mesh, 1);
371
372 BFT_MALLOC(face_cell_g, n_faces*2, cs_gnum_t);
373
374 for (i = 0; i < n_i_faces; i++) {
375 face_cell_g[i*2] = cell_gnum[mesh->i_face_cells[i][0]];
376 face_cell_g[i*2 + 1] = cell_gnum[mesh->i_face_cells[i][1]];
377 }
378 for (i = 0, j = n_i_faces; i < n_b_faces; i++, j++) {
379 face_cell_g[j*2] = cell_gnum[mesh->b_face_cells[i]];
380 face_cell_g[j*2 + 1] = 0;
381 }
382
383 BFT_FREE(cell_gnum);
384
385 if (transfer == true) {
386 BFT_FREE(mesh->global_cell_num);
387 BFT_FREE(mesh->i_face_cells);
388 BFT_FREE(mesh->b_face_cells);
389 }
390
391 /* Distribute to blocks and write */
392
393 BFT_MALLOC(mb->face_cells,
394 (mb->face_bi.gnum_range[1] - mb->face_bi.gnum_range[0]) * 2,
395 cs_gnum_t);
396
397 cs_part_to_block_copy_array(d,
398 gnum_type,
399 2,
400 face_cell_g,
401 mb->face_cells);
402
403 BFT_FREE(face_cell_g);
404
405 if (pp_out != NULL) {
406 if (transfer == true)
407 cs_io_write_block("face_cells",
408 mb->n_g_faces,
409 mb->face_bi.gnum_range[0],
410 mb->face_bi.gnum_range[1],
411 2, /* location_id, */
412 0, /* index id */
413 2, /* n_location_vals */
414 gnum_type,
415 mb->face_cells,
416 pp_out);
417 else
418 cs_io_write_block_buffer("face_cells",
419 mb->n_g_faces,
420 mb->face_bi.gnum_range[0],
421 mb->face_bi.gnum_range[1],
422 2, /* location_id, */
423 0, /* index id */
424 2, /* n_location_vals */
425 gnum_type,
426 mb->face_cells,
427 pp_out);
428 }
429
430 if (transfer == false)
431 BFT_FREE(mb->face_cells);
432
433 /* Now write pre-distributed blocks for cell group classes if required */
434
435 if (pp_out != NULL) {
436 if (transfer == true)
437 cs_io_write_block("cell_group_class_id",
438 mesh->n_g_cells,
439 mb->cell_bi.gnum_range[0],
440 mb->cell_bi.gnum_range[1],
441 1, /* location_id, */
442 0, /* index id */
443 1, /* n_location_vals */
444 int_type,
445 mb->cell_gc_id,
446 pp_out);
447 else
448 cs_io_write_block_buffer("cell_group_class_id",
449 mesh->n_g_cells,
450 mb->cell_bi.gnum_range[0],
451 mb->cell_bi.gnum_range[1],
452 1, /* location_id, */
453 0, /* index id */
454 1, /* n_location_vals */
455 int_type,
456 mb->cell_gc_id,
457 pp_out);
458 }
459
460 if (transfer == false)
461 BFT_FREE(mb->cell_gc_id);
462
463 /* Face group classes */
464 /*--------------------*/
465
466 BFT_MALLOC(mb->face_gc_id,
467 (mb->face_bi.gnum_range[1] - mb->face_bi.gnum_range[0]),
468 int);
469
470 BFT_MALLOC(face_gc_id, n_faces, int);
471
472 for (i = 0; i < n_i_faces; i++)
473 face_gc_id[i] = mesh->i_face_family[i];
474 for (i = 0, j = n_i_faces; i < n_b_faces; i++, j++)
475 face_gc_id[j] = mesh->b_face_family[i];
476
477 if (transfer == true) {
478 BFT_FREE(mesh->i_face_family);
479 BFT_FREE(mesh->b_face_family);
480 }
481
482 /* Distribute to blocks, write if required */
483
484 cs_part_to_block_copy_array(d,
485 int_type,
486 1,
487 face_gc_id,
488 mb->face_gc_id);
489
490 BFT_FREE(face_gc_id);
491
492 if (pp_out != NULL) {
493 if (transfer == true)
494 cs_io_write_block("face_group_class_id",
495 mb->n_g_faces,
496 mb->face_bi.gnum_range[0],
497 mb->face_bi.gnum_range[1],
498 2, /* location_id, */
499 0, /* index id */
500 1, /* n_location_vals */
501 int_type,
502 mb->face_gc_id,
503 pp_out);
504 else
505 cs_io_write_block_buffer("face_group_class_id",
506 mb->n_g_faces,
507 mb->face_bi.gnum_range[0],
508 mb->face_bi.gnum_range[1],
509 2, /* location_id, */
510 0, /* index id */
511 1, /* n_location_vals */
512 int_type,
513 mb->face_gc_id,
514 pp_out);
515 }
516
517 if (transfer == false)
518 BFT_FREE(mb->face_gc_id);
519
520 /* Face refinement generation */
521 /*----------------------------*/
522
523 if (mb->have_face_r_gen) {
524
525 char *face_r_gen = NULL;
526
527 BFT_MALLOC(mb->face_r_gen,
528 (mb->face_bi.gnum_range[1] - mb->face_bi.gnum_range[0]),
529 char);
530
531 BFT_MALLOC(face_r_gen, n_faces, char);
532
533 for (i = 0; i < n_i_faces; i++)
534 face_r_gen[i] = mesh->i_face_r_gen[i];
535 for (i = 0, j = n_i_faces; i < n_b_faces; i++, j++)
536 face_r_gen[j] = 0;
537
538 /* Distribute to blocks, write if required */
539
540 cs_part_to_block_copy_array(d,
541 CS_CHAR,
542 1,
543 face_r_gen,
544 mb->face_r_gen);
545
546 BFT_FREE(face_r_gen);
547
548 if (pp_out != NULL) {
549 if (transfer == true)
550 cs_io_write_block("face_refinement_generation",
551 mb->n_g_faces,
552 mb->face_bi.gnum_range[0],
553 mb->face_bi.gnum_range[1],
554 2, /* location_id, */
555 0, /* index id */
556 1, /* n_location_vals */
557 CS_CHAR,
558 mb->face_r_gen,
559 pp_out);
560 else
561 cs_io_write_block_buffer("face_refinement_generation",
562 mb->n_g_faces,
563 mb->face_bi.gnum_range[0],
564 mb->face_bi.gnum_range[1],
565 2, /* location_id, */
566 0, /* index id */
567 1, /* n_location_vals */
568 CS_CHAR,
569 mb->face_r_gen,
570 pp_out);
571
572 }
573 }
574
575 if (transfer == true)
576 BFT_FREE(mesh->i_face_r_gen);
577
578 if (transfer == false)
579 BFT_FREE(mb->face_r_gen);
580
581 /* Face -> vertex connectivity */
582 /*-----------------------------*/
583
584 _mesh_to_builder_face_vertices_g(mesh, mb, d);
585
586 if (transfer == true) {
587 BFT_FREE(mesh->i_face_vtx_idx);
588 BFT_FREE(mesh->i_face_vtx_lst);
589 BFT_FREE(mesh->b_face_vtx_idx);
590 BFT_FREE(mesh->b_face_vtx_lst);
591 }
592
593 if (pp_out != NULL)
594 _write_face_vertices_g(mesh, mb, transfer, pp_out);
595
596 if (transfer == false) {
597 BFT_FREE(mb->face_vertices_idx);
598 BFT_FREE(mb->face_vertices);
599 }
600
601 /* Free face part to block distribution structures */
602
603 cs_part_to_block_destroy(&d);
604
605 /* Vertex coordinates */
606 /*--------------------*/
607
608 d = cs_part_to_block_create_by_gnum(cs_glob_mpi_comm,
609 mb->vertex_bi,
610 mesh->n_vertices,
611 mesh->global_vtx_num);
612
613 BFT_MALLOC(mb->vertex_coords,
614 (mb->vertex_bi.gnum_range[1] - mb->vertex_bi.gnum_range[0]) * 3,
615 cs_real_t);
616
617 cs_part_to_block_copy_array(d,
618 real_type,
619 3,
620 mesh->vtx_coord,
621 mb->vertex_coords);
622
623 if (transfer == true)
624 BFT_FREE(mesh->vtx_coord);
625
626 if (pp_out != NULL) {
627 if (transfer == true)
628 cs_io_write_block("vertex_coords",
629 mesh->n_g_vertices,
630 mb->vertex_bi.gnum_range[0],
631 mb->vertex_bi.gnum_range[1],
632 3, /* location_id, */
633 0, /* index id */
634 3, /* n_location_vals */
635 real_type,
636 mb->vertex_coords,
637 pp_out);
638 else
639 cs_io_write_block_buffer("vertex_coords",
640 mesh->n_g_vertices,
641 mb->vertex_bi.gnum_range[0],
642 mb->vertex_bi.gnum_range[1],
643 3, /* location_id, */
644 0, /* index id */
645 3, /* n_location_vals */
646 real_type,
647 mb->vertex_coords,
648 pp_out);
649 }
650
651 if (transfer == true)
652 BFT_FREE(mesh->global_vtx_num);
653 else
654 BFT_FREE(mb->vertex_coords);
655
656 if (transfer == true)
657 BFT_FREE(mesh->global_vtx_num);
658 else
659 BFT_FREE(mb->vertex_coords);
660
661 cs_part_to_block_destroy(&d);
662 }
663
664 #endif /* defined(HAVE_MPI) */
665
666 /*----------------------------------------------------------------------------
667 * Transfer or save mesh data to builder in serial mode.
668 *
669 * parameters:
670 * mesh <-> pointer to mesh structure
671 * mb <-> mesh builder
672 * transfer <-- if true, data is transferred from mesh to builder;
673 * if false, builder fields are only used as a temporary
674 * arrays.
675 * pp_out <-> optional output file, or NULL
676 *----------------------------------------------------------------------------*/
677
678 static void
_mesh_to_builder_l(cs_mesh_t * mesh,cs_mesh_builder_t * mb,bool transfer,cs_io_t * pp_out)679 _mesh_to_builder_l(cs_mesh_t *mesh,
680 cs_mesh_builder_t *mb,
681 bool transfer,
682 cs_io_t *pp_out)
683 {
684 cs_lnum_t i, j, k, l, n_face_vertices, cell_id_0, cell_id_1;
685 cs_gnum_t g_vtx_connect_size;
686
687 cs_lnum_t *order = NULL, *i_order = NULL, *b_order = NULL;
688
689 const cs_lnum_t n_i_faces = mesh->n_i_faces;
690 const cs_lnum_t n_b_faces = mesh->n_b_faces;
691 const cs_lnum_t n_faces = mesh->n_i_faces + mesh->n_b_faces;
692
693 const cs_datatype_t real_type
694 = (sizeof(cs_real_t) == 8) ? CS_DOUBLE : CS_FLOAT;
695 const cs_datatype_t gnum_type
696 = (sizeof(cs_gnum_t) == 8) ? CS_UINT64 : CS_UINT32;
697 const cs_datatype_t int_type
698 = (sizeof(int) == 8) ? CS_INT64 : CS_INT32;
699
700 /* Face ordering */
701
702 i_order = cs_order_gnum(NULL, mesh->global_i_face_num, mesh->n_i_faces);
703 if (mesh->n_b_faces > 0)
704 b_order = cs_order_gnum(NULL, mesh->global_b_face_num, mesh->n_b_faces);
705
706 if (transfer == true) {
707 BFT_FREE(mesh->global_i_face_num);
708 BFT_FREE(mesh->global_b_face_num);
709 }
710
711 /* Face -> cell connectivity (excluding periodic values )*/
712
713 BFT_MALLOC(mb->face_cells, n_faces*2, cs_gnum_t);
714
715 if (mesh->global_cell_num != NULL) {
716
717 for (i = 0; i < n_i_faces; i++) {
718 k = i_order[i];
719 cell_id_0 = mesh->i_face_cells[k][0];
720 cell_id_1 = mesh->i_face_cells[k][1];
721 if (cell_id_0 < mesh->n_cells)
722 mb->face_cells[i*2] = mesh->global_cell_num[cell_id_0];
723 else
724 mb->face_cells[i*2] = 0;
725 if (cell_id_1 < mesh->n_cells)
726 mb->face_cells[i*2 + 1] = mesh->global_cell_num[cell_id_1];
727 else
728 mb->face_cells[i*2 + 1] = 0;
729 }
730 for (i = 0, j = n_i_faces; i < n_b_faces; i++, j++) {
731 k = b_order[i];
732 mb->face_cells[j*2] = mesh->global_cell_num[mesh->b_face_cells[k]];
733 mb->face_cells[j*2 + 1] = 0;
734 }
735
736 if (transfer == true)
737 BFT_FREE(mesh->global_cell_num);
738
739 }
740 else { /* if (mesh->global_cell_num == NULL) */
741
742 for (i = 0; i < n_i_faces; i++) {
743 k = i_order[i];
744 cell_id_0 = mesh->i_face_cells[k][0];
745 cell_id_1 = mesh->i_face_cells[k][1];
746 if (cell_id_0 < mesh->n_cells)
747 mb->face_cells[i*2] = cell_id_0 + 1;
748 else
749 mb->face_cells[i*2] = 0;
750 if (cell_id_1 < mesh->n_cells)
751 mb->face_cells[i*2 + 1] = cell_id_1 + 1;
752 else
753 mb->face_cells[i*2 + 1] = 0;
754 }
755 for (i = 0, j = n_i_faces; i < n_b_faces; i++, j++) {
756 k = b_order[i];
757 mb->face_cells[j*2] = mesh->b_face_cells[k] + 1;
758 mb->face_cells[j*2 + 1] = 0;
759 }
760
761 }
762
763 if (transfer == true) {
764 BFT_FREE(mesh->i_face_cells);
765 BFT_FREE(mesh->b_face_cells);
766 }
767
768 if (pp_out != NULL) {
769 if (transfer == true)
770 cs_io_write_block("face_cells",
771 mb->n_g_faces,
772 1,
773 mb->n_g_faces + 1,
774 2, /* location_id, */
775 0, /* index id */
776 2, /* n_location_vals */
777 gnum_type,
778 mb->face_cells,
779 pp_out);
780 else
781 cs_io_write_block_buffer("face_cells",
782 mb->n_g_faces,
783 1,
784 mb->n_g_faces + 1,
785 2, /* location_id, */
786 0, /* index id */
787 2, /* n_location_vals */
788 gnum_type,
789 mb->face_cells,
790 pp_out);
791 }
792
793 if (transfer == false)
794 BFT_FREE(mb->face_cells);
795
796 /* Cell group classes */
797
798 BFT_MALLOC(mb->cell_gc_id, mesh->n_cells, int);
799
800 order = cs_order_gnum(NULL, mesh->global_cell_num, mesh->n_cells);
801 for (i = 0; i < mesh->n_cells; i++)
802 mb->cell_gc_id[i] = mesh->cell_family[order[i]];
803 BFT_FREE(order);
804
805 if (transfer == true)
806 BFT_FREE(mesh->cell_family);
807
808 if (pp_out != NULL) {
809 if (transfer == true)
810 cs_io_write_block("cell_group_class_id",
811 mesh->n_g_cells,
812 1,
813 mesh->n_g_cells + 1,
814 1, /* location_id, */
815 0, /* index id */
816 1, /* n_location_vals */
817 int_type,
818 mb->cell_gc_id,
819 pp_out);
820 else
821 cs_io_write_block_buffer("cell_group_class_id",
822 mesh->n_g_cells,
823 1,
824 mesh->n_g_cells + 1,
825 1, /* location_id, */
826 0, /* index id */
827 1, /* n_location_vals */
828 int_type,
829 mb->cell_gc_id,
830 pp_out);
831 }
832
833 if (transfer == false)
834 BFT_FREE(mb->cell_gc_id);
835
836 /* Face group classes */
837
838 BFT_MALLOC(mb->face_gc_id, n_faces, int);
839
840 for (i = 0; i < n_i_faces; i++)
841 mb->face_gc_id[i] = mesh->i_face_family[i_order[i]];
842 for (i = 0, j = n_i_faces; i < n_b_faces; i++, j++)
843 mb->face_gc_id[j] = mesh->b_face_family[b_order[i]];
844
845 if (transfer == true) {
846 BFT_FREE(mesh->i_face_family);
847 BFT_FREE(mesh->b_face_family);
848 }
849
850 if (pp_out != NULL) {
851 if (transfer == true)
852 cs_io_write_block("face_group_class_id",
853 mb->n_g_faces,
854 1,
855 mb->n_g_faces + 1,
856 2, /* location_id, */
857 0, /* index id */
858 1, /* n_location_vals */
859 int_type,
860 mb->face_gc_id,
861 pp_out);
862 else
863 cs_io_write_block_buffer("face_group_class_id",
864 mb->n_g_faces,
865 1,
866 mb->n_g_faces + 1,
867 2, /* location_id, */
868 0, /* index id */
869 1, /* n_location_vals */
870 int_type,
871 mb->face_gc_id,
872 pp_out);
873 }
874
875 if (transfer == false)
876 BFT_FREE(mb->face_gc_id);
877
878 /* Face refinement generation */
879
880 if (mb->have_face_r_gen) {
881
882 BFT_MALLOC(mb->face_r_gen, n_faces, char);
883
884 for (i = 0; i < n_i_faces; i++)
885 mb->face_r_gen[i] = mesh->i_face_r_gen[i_order[i]];
886 for (i = 0, j = n_i_faces; i < n_b_faces; i++, j++)
887 mb->face_r_gen[j] = 0;
888
889 if (transfer == true)
890 BFT_FREE(mesh->i_face_r_gen);
891
892 if (pp_out != NULL) {
893 if (transfer == true)
894 cs_io_write_block("face_refinement_generation",
895 mb->n_g_faces,
896 1,
897 mb->n_g_faces + 1,
898 2, /* location_id, */
899 0, /* index id */
900 1, /* n_location_vals */
901 CS_CHAR,
902 mb->face_r_gen,
903 pp_out);
904 else
905 cs_io_write_block_buffer("face_refinement_generation",
906 mb->n_g_faces,
907 1,
908 mb->n_g_faces + 1,
909 2, /* location_id, */
910 0, /* index id */
911 1, /* n_location_vals */
912 CS_CHAR,
913 mb->face_r_gen,
914 pp_out);
915 }
916
917 if (transfer == false)
918 BFT_FREE(mb->face_r_gen);
919 }
920
921 /* Face -> vertex connectivity */
922 /*-----------------------------*/
923
924 if (pp_out != NULL) {
925
926 cs_gnum_t *face_vtx_idx_g = NULL;
927
928 BFT_MALLOC(face_vtx_idx_g, n_faces + 1, cs_gnum_t);
929
930 face_vtx_idx_g[0] = 1;
931 for (i = 0; i < n_i_faces; i++) {
932 k = i_order[i];
933 n_face_vertices = mesh->i_face_vtx_idx[k+1] - mesh->i_face_vtx_idx[k];
934 face_vtx_idx_g[i+1] = face_vtx_idx_g[i] + n_face_vertices;
935 }
936 for (i = 0, j = n_i_faces; i < n_b_faces; i++, j++) {
937 k = b_order[i];
938 n_face_vertices = mesh->b_face_vtx_idx[k+1] - mesh->b_face_vtx_idx[k];
939 face_vtx_idx_g[j+1] = face_vtx_idx_g[j] + n_face_vertices;
940 }
941
942 /* Now write buffer */
943
944 cs_io_write_block_buffer("face_vertices_index",
945 mb->n_g_faces + 1,
946 1,
947 n_faces + 2,
948 2, /* location_id, */
949 1, /* index id */
950 1, /* n_location_vals */
951 gnum_type,
952 face_vtx_idx_g,
953 pp_out);
954
955 BFT_FREE(face_vtx_idx_g);
956
957 }
958
959 if (transfer == true) {
960
961 BFT_MALLOC(mb->face_vertices_idx, n_faces + 1, cs_lnum_t);
962
963 mb->face_vertices_idx[0] = 0;
964 for (i = 0; i < n_i_faces; i++) {
965 k = i_order[i];
966 n_face_vertices = mesh->i_face_vtx_idx[k+1] - mesh->i_face_vtx_idx[k];
967 mb->face_vertices_idx[i+1] = mb->face_vertices_idx[i] + n_face_vertices;
968 }
969 for (i = 0, j = n_i_faces; i < n_b_faces; i++, j++) {
970 k = b_order[i];
971 n_face_vertices = mesh->b_face_vtx_idx[k+1] - mesh->b_face_vtx_idx[k];
972 mb->face_vertices_idx[j+1] = mb->face_vertices_idx[j] + n_face_vertices;
973 }
974
975 }
976
977 /* Build connectivity */
978
979 g_vtx_connect_size = ( mesh->i_face_vtx_connect_size
980 + mesh->b_face_vtx_connect_size);
981
982 BFT_MALLOC(mb->face_vertices, g_vtx_connect_size, cs_gnum_t);
983
984 if (mesh->global_vtx_num != NULL) {
985 l = 0;
986 for (i = 0; i < n_i_faces; i++) {
987 k = i_order[i];
988 for (j = mesh->i_face_vtx_idx[k]; j < mesh->i_face_vtx_idx[k+1]; j++)
989 mb->face_vertices[l++]
990 = mesh->global_vtx_num[mesh->i_face_vtx_lst[j - 1] - 1];
991 }
992 for (i = 0; i < n_b_faces; i++) {
993 k = b_order[i];
994 for (j = mesh->b_face_vtx_idx[k]; j < mesh->b_face_vtx_idx[k+1]; j++)
995 mb->face_vertices[l++]
996 = mesh->global_vtx_num[mesh->b_face_vtx_lst[j]];
997 }
998 }
999 else { /* if (mesh->global_vtx_num == NULL) */
1000 l = 0;
1001 for (i = 0; i < n_i_faces; i++) {
1002 k = i_order[i];
1003 for (j = mesh->i_face_vtx_idx[k]; j < mesh->i_face_vtx_idx[k+1]; j++)
1004 mb->face_vertices[l++] = mesh->i_face_vtx_lst[j] + 1;
1005 }
1006 for (i = 0; i < n_b_faces; i++) {
1007 k = b_order[i];
1008 for (j = mesh->b_face_vtx_idx[k]; j < mesh->b_face_vtx_idx[k+1]; j++)
1009 mb->face_vertices[l++] = mesh->b_face_vtx_lst[j] + 1;
1010 }
1011 }
1012
1013 BFT_FREE(b_order);
1014 BFT_FREE(i_order);
1015
1016 if (transfer == true) {
1017 BFT_FREE(mesh->i_face_vtx_idx);
1018 BFT_FREE(mesh->i_face_vtx_lst);
1019 BFT_FREE(mesh->b_face_vtx_idx);
1020 BFT_FREE(mesh->b_face_vtx_lst);
1021 }
1022
1023 if (pp_out != NULL) {
1024 if (transfer == true)
1025 cs_io_write_block("face_vertices",
1026 g_vtx_connect_size,
1027 1,
1028 g_vtx_connect_size + 1,
1029 2, /* location_id, */
1030 1, /* index id */
1031 1, /* n_location_vals */
1032 gnum_type,
1033 mb->face_vertices,
1034 pp_out);
1035 else
1036 cs_io_write_block_buffer("face_vertices",
1037 g_vtx_connect_size,
1038 1,
1039 g_vtx_connect_size + 1,
1040 2, /* location_id, */
1041 1, /* index id */
1042 1, /* n_location_vals */
1043 gnum_type,
1044 mb->face_vertices,
1045 pp_out);
1046 }
1047
1048 if (transfer == false)
1049 BFT_FREE(mb->face_vertices);
1050
1051 /* Vertex coordinates */
1052 /*--------------------*/
1053
1054 BFT_MALLOC(mb->vertex_coords, mesh->n_vertices*3, cs_real_t);
1055
1056 order = cs_order_gnum(NULL, mesh->global_vtx_num, mesh->n_vertices);
1057 for (i = 0; i < mesh->n_vertices; i++) {
1058 j = order[i];
1059 for (k = 0; k <3; k++)
1060 mb->vertex_coords[i*3 + k] = mesh->vtx_coord[j*3 + k];
1061 }
1062 BFT_FREE(order);
1063
1064 if (transfer == true)
1065 BFT_FREE(mesh->vtx_coord);
1066
1067 if (pp_out != NULL) {
1068 if (transfer == true)
1069 cs_io_write_block("vertex_coords",
1070 mesh->n_g_vertices,
1071 1,
1072 mesh->n_vertices + 1,
1073 3, /* location_id, */
1074 0, /* index id */
1075 3, /* n_location_vals */
1076 real_type,
1077 mb->vertex_coords,
1078 pp_out);
1079 else
1080 cs_io_write_block_buffer("vertex_coords",
1081 mesh->n_g_vertices,
1082 1,
1083 mesh->n_vertices + 1,
1084 3, /* location_id, */
1085 0, /* index id */
1086 3, /* n_location_vals */
1087 real_type,
1088 mb->vertex_coords,
1089 pp_out);
1090 }
1091
1092 if (transfer == true)
1093 BFT_FREE(mesh->global_vtx_num);
1094 else
1095 BFT_FREE(mb->vertex_coords);
1096 }
1097
1098 /*----------------------------------------------------------------------------
1099 * Save a mesh as preprocessor data.
1100 *
1101 * parameters:
1102 * mesh <-- pointer to mesh structure
1103 * mb <-- pointer to optional mesh builder structure, or NULL
1104 * pp_out <-> output file
1105 *----------------------------------------------------------------------------*/
1106
1107 static void
_write_dimensions(cs_mesh_t * mesh,cs_mesh_builder_t * mb,cs_io_t * pp_out)1108 _write_dimensions(cs_mesh_t *mesh,
1109 cs_mesh_builder_t *mb,
1110 cs_io_t *pp_out)
1111 {
1112 cs_lnum_t i;
1113
1114 const cs_gnum_t n_g_faces = mesh->n_g_i_faces + mesh->n_g_b_faces;
1115 const cs_datatype_t gnum_type
1116 = (sizeof(cs_gnum_t) == 8) ? CS_UINT64 : CS_UINT32;
1117 const cs_datatype_t lnum_type
1118 = (sizeof(cs_lnum_t) == 8) ? CS_INT64 : CS_INT32;
1119
1120 /* Write headers */
1121 /*---------------*/
1122
1123 cs_io_write_global("start_block:dimensions", 0, 0, 0, 0, CS_DATATYPE_NULL,
1124 NULL, pp_out);
1125
1126 cs_io_write_global("n_cells", 1, 1, 0, 1, gnum_type,
1127 &(mesh->n_g_cells), pp_out);
1128
1129 cs_io_write_global("n_faces", 1, 2, 0, 1, gnum_type,
1130 &n_g_faces, pp_out);
1131
1132 cs_io_write_global("n_vertices", 1, 3, 0, 1, gnum_type,
1133 &(mesh->n_g_vertices), pp_out);
1134
1135 cs_io_write_global("face_vertices_size", 1, 0, 0, 1, gnum_type,
1136 &(mb->n_g_face_connect_size), pp_out);
1137
1138 cs_io_write_global("n_group_classes", 1, 0, 0, 1, lnum_type,
1139 &(mesh->n_families), pp_out);
1140
1141 cs_io_write_global("n_group_class_props_max", 1, 0, 0, 1, lnum_type,
1142 &(mesh->n_max_family_items), pp_out);
1143
1144 if (mesh->n_groups > 0) {
1145
1146 cs_io_write_global("n_groups", 1, 0, 0, 1, lnum_type,
1147 &(mesh->n_groups), pp_out);
1148
1149 cs_lnum_t *_group_idx;
1150 BFT_MALLOC(_group_idx, mesh->n_groups + 1, cs_lnum_t);
1151 for (i = 0; i < mesh->n_groups + 1; i++)
1152 _group_idx[i] = mesh->group_idx[i] + 1;
1153 cs_io_write_global("group_name_index",
1154 mesh->n_groups + 1, 0, 0, 1, lnum_type,
1155 _group_idx, pp_out);
1156 BFT_FREE(_group_idx);
1157
1158 cs_io_write_global("group_name",
1159 mesh->group_idx[mesh->n_groups], 0, 0, 1, CS_CHAR,
1160 mesh->group, pp_out);
1161 }
1162
1163 cs_io_write_global("group_class_properties",
1164 (mesh->n_families * mesh->n_max_family_items), 0, 0, 1,
1165 CS_INT_TYPE,
1166 mesh->family_item, pp_out);
1167
1168 if (mesh->n_init_perio > 0) {
1169
1170 int n_rot_perio = 0;
1171
1172 /* Count rotation periodicities */
1173
1174 for (i = 0; i < mesh->n_init_perio; i++) {
1175 if ( fvm_periodicity_get_type(mesh->periodicity, i*2)
1176 >= FVM_PERIODICITY_ROTATION)
1177 n_rot_perio += 1;
1178 }
1179
1180 /* Output periodicity information */
1181
1182 cs_io_write_global("n_periodic_directions", 1, 0, 0, 1, lnum_type,
1183 &(mesh->n_init_perio), pp_out);
1184
1185 cs_io_write_global("n_periodic_rotations", 1, 0, 0, 1, lnum_type,
1186 &n_rot_perio, pp_out);
1187
1188 for (i = 0; i < mesh->n_init_perio; i++)
1189 mb->n_g_per_face_couples[i] *= 2;
1190
1191 cs_io_write_global("n_periodic_faces", mesh->n_init_perio, 0, 0, 1,
1192 gnum_type,
1193 mb->n_g_per_face_couples, pp_out);
1194
1195 for (i = 0; i < mesh->n_init_perio; i++)
1196 mb->n_g_per_face_couples[i] /= 2;
1197
1198 }
1199
1200 cs_io_write_global("end_block:dimensions", 0, 0, 0, 0, CS_DATATYPE_NULL,
1201 NULL, pp_out);
1202 }
1203
1204 /*----------------------------------------------------------------------------
1205 * Write mesh periodicity metadata.
1206 *
1207 * parameters:
1208 * mesh <-- pointer to mesh structure
1209 * perio_num <-- periodicity num
1210 * pp_out <-> output file
1211 *----------------------------------------------------------------------------*/
1212
1213 static void
_write_mesh_perio_metadata(const cs_mesh_t * mesh,int perio_num,cs_io_t * pp_out)1214 _write_mesh_perio_metadata(const cs_mesh_t *mesh,
1215 int perio_num,
1216 cs_io_t *pp_out)
1217 {
1218 char section_name[32];
1219
1220 const cs_datatype_t lnum_type
1221 = (sizeof(cs_lnum_t) == 8) ? CS_INT64 : CS_INT32;
1222
1223 double matrix[3][4];
1224 cs_lnum_t perio_type = 0;
1225
1226 const int tr_id = (perio_num - 1)*2;
1227 const fvm_periodicity_t *perio = mesh->periodicity;
1228
1229 assert(perio != NULL);
1230
1231 /* Get periodicity type and matrix */
1232
1233 perio_type = fvm_periodicity_get_type(perio, tr_id);
1234
1235 fvm_periodicity_get_matrix(perio, tr_id, matrix);
1236
1237 /* Write global data */
1238
1239 sprintf(section_name, "periodicity_type_%02d", perio_num);
1240
1241 cs_io_write_global(section_name, 1, 0, 0, 1, lnum_type,
1242 &perio_type, pp_out);
1243
1244 sprintf(section_name, "periodicity_matrix_%02d", perio_num);
1245
1246 cs_io_write_global(section_name, 12, 0, 0, 1, CS_DOUBLE,
1247 matrix, pp_out);
1248 }
1249
1250 #if defined(HAVE_MPI)
1251
1252 /*----------------------------------------------------------------------------
1253 * Write mesh periodicity data in parallel.
1254 *
1255 * parameters:
1256 * perio_num <-- periodicity number
1257 * n_perio_couples <-- number of periodic face couples for this periodicity
1258 * perio_couples <-> periodic face couples for this periodicity
1259 * min_rank_step <-- minimum rank step between blocks
1260 * transfer <-- if true, mesh transferred to builder;
1261 * if false, builder is a temporary copy
1262 * pp_out <-> output file
1263 *----------------------------------------------------------------------------*/
1264
1265 static void
_write_mesh_perio_data_g(int perio_num,cs_lnum_t n_perio_couples,cs_gnum_t perio_couples[],int min_rank_step,bool transfer,cs_io_t * pp_out)1266 _write_mesh_perio_data_g(int perio_num,
1267 cs_lnum_t n_perio_couples,
1268 cs_gnum_t perio_couples[],
1269 int min_rank_step,
1270 bool transfer,
1271 cs_io_t *pp_out)
1272 {
1273 char section_name[32];
1274 cs_block_dist_info_t bi;
1275
1276 cs_gnum_t n_g_couples = 0;
1277 cs_gnum_t *_perio_couples = NULL;
1278 const cs_gnum_t *couple_g_num = NULL;
1279 cs_part_to_block_t *d = NULL;
1280
1281 const cs_datatype_t gnum_type
1282 = (sizeof(cs_gnum_t) == 8) ? CS_UINT64 : CS_UINT32;
1283
1284 /* Create global couple numbering */
1285
1286 fvm_io_num_t *c_io_num = fvm_io_num_create_from_adj_s(NULL,
1287 perio_couples,
1288 n_perio_couples,
1289 2);
1290
1291 n_g_couples = fvm_io_num_get_global_count(c_io_num);
1292 couple_g_num = fvm_io_num_get_global_num(c_io_num);
1293
1294 /* Create associated block info and distribution */
1295
1296 bi = cs_block_dist_compute_sizes(cs_glob_rank_id,
1297 cs_glob_n_ranks,
1298 min_rank_step,
1299 0,
1300 n_g_couples);
1301
1302 d = cs_part_to_block_create_by_gnum(cs_glob_mpi_comm,
1303 bi,
1304 n_perio_couples,
1305 couple_g_num);
1306
1307 BFT_MALLOC(_perio_couples,
1308 (bi.gnum_range[1] - bi.gnum_range[0]) * 2,
1309 cs_gnum_t);
1310
1311 cs_part_to_block_copy_array(d,
1312 gnum_type,
1313 2,
1314 perio_couples,
1315 _perio_couples);
1316
1317 cs_part_to_block_destroy(&d);
1318
1319 c_io_num = fvm_io_num_destroy(c_io_num);
1320
1321 /* Write face couples */
1322
1323 sprintf(section_name, "periodicity_faces_%02d", perio_num);
1324
1325 if (transfer == true)
1326 cs_io_write_block(section_name,
1327 n_g_couples,
1328 bi.gnum_range[0],
1329 bi.gnum_range[1],
1330 0, /* location_id, */
1331 0, /* index id */
1332 2, /* n_location_vals */
1333 gnum_type,
1334 _perio_couples,
1335 pp_out);
1336 else
1337 cs_io_write_block_buffer(section_name,
1338 n_g_couples,
1339 bi.gnum_range[0],
1340 bi.gnum_range[1],
1341 0, /* location_id, */
1342 0, /* index id */
1343 2, /* n_location_vals */
1344 gnum_type,
1345 _perio_couples,
1346 pp_out);
1347
1348 BFT_FREE(_perio_couples);
1349 }
1350
1351 #endif /* defined(HAVE_MPI) */
1352
1353 /*----------------------------------------------------------------------------
1354 * Write mesh periodicity data in single-processor mode.
1355 *
1356 * parameters:
1357 * perio_num <-- periodicity number
1358 * n_perio_couples <-- number of periodic face couples for this periodicity
1359 * perio_couples <-> periodic face couples for this periodicity
1360 * transfer <-- if true, mesh transferred to builder;
1361 * if false, builder is a temporary copy
1362 * pp_out <-> output file
1363 *----------------------------------------------------------------------------*/
1364
1365 static void
_write_mesh_perio_data_l(int perio_num,cs_lnum_t n_perio_couples,cs_gnum_t perio_couples[],bool transfer,cs_io_t * pp_out)1366 _write_mesh_perio_data_l(int perio_num,
1367 cs_lnum_t n_perio_couples,
1368 cs_gnum_t perio_couples[],
1369 bool transfer,
1370 cs_io_t *pp_out)
1371 {
1372 char section_name[32];
1373
1374 const cs_datatype_t gnum_type
1375 = (sizeof(cs_gnum_t) == 8) ? CS_UINT64 : CS_UINT32;
1376
1377 /* Write face couples */
1378
1379 sprintf(section_name, "periodicity_faces_%02d", perio_num);
1380
1381 if (transfer == true)
1382 cs_io_write_block(section_name,
1383 n_perio_couples,
1384 1,
1385 n_perio_couples + 1,
1386 0, /* location_id, */
1387 0, /* index id */
1388 2, /* n_location_vals */
1389 gnum_type,
1390 perio_couples,
1391 pp_out);
1392 else
1393 cs_io_write_block_buffer(section_name,
1394 n_perio_couples,
1395 1,
1396 n_perio_couples + 1,
1397 0, /* location_id, */
1398 0, /* index id */
1399 2, /* n_location_vals */
1400 gnum_type,
1401 perio_couples,
1402 pp_out);
1403 }
1404
1405 /*! (DOXYGEN_SHOULD_SKIP_THIS) \endcond */
1406
1407 /*=============================================================================
1408 * Public function definitions
1409 *============================================================================*/
1410
1411 /*----------------------------------------------------------------------------*/
1412 /*!
1413 * \brief Transfer mesh to mesh builder structure.
1414 *
1415 * As the dataflow is very similar, but may be done array-by array to minimize
1416 * memory overhead, this function also handles a part of the output
1417 * to file needed to save a mesh file.
1418 *
1419 * \param[in, out] mesh pointer to mesh structure
1420 * \param[in, out] mb pointer to mesh builder structure
1421 * \param[in] transfer if true, data is transferred from mesh to builder;
1422 * if false, builder fields are only used as a
1423 * temporary arrays.
1424 * \param[in, out] pp_out optional output file, or NULL
1425 */
1426 /*----------------------------------------------------------------------------*/
1427
1428 void
cs_mesh_to_builder(cs_mesh_t * mesh,cs_mesh_builder_t * mb,bool transfer,cs_io_t * pp_out)1429 cs_mesh_to_builder(cs_mesh_t *mesh,
1430 cs_mesh_builder_t *mb,
1431 bool transfer,
1432 cs_io_t *pp_out)
1433 {
1434 int i;
1435 cs_gnum_t g_i_face_vertices_size = 0, g_b_face_vertices_size = 0;
1436
1437 /* Clear rebuildable mesh data in case of transfer;
1438 halos are kept at this stage, as they are required for this algorithm */
1439
1440 if (transfer)
1441 cs_mesh_free_rebuildable(mesh, false);
1442
1443 /* Clear previous builder data if present (periodicity done separately) */
1444
1445 mb->have_cell_rank = false;
1446
1447 BFT_FREE(mb->face_cells);
1448 BFT_FREE(mb->face_vertices_idx);
1449 BFT_FREE(mb->face_vertices);
1450 BFT_FREE(mb->cell_gc_id);
1451 BFT_FREE(mb->face_gc_id);
1452 BFT_FREE(mb->face_r_gen);
1453 BFT_FREE(mb->vertex_coords);
1454
1455 /* Precompute some sizes */
1456
1457 mb->n_g_faces = mesh->n_g_i_faces + mesh->n_g_b_faces;
1458
1459 cs_mesh_builder_define_block_dist(mb,
1460 cs_glob_rank_id,
1461 cs_glob_n_ranks,
1462 mb->min_rank_step,
1463 0,
1464 mesh->n_g_cells,
1465 mb->n_g_faces,
1466 mesh->n_g_vertices);
1467
1468 cs_mesh_g_face_vertices_sizes(mesh,
1469 &g_i_face_vertices_size,
1470 &g_b_face_vertices_size);
1471
1472 mb->n_g_face_connect_size = g_i_face_vertices_size + g_b_face_vertices_size;
1473
1474 /* Get refinement info if needed */
1475
1476 int r_flag = 0;
1477 for (cs_lnum_t j = 0; j < mesh->n_i_faces; j++) {
1478 if (mesh->i_face_r_gen[j] != 0) {
1479 r_flag = 1;
1480 break;
1481 }
1482 }
1483
1484 #if defined(HAVE_MPI)
1485 if (cs_glob_n_ranks > 1) {
1486 int _r_flag = r_flag;
1487 MPI_Allreduce(&_r_flag, &r_flag, 1,
1488 MPI_INT, MPI_MAX, cs_glob_mpi_comm);
1489 }
1490 #endif
1491
1492 mb->have_face_r_gen = (r_flag) ? true : false;
1493
1494 /* Get periodic faces information if required */
1495
1496 cs_mesh_to_builder_perio_faces(mesh, mb);
1497
1498 /* Write metadata if output is required */
1499
1500 if (pp_out != NULL)
1501 _write_dimensions(mesh, mb, pp_out);
1502
1503 /* Main mesh data */
1504
1505 if (pp_out != NULL)
1506 cs_io_write_global("start_block:data", 0, 0, 0, 0, CS_DATATYPE_NULL,
1507 NULL, pp_out);
1508
1509 #if defined(HAVE_MPI)
1510
1511 if (cs_glob_n_ranks > 1)
1512 _mesh_to_builder_g(mesh, mb, transfer, pp_out);
1513
1514 #endif
1515
1516 if (cs_glob_n_ranks == 1)
1517 _mesh_to_builder_l(mesh, mb, transfer, pp_out);
1518
1519 /* Finish clearing rebuildable mesh data in case of transfer */
1520
1521 if (transfer) {
1522 if (mesh->vtx_interfaces != NULL)
1523 cs_interface_set_destroy(&(mesh->vtx_interfaces));
1524 if (mesh->halo != NULL)
1525 cs_halo_destroy(&(mesh->halo));
1526 }
1527
1528 /* Periodicity data */
1529
1530 if (pp_out != NULL) {
1531
1532 for (i = 0; i < mesh->n_init_perio; i++) {
1533
1534 _write_mesh_perio_metadata(mesh, i+1, pp_out);
1535
1536 #if defined(HAVE_MPI)
1537 if (cs_glob_n_ranks > 1)
1538 _write_mesh_perio_data_g(i+1,
1539 mb->n_per_face_couples[i],
1540 mb->per_face_couples[i],
1541 mb->min_rank_step,
1542 transfer,
1543 pp_out);
1544 #endif
1545
1546 if (cs_glob_n_ranks == 1)
1547 _write_mesh_perio_data_l(i+1,
1548 mb->n_per_face_couples[i],
1549 mb->per_face_couples[i],
1550 transfer,
1551 pp_out);
1552
1553 if (transfer == false)
1554 BFT_FREE(mb->per_face_couples[i]);
1555
1556 }
1557
1558 cs_io_write_global("end_block:data", 0, 0, 0, 0, CS_DATATYPE_NULL,
1559 NULL, pp_out);
1560
1561 }
1562 }
1563
1564 /*----------------------------------------------------------------------------*/
1565 /*!
1566 * \brief Transfer mesh partitioning info to mesh builder structure.
1567 *
1568 * \param[in] mesh pointer to mesh structure
1569 * \param[in, out] mb pointer to mesh builder structure
1570 */
1571 /*----------------------------------------------------------------------------*/
1572
1573 void
cs_mesh_to_builder_partition(const cs_mesh_t * mesh,cs_mesh_builder_t * mb)1574 cs_mesh_to_builder_partition(const cs_mesh_t *mesh,
1575 cs_mesh_builder_t *mb)
1576 {
1577 #if defined(HAVE_MPI)
1578
1579 if (cs_glob_n_ranks > 1) {
1580
1581 const cs_datatype_t int_type
1582 = (sizeof(int) == 8) ? CS_INT64 : CS_INT32;
1583
1584 /* Distribute cell group class info to blocks */
1585 /*---------------------------------------------*/
1586
1587 mb->cell_bi = cs_block_dist_compute_sizes(cs_glob_rank_id,
1588 cs_glob_n_ranks,
1589 mb->min_rank_step,
1590 0,
1591 mesh->n_g_cells);
1592
1593 mb->have_cell_rank = true;
1594 BFT_REALLOC(mb->cell_rank,
1595 (mb->cell_bi.gnum_range[1] - mb->cell_bi.gnum_range[0]),
1596 int);
1597
1598 int *cell_rank;
1599 BFT_MALLOC(cell_rank, mesh->n_cells, int);
1600 for (cs_lnum_t i = 0; i < mesh->n_cells; i++)
1601 cell_rank[i] = cs_glob_rank_id;
1602
1603 cs_part_to_block_t *d
1604 = cs_part_to_block_create_by_gnum(cs_glob_mpi_comm,
1605 mb->cell_bi,
1606 mesh->n_cells,
1607 mesh->global_cell_num);
1608
1609 cs_part_to_block_copy_array(d,
1610 int_type,
1611 1,
1612 cell_rank,
1613 mb->cell_rank);
1614
1615 cs_part_to_block_destroy(&d);
1616
1617 BFT_FREE(cell_rank);
1618
1619 }
1620
1621 #endif
1622 }
1623
1624 /*----------------------------------------------------------------------------*/
1625 /*!
1626 * \brief Reconstruct periodic faces info from mesh to builder.
1627 *
1628 * \param[in] mesh pointer to mesh structure
1629 * \param[in, out] mb pointer to mesh builder structure
1630 */
1631 /*----------------------------------------------------------------------------*/
1632
1633 void
cs_mesh_to_builder_perio_faces(const cs_mesh_t * mesh,cs_mesh_builder_t * mb)1634 cs_mesh_to_builder_perio_faces(const cs_mesh_t *mesh,
1635 cs_mesh_builder_t *mb)
1636 {
1637 cs_lnum_t i;
1638
1639 /* Get periodic faces information if required */
1640
1641 mb->n_perio = mesh->n_init_perio;
1642
1643 if (mesh->n_init_perio < 1)
1644 return;
1645
1646 cs_mesh_get_perio_faces(mesh,
1647 &(mb->n_per_face_couples),
1648 &(mb->per_face_couples));
1649
1650 /* In parallel, each periodic couple returned by
1651 cs_mesh_get_perio_faces() should appear on one rank only,
1652 so the global number of couples is simply the sum over all
1653 ranks of the local counts. */
1654
1655 BFT_MALLOC(mb->n_g_per_face_couples, mesh->n_init_perio, cs_gnum_t);
1656
1657 #if defined(HAVE_MPI)
1658
1659 if (cs_glob_n_ranks > 1) {
1660
1661 cs_gnum_t *_n_l_perio_faces = NULL;
1662
1663 BFT_MALLOC(_n_l_perio_faces, mesh->n_init_perio, cs_gnum_t);
1664
1665 for (i = 0; i < mesh->n_init_perio; i++)
1666 _n_l_perio_faces[i] = mb->n_per_face_couples[i];
1667
1668 MPI_Allreduce(_n_l_perio_faces, mb->n_g_per_face_couples, mesh->n_init_perio,
1669 CS_MPI_GNUM, MPI_SUM, cs_glob_mpi_comm);
1670
1671 BFT_FREE(_n_l_perio_faces);
1672 }
1673
1674 #endif
1675
1676 if (cs_glob_n_ranks == 1) {
1677 for (i = 0; i < mesh->n_init_perio; i++)
1678 mb->n_g_per_face_couples[i] = mb->n_per_face_couples[i];
1679 }
1680 }
1681
1682 /*----------------------------------------------------------------------------*/
1683
1684 END_C_DECLS
1685