1 /*============================================================================
2 * Write a nodal representation associated with a mesh and associated
3 * variables to EnSight Gold files
4 *============================================================================*/
5
6 /*
7 This file is part of Code_Saturne, a general-purpose CFD tool.
8
9 Copyright (C) 1998-2021 EDF S.A.
10
11 This program is free software; you can redistribute it and/or modify it under
12 the terms of the GNU General Public License as published by the Free Software
13 Foundation; either version 2 of the License, or (at your option) any later
14 version.
15
16 This program is distributed in the hope that it will be useful, but WITHOUT
17 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
18 FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
19 details.
20
21 You should have received a copy of the GNU General Public License along with
22 this program; if not, write to the Free Software Foundation, Inc., 51 Franklin
23 Street, Fifth Floor, Boston, MA 02110-1301, USA.
24 */
25
26 /*----------------------------------------------------------------------------*/
27
28 #include "cs_defs.h"
29
30 /*----------------------------------------------------------------------------
31 * Standard C library headers
32 *----------------------------------------------------------------------------*/
33
34 #include <assert.h>
35 #include <errno.h>
36 #include <stdio.h>
37 #include <stdlib.h>
38 #include <string.h>
39
40 /*----------------------------------------------------------------------------
41 * Local headers
42 *----------------------------------------------------------------------------*/
43
44 #include "bft_error.h"
45 #include "bft_mem.h"
46
47 #include "fvm_defs.h"
48 #include "fvm_io_num.h"
49 #include "fvm_nodal.h"
50 #include "fvm_nodal_priv.h"
51 #include "fvm_to_ensight_case.h"
52 #include "fvm_writer_helper.h"
53 #include "fvm_writer_priv.h"
54
55 #include "cs_block_dist.h"
56 #include "cs_file.h"
57 #include "cs_parall.h"
58 #include "cs_part_to_block.h"
59
60 /*----------------------------------------------------------------------------
61 * Header for the current file
62 *----------------------------------------------------------------------------*/
63
64 #include "fvm_to_ensight.h"
65
66 /*----------------------------------------------------------------------------*/
67
68 BEGIN_C_DECLS
69
70 /*! \cond DOXYGEN_SHOULD_SKIP_THIS */
71
72 /*============================================================================
73 * Local Type Definitions
74 *============================================================================*/
75
76 /*----------------------------------------------------------------------------
77 * EnSight Gold writer structure
78 *----------------------------------------------------------------------------*/
79
80 typedef struct {
81
82 char *name; /* Writer name */
83
84 int rank; /* Rank of current process in communicator */
85 int n_ranks; /* Number of processes in communicator */
86
87 bool text_mode; /* true if using text output */
88 bool swap_endian; /* true if binary file endianness must
89 be changed */
90
91 bool discard_polygons; /* Option to discard polygonal elements */
92 bool discard_polyhedra; /* Option to discard polyhedral elements */
93
94 bool divide_polygons; /* Option to tesselate polygonal elements */
95 bool divide_polyhedra; /* Option to tesselate polyhedral elements */
96
97 fvm_to_ensight_case_t *case_info; /* Associated case structure */
98
99 #if defined(HAVE_MPI)
100 int min_rank_step; /* Minimum rank step */
101 int min_block_size; /* Minimum block buffer size */
102 MPI_Comm block_comm; /* Associated MPI block communicator */
103 MPI_Comm comm; /* Associated MPI communicator */
104 #endif
105
106 } fvm_to_ensight_writer_t;
107
108 /*----------------------------------------------------------------------------
109 * Indirect file structure to handle both text and binary files
110 *----------------------------------------------------------------------------*/
111
112 typedef struct {
113
114 FILE *tf; /* Text file handing structure */
115 cs_file_t *bf; /* Binary file handling structure */
116
117 } _ensight_file_t;
118
119 /*----------------------------------------------------------------------------
120 * Context structure for fvm_writer_field_helper_output_* functions.
121 *----------------------------------------------------------------------------*/
122
123 typedef struct {
124
125 fvm_to_ensight_writer_t *writer; /* Pointer to writer structure */
126 _ensight_file_t *file; /* Pointer to file handler structure */
127
128 } _ensight_context_t;
129
130 /*============================================================================
131 * Static global variables
132 *============================================================================*/
133
134 static const char *_ensight_type_name[FVM_N_ELEMENT_TYPES] = {"bar2",
135 "tria3",
136 "quad4",
137 "nsided",
138 "tetra4",
139 "pyramid5",
140 "penta6",
141 "hexa8",
142 "nfaced"};
143
144 /* for symetric tensors, Code_Saturne assumes xx, yy, zz, xy, yz, xy,
145 xhere EnSight assumes xx, yy, zz, xy, xy, yz, so permutation is required */
146
147 static const int _ensight_c_order_6[6] = {0, 1, 2, 3, 5, 4};
148
149 /*============================================================================
150 * Private function definitions
151 *============================================================================*/
152
153 /*----------------------------------------------------------------------------
154 * Open an EnSight Gold geometry or variable file
155 *
156 * parameters:
157 * this_writer <-- pointer to Ensight Gold writer structure.
158 * filename <-- name of file to open.
159 * apend <-- if true, append to file instead of overwriting
160 *----------------------------------------------------------------------------*/
161
162 static _ensight_file_t
_open_ensight_file(const fvm_to_ensight_writer_t * this_writer,const char * filename,bool append)163 _open_ensight_file(const fvm_to_ensight_writer_t *this_writer,
164 const char *filename,
165 bool append)
166 {
167 _ensight_file_t f = {NULL, NULL};
168
169 if (this_writer->text_mode == true) {
170 if (this_writer->rank == 0) {
171 if (append)
172 f.tf = fopen(filename, "a");
173 else
174 f.tf = fopen(filename, "w");
175 if (f.tf == NULL)
176 bft_error(__FILE__, __LINE__, 0,
177 _("Error opening file \"%s\":\n\n"
178 " %s"), filename, strerror(errno));
179 }
180 }
181 else {
182
183 cs_file_mode_t mode = append ? CS_FILE_MODE_APPEND : CS_FILE_MODE_WRITE;
184 cs_file_access_t method;
185
186 #if defined(HAVE_MPI)
187
188 MPI_Info hints;
189 cs_file_get_default_access(CS_FILE_MODE_WRITE, &method, &hints);
190 f.bf = cs_file_open(filename,
191 mode,
192 method,
193 hints,
194 this_writer->block_comm,
195 this_writer->comm);
196
197 #else
198
199 cs_file_get_default_access(CS_FILE_MODE_WRITE, &method);
200 f.bf = cs_file_open(filename, mode, method);
201
202 #endif
203
204 if (this_writer->swap_endian == true)
205 cs_file_set_swap_endian(f.bf, 1);
206 }
207
208 return f;
209 }
210
211 /*----------------------------------------------------------------------------
212 * close an EnSight Gold geometry or variable file
213 *
214 * parameters:
215 * f <-- pointer to file handler structure.
216 *----------------------------------------------------------------------------*/
217
218 static void
_free_ensight_file(_ensight_file_t * f)219 _free_ensight_file(_ensight_file_t *f)
220 {
221 if (f->tf != NULL) {
222 if (fclose(f->tf) != 0)
223 bft_error(__FILE__, __LINE__, 0,
224 _("Error closing EnSight output file (text mode):\n\n"
225 " %s"), strerror(errno));
226 f->tf = NULL;
227 }
228
229 else if (f->bf != NULL)
230 f->bf = cs_file_free(f->bf);
231 }
232
233 /*----------------------------------------------------------------------------
234 * Write string to a text or C binary EnSight Gold file
235 *
236 * parameters:
237 * f <-- file to write to
238 * s <-- string to write
239 *----------------------------------------------------------------------------*/
240
241 static void
_write_string(_ensight_file_t f,const char * s)242 _write_string(_ensight_file_t f,
243 const char *s)
244 {
245 size_t i;
246 char buf[82];
247
248 if (f.tf != NULL) {
249 strncpy(buf, s, 80);
250 buf[80] = '\0';
251 fprintf(f.tf, "%s\n", buf);
252 }
253
254 else if (f.bf != NULL) {
255 strncpy(buf, s, 80);
256 buf[80] = '\0';
257 for (i = strlen(buf); i < 80; i++)
258 buf[i] = '\0';
259 cs_file_write_global(f.bf, buf, 1, 80);
260 }
261 }
262
263 /*----------------------------------------------------------------------------
264 * Write integer to a text or C binary EnSight Gold file
265 *
266 * parameters:
267 * f <-- file to write to
268 * n <-- integer value to write
269 *----------------------------------------------------------------------------*/
270
271 inline static void
_write_int(_ensight_file_t f,int32_t n)272 _write_int(_ensight_file_t f,
273 int32_t n)
274 {
275 if (f.tf != NULL)
276 fprintf(f.tf, "%10d\n", (int)n);
277
278 else if (f.bf != NULL) {
279 int _n = n;
280 cs_file_write_global(f.bf, &_n, sizeof(int32_t), 1);
281 }
282 }
283
284 /*----------------------------------------------------------------------------
285 * Write headers to an EnSight Gold geometry file.
286 *
287 * parameters:
288 * this_writer <-- pointer to Ensight Gold writer structure.
289 * f <-- pointer to file to initialize.
290 *----------------------------------------------------------------------------*/
291
292 static void
_write_geom_headers(fvm_to_ensight_writer_t * this_writer,_ensight_file_t f)293 _write_geom_headers(fvm_to_ensight_writer_t *this_writer,
294 _ensight_file_t f)
295 {
296 if (f.bf != NULL)
297 _write_string(f, "C Binary");
298
299 /* 1st description line */
300 {
301 char buf[81] = "";
302 if (this_writer->name != NULL)
303 strncpy(buf, this_writer->name, 80);
304 buf[80] = '\0';
305 _write_string(f, buf);
306 }
307 /* 2nd description line */
308 _write_string(f, "Output by Code_Saturne version "VERSION);
309 _write_string(f, "node id assign");
310 _write_string(f, "element id assign");
311 }
312
313 #if defined(HAVE_MPI)
314
315 /*----------------------------------------------------------------------------
316 * Write block of a vector of floats to an EnSight Gold file.
317 *
318 * This variant is called in parallel mode, where the values are already
319 * in a temporary block buffer and may be discarded.
320 *
321 * parameters:
322 * num_start <-- global number of first element for this block
323 * num_end <-- global number of past the last element for this block
324 * values <-- pointer to values block array
325 * comm <-- associated MPI communicator
326 * f <-- file to write to
327 *----------------------------------------------------------------------------*/
328
329 static void
_write_block_floats_g(cs_gnum_t num_start,cs_gnum_t num_end,float values[],MPI_Comm comm,_ensight_file_t f)330 _write_block_floats_g(cs_gnum_t num_start,
331 cs_gnum_t num_end,
332 float values[],
333 MPI_Comm comm,
334 _ensight_file_t f)
335 {
336 size_t i;
337 cs_gnum_t j;
338
339 /* In Binary mode, all ranks have a file structure,
340 we may use use a collective call */
341
342 if (f.bf != NULL)
343 cs_file_write_block_buffer(f.bf,
344 values,
345 sizeof(float),
346 1,
347 num_start,
348 num_end);
349
350 /* If all ranks do not have a binary file structure pointer, then
351 we are using a text file, open only on rank 0 */
352
353 else {
354
355 float *_values = NULL;
356 cs_file_serializer_t *s
357 = cs_file_serializer_create(sizeof(float),
358 1,
359 num_start,
360 num_end,
361 0,
362 values,
363 comm);
364
365 do {
366
367 cs_gnum_t range[2] = {num_start, num_end};
368
369 _values = cs_file_serializer_advance(s, range);
370
371 if (_values != NULL) { /* only possible on rank 0 */
372 assert(f.tf != NULL);
373 for (i = 0, j = range[0]; j < range[1]; i++, j++)
374 fprintf(f.tf, "%12.5e\n", _values[i]);
375 }
376
377 } while (_values != NULL);
378
379 cs_file_serializer_destroy(&s);
380 }
381
382 }
383
384 #endif /* defined(HAVE_MPI) */
385
386 /*----------------------------------------------------------------------------
387 * Write block of a vector of floats to an EnSight Gold file
388 *
389 * parameters:
390 * n_values <-- number of values to write
391 * num_end <-- global number of past the last element for this block
392 * values <-- pointer to values block array
393 * f <-- file to write to
394 *----------------------------------------------------------------------------*/
395
396 static void
_write_block_floats_l(size_t n_values,const float values[],_ensight_file_t f)397 _write_block_floats_l(size_t n_values,
398 const float values[],
399 _ensight_file_t f)
400 {
401 size_t i;
402
403 if (f.tf != NULL) {
404 for (i = 0; i < n_values; i++)
405 fprintf(f.tf, "%12.5e\n", values[i]);
406 }
407 else if (f.bf != NULL)
408 cs_file_write_global(f.bf, values, sizeof(float), n_values);
409 }
410
411 #if defined(HAVE_MPI)
412
413 /*----------------------------------------------------------------------------
414 * Write vertex coordinates to an EnSight Gold file in parallel mode
415 *
416 * parameters:
417 * this_writer <-- pointer to associated writer
418 * mesh <-- pointer to nodal mesh structure
419 * f <-- associated file handle
420 *----------------------------------------------------------------------------*/
421
422 static void
_export_vertex_coords_g(const fvm_to_ensight_writer_t * this_writer,const fvm_nodal_t * mesh,_ensight_file_t f)423 _export_vertex_coords_g(const fvm_to_ensight_writer_t *this_writer,
424 const fvm_nodal_t *mesh,
425 _ensight_file_t f)
426 {
427 cs_lnum_t i, j;
428 size_t stride;
429 cs_block_dist_info_t bi;
430
431 cs_gnum_t n_g_extra_vertices = 0, n_g_vertices_tot = 0;
432 cs_lnum_t n_extra_vertices = 0, n_vertices_tot = 0;
433 cs_coord_t *extra_vertex_coords = NULL;
434 float *part_coords = NULL, *block_coords = NULL;
435
436 cs_part_to_block_t *d = NULL;
437 size_t block_buf_size = 0;
438
439 const cs_coord_t *vertex_coords = mesh->vertex_coords;
440 const cs_lnum_t *parent_vertex_num = mesh->parent_vertex_num;
441 const cs_lnum_t n_vertices
442 = fvm_io_num_get_local_count(mesh->global_vertex_num);
443 cs_gnum_t n_g_vertices
444 = fvm_io_num_get_global_count(mesh->global_vertex_num);
445
446 /* Check for extra vertex coordinates */
447
448 fvm_writer_count_extra_vertices(mesh,
449 this_writer->divide_polyhedra,
450 &n_g_extra_vertices,
451 &n_extra_vertices);
452
453 n_vertices_tot = n_vertices + n_extra_vertices;
454 n_g_vertices_tot = n_g_vertices + n_g_extra_vertices;
455
456 /* Initialize distribution info */
457
458 fvm_writer_vertex_part_to_block_create(this_writer->min_rank_step,
459 this_writer->min_block_size,
460 n_g_extra_vertices,
461 n_extra_vertices,
462 mesh,
463 &bi,
464 &d,
465 this_writer->comm);
466
467 /* Compute extra vertex coordinates if present */
468
469 extra_vertex_coords = fvm_writer_extra_vertex_coords(mesh, n_extra_vertices);
470
471 /* Build arrays */
472
473 block_buf_size = (bi.gnum_range[1] - bi.gnum_range[0]);
474 BFT_MALLOC(block_coords, block_buf_size, float);
475 BFT_MALLOC(part_coords, n_vertices_tot, float);
476
477 /* Vertex coordinates */
478 /*--------------------*/
479
480 stride = (size_t)(mesh->dim);
481
482 _write_string(f, "coordinates");
483 _write_int(f, (int)(n_g_vertices_tot));
484
485 /* Loop on dimension (de-interlace coordinates, always 3D for EnSight) */
486
487 for (j = 0; j < 3; j++) {
488
489 if (j < mesh->dim) {
490
491 if (parent_vertex_num != NULL) {
492 for (i = 0; i < n_vertices; i++)
493 part_coords[i] = vertex_coords[(parent_vertex_num[i]-1)*stride + j];
494 }
495 else {
496 for (i = 0; i < n_vertices; i++)
497 part_coords[i] = vertex_coords[i*stride + j];
498 }
499
500 for (i = 0; i < n_extra_vertices; i++)
501 part_coords[n_vertices + i] = extra_vertex_coords[(i*stride) + j];
502 }
503 else {
504 for (i = 0; i < n_vertices_tot; i++)
505 part_coords[i] = 0.0;
506 }
507
508 cs_part_to_block_copy_array(d,
509 CS_FLOAT,
510 1,
511 part_coords,
512 block_coords);
513
514 _write_block_floats_g(bi.gnum_range[0],
515 bi.gnum_range[1],
516 block_coords,
517 this_writer->comm,
518 f);
519
520 } /* end of loop on spatial dimension */
521
522 cs_part_to_block_destroy(&d);
523
524 BFT_FREE(block_coords);
525 BFT_FREE(part_coords);
526 if (extra_vertex_coords != NULL)
527 BFT_FREE(extra_vertex_coords);
528 }
529
530 #endif /* defined(HAVE_MPI) */
531
532 /*----------------------------------------------------------------------------
533 * Write vertex coordinates to an EnSight Gold file in serial mode
534 *
535 * parameters:
536 * this_writer <-- pointer to associated writer
537 * mesh <-- pointer to nodal mesh structure
538 * f <-- associated file handle
539 *----------------------------------------------------------------------------*/
540
541 static void
_export_vertex_coords_l(const fvm_to_ensight_writer_t * this_writer,const fvm_nodal_t * mesh,_ensight_file_t f)542 _export_vertex_coords_l(const fvm_to_ensight_writer_t *this_writer,
543 const fvm_nodal_t *mesh,
544 _ensight_file_t f)
545 {
546 cs_lnum_t i, j;
547 cs_lnum_t n_extra_vertices = 0;
548 cs_coord_t *extra_vertex_coords = NULL;
549 float *coords_tmp = NULL;
550
551 const cs_lnum_t n_vertices = mesh->n_vertices;
552 const cs_coord_t *vertex_coords = mesh->vertex_coords;
553 const cs_lnum_t *parent_vertex_num = mesh->parent_vertex_num;
554
555 const size_t stride = (size_t)(mesh->dim);
556
557 /* Compute extra vertex coordinates if present */
558
559 fvm_writer_count_extra_vertices(mesh,
560 this_writer->divide_polyhedra,
561 NULL,
562 &n_extra_vertices);
563
564 extra_vertex_coords = fvm_writer_extra_vertex_coords(mesh,
565 n_extra_vertices);
566
567 /* Vertex coordinates */
568 /*--------------------*/
569
570 BFT_MALLOC(coords_tmp, CS_MAX(n_vertices, n_extra_vertices), float);
571
572 _write_string(f, "coordinates");
573 _write_int(f, n_vertices + n_extra_vertices);
574
575 /* Loop on dimension (de-interlace coordinates, always 3D for EnSight) */
576
577 for (j = 0; j < 3; j++) {
578
579 /* First, handle regular vertices */
580
581 if (j < mesh->dim) {
582 if (parent_vertex_num != NULL) {
583 for (i = 0; i < n_vertices; i++) {
584 assert(parent_vertex_num[i] != 0);
585 coords_tmp[i]
586 = (float)(vertex_coords[(parent_vertex_num[i]-1)*stride + j]);
587 }
588 }
589 else {
590 for (i = 0; i < n_vertices; i++)
591 coords_tmp[i] = (float)(vertex_coords[i*stride + j]);
592 }
593 }
594 else {
595 for (i = 0; i < (n_vertices); i++)
596 coords_tmp[i] = 0.0;
597 }
598
599 _write_block_floats_l(n_vertices,
600 coords_tmp,
601 f);
602
603 /* Handle extra vertices (only occur with polyhedra tesselations in 3d) */
604
605 for (i = 0; i < n_extra_vertices; i++)
606 coords_tmp[i] = (float)(extra_vertex_coords[i*stride + j]);
607
608 if (n_extra_vertices > 0)
609 _write_block_floats_l(n_extra_vertices,
610 coords_tmp,
611 f);
612
613 } /* end of loop on mesh dimension */
614
615 BFT_FREE(coords_tmp);
616
617 if (extra_vertex_coords != NULL)
618 BFT_FREE(extra_vertex_coords);
619 }
620
621 #if defined(HAVE_MPI)
622
623 /*----------------------------------------------------------------------------
624 * Write strided connectivity block to an EnSight Gold file in text mode
625 *
626 * parameters:
627 * stride <-- number of vertices per element type
628 * n_elems <-- number of elements
629 * connect <-- connectivity array
630 * tf <-- file to write to
631 *----------------------------------------------------------------------------*/
632
633 static void
_write_connect_block_gt(int stride,cs_lnum_t n_elems,const int32_t connect[],FILE * tf)634 _write_connect_block_gt(int stride,
635 cs_lnum_t n_elems,
636 const int32_t connect[],
637 FILE *tf)
638 {
639 cs_lnum_t i;
640
641 assert(tf != NULL);
642
643 switch(stride) {
644
645 case 1: /* length */
646 for (i = 0; i < n_elems; i++)
647 fprintf(tf, "%10d\n",
648 (int)connect[i]);
649 break;
650
651 case 2: /* edge */
652 for (i = 0; i < n_elems; i++)
653 fprintf(tf, "%10d%10d\n",
654 (int)connect[i*2],
655 (int)connect[i*2+1]);
656 break;
657
658 case 3: /* FVM_FACE_TRIA */
659 for (i = 0; i < n_elems; i++)
660 fprintf(tf, "%10d%10d%10d\n",
661 (int)connect[i*3],
662 (int)connect[i*3+1],
663 (int)connect[i*3+2]);
664 break;
665
666 case 4: /* FVM_FACE_QUAD or FVM_CELL_TETRA */
667 for (i = 0; i < n_elems; i++)
668 fprintf(tf, "%10d%10d%10d%10d\n",
669 (int)connect[i*4],
670 (int)connect[i*4+1],
671 (int)connect[i*4+2],
672 (int)connect[i*4+3]);
673 break;
674
675 case 5: /* FVM_CELL_PYRAM */
676 for (i = 0; i < n_elems; i++)
677 fprintf(tf, "%10d%10d%10d%10d%10d\n",
678 (int)connect[i*5],
679 (int)connect[i*5+1],
680 (int)connect[i*5+2],
681 (int)connect[i*5+3],
682 (int)connect[i*5+4]);
683 break;
684
685 case 6: /* FVM_CELL_PRISM */
686 for (i = 0; i < n_elems; i++)
687 fprintf(tf, "%10d%10d%10d%10d%10d%10d\n",
688 (int)connect[i*6],
689 (int)connect[i*6+1],
690 (int)connect[i*6+2],
691 (int)connect[i*6+3],
692 (int)connect[i*6+4],
693 (int)connect[i*6+5]);
694 break;
695
696 case 8: /* FVM_CELL_HEXA */
697 for (i = 0; i < n_elems; i++)
698 fprintf(tf,
699 "%10d%10d%10d%10d%10d%10d%10d%10d\n",
700 (int)connect[i*8],
701 (int)connect[i*8+1],
702 (int)connect[i*8+2],
703 (int)connect[i*8+3],
704 (int)connect[i*8+4],
705 (int)connect[i*8+5],
706 (int)connect[i*8+6],
707 (int)connect[i*8+7]);
708 break;
709
710 default:
711 assert(0);
712 }
713 }
714
715 /*----------------------------------------------------------------------------
716 * Write strided global connectivity block to an EnSight Gold file
717 *
718 * parameters:
719 * stride <-- number of vertices per element type
720 * num_start <-- global number of first element for this block
721 * num_end <-- global number of past last element for this block
722 * block_connect <-> global connectivity block array
723 * comm <-- associated MPI communicator
724 * f <-- associated file handle
725 *----------------------------------------------------------------------------*/
726
727 static void
_write_block_connect_g(int stride,cs_gnum_t num_start,cs_gnum_t num_end,int32_t block_connect[],MPI_Comm comm,_ensight_file_t f)728 _write_block_connect_g(int stride,
729 cs_gnum_t num_start,
730 cs_gnum_t num_end,
731 int32_t block_connect[],
732 MPI_Comm comm,
733 _ensight_file_t f)
734 {
735 /* In Binary mode, all ranks have a file structure,
736 we may use use a collective call */
737
738 if (f.bf != NULL)
739 cs_file_write_block_buffer(f.bf,
740 block_connect,
741 sizeof(int32_t),
742 stride,
743 num_start,
744 num_end);
745
746 /* If all ranks do not have a binary file structure pointer, then
747 we are using a text file, open only on rank 0 */
748
749 else {
750
751 int32_t *_block_connect = NULL;
752
753 cs_file_serializer_t *s = cs_file_serializer_create(sizeof(int32_t),
754 stride,
755 num_start,
756 num_end,
757 0,
758 block_connect,
759 comm);
760
761 do {
762 cs_gnum_t range[2] = {num_start, num_end};
763
764 _block_connect = cs_file_serializer_advance(s, range);
765
766 if (_block_connect != NULL) /* only possible on rank 0 */
767 _write_connect_block_gt(stride,
768 (range[1] - range[0]),
769 _block_connect,
770 f.tf);
771
772 } while (_block_connect != NULL);
773
774 cs_file_serializer_destroy(&s);
775 }
776 }
777
778 #endif /* defined(HAVE_MPI) */
779
780 /*----------------------------------------------------------------------------
781 * Write strided local connectivity to an EnSight Gold file
782 *
783 * parameters:
784 * stride <-- number of vertices per element type
785 * n_elems <-- number of elements
786 * connect <-- connectivity array
787 * f <-- file to write to
788 *----------------------------------------------------------------------------*/
789
790 static void
_write_connect_l(int stride,cs_lnum_t n_elems,const cs_lnum_t connect[],_ensight_file_t f)791 _write_connect_l(int stride,
792 cs_lnum_t n_elems,
793 const cs_lnum_t connect[],
794 _ensight_file_t f)
795 {
796 cs_lnum_t i;
797
798 if (f.tf != NULL) { /* Text mode */
799
800 switch(stride) {
801
802 case 2: /* edge */
803 for (i = 0; i < n_elems; i++)
804 fprintf(f.tf, "%10d%10d\n",
805 (int)connect[i*2],
806 (int)connect[i*2+1]);
807 break;
808
809 case 3: /* FVM_FACE_TRIA */
810 for (i = 0; i < n_elems; i++)
811 fprintf(f.tf, "%10d%10d%10d\n",
812 (int)connect[i*3],
813 (int)connect[i*3+1],
814 (int)connect[i*3+2]);
815 break;
816
817 case 4: /* FVM_FACE_QUAD or FVM_CELL_TETRA */
818 for (i = 0; i < n_elems; i++)
819 fprintf(f.tf, "%10d%10d%10d%10d\n",
820 (int)connect[i*4],
821 (int)connect[i*4+1],
822 (int)connect[i*4+2],
823 (int)connect[i*4+3]);
824 break;
825
826 case 5: /* FVM_CELL_PYRAM */
827 for (i = 0; i < n_elems; i++)
828 fprintf(f.tf, "%10d%10d%10d%10d%10d\n",
829 (int)connect[i*5],
830 (int)connect[i*5+1],
831 (int)connect[i*5+2],
832 (int)connect[i*5+3],
833 (int)connect[i*5+4]);
834 break;
835
836 case 6: /* FVM_CELL_PRISM */
837 for (i = 0; i < n_elems; i++)
838 fprintf(f.tf, "%10d%10d%10d%10d%10d%10d\n",
839 (int)connect[i*6],
840 (int)connect[i*6+1],
841 (int)connect[i*6+2],
842 (int)connect[i*6+3],
843 (int)connect[i*6+4],
844 (int)connect[i*6+5]);
845 break;
846
847 case 8: /* FVM_CELL_HEXA */
848 for (i = 0; i < n_elems; i++)
849 fprintf(f.tf,
850 "%10d%10d%10d%10d%10d%10d%10d%10d\n",
851 (int)connect[i*8],
852 (int)connect[i*8+1],
853 (int)connect[i*8+2],
854 (int)connect[i*8+3],
855 (int)connect[i*8+4],
856 (int)connect[i*8+5],
857 (int)connect[i*8+6],
858 (int)connect[i*8+7]);
859 break;
860
861 default:
862 assert(0);
863 }
864
865 }
866 else if (f.bf != NULL) { /* Binary mode */
867
868 size_t j;
869 size_t k = 0;
870 int32_t *buffer = NULL;
871 const size_t n_values = n_elems*stride;
872 const size_t buffer_size = n_values > 64 ? (n_values / 8) : n_values;
873
874 BFT_MALLOC(buffer, buffer_size, int32_t);
875
876 while (k < n_values) {
877 for (j = 0; j < buffer_size && k < n_values; j++)
878 buffer[j] = (int)(connect[k++]);
879 cs_file_write_global(f.bf, buffer, sizeof(int32_t), j);
880 }
881
882 BFT_FREE(buffer);
883 }
884
885 }
886
887 #if defined(HAVE_MPI)
888
889 /*----------------------------------------------------------------------------
890 * Output function for field values.
891 *
892 * This function is passed to fvm_writer_field_helper_output_* functions.
893 *
894 * parameters:
895 * context <-> pointer to writer and field context
896 * datatype <-- output datatype
897 * dimension <-- output field dimension
898 * component_id <-- output component id (if non-interleaved)
899 * block_start <-- start global number of element for current block
900 * block_end <-- past-the-end global number of element for current block
901 * buffer <-> associated output buffer
902 *----------------------------------------------------------------------------*/
903
904 static void
_field_output_g(void * context,cs_datatype_t datatype,int dimension,int component_id,cs_gnum_t block_start,cs_gnum_t block_end,void * buffer)905 _field_output_g(void *context,
906 cs_datatype_t datatype,
907 int dimension,
908 int component_id,
909 cs_gnum_t block_start,
910 cs_gnum_t block_end,
911 void *buffer)
912 {
913 CS_NO_WARN_IF_UNUSED(datatype);
914 CS_UNUSED(dimension);
915 CS_UNUSED(component_id);
916
917 _ensight_context_t *c = context;
918
919 fvm_to_ensight_writer_t *w = c->writer;
920 _ensight_file_t *f = c->file;
921
922 assert(datatype == CS_FLOAT);
923
924 _write_block_floats_g(block_start,
925 block_end,
926 buffer,
927 w->comm,
928 *f);
929 }
930
931 /*----------------------------------------------------------------------------
932 * Write "trivial" point elements to an EnSight Gold file in parallel mode
933 *
934 * parameters:
935 * w <-- pointer to writer structure
936 * mesh <-- pointer to nodal mesh structure
937 * f <-- file to write to
938 *----------------------------------------------------------------------------*/
939
940 static void
_export_point_elements_g(const fvm_to_ensight_writer_t * w,const fvm_nodal_t * mesh,_ensight_file_t f)941 _export_point_elements_g(const fvm_to_ensight_writer_t *w,
942 const fvm_nodal_t *mesh,
943 _ensight_file_t f)
944 {
945 const cs_gnum_t n_g_vertices
946 = fvm_io_num_get_global_count(mesh->global_vertex_num);
947
948 _write_string(f, "point");
949 _write_int(f, (int)n_g_vertices);
950
951 if (f.tf != NULL) { /* Text mode, rank 0 only */
952
953 cs_gnum_t i;
954 int32_t j = 1;
955
956 for (i = 0; i < n_g_vertices; i++)
957 fprintf(f.tf, "%10d\n", j++);
958
959 }
960 else if (f.bf != NULL) { /* Binary mode */
961
962 cs_lnum_t i;
963 cs_gnum_t j;
964 cs_block_dist_info_t bi;
965
966 size_t min_block_size = w->min_block_size / sizeof(float);
967 int32_t *connect = NULL;
968
969 bi = cs_block_dist_compute_sizes(w->rank,
970 w->n_ranks,
971 w->min_rank_step,
972 min_block_size,
973 n_g_vertices);
974
975 BFT_MALLOC(connect, bi.gnum_range[1] - bi.gnum_range[0], int32_t);
976
977 for (i = 0, j = bi.gnum_range[0]; j < bi.gnum_range[1]; i++, j++)
978 connect[i] = j;
979
980 _write_block_connect_g(1,
981 bi.gnum_range[0],
982 bi.gnum_range[1],
983 connect,
984 w->comm,
985 f);
986
987 BFT_FREE(connect);
988 }
989 }
990
991 #endif /* defined(HAVE_MPI) */
992
993 /*----------------------------------------------------------------------------
994 * Write "trivial" point elements to an EnSight Gold file in serial mode
995 *
996 * parameters:
997 * mesh <-- pointer to nodal mesh structure
998 * f <-- file to write to
999 *----------------------------------------------------------------------------*/
1000
1001 static void
_export_point_elements_l(const fvm_nodal_t * mesh,_ensight_file_t f)1002 _export_point_elements_l(const fvm_nodal_t *mesh,
1003 _ensight_file_t f)
1004 {
1005 const cs_lnum_t n_vertices = mesh->n_vertices;
1006
1007 _write_string(f, "point");
1008 _write_int(f, (int)n_vertices);
1009
1010 if (n_vertices == 0)
1011 return;
1012
1013 if (f.tf != NULL) { /* Text mode */
1014 int i;
1015 for (i = 0; i < n_vertices; i++)
1016 fprintf(f.tf, "%10d\n", i+1);
1017 }
1018
1019 else if (f.bf != NULL) { /* Binary mode */
1020
1021 int32_t k, j_end;
1022 int32_t j = 1;
1023 int32_t *buf = NULL;
1024 const int32_t bufsize = n_vertices > 64 ? (n_vertices / 8) : n_vertices;
1025
1026 BFT_MALLOC(buf, bufsize, int32_t);
1027
1028 j_end = n_vertices + 1;
1029 while (j < j_end) {
1030 for (k = 0; j < j_end && k < bufsize; k++)
1031 buf[k] = j++;
1032 cs_file_write_global(f.bf, buf, sizeof(int32_t), k);
1033 }
1034
1035 BFT_FREE(buf);
1036 }
1037 }
1038
1039 #if defined(HAVE_MPI)
1040
1041 /*----------------------------------------------------------------------------
1042 * Write indexed element lengths from a nodal mesh to an EnSight Gold
1043 * file in parallel mode
1044 *
1045 * parameters:
1046 * w <-- pointer to writer structure
1047 * global_element_num <-- global element numbering
1048 * vertex_index <-- pointer to element -> vertex index
1049 * n_ranks <-- number of processes in communicator
1050 * f <-- associated file handle
1051 *----------------------------------------------------------------------------*/
1052
1053 static void
_write_lengths_g(const fvm_to_ensight_writer_t * w,const fvm_io_num_t * global_element_num,const cs_lnum_t vertex_index[],_ensight_file_t f)1054 _write_lengths_g(const fvm_to_ensight_writer_t *w,
1055 const fvm_io_num_t *global_element_num,
1056 const cs_lnum_t vertex_index[],
1057 _ensight_file_t f)
1058 {
1059 cs_lnum_t i;
1060 cs_block_dist_info_t bi;
1061
1062 int32_t *part_lengths = NULL;
1063 int32_t *block_lengths = NULL;
1064
1065 cs_part_to_block_t *d = NULL;
1066
1067 const size_t min_block_size = w->min_block_size / sizeof(int32_t);
1068 const cs_lnum_t n_elements
1069 = fvm_io_num_get_local_count(global_element_num);
1070 const cs_lnum_t n_g_elements
1071 = fvm_io_num_get_global_count(global_element_num);
1072 const cs_gnum_t *g_num
1073 = fvm_io_num_get_global_num(global_element_num);
1074
1075 /* Allocate block buffer */
1076
1077 bi = cs_block_dist_compute_sizes(w->rank,
1078 w->n_ranks,
1079 w->min_rank_step,
1080 min_block_size,
1081 n_g_elements);
1082
1083 /* Build distribution structures */
1084
1085 BFT_MALLOC(block_lengths, bi.gnum_range[1] - bi.gnum_range[0], int);
1086 BFT_MALLOC(part_lengths, n_elements, int32_t);
1087
1088 for (i = 0; i < n_elements; i++)
1089 part_lengths[i] = vertex_index[i+1] - vertex_index[i];
1090
1091 d = cs_part_to_block_create_by_gnum(w->comm, bi, n_elements, g_num);
1092
1093 cs_part_to_block_copy_array(d,
1094 CS_INT32,
1095 1,
1096 part_lengths,
1097 block_lengths);
1098
1099 cs_part_to_block_destroy(&d);
1100 BFT_FREE(part_lengths);
1101
1102 /* Write to file */
1103
1104 _write_block_connect_g(1,
1105 bi.gnum_range[0],
1106 bi.gnum_range[1],
1107 block_lengths,
1108 w->comm,
1109 f);
1110
1111 BFT_FREE(block_lengths);
1112 }
1113
1114 /*----------------------------------------------------------------------------
1115 * Write block-distributed indexed element (polygons or polyhedra)
1116 * cell -> vertex connectivity to an EnSight Gold file in parallel mode.
1117 *
1118 * In text mode, zeroes may be used in place of extra vertex numbers
1119 * to indicate extra newlines between face -> vertex definitions.
1120 *
1121 * parameters:
1122 * num_start <-- global number of first element for this block
1123 * num_end <-- global number of past last element for this block
1124 * block_index <-- global connectivity block array
1125 * block_connect <-> global connectivity block array
1126 * comm <-- associated MPI communicator
1127 * f <-- associated file handle
1128 *----------------------------------------------------------------------------*/
1129
1130 static void
_write_block_indexed(cs_gnum_t num_start,cs_gnum_t num_end,const cs_lnum_t block_index[],int32_t block_connect[],MPI_Comm comm,_ensight_file_t f)1131 _write_block_indexed(cs_gnum_t num_start,
1132 cs_gnum_t num_end,
1133 const cs_lnum_t block_index[],
1134 int32_t block_connect[],
1135 MPI_Comm comm,
1136 _ensight_file_t f)
1137 {
1138 cs_gnum_t block_size = 0, block_start = 0, block_end = 0;
1139
1140 /* Prepare write to file */
1141
1142 block_size = block_index[num_end - num_start];
1143
1144 MPI_Scan(&block_size, &block_end, 1, CS_MPI_GNUM, MPI_SUM, comm);
1145 block_end += 1;
1146 block_start = block_end - block_size;
1147
1148 /* In Binary mode, all ranks have a file structure,
1149 we may use use a collective call */
1150
1151 if (f.bf != NULL)
1152 cs_file_write_block_buffer(f.bf,
1153 block_connect,
1154 sizeof(int32_t),
1155 1,
1156 block_start,
1157 block_end);
1158
1159 /* If all ranks do not have a binary file structure pointer, then
1160 we are using a text file, open only on rank 0 */
1161
1162 else {
1163 cs_lnum_t i;
1164 int32_t *_block_vtx_num = NULL;
1165 cs_file_serializer_t *s = cs_file_serializer_create(sizeof(int32_t),
1166 1,
1167 block_start,
1168 block_end,
1169 0,
1170 block_connect,
1171 comm);
1172
1173 do {
1174 cs_gnum_t j;
1175 cs_gnum_t range[2] = {block_start, block_end};
1176 _block_vtx_num = cs_file_serializer_advance(s, range);
1177 if (_block_vtx_num != NULL) { /* only possible on rank 0 */
1178 assert(f.tf != NULL);
1179 for (i = 0, j = range[0]; j < range[1]; i++, j++) {
1180 if (_block_vtx_num[i] != 0)
1181 fprintf(f.tf, "%10d", _block_vtx_num[i]);
1182 else
1183 fprintf(f.tf, "\n");
1184 }
1185 }
1186 } while (_block_vtx_num != NULL);
1187
1188 cs_file_serializer_destroy(&s);
1189 }
1190 }
1191
1192 /*----------------------------------------------------------------------------
1193 * Write indexed element (polygons or polyhedra) cell -> vertex connectivity
1194 * to an EnSight Gold file in parallel mode.
1195 *
1196 * In text mode, zeroes may be used in place of extra vertex numbers
1197 * to indicate extra newlines between face -> vertex definitions.
1198 *
1199 * parameters:
1200 * w <-- pointer to writer structure
1201 * global_vertex_num <-- vertex global numbering
1202 * global_element_num <-- global element numbering
1203 * vertex_index <-- element -> vertex index
1204 * vertex_num <-- element -> vertex number
1205 * f <-- associated file handle
1206 *----------------------------------------------------------------------------*/
1207
1208 static void
_write_indexed_connect_g(const fvm_to_ensight_writer_t * w,const fvm_io_num_t * global_element_num,const cs_lnum_t vertex_index[],const int32_t vertex_num[],_ensight_file_t f)1209 _write_indexed_connect_g(const fvm_to_ensight_writer_t *w,
1210 const fvm_io_num_t *global_element_num,
1211 const cs_lnum_t vertex_index[],
1212 const int32_t vertex_num[],
1213 _ensight_file_t f)
1214 {
1215 cs_block_dist_info_t bi;
1216
1217 cs_gnum_t loc_size = 0, tot_size = 0, block_size = 0;
1218 cs_part_to_block_t *d = NULL;
1219 cs_lnum_t *block_index = NULL;
1220 int32_t *block_vtx_num = NULL;
1221 size_t min_block_size = w->min_block_size / sizeof(int32_t);
1222
1223 const cs_gnum_t n_g_elements
1224 = fvm_io_num_get_global_count(global_element_num);
1225 const cs_lnum_t n_elements
1226 = fvm_io_num_get_local_count(global_element_num);
1227 const cs_gnum_t *g_elt_num
1228 = fvm_io_num_get_global_num(global_element_num);
1229
1230 /* Adjust min block size based on minimum element size */
1231
1232 loc_size = vertex_index[n_elements];
1233 MPI_Allreduce(&loc_size, &tot_size, 1, CS_MPI_GNUM, MPI_SUM, w->comm);
1234
1235 min_block_size /= (tot_size / n_g_elements);
1236
1237 /* Allocate memory for additionnal indexes */
1238
1239 bi = cs_block_dist_compute_sizes(w->rank,
1240 w->n_ranks,
1241 w->min_rank_step,
1242 min_block_size,
1243 n_g_elements);
1244
1245 BFT_MALLOC(block_index, bi.gnum_range[1] - bi.gnum_range[0] + 1, cs_lnum_t);
1246
1247 d = cs_part_to_block_create_by_gnum(w->comm, bi, n_elements, g_elt_num);
1248
1249 cs_part_to_block_copy_index(d,
1250 vertex_index,
1251 block_index);
1252
1253 block_size = block_index[bi.gnum_range[1] - bi.gnum_range[0]];
1254
1255 BFT_MALLOC(block_vtx_num, block_size, int32_t);
1256
1257 cs_part_to_block_copy_indexed(d,
1258 CS_INT32,
1259 vertex_index,
1260 vertex_num,
1261 block_index,
1262 block_vtx_num);
1263
1264 /* Write to file */
1265
1266 _write_block_indexed(bi.gnum_range[0],
1267 bi.gnum_range[1],
1268 block_index,
1269 block_vtx_num,
1270 w->comm,
1271 f);
1272
1273 /* Free memory */
1274
1275 BFT_FREE(block_vtx_num);
1276 cs_part_to_block_destroy(&d);
1277 BFT_FREE(block_index);
1278 }
1279
1280 /*----------------------------------------------------------------------------
1281 * Write polyhedra from a nodal mesh to an EnSight Gold file in parallel mode
1282 *
1283 * parameters:
1284 * w <-- pointer to writer structure
1285 * export_section <-- pointer to EnSight section helper structure
1286 * global_vertex_num <-- pointer to vertex global numbering
1287 * f <-- associated file handle
1288 *
1289 * returns:
1290 * pointer to next EnSight section helper structure in list
1291 *----------------------------------------------------------------------------*/
1292
1293 static const fvm_writer_section_t *
_export_nodal_polyhedra_g(const fvm_to_ensight_writer_t * w,const fvm_writer_section_t * export_section,const fvm_io_num_t * global_vertex_num,_ensight_file_t f)1294 _export_nodal_polyhedra_g(const fvm_to_ensight_writer_t *w,
1295 const fvm_writer_section_t *export_section,
1296 const fvm_io_num_t *global_vertex_num,
1297 _ensight_file_t f)
1298 {
1299 cs_lnum_t i, j, k, l, face_id;
1300
1301 cs_lnum_t face_length, cell_length;
1302 cs_block_dist_info_t bi;
1303
1304 cs_part_to_block_t *d = NULL;
1305 const fvm_writer_section_t *current_section;
1306
1307 /* Export number of faces per polyhedron */
1308 /*---------------------------------------*/
1309
1310 current_section = export_section;
1311
1312 do { /* loop on sections which should be appended */
1313
1314 const fvm_nodal_section_t *section = current_section->section;
1315
1316 _write_lengths_g(w,
1317 section->global_element_num,
1318 section->face_index,
1319 f);
1320
1321 current_section = current_section->next;
1322
1323 } while ( current_section != NULL
1324 && current_section->continues_previous == true);
1325
1326 /* Export number of vertices per face per polyhedron */
1327 /*---------------------------------------------------*/
1328
1329 current_section = export_section;
1330
1331 do { /* loop on sections which should be appended */
1332
1333 cs_gnum_t block_size = 0, block_start = 0, block_end = 0;
1334 cs_lnum_t *block_index = NULL;
1335
1336 size_t min_block_size = w->min_block_size / sizeof(int32_t);
1337 int32_t *part_face_len = NULL, *block_face_len = NULL;
1338
1339 const fvm_nodal_section_t *section = current_section->section;
1340 const cs_lnum_t n_elements
1341 = fvm_io_num_get_local_count(section->global_element_num);
1342 const cs_gnum_t n_g_elements
1343 = fvm_io_num_get_global_count(section->global_element_num);
1344 const cs_gnum_t *g_elt_num
1345 = fvm_io_num_get_global_num(section->global_element_num);
1346
1347 /* Build local polyhedron face lengths information */
1348
1349 BFT_MALLOC(part_face_len,
1350 section->face_index[section->n_elements],
1351 int32_t);
1352
1353 k = 0;
1354 for (i = 0; i < section->n_elements; i++) {
1355 for (j = section->face_index[i]; j < section->face_index[i+1]; j++) {
1356 face_id = CS_ABS(section->face_num[j]) - 1;
1357 face_length = ( section->vertex_index[face_id+1]
1358 - section->vertex_index[face_id]);
1359 part_face_len[k++] = face_length;
1360 }
1361 }
1362 assert(k == section->face_index[section->n_elements]);
1363
1364 /* Prepare distribution structures */
1365
1366 bi = cs_block_dist_compute_sizes(w->rank,
1367 w->n_ranks,
1368 w->min_rank_step,
1369 min_block_size,
1370 n_g_elements);
1371
1372 d = cs_part_to_block_create_by_gnum(w->comm,
1373 bi,
1374 n_elements,
1375 g_elt_num);
1376
1377 BFT_MALLOC(block_index, bi.gnum_range[1] - bi.gnum_range[0] + 1, cs_lnum_t);
1378
1379 cs_part_to_block_copy_index(d,
1380 section->face_index,
1381 block_index);
1382
1383 block_size = block_index[bi.gnum_range[1] - bi.gnum_range[0]];
1384
1385 BFT_MALLOC(block_face_len, block_size, int32_t);
1386
1387 cs_part_to_block_copy_indexed(d,
1388 CS_INT32,
1389 section->face_index,
1390 part_face_len,
1391 block_index,
1392 block_face_len);
1393
1394 MPI_Scan(&block_size, &block_end, 1, CS_MPI_GNUM, MPI_SUM, w->comm);
1395 block_end += 1;
1396 block_start = block_end - block_size;
1397
1398 _write_block_connect_g(1,
1399 block_start,
1400 block_end,
1401 block_face_len,
1402 w->comm,
1403 f);
1404
1405 BFT_FREE(block_face_len);
1406
1407 cs_part_to_block_destroy(&d);
1408
1409 BFT_FREE(block_index);
1410 BFT_FREE(part_face_len);
1411
1412 current_section = current_section->next;
1413
1414 } while ( current_section != NULL
1415 && current_section->continues_previous == true);
1416
1417 /* Export cell->vertex connectivity by blocks */
1418 /*--------------------------------------------*/
1419
1420 current_section = export_section;
1421
1422 do { /* loop on sections which should be appended */
1423
1424 cs_lnum_t *part_vtx_idx = NULL;
1425 int32_t *part_vtx_num = NULL;
1426
1427 const fvm_nodal_section_t *section = current_section->section;
1428 const cs_gnum_t *g_vtx_num
1429 = fvm_io_num_get_global_num(global_vertex_num);
1430
1431 BFT_MALLOC(part_vtx_idx, section->n_elements + 1, cs_lnum_t);
1432
1433 l = 0;
1434
1435 if (f.bf != NULL) { /* In binary mode, build cell -> vertex connectivity */
1436
1437 part_vtx_idx[0] = 0;
1438 for (i = 0; i < section->n_elements; i++) {
1439 cell_length = 0;
1440 for (j = section->face_index[i]; j < section->face_index[i+1]; j++) {
1441 face_id = CS_ABS(section->face_num[j]) - 1;
1442 face_length = ( section->vertex_index[face_id+1]
1443 - section->vertex_index[face_id]);
1444 cell_length += face_length;
1445 }
1446 part_vtx_idx[i+1] = part_vtx_idx[i] + cell_length;
1447 }
1448
1449 }
1450
1451 /* In text mode, add zeroes to cell vertex connectivity to mark face
1452 bounds (so as to add newlines) */
1453
1454 else { /* we are in text mode if f.bf == NULL */
1455
1456 part_vtx_idx[0] = 0;
1457 for (i = 0; i < section->n_elements; i++) {
1458 cell_length = 0;
1459 for (j = section->face_index[i]; j < section->face_index[i+1]; j++) {
1460 face_id = CS_ABS(section->face_num[j]) - 1;
1461 face_length = ( section->vertex_index[face_id+1]
1462 - section->vertex_index[face_id]);
1463 cell_length += face_length + 1;
1464 }
1465 part_vtx_idx[i+1] = part_vtx_idx[i] + cell_length;
1466 }
1467
1468 }
1469
1470 BFT_MALLOC(part_vtx_num, part_vtx_idx[section->n_elements], int32_t);
1471
1472 l = 0;
1473
1474 for (i = 0; i < section->n_elements; i++) {
1475 for (j = section->face_index[i]; j < section->face_index[i+1]; j++) {
1476 if (section->face_num[j] > 0) {
1477 face_id = section->face_num[j] - 1;
1478 for (k = section->vertex_index[face_id];
1479 k < section->vertex_index[face_id+1];
1480 k++)
1481 part_vtx_num[l++] = g_vtx_num[section->vertex_num[k] - 1];
1482 }
1483 else {
1484 face_id = -section->face_num[j] - 1;
1485 k = section->vertex_index[face_id];
1486 part_vtx_num[l++] = g_vtx_num[section->vertex_num[k] - 1];
1487 for (k = section->vertex_index[face_id+1] - 1;
1488 k > section->vertex_index[face_id];
1489 k--)
1490 part_vtx_num[l++] = g_vtx_num[section->vertex_num[k] - 1];
1491 }
1492 if (f.bf == NULL)
1493 part_vtx_num[l++] = 0; /* mark face limits in text mode */
1494 }
1495
1496 }
1497
1498 /* Now distribute and write cells -> vertices connectivity */
1499
1500 _write_indexed_connect_g(w,
1501 section->global_element_num,
1502 part_vtx_idx,
1503 part_vtx_num,
1504 f);
1505
1506 BFT_FREE(part_vtx_num);
1507 BFT_FREE(part_vtx_idx);
1508
1509 current_section = current_section->next;
1510
1511 } while ( current_section != NULL
1512 && current_section->continues_previous == true);
1513
1514 return current_section;
1515 }
1516
1517 #endif /* defined(HAVE_MPI) */
1518
1519 /*----------------------------------------------------------------------------
1520 * Write polyhedra from a nodal mesh to an EnSight Gold file in serial mode
1521 *
1522 * parameters:
1523 * export_section <-- pointer to EnSight section helper structure
1524 * f <-- associated file handle
1525 *
1526 * returns:
1527 * pointer to next EnSight section helper structure in list
1528 *----------------------------------------------------------------------------*/
1529
1530 static const fvm_writer_section_t *
_export_nodal_polyhedra_l(const fvm_writer_section_t * export_section,_ensight_file_t f)1531 _export_nodal_polyhedra_l(const fvm_writer_section_t *export_section,
1532 _ensight_file_t f)
1533 {
1534 int face_sgn;
1535 cs_lnum_t i, j, k, l;
1536 size_t i_buf;
1537
1538 cs_lnum_t face_length, face_id;
1539
1540 size_t buffer_size = 0;
1541 int32_t *buffer = NULL;
1542
1543 const fvm_writer_section_t *current_section;
1544
1545 /* Write number of faces per cell */
1546 /*--------------------------------*/
1547
1548 current_section = export_section;
1549
1550 do { /* loop on sections which should be appended */
1551
1552 const fvm_nodal_section_t *section = current_section->section;
1553
1554 if (f.tf != NULL) { /* Text mode */
1555 for (i = 0; i < section->n_elements; i++)
1556 fprintf(f.tf, "%10d\n",
1557 (int)( section->face_index[i+1]
1558 - section->face_index[i]));
1559 }
1560 else { /* binary mode */
1561
1562 /* First, allocate a buffer large enough so that the number of
1563 writes is limited, small enough so that the memory overhead is
1564 minimal; polyhedral connectivity is at least 4 faces x 3 vertices
1565 per cell, usually quite a bit more, so this is 1/3 of the minimum */
1566
1567 if (buffer_size < (size_t)section->n_elements * 4) {
1568 buffer_size = section->n_elements * 4;
1569 BFT_REALLOC(buffer, buffer_size, int32_t);
1570 }
1571
1572 /* Now fill buffer and write */
1573
1574 for (i = 0, i_buf = 0; i < section->n_elements; i++) {
1575 if (i_buf == buffer_size) {
1576 cs_file_write_global(f.bf, buffer, sizeof(int32_t), i_buf);
1577 i_buf = 0;
1578 }
1579 buffer[i_buf++] = (int)( section->face_index[i+1]
1580 - section->face_index[i]);
1581 }
1582 if (i_buf > 0)
1583 cs_file_write_global(f.bf, buffer, sizeof(int32_t), i_buf);
1584
1585 }
1586
1587 current_section = current_section->next;
1588
1589 } while ( current_section != NULL
1590 && current_section->continues_previous == true);
1591
1592 /* Write number of vertices/face */
1593 /*-------------------------------*/
1594
1595 current_section = export_section;
1596
1597 do { /* loop on sections which should be appended */
1598
1599 const fvm_nodal_section_t *section = current_section->section;
1600
1601 for (i = 0, i_buf = 0; i < section->n_elements; i++) {
1602
1603 /* Loop on cell faces */
1604
1605 for (j = section->face_index[i];
1606 j < section->face_index[i+1];
1607 j++) {
1608
1609 if (section->face_num[j] > 0)
1610 face_id = section->face_num[j] - 1;
1611 else
1612 face_id = -section->face_num[j] - 1;
1613
1614 face_length = ( section->vertex_index[face_id+1]
1615 - section->vertex_index[face_id]);
1616
1617 if (f.tf != NULL)
1618 fprintf(f.tf, "%10d\n",
1619 (int)face_length);
1620 else {
1621 if (i_buf == buffer_size) {
1622 cs_file_write_global(f.bf, buffer, sizeof(int32_t), i_buf);
1623 i_buf = 0;
1624 }
1625 buffer[i_buf++] = (int)face_length;
1626 }
1627
1628 }
1629
1630 }
1631
1632 if (f.bf != NULL && i_buf > 0)
1633 cs_file_write_global(f.bf, buffer, sizeof(int32_t), i_buf);
1634
1635 current_section = current_section->next;
1636
1637 } while ( current_section != NULL
1638 && current_section->continues_previous == true);
1639
1640 /* Write cell/vertex connectivity */
1641 /*--------------------------------*/
1642
1643 current_section = export_section;
1644
1645 do { /* loop on sections which should be appended */
1646
1647 const fvm_nodal_section_t *section = current_section->section;
1648
1649 for (i = 0, i_buf = 0; i < section->n_elements; i++) {
1650
1651 /* Loop on cell faces */
1652
1653 for (j = section->face_index[i];
1654 j < section->face_index[i+1];
1655 j++) {
1656
1657 /* Print face vertex numbers */
1658
1659 if (section->face_num[j] > 0) {
1660 face_id = section->face_num[j] - 1;
1661 face_sgn = 1;
1662 }
1663 else {
1664 face_id = -section->face_num[j] - 1;
1665 face_sgn = -1;
1666 }
1667
1668 face_length = ( section->vertex_index[face_id+1]
1669 - section->vertex_index[face_id]);
1670
1671 if (f.tf != NULL) { /* text mode */
1672 for (k = 0; k < face_length; k++) {
1673 l = section->vertex_index[face_id]
1674 + (face_length + (k*face_sgn))%face_length;
1675 fprintf(f.tf, "%10d", (int)section->vertex_num[l]);
1676 }
1677 fprintf(f.tf, "\n");
1678 }
1679 else { /* binary mode */
1680 for (k = 0; k < face_length; k++) {
1681 l = section->vertex_index[face_id]
1682 + (face_length + (k*face_sgn))%face_length;
1683 if (i_buf == buffer_size) {
1684 cs_file_write_global(f.bf, buffer, sizeof(int32_t), i_buf);
1685 i_buf = 0;
1686 }
1687 buffer[i_buf++] = (int)section->vertex_num[l];
1688 }
1689 }
1690
1691 } /* End of loop on cell faces */
1692
1693 } /* End of loop on polyhedral cells */
1694
1695 if (f.bf != NULL && i_buf > 0)
1696 cs_file_write_global(f.bf, buffer, sizeof(int32_t), i_buf);
1697
1698 current_section = current_section->next;
1699
1700 } while ( current_section != NULL
1701 && current_section->continues_previous == true);
1702
1703 if (buffer != NULL)
1704 BFT_FREE(buffer);
1705
1706 return current_section;
1707 }
1708
1709 #if defined(HAVE_MPI)
1710
1711 /*----------------------------------------------------------------------------
1712 * Write polygons from a nodal mesh to an EnSight Gold file in parallel mode
1713 *
1714 * parameters:
1715 * w <-- pointer to writer structure
1716 * export_section <-- pointer to EnSight section helper structure
1717 * global_vertex_num <-- pointer to vertex global numbering
1718 * f <-- associated file handle
1719 *
1720 * returns:
1721 * pointer to next EnSight section helper structure in list
1722 *----------------------------------------------------------------------------*/
1723
1724 static const fvm_writer_section_t *
_export_nodal_polygons_g(const fvm_to_ensight_writer_t * w,const fvm_writer_section_t * export_section,const fvm_io_num_t * global_vertex_num,_ensight_file_t f)1725 _export_nodal_polygons_g(const fvm_to_ensight_writer_t *w,
1726 const fvm_writer_section_t *export_section,
1727 const fvm_io_num_t *global_vertex_num,
1728 _ensight_file_t f)
1729 {
1730 const fvm_writer_section_t *current_section;
1731
1732 /* Export number of vertices per polygon by blocks */
1733 /*-------------------------------------------------*/
1734
1735 current_section = export_section;
1736
1737 do { /* loop on sections which should be appended */
1738
1739 const fvm_nodal_section_t *section = current_section->section;
1740
1741 _write_lengths_g(w,
1742 section->global_element_num,
1743 section->vertex_index,
1744 f);
1745
1746 current_section = current_section->next;
1747
1748 } while ( current_section != NULL
1749 && current_section->continues_previous == true);
1750
1751 /* Export face->vertex connectivity by blocks */
1752 /*--------------------------------------------*/
1753
1754 current_section = export_section;
1755
1756 do { /* loop on sections which should be appended */
1757
1758 cs_lnum_t i, j, k;
1759 cs_lnum_t *_part_vtx_idx = NULL;
1760 const cs_lnum_t *part_vtx_idx = NULL;
1761 int32_t *part_vtx_num = NULL;
1762
1763 const fvm_nodal_section_t *section = current_section->section;
1764 const cs_gnum_t *g_vtx_num
1765 = fvm_io_num_get_global_num(global_vertex_num);
1766
1767 if (f.bf != NULL) /* In binary mode, use existing index */
1768 part_vtx_idx = section->vertex_index;
1769
1770 /* In text mode, add zeroes to cell vertex connectivity to mark face
1771 bounds (so as to add newlines) */
1772
1773 else { /* we are in text mode if f.bf == NULL */
1774
1775 BFT_MALLOC(_part_vtx_idx, section->n_elements + 1, cs_lnum_t);
1776
1777 _part_vtx_idx[0] = 0;
1778 for (i = 0; i < section->n_elements; i++)
1779 _part_vtx_idx[i+1] = _part_vtx_idx[i] + ( section->vertex_index[i+1]
1780 - section->vertex_index[i]) + 1;
1781
1782 part_vtx_idx = _part_vtx_idx;
1783 }
1784
1785 /* Build connectivity array */
1786
1787 BFT_MALLOC(part_vtx_num, part_vtx_idx[section->n_elements], int32_t);
1788
1789 if (f.bf != NULL) { /* binary mode */
1790 for (i = 0, k = 0; i < section->n_elements; i++) {
1791 for (j = section->vertex_index[i];
1792 j < section->vertex_index[i+1];
1793 j++)
1794 part_vtx_num[k++] = g_vtx_num[section->vertex_num[j] - 1];
1795 }
1796 }
1797
1798 else { /* text mode */
1799 for (i = 0, k = 0; i < section->n_elements; i++) {
1800 for (j = section->vertex_index[i];
1801 j < section->vertex_index[i+1];
1802 j++)
1803 part_vtx_num[k++] = g_vtx_num[section->vertex_num[j] - 1];
1804 part_vtx_num[k++] = 0; /* mark face bounds in text mode */
1805 }
1806 }
1807
1808 /* Now distribute and write cell -> vertices connectivity */
1809
1810 _write_indexed_connect_g(w,
1811 section->global_element_num,
1812 part_vtx_idx,
1813 part_vtx_num,
1814 f);
1815
1816 BFT_FREE(part_vtx_num);
1817 if (_part_vtx_idx != NULL)
1818 BFT_FREE(_part_vtx_idx);
1819
1820 current_section = current_section->next;
1821
1822 } while ( current_section != NULL
1823 && current_section->continues_previous == true);
1824
1825 return current_section;
1826 }
1827
1828 #endif /* defined(HAVE_MPI) */
1829
1830 /*----------------------------------------------------------------------------
1831 * Write polygons from a nodal mesh to a text file in serial mode
1832 *
1833 * parameters:
1834 * export_section <-- pointer to EnSight section helper structure
1835 * f <-- associated file handle
1836 *
1837 * returns:
1838 * pointer to next EnSight section helper structure in list
1839 *----------------------------------------------------------------------------*/
1840
1841 static const fvm_writer_section_t *
_export_nodal_polygons_l(const fvm_writer_section_t * export_section,_ensight_file_t f)1842 _export_nodal_polygons_l(const fvm_writer_section_t *export_section,
1843 _ensight_file_t f)
1844
1845
1846 {
1847 cs_lnum_t i, j;
1848 size_t i_buf;
1849
1850 size_t buffer_size = 0;
1851 int32_t *buffer = NULL;
1852
1853 const fvm_writer_section_t *current_section = NULL;
1854
1855 /* Print face connectivity directly, without using extra buffers */
1856
1857 /* First loop on all polygonal faces, to write number of vertices */
1858 /*----------------------------------------------------------------*/
1859
1860 current_section = export_section;
1861
1862 do { /* Loop on sections that should be grouped */
1863
1864 const fvm_nodal_section_t *section = current_section->section;
1865
1866 if (f.tf != NULL) { /* Text mode */
1867 for (i = 0; i < section->n_elements; i++)
1868 fprintf(f.tf, "%10d\n", (int)( section->vertex_index[i+1]
1869 - section->vertex_index[i]));
1870 }
1871 else { /* binary mode */
1872
1873 /* First, allocate a buffer large enough so that the number of
1874 writes is limited, small enough so that the memory overhead is
1875 minimal; polygonal connectivity is at least 3 vertices per face,
1876 usually 5 or more, so this is 1/3 of the minimum */
1877
1878 if (buffer_size < (size_t)section->n_elements) {
1879 buffer_size = section->n_elements;
1880 BFT_REALLOC(buffer, buffer_size, int32_t);
1881 }
1882
1883 /* Now fill buffer and write */
1884
1885 for (i = 0, i_buf = 0; i < section->n_elements; i++) {
1886 if (i_buf == buffer_size) {
1887 cs_file_write_global(f.bf, buffer, sizeof(int32_t), i_buf);
1888 i_buf = 0;
1889 }
1890 buffer[i_buf++] = (int)( section->vertex_index[i+1]
1891 - section->vertex_index[i]);
1892 }
1893 if (i_buf > 0)
1894 cs_file_write_global(f.bf, buffer, sizeof(int32_t), i_buf);
1895 }
1896
1897 current_section = current_section->next;
1898
1899 } while ( current_section != NULL
1900 && current_section->continues_previous == true);
1901
1902 /* Loop on all polygonal faces */
1903 /*-----------------------------*/
1904
1905 current_section = export_section;
1906
1907 do { /* Loop on sections that should be grouped */
1908
1909 const fvm_nodal_section_t *section = current_section->section;
1910
1911 for (i = 0, i_buf = 0; i < section->n_elements; i++) {
1912
1913 /* Print face vertex numbers */
1914
1915 if (f.tf != NULL) { /* text mode */
1916 for (j = section->vertex_index[i];
1917 j < section->vertex_index[i+1];
1918 j++)
1919 fprintf(f.tf, "%10d", (int)section->vertex_num[j]);
1920 fprintf(f.tf, "\n");
1921 }
1922 else { /* binary mode */
1923 for (j = section->vertex_index[i];
1924 j < section->vertex_index[i+1];
1925 j++) {
1926 if (i_buf == buffer_size) {
1927 cs_file_write_global(f.bf, buffer, sizeof(int32_t), i_buf);
1928 i_buf = 0;
1929 }
1930 buffer[i_buf++] = (int)section->vertex_num[j];
1931 }
1932 }
1933
1934 } /* End of loop on polygonal faces */
1935
1936 if (f.bf != NULL && i_buf > 0)
1937 cs_file_write_global(f.bf, buffer, sizeof(int32_t), i_buf);
1938
1939 current_section = current_section->next;
1940
1941 } while ( current_section != NULL
1942 && current_section->continues_previous == true);
1943
1944 if (buffer != NULL)
1945 BFT_FREE(buffer);
1946
1947 return current_section;
1948 }
1949
1950 #if defined(HAVE_MPI)
1951
1952 /*----------------------------------------------------------------------------
1953 * Write tesselated element cell -> vertex connectivity to an EnSight Gold
1954 * file in parallel mode.
1955 *
1956 * parameters:
1957 * w <-- pointer to writer structure
1958 * global_vertex_num <-- vertex global numbering
1959 * global_element_num <-- global element numbering
1960 * tesselation <-- element tesselation description
1961 * type <-- tesselated sub-element type
1962 * extra_vertex_base <-- starting number for added vertices
1963 * f <-- associated file handle
1964 *----------------------------------------------------------------------------*/
1965
1966 static void
_write_tesselated_connect_g(const fvm_to_ensight_writer_t * w,const fvm_io_num_t * global_vertex_num,const fvm_io_num_t * global_element_num,const fvm_tesselation_t * tesselation,fvm_element_t type,const cs_gnum_t extra_vertex_base,_ensight_file_t f)1967 _write_tesselated_connect_g(const fvm_to_ensight_writer_t *w,
1968 const fvm_io_num_t *global_vertex_num,
1969 const fvm_io_num_t *global_element_num,
1970 const fvm_tesselation_t *tesselation,
1971 fvm_element_t type,
1972 const cs_gnum_t extra_vertex_base,
1973 _ensight_file_t f)
1974 {
1975 cs_block_dist_info_t bi;
1976
1977 cs_lnum_t part_size = 0;
1978
1979 cs_gnum_t n_g_sub_elements = 0;
1980 cs_gnum_t block_size = 0, block_start = 0, block_end = 0;
1981
1982 cs_part_to_block_t *d = NULL;
1983 cs_lnum_t *part_index, *block_index = NULL;
1984 int32_t *part_vtx_num = NULL, *block_vtx_num = NULL;
1985 cs_gnum_t *part_vtx_gnum = NULL;
1986
1987 size_t min_block_size = w->min_block_size / sizeof(int32_t);
1988
1989 const int stride = fvm_nodal_n_vertices_element[type];
1990 const cs_lnum_t n_elements = fvm_tesselation_n_elements(tesselation);
1991 const cs_gnum_t n_g_elements
1992 = fvm_io_num_get_global_count(global_element_num);
1993 const cs_lnum_t n_sub_elements
1994 = fvm_tesselation_n_sub_elements(tesselation, type);
1995 const cs_lnum_t *sub_element_idx
1996 = fvm_tesselation_sub_elt_index(tesselation, type);
1997 const cs_gnum_t *g_elt_num
1998 = fvm_io_num_get_global_num(global_element_num);
1999
2000 /* Adjust min block size based on mean number of sub-elements */
2001
2002 fvm_tesselation_get_global_size(tesselation,
2003 type,
2004 &n_g_sub_elements,
2005 NULL);
2006
2007 min_block_size /= ((n_g_sub_elements*1.)/n_g_elements) * stride;
2008
2009 /* Decode connectivity */
2010
2011 part_size = n_sub_elements * stride;
2012 assert(sub_element_idx[n_elements]*stride == part_size);
2013
2014 if (n_elements > 0) {
2015 BFT_MALLOC(part_vtx_num, part_size, int32_t);
2016 BFT_MALLOC(part_vtx_gnum, part_size, cs_gnum_t);
2017 }
2018
2019 fvm_tesselation_decode_g(tesselation,
2020 type,
2021 global_vertex_num,
2022 extra_vertex_base,
2023 part_vtx_gnum);
2024
2025 /* Convert to write type */
2026
2027 if (n_elements > 0) {
2028 for (cs_lnum_t i = 0; i < part_size; i++)
2029 part_vtx_num[i] = part_vtx_gnum[i];
2030 BFT_FREE(part_vtx_gnum);
2031 }
2032
2033 /* Allocate memory for additionnal indexes and decoded connectivity */
2034
2035 bi = cs_block_dist_compute_sizes(w->rank,
2036 w->n_ranks,
2037 w->min_rank_step,
2038 min_block_size,
2039 n_g_elements);
2040
2041 BFT_MALLOC(block_index, bi.gnum_range[1] - bi.gnum_range[0] + 1, cs_lnum_t);
2042 BFT_MALLOC(part_index, n_elements + 1, cs_lnum_t);
2043
2044 d = cs_part_to_block_create_by_gnum(w->comm, bi, n_elements, g_elt_num);
2045
2046 part_index[0] = 0;
2047 for (cs_lnum_t i = 0; i < n_elements; i++) {
2048 part_index[i+1] = part_index[i] + ( sub_element_idx[i+1]
2049 - sub_element_idx[i]) * stride;
2050 }
2051
2052 /* Copy index */
2053
2054 cs_part_to_block_copy_index(d,
2055 part_index,
2056 block_index);
2057
2058 block_size = (block_index[bi.gnum_range[1] - bi.gnum_range[0]]);
2059
2060 /* Copy connectivity */
2061
2062 BFT_MALLOC(block_vtx_num, block_size, int32_t);
2063
2064 cs_part_to_block_copy_indexed(d,
2065 CS_INT32,
2066 part_index,
2067 part_vtx_num,
2068 block_index,
2069 block_vtx_num);
2070
2071 cs_part_to_block_destroy(&d);
2072
2073 BFT_FREE(part_vtx_num);
2074 BFT_FREE(part_index);
2075 BFT_FREE(block_index);
2076
2077 /* Write to file */
2078
2079 block_size /= stride;
2080
2081 MPI_Scan(&block_size, &block_end, 1, CS_MPI_GNUM, MPI_SUM, w->comm);
2082 block_end += 1;
2083 block_start = block_end - block_size;
2084
2085 _write_block_connect_g(stride,
2086 block_start,
2087 block_end,
2088 block_vtx_num,
2089 w->comm,
2090 f);
2091
2092 /* Free remaining memory */
2093
2094 BFT_FREE(block_vtx_num);
2095 }
2096
2097 /*----------------------------------------------------------------------------
2098 * Write tesselated element connectivity from a nodal mesh to an EnSight Gold
2099 * file in parallel mode
2100 *
2101 * parameters:
2102 * w <-- pointer to writer structure
2103 * export_section <-- pointer to EnSight section helper structure
2104 * global_vertex_num <-- pointer to vertex global numbering
2105 * f <-- associated file handle
2106 *
2107 * returns:
2108 * pointer to next EnSight section helper structure in list
2109 *----------------------------------------------------------------------------*/
2110
2111 static const fvm_writer_section_t *
_export_nodal_tesselated_g(const fvm_to_ensight_writer_t * w,const fvm_writer_section_t * export_section,const fvm_io_num_t * global_vertex_num,_ensight_file_t f)2112 _export_nodal_tesselated_g(const fvm_to_ensight_writer_t *w,
2113 const fvm_writer_section_t *export_section,
2114 const fvm_io_num_t *global_vertex_num,
2115 _ensight_file_t f)
2116 {
2117 const fvm_writer_section_t *current_section;
2118
2119 /* Export face->vertex connectivity by blocks */
2120 /*--------------------------------------------*/
2121
2122 current_section = export_section;
2123
2124 do { /* loop on sections which should be appended */
2125
2126 const fvm_nodal_section_t *section = current_section->section;
2127
2128 _write_tesselated_connect_g(w,
2129 global_vertex_num,
2130 section->global_element_num,
2131 section->tesselation,
2132 current_section->type,
2133 current_section->extra_vertex_base,
2134 f);
2135
2136 current_section = current_section->next;
2137
2138 } while ( current_section != NULL
2139 && current_section->continues_previous == true
2140 && ( current_section->section->type
2141 == export_section->section->type));
2142
2143 return current_section;
2144 }
2145
2146 #endif /* defined(HAVE_MPI) */
2147
2148 /*----------------------------------------------------------------------------
2149 * Write tesselated element connectivity from a nodal mesh to an EnSight Gold
2150 * file in parallel mode
2151 *
2152 * parameters:
2153 * export_section <-- pointer to EnSight section helper structure
2154 * f <-- associated file handle
2155 *
2156 * returns:
2157 * pointer to next EnSight section helper structure in list
2158 *----------------------------------------------------------------------------*/
2159
2160 static const fvm_writer_section_t *
_export_nodal_tesselated_l(const fvm_writer_section_t * export_section,_ensight_file_t f)2161 _export_nodal_tesselated_l(const fvm_writer_section_t *export_section,
2162 _ensight_file_t f)
2163 {
2164 const fvm_writer_section_t *current_section;
2165
2166 current_section = export_section;
2167
2168 do { /* loop on sections which should be appended */
2169
2170 const fvm_nodal_section_t *section = current_section->section;
2171
2172 cs_lnum_t start_id, end_id;
2173 cs_lnum_t n_sub_elements_max;
2174 cs_lnum_t n_buffer_elements_max = section->n_elements;
2175 cs_lnum_t *vertex_num = NULL;
2176
2177 const cs_lnum_t *sub_element_idx
2178 = fvm_tesselation_sub_elt_index(section->tesselation,
2179 export_section->type);
2180
2181 fvm_tesselation_get_global_size(section->tesselation,
2182 export_section->type,
2183 NULL,
2184 &n_sub_elements_max);
2185 if (n_sub_elements_max > n_buffer_elements_max)
2186 n_buffer_elements_max = n_sub_elements_max;
2187
2188 BFT_MALLOC(vertex_num,
2189 ( n_buffer_elements_max
2190 * fvm_nodal_n_vertices_element[export_section->type]),
2191 cs_lnum_t);
2192
2193 for (start_id = 0;
2194 start_id < section->n_elements;
2195 start_id = end_id) {
2196
2197 end_id
2198 = fvm_tesselation_decode(section->tesselation,
2199 current_section->type,
2200 start_id,
2201 n_buffer_elements_max,
2202 export_section->extra_vertex_base,
2203 vertex_num);
2204
2205 _write_connect_l(fvm_nodal_n_vertices_element[export_section->type],
2206 ( sub_element_idx[end_id]
2207 - sub_element_idx[start_id]),
2208 vertex_num,
2209 f);
2210
2211 }
2212
2213 BFT_FREE(vertex_num);
2214
2215 current_section = current_section->next;
2216
2217 } while ( current_section != NULL
2218 && current_section->continues_previous == true
2219 && ( current_section->section->type
2220 == export_section->section->type));
2221
2222 return current_section;
2223 }
2224
2225 #if defined(HAVE_MPI)
2226
2227 /*----------------------------------------------------------------------------
2228 * Write strided elements from a nodal mesh to an EnSight Gold file in
2229 * parallel mode
2230 *
2231 * parameters:
2232 * w <-- pointer to writer structure
2233 * export_section <-- pointer to EnSight section helper structure
2234 * global_vertex_num <-- pointer to vertex global numbering
2235 * f <-- associated file handle
2236 *
2237 * returns:
2238 * pointer to next EnSight section helper structure in list
2239 *----------------------------------------------------------------------------*/
2240
2241 static const fvm_writer_section_t *
_export_nodal_strided_g(const fvm_to_ensight_writer_t * w,const fvm_writer_section_t * export_section,const fvm_io_num_t * global_vertex_num,_ensight_file_t f)2242 _export_nodal_strided_g(const fvm_to_ensight_writer_t *w,
2243 const fvm_writer_section_t *export_section,
2244 const fvm_io_num_t *global_vertex_num,
2245 _ensight_file_t f)
2246 {
2247 cs_lnum_t i, j;
2248
2249 const fvm_writer_section_t *current_section;
2250
2251 /* Export vertex connectivity */
2252
2253 current_section = export_section;
2254
2255 do { /* loop on sections which should be appended */
2256
2257 cs_block_dist_info_t bi;
2258
2259 cs_lnum_t block_size = 0;
2260 cs_part_to_block_t *d = NULL;
2261 int32_t *part_vtx_num = NULL, *block_vtx_num = NULL;
2262
2263 const fvm_nodal_section_t *section = current_section->section;
2264 const int stride = fvm_nodal_n_vertices_element[section->type];
2265
2266 const size_t min_block_size
2267 = w->min_block_size / (sizeof(int32_t) * stride);
2268
2269 const cs_lnum_t n_elements
2270 = fvm_io_num_get_local_count(section->global_element_num);
2271 const cs_gnum_t n_g_elements
2272 = fvm_io_num_get_global_count(section->global_element_num);
2273 const cs_gnum_t *g_elt_num
2274 = fvm_io_num_get_global_num(section->global_element_num);
2275 const cs_gnum_t *g_vtx_num
2276 = fvm_io_num_get_global_num(global_vertex_num);
2277
2278 /* Prepare distribution structures */
2279
2280 bi = cs_block_dist_compute_sizes(w->rank,
2281 w->n_ranks,
2282 w->min_rank_step,
2283 min_block_size,
2284 n_g_elements);
2285
2286 d = cs_part_to_block_create_by_gnum(w->comm,
2287 bi,
2288 n_elements,
2289 g_elt_num);
2290
2291 /* Build connectivity */
2292
2293 block_size = bi.gnum_range[1] - bi.gnum_range[0];
2294
2295 BFT_MALLOC(block_vtx_num, block_size*stride, int32_t);
2296 BFT_MALLOC(part_vtx_num, n_elements*stride, int32_t);
2297
2298 for (i = 0; i < n_elements; i++) {
2299 for (j = 0; j < stride; j++) {
2300 part_vtx_num[i*stride + j]
2301 = g_vtx_num[section->vertex_num[i*stride + j] - 1];
2302 }
2303 }
2304
2305 cs_part_to_block_copy_array(d,
2306 CS_INT32,
2307 stride,
2308 part_vtx_num,
2309 block_vtx_num);
2310
2311 BFT_FREE(part_vtx_num);
2312
2313 _write_block_connect_g(stride,
2314 bi.gnum_range[0],
2315 bi.gnum_range[1],
2316 block_vtx_num,
2317 w->comm,
2318 f);
2319
2320 BFT_FREE(block_vtx_num);
2321
2322 cs_part_to_block_destroy(&d);
2323
2324 current_section = current_section->next;
2325
2326 } while ( current_section != NULL
2327 && current_section->continues_previous == true
2328 && ( current_section->section->type
2329 == export_section->section->type));
2330
2331 return current_section;
2332 }
2333
2334 #endif /* defined(HAVE_MPI) */
2335
2336 /*----------------------------------------------------------------------------
2337 * Write field values associated with nodal values of a nodal mesh to
2338 * an EnSight Gold file in serial mode.
2339 *
2340 * Output fields ar either scalar or 3d vectors or scalars, and are
2341 * non interlaced. Input arrays may be less than 2d, in which case the z
2342 * values are set to 0, and may be interlaced or not.
2343 *
2344 * parameters:
2345 * n_entities <-- number of entities
2346 * input_dim <-- input field dimension
2347 * output_dim <-- output field dimension
2348 * interlace <-- indicates if field in memory is interlaced
2349 * n_parent_lists <-- indicates if field values are to be obtained
2350 * directly through the local entity index (when 0) or
2351 * through the parent entity numbers (when 1 or more)
2352 * parent_num_shift <-- parent list to common number index shifts;
2353 * size: n_parent_lists
2354 * datatype <-- input data type (output is real)
2355 * field_values <-- array of associated field value arrays
2356 * f <-- associated file handle
2357 *----------------------------------------------------------------------------*/
2358
2359 static void
_export_field_values_nl(const fvm_nodal_t * mesh,fvm_writer_field_helper_t * helper,int input_dim,cs_interlace_t interlace,int n_parent_lists,const cs_lnum_t parent_num_shift[],cs_datatype_t datatype,const void * const field_values[],_ensight_file_t f)2360 _export_field_values_nl(const fvm_nodal_t *mesh,
2361 fvm_writer_field_helper_t *helper,
2362 int input_dim,
2363 cs_interlace_t interlace,
2364 int n_parent_lists,
2365 const cs_lnum_t parent_num_shift[],
2366 cs_datatype_t datatype,
2367 const void *const field_values[],
2368 _ensight_file_t f)
2369 {
2370 int i;
2371 size_t output_size;
2372 float *output_buffer;
2373
2374 int output_dim = fvm_writer_field_helper_field_dim(helper);
2375
2376 const size_t output_buffer_size
2377 = mesh->n_vertices > 16 ? (mesh->n_vertices / 4) : mesh->n_vertices;
2378
2379 BFT_MALLOC(output_buffer, output_buffer_size, float);
2380
2381 for (i = 0; i < output_dim; i++) {
2382
2383 const int i_in = (input_dim == 6) ? _ensight_c_order_6[i] : i;
2384
2385 while (fvm_writer_field_helper_step_nl(helper,
2386 mesh,
2387 input_dim,
2388 i_in,
2389 interlace,
2390 n_parent_lists,
2391 parent_num_shift,
2392 datatype,
2393 field_values,
2394 output_buffer,
2395 output_buffer_size,
2396 &output_size) == 0) {
2397
2398 _write_block_floats_l(output_size,
2399 output_buffer,
2400 f);
2401
2402 }
2403 }
2404
2405 BFT_FREE(output_buffer);
2406 }
2407
2408 /*----------------------------------------------------------------------------
2409 * Write field values associated with element values of a nodal mesh to
2410 * an EnSight Gold file.
2411 *
2412 * Output fields ar either scalar or 3d vectors or scalars, and are
2413 * non interlaced. Input arrays may be less than 2d, in which case the z
2414 * values are set to 0, and may be interlaced or not.
2415 *
2416 * parameters:
2417 * export_section <-- pointer to EnSight section helper structure
2418 * helper <-- pointer to general writer helper structure
2419 * input_dim <-- input field dimension
2420 * interlace <-- indicates if field in memory is interlaced
2421 * n_parent_lists <-- indicates if field values are to be obtained
2422 * directly through the local entity index (when 0) or
2423 * through the parent entity numbers (when 1 or more)
2424 * parent_num_shift <-- parent list to common number index shifts;
2425 * size: n_parent_lists
2426 * datatype <-- indicates the data type of (source) field values
2427 * field_values <-- array of associated field value arrays
2428 * f <-- associated file handle
2429 *
2430 * returns:
2431 * pointer to next EnSight section helper structure in list
2432 *----------------------------------------------------------------------------*/
2433
2434 static const fvm_writer_section_t *
_export_field_values_el(const fvm_writer_section_t * export_section,fvm_writer_field_helper_t * helper,int input_dim,cs_interlace_t interlace,int n_parent_lists,const cs_lnum_t parent_num_shift[],cs_datatype_t datatype,const void * const field_values[],_ensight_file_t f)2435 _export_field_values_el(const fvm_writer_section_t *export_section,
2436 fvm_writer_field_helper_t *helper,
2437 int input_dim,
2438 cs_interlace_t interlace,
2439 int n_parent_lists,
2440 const cs_lnum_t parent_num_shift[],
2441 cs_datatype_t datatype,
2442 const void *const field_values[],
2443 _ensight_file_t f)
2444 {
2445 int i;
2446 size_t input_size = 0, output_size = 0;
2447 size_t min_output_buffer_size = 0, output_buffer_size = 0;
2448 float *output_buffer = NULL;
2449
2450 const fvm_writer_section_t *current_section = NULL;
2451
2452 int output_dim = fvm_writer_field_helper_field_dim(helper);
2453
2454 /* Blocking for arbitrary buffer size, but should be small enough
2455 to add little additional memory requirement (in proportion), large
2456 enough to limit number of write calls. */
2457
2458 fvm_writer_field_helper_get_size(helper,
2459 &input_size,
2460 &output_size,
2461 &min_output_buffer_size);
2462
2463 output_buffer_size = input_size / 4;
2464 output_buffer_size = CS_MAX(output_buffer_size, min_output_buffer_size);
2465 output_buffer_size = CS_MAX(output_buffer_size, 128);
2466 output_buffer_size = CS_MIN(output_buffer_size, output_size);
2467
2468 BFT_MALLOC(output_buffer, output_buffer_size, float);
2469
2470 /* Loop on dimension (de-interlace vectors, always 3D for EnSight) */
2471
2472 for (i = 0; i < output_dim; i++) {
2473
2474 bool loop_on_sections = true;
2475
2476 const int i_in = (input_dim == 6) ? _ensight_c_order_6[i] : i;
2477
2478 current_section = export_section;
2479
2480 while (loop_on_sections == true) {
2481
2482 while (fvm_writer_field_helper_step_el(helper,
2483 current_section,
2484 input_dim,
2485 i_in,
2486 interlace,
2487 n_parent_lists,
2488 parent_num_shift,
2489 datatype,
2490 field_values,
2491 output_buffer,
2492 output_buffer_size,
2493 &output_size) == 0) {
2494
2495 _write_block_floats_l(output_size,
2496 output_buffer,
2497 f);
2498
2499 }
2500
2501 current_section = current_section->next;
2502
2503 if ( current_section == NULL
2504 || current_section->continues_previous == false)
2505 loop_on_sections = false;
2506
2507 } /* while (loop on sections) */
2508
2509 } /* end of loop on spatial dimension */
2510
2511 BFT_FREE(output_buffer);
2512
2513 return current_section;
2514 }
2515
2516 /*! (DOXYGEN_SHOULD_SKIP_THIS) \endcond */
2517
2518 /*============================================================================
2519 * Public function definitions
2520 *============================================================================*/
2521
2522 /*----------------------------------------------------------------------------
2523 * Initialize FVM to EnSight Gold file writer.
2524 *
2525 * Options are:
2526 * text output text files
2527 * binary output binary files (default)
2528 * big_endian force binary files to big-endian
2529 * discard_polygons do not output polygons or related values
2530 * discard_polyhedra do not output polyhedra or related values
2531 * divide_polygons tesselate polygons with triangles
2532 * divide_polyhedra tesselate polyhedra with tetrahedra and pyramids
2533 * (adding a vertex near each polyhedron's center)
2534 *
2535 * parameters:
2536 * name <-- base output case name.
2537 * options <-- whitespace separated, lowercase options list
2538 * time_dependecy <-- indicates if and how meshes will change with time
2539 * comm <-- associated MPI communicator.
2540 *
2541 * returns:
2542 * pointer to opaque EnSight Gold writer structure.
2543 *----------------------------------------------------------------------------*/
2544
2545 #if defined(HAVE_MPI)
2546 void *
fvm_to_ensight_init_writer(const char * name,const char * path,const char * options,fvm_writer_time_dep_t time_dependency,MPI_Comm comm)2547 fvm_to_ensight_init_writer(const char *name,
2548 const char *path,
2549 const char *options,
2550 fvm_writer_time_dep_t time_dependency,
2551 MPI_Comm comm)
2552 #else
2553 void *
2554 fvm_to_ensight_init_writer(const char *name,
2555 const char *path,
2556 const char *options,
2557 fvm_writer_time_dep_t time_dependency)
2558 #endif
2559 {
2560 fvm_to_ensight_writer_t *this_writer = NULL;
2561
2562 /* Initialize writer */
2563
2564 BFT_MALLOC(this_writer, 1, fvm_to_ensight_writer_t);
2565
2566 BFT_MALLOC(this_writer->name, strlen(name) + 1, char);
2567 strcpy(this_writer->name, name);
2568
2569 this_writer->text_mode = false;
2570 this_writer->swap_endian = false;
2571 this_writer->discard_polygons = false;
2572 this_writer->discard_polyhedra = false;
2573 this_writer->divide_polygons = false;
2574 this_writer->divide_polyhedra = false;
2575
2576 this_writer->rank = 0;
2577 this_writer->n_ranks = 1;
2578
2579 #if defined(HAVE_MPI)
2580 {
2581 int mpi_flag, rank, n_ranks;
2582 MPI_Comm w_block_comm, w_comm;
2583 this_writer->min_rank_step = 1;
2584 this_writer->min_block_size = 0;
2585 this_writer->block_comm = MPI_COMM_NULL;
2586 this_writer->comm = MPI_COMM_NULL;
2587 MPI_Initialized(&mpi_flag);
2588 if (mpi_flag && comm != MPI_COMM_NULL) {
2589 size_t min_block_size = cs_parall_get_min_coll_buf_size();
2590 this_writer->comm = comm;
2591 MPI_Comm_rank(this_writer->comm, &rank);
2592 MPI_Comm_size(this_writer->comm, &n_ranks);
2593 this_writer->rank = rank;
2594 this_writer->n_ranks = n_ranks;
2595 cs_file_get_default_comm(NULL, &w_block_comm, &w_comm);
2596 if (comm == w_comm) {
2597 this_writer->min_block_size = min_block_size;
2598 this_writer->block_comm = w_block_comm;
2599 }
2600 this_writer->comm = comm;
2601 }
2602 }
2603 #endif /* defined(HAVE_MPI) */
2604
2605 /* Parse options */
2606
2607 if (options != NULL) {
2608
2609 int i1, i2, l_opt;
2610 int l_tot = strlen(options);
2611
2612 i1 = 0; i2 = 0;
2613 while (i1 < l_tot) {
2614
2615 for (i2 = i1; i2 < l_tot && options[i2] != ' '; i2++);
2616 l_opt = i2 - i1;
2617
2618 if ((l_opt == 4) && (strncmp(options + i1, "text", l_opt) == 0))
2619 this_writer->text_mode = true;
2620 else if ((l_opt == 6) && (strncmp(options + i1, "binary", l_opt) == 0))
2621 this_writer->text_mode = false;
2622
2623 else if ( (l_opt == 10)
2624 && (strncmp(options + i1, "big_endian", l_opt) == 0)) {
2625 int int_endian = 0;
2626 this_writer->text_mode = false;
2627 /* Check if system is "big-endian" or "little-endian" */
2628 *((char *)(&int_endian)) = '\1';
2629 if (int_endian == 1)
2630 this_writer->swap_endian = 1;
2631 }
2632
2633 else if ( (l_opt == 16)
2634 && (strncmp(options + i1, "discard_polygons", l_opt) == 0))
2635 this_writer->discard_polygons = true;
2636 else if ( (l_opt == 17)
2637 && (strncmp(options + i1, "discard_polyhedra", l_opt) == 0))
2638 this_writer->discard_polyhedra = true;
2639
2640 else if ( (l_opt == 15)
2641 && (strncmp(options + i1, "divide_polygons", l_opt) == 0))
2642 this_writer->divide_polygons = true;
2643 else if ( (l_opt == 16)
2644 && (strncmp(options + i1, "divide_polyhedra", l_opt) == 0))
2645 this_writer->divide_polyhedra = true;
2646
2647 for (i1 = i2 + 1; i1 < l_tot && options[i1] == ' '; i1++);
2648
2649 }
2650
2651 }
2652
2653 this_writer->case_info = fvm_to_ensight_case_create(name,
2654 path,
2655 time_dependency);
2656
2657 /* Return writer */
2658
2659 return this_writer;
2660 }
2661
2662 /*----------------------------------------------------------------------------
2663 * Finalize FVM to EnSight Gold file writer.
2664 *
2665 * parameters:
2666 * this_writer_p <-- pointer to opaque Ensight Gold writer structure.
2667 *
2668 * returns:
2669 * NULL pointer
2670 *----------------------------------------------------------------------------*/
2671
2672 void *
fvm_to_ensight_finalize_writer(void * this_writer_p)2673 fvm_to_ensight_finalize_writer(void *this_writer_p)
2674 {
2675 fvm_to_ensight_writer_t *this_writer
2676 = (fvm_to_ensight_writer_t *)this_writer_p;
2677
2678 BFT_FREE(this_writer->name);
2679
2680 fvm_to_ensight_case_destroy(this_writer->case_info);
2681
2682 BFT_FREE(this_writer);
2683
2684 return NULL;
2685 }
2686
2687 /*----------------------------------------------------------------------------
2688 * Associate new time step with an EnSight geometry.
2689 *
2690 * parameters:
2691 * this_writer_p <-- pointer to associated writer
2692 * time_step <-- time step number
2693 * time_value <-- time_value number
2694 *----------------------------------------------------------------------------*/
2695
2696 void
fvm_to_ensight_set_mesh_time(void * this_writer_p,const int time_step,const double time_value)2697 fvm_to_ensight_set_mesh_time(void *this_writer_p,
2698 const int time_step,
2699 const double time_value)
2700 {
2701 fvm_to_ensight_writer_t *this_writer
2702 = (fvm_to_ensight_writer_t *)this_writer_p;
2703
2704 fvm_to_ensight_case_set_geom_time(this_writer->case_info,
2705 time_step,
2706 time_value);
2707 }
2708
2709 /*----------------------------------------------------------------------------
2710 * Indicate if a elements of a given type in a mesh associated to a given
2711 * EnSight Gold file writer need to be tesselated.
2712 *
2713 * parameters:
2714 * this_writer_p <-- pointer to associated writer
2715 * mesh <-- pointer to nodal mesh structure that should be written
2716 * element_type <-- element type we are interested in
2717 *
2718 * returns:
2719 * 1 if tesselation of the given element type is needed, 0 otherwise
2720 *----------------------------------------------------------------------------*/
2721
2722 int
fvm_to_ensight_needs_tesselation(void * this_writer_p,const fvm_nodal_t * mesh,fvm_element_t element_type)2723 fvm_to_ensight_needs_tesselation(void *this_writer_p,
2724 const fvm_nodal_t *mesh,
2725 fvm_element_t element_type)
2726 {
2727 int i;
2728 int retval = 0;
2729 fvm_to_ensight_writer_t *this_writer
2730 = (fvm_to_ensight_writer_t *)this_writer_p;
2731
2732 const int export_dim = fvm_nodal_get_max_entity_dim(mesh);
2733
2734 /* Initial count and allocation */
2735
2736 if ( ( element_type == FVM_FACE_POLY
2737 && this_writer->divide_polygons == true)
2738 || ( element_type == FVM_CELL_POLY
2739 && this_writer->divide_polyhedra == true)) {
2740
2741 for (i = 0; i < mesh->n_sections; i++) {
2742
2743 const fvm_nodal_section_t *const section = mesh->sections[i];
2744
2745 /* Output if entity dimension equal to highest in mesh
2746 (i.e. no output of faces if cells present, or edges
2747 if cells or faces) */
2748
2749 if (section->entity_dim == export_dim) {
2750 if (section->type == element_type)
2751 retval = 1;
2752 }
2753
2754 }
2755
2756 }
2757
2758 return retval;
2759 }
2760
2761 /*----------------------------------------------------------------------------
2762 * Write nodal mesh to a an EnSight Gold file
2763 *
2764 * parameters:
2765 * this_writer_p <-- pointer to associated writer
2766 * mesh <-- pointer to nodal mesh structure that should be written
2767 *----------------------------------------------------------------------------*/
2768
2769 void
fvm_to_ensight_export_nodal(void * this_writer_p,const fvm_nodal_t * mesh)2770 fvm_to_ensight_export_nodal(void *this_writer_p,
2771 const fvm_nodal_t *mesh)
2772 {
2773 int part_num;
2774
2775 const fvm_writer_section_t *export_section = NULL;
2776 fvm_writer_section_t *export_list = NULL;
2777 fvm_to_ensight_writer_t *this_writer
2778 = (fvm_to_ensight_writer_t *)this_writer_p;
2779 _ensight_file_t f = {NULL, NULL};
2780
2781 const int rank = this_writer->rank;
2782 const int n_ranks = this_writer->n_ranks;
2783
2784 /* Initialization */
2785 /*----------------*/
2786
2787 fvm_to_ensight_case_file_info_t file_info;
2788
2789 /* Get part number */
2790
2791 part_num = fvm_to_ensight_case_get_part_num(this_writer->case_info,
2792 mesh->name);
2793 if (part_num == 0)
2794 part_num = fvm_to_ensight_case_add_part(this_writer->case_info,
2795 mesh->name);
2796
2797 /* Open geometry file in append mode */
2798
2799 file_info = fvm_to_ensight_case_get_geom_file(this_writer->case_info);
2800
2801 f = _open_ensight_file(this_writer,
2802 file_info.name,
2803 file_info.queried);
2804
2805 if (file_info.queried == false)
2806 _write_geom_headers(this_writer, f);
2807
2808 /* Part header */
2809
2810 _write_string(f, "part");
2811 _write_int(f, part_num);
2812 if (mesh->name != NULL)
2813 _write_string(f, mesh->name);
2814 else
2815 _write_string(f, "unnamed");
2816
2817 /* Build list of sections that are used here, in order of output */
2818
2819 int export_dim = fvm_nodal_get_max_entity_dim(mesh);
2820
2821 export_list = fvm_writer_export_list(mesh,
2822 export_dim,
2823 export_dim,
2824 -1,
2825 true,
2826 false,
2827 this_writer->discard_polygons,
2828 this_writer->discard_polyhedra,
2829 this_writer->divide_polygons,
2830 this_writer->divide_polyhedra);
2831
2832 /* Vertex coordinates */
2833 /*--------------------*/
2834
2835 #if defined(HAVE_MPI)
2836 if (n_ranks > 1)
2837 _export_vertex_coords_g(this_writer, mesh, f);
2838 #endif
2839
2840 if (n_ranks == 1)
2841 _export_vertex_coords_l(this_writer, mesh, f);
2842
2843 /* If no sections are present (i.e. we may only have vertices),
2844 add "point" elements */
2845
2846 if (export_list == NULL) {
2847
2848 #if defined(HAVE_MPI)
2849 if (n_ranks > 1)
2850 _export_point_elements_g(this_writer, mesh, f);
2851 #endif
2852 if (n_ranks == 1)
2853 _export_point_elements_l(mesh, f);
2854
2855 }
2856
2857 /* Element connectivity */
2858 /*----------------------*/
2859
2860 export_section = export_list;
2861
2862 while (export_section != NULL) {
2863
2864 const fvm_nodal_section_t *section = export_section->section;
2865
2866 /* Print header if start of corresponding EnSight section */
2867
2868 if (export_section->continues_previous == false) {
2869
2870 cs_gnum_t n_g_elements = 0;
2871 const fvm_writer_section_t *next_section = export_section;
2872
2873 do {
2874
2875 if (next_section->section->type == export_section->type)
2876 n_g_elements += fvm_nodal_section_n_g_elements(next_section->section);
2877
2878 else {
2879 cs_gnum_t n_g_sub_elements = 0;
2880 fvm_tesselation_get_global_size(next_section->section->tesselation,
2881 next_section->type,
2882 &n_g_sub_elements,
2883 NULL);
2884 n_g_elements += n_g_sub_elements;
2885 }
2886
2887 next_section = next_section->next;
2888
2889 } while (next_section != NULL && next_section->continues_previous == true);
2890
2891 _write_string(f, _ensight_type_name[export_section->type]);
2892 _write_int(f, n_g_elements);
2893 }
2894
2895 /* Output for strided (regular) element types */
2896 /*--------------------------------------------*/
2897
2898 if (section->stride > 0) {
2899
2900 #if defined(HAVE_MPI)
2901
2902 if (n_ranks > 1)
2903 export_section = _export_nodal_strided_g(this_writer,
2904 export_section,
2905 mesh->global_vertex_num,
2906 f);
2907
2908 #endif /* defined(HAVE_MPI) */
2909
2910 if (n_ranks == 1) { /* start of output in serial mode */
2911
2912 _write_connect_l(section->stride,
2913 section->n_elements,
2914 section->vertex_num,
2915 f);
2916
2917 export_section = export_section->next;
2918
2919 }
2920
2921 } /* end of output for strided element types */
2922
2923 /* Output for tesselated polygons or polyhedra */
2924 /*---------------------------------------------*/
2925
2926 else if (export_section->type != section->type) {
2927
2928 #if defined(HAVE_MPI)
2929
2930 /* output in parallel mode */
2931
2932 if (n_ranks > 1)
2933 export_section = _export_nodal_tesselated_g(this_writer,
2934 export_section,
2935 mesh->global_vertex_num,
2936 f);
2937 #endif /* defined(HAVE_MPI) */
2938
2939 if (n_ranks == 1)
2940 export_section = _export_nodal_tesselated_l(export_section,
2941 f);
2942
2943 }
2944
2945 /* Output for polygons */
2946 /*---------------------*/
2947
2948 else if (export_section->type == FVM_FACE_POLY) {
2949 #if defined(HAVE_MPI)
2950
2951 /* output in parallel mode */
2952
2953 if (n_ranks > 1)
2954 export_section = _export_nodal_polygons_g(this_writer,
2955 export_section,
2956 mesh->global_vertex_num,
2957 f);
2958 #endif /* defined(HAVE_MPI) */
2959
2960 if (n_ranks == 1)
2961 export_section = _export_nodal_polygons_l(export_section,
2962 f);
2963
2964 }
2965
2966 /* Output for polyhedra */
2967 /*----------------------*/
2968
2969 else if (export_section->type == FVM_CELL_POLY) {
2970
2971 #if defined(HAVE_MPI)
2972
2973 /* output in parallel mode */
2974
2975 if (n_ranks > 1)
2976 export_section =_export_nodal_polyhedra_g(this_writer,
2977 export_section,
2978 mesh->global_vertex_num,
2979 f);
2980
2981 #endif /* defined(HAVE_MPI) */
2982
2983 if (n_ranks == 1)
2984 export_section = _export_nodal_polyhedra_l(export_section,
2985 f);
2986
2987 }
2988
2989 } /* End of loop on sections */
2990
2991 BFT_FREE(export_list);
2992
2993 /* Close geometry file and update case file */
2994 /*------------------------------------------*/
2995
2996 _free_ensight_file(&f);
2997
2998 fvm_to_ensight_case_write_case(this_writer->case_info, rank);
2999 }
3000
3001 /*----------------------------------------------------------------------------
3002 * Write field associated with a nodal mesh to an EnSight Gold file.
3003 *
3004 * Assigning a negative value to the time step indicates a time-independent
3005 * field (in which case the time_value argument is unused).
3006 *
3007 * parameters:
3008 * this_writer_p <-- pointer to associated writer
3009 * mesh <-- pointer to associated nodal mesh structure
3010 * name <-- variable name
3011 * location <-- variable definition location (nodes or elements)
3012 * dimension <-- variable dimension (0: constant, 1: scalar,
3013 * 3: vector, 6: sym. tensor, 9: asym. tensor)
3014 * interlace <-- indicates if variable in memory is interlaced
3015 * n_parent_lists <-- indicates if variable values are to be obtained
3016 * directly through the local entity index (when 0) or
3017 * through the parent entity numbers (when 1 or more)
3018 * parent_num_shift <-- parent number to value array index shifts;
3019 * size: n_parent_lists
3020 * datatype <-- indicates the data type of (source) field values
3021 * time_step <-- number of the current time step
3022 * time_value <-- associated time value
3023 * field_values <-- array of associated field value arrays
3024 *----------------------------------------------------------------------------*/
3025
3026 void
fvm_to_ensight_export_field(void * this_writer_p,const fvm_nodal_t * mesh,const char * name,fvm_writer_var_loc_t location,int dimension,cs_interlace_t interlace,int n_parent_lists,const cs_lnum_t parent_num_shift[],cs_datatype_t datatype,int time_step,double time_value,const void * const field_values[])3027 fvm_to_ensight_export_field(void *this_writer_p,
3028 const fvm_nodal_t *mesh,
3029 const char *name,
3030 fvm_writer_var_loc_t location,
3031 int dimension,
3032 cs_interlace_t interlace,
3033 int n_parent_lists,
3034 const cs_lnum_t parent_num_shift[],
3035 cs_datatype_t datatype,
3036 int time_step,
3037 double time_value,
3038 const void *const field_values[])
3039 {
3040 int output_dim, part_num;
3041 fvm_to_ensight_case_file_info_t file_info;
3042
3043 const fvm_writer_section_t *export_section = NULL;
3044 fvm_writer_field_helper_t *helper = NULL;
3045 fvm_writer_section_t *export_list = NULL;
3046 fvm_to_ensight_writer_t *w = (fvm_to_ensight_writer_t *)this_writer_p;
3047 _ensight_file_t f = {NULL, NULL};
3048
3049 const int rank = w->rank;
3050 const int n_ranks = w->n_ranks;
3051
3052 /* Initialization */
3053 /*----------------*/
3054
3055 /* Dimension */
3056
3057 output_dim = dimension;
3058 if (dimension == 2)
3059 output_dim = 3;
3060 else if (dimension > 3 && dimension != 6 && dimension != 9)
3061 bft_error(__FILE__, __LINE__, 0,
3062 _("Data of dimension %d not handled"), dimension);
3063
3064 const int *comp_order = (dimension == 6) ? _ensight_c_order_6 : NULL;
3065
3066 /* Get part number */
3067
3068 part_num = fvm_to_ensight_case_get_part_num(w->case_info,
3069 mesh->name);
3070 if (part_num == 0)
3071 part_num = fvm_to_ensight_case_add_part(w->case_info,
3072 mesh->name);
3073
3074 /* Open variable file */
3075
3076 file_info = fvm_to_ensight_case_get_var_file(w->case_info,
3077 name,
3078 output_dim,
3079 location,
3080 time_step,
3081 time_value);
3082
3083 f = _open_ensight_file(w, file_info.name, file_info.queried);
3084
3085 if (file_info.queried == false) {
3086
3087 char buf[81] = "";
3088
3089 /* New files start with description line */
3090 #if HAVE_SNPRINTF
3091 if (time_step > -1)
3092 snprintf(buf, 80, "%s (time values: %d, %g)",
3093 name, time_step, time_value);
3094 else
3095 strncpy(buf, name, 80);
3096 #else
3097 strncpy(buf, name, 80);
3098 #endif
3099 buf[80] = '\0';
3100 _write_string(f, buf);
3101 }
3102
3103 /* Initialize writer helper */
3104 /*--------------------------*/
3105
3106 /* Build list of sections that are used here, in order of output */
3107
3108 int export_dim = fvm_nodal_get_max_entity_dim(mesh);
3109
3110 export_list = fvm_writer_export_list(mesh,
3111 export_dim,
3112 export_dim,
3113 -1,
3114 true,
3115 false,
3116 w->discard_polygons,
3117 w->discard_polyhedra,
3118 w->divide_polygons,
3119 w->divide_polyhedra);
3120
3121 helper = fvm_writer_field_helper_create(mesh,
3122 export_list,
3123 output_dim,
3124 CS_NO_INTERLACE,
3125 CS_FLOAT,
3126 location);
3127
3128 #if defined(HAVE_MPI)
3129
3130 if (n_ranks > 1)
3131 fvm_writer_field_helper_init_g(helper,
3132 w->min_rank_step,
3133 w->min_block_size,
3134 w->comm);
3135
3136 #endif
3137
3138 /* Part header */
3139
3140 _write_string(f, "part");
3141 _write_int(f, part_num);
3142
3143 /* Per node variable */
3144 /*-------------------*/
3145
3146 if (location == FVM_WRITER_PER_NODE) {
3147
3148 _write_string(f, "coordinates");
3149
3150 #if defined(HAVE_MPI)
3151
3152 if (n_ranks > 1) {
3153
3154 _ensight_context_t c;
3155 c.writer = w;
3156 c.file = &f;
3157
3158 fvm_writer_field_helper_output_n(helper,
3159 &c,
3160 mesh,
3161 dimension,
3162 interlace,
3163 comp_order,
3164 n_parent_lists,
3165 parent_num_shift,
3166 datatype,
3167 field_values,
3168 _field_output_g);
3169
3170 }
3171
3172 #endif /* defined(HAVE_MPI) */
3173
3174 if (n_ranks == 1)
3175 _export_field_values_nl(mesh,
3176 helper,
3177 dimension,
3178 interlace,
3179 n_parent_lists,
3180 parent_num_shift,
3181 datatype,
3182 field_values,
3183 f);
3184 }
3185
3186 /* Per element variable */
3187 /*----------------------*/
3188
3189 else if (location == FVM_WRITER_PER_ELEMENT) {
3190
3191 export_section = export_list;
3192
3193 while (export_section != NULL) {
3194
3195 /* Print header if start of corresponding EnSight section */
3196
3197 if (export_section->continues_previous == false)
3198 _write_string(f, _ensight_type_name[export_section->type]);
3199
3200 /* Output per grouped sections */
3201
3202 #if defined(HAVE_MPI)
3203
3204 if (n_ranks > 1) {
3205
3206 _ensight_context_t c;
3207 c.writer = w;
3208 c.file = &f;
3209
3210 export_section = fvm_writer_field_helper_output_e(helper,
3211 &c,
3212 export_section,
3213 dimension,
3214 interlace,
3215 comp_order,
3216 n_parent_lists,
3217 parent_num_shift,
3218 datatype,
3219 field_values,
3220 _field_output_g);
3221
3222 }
3223
3224 #endif /* defined(HAVE_MPI) */
3225
3226 if (n_ranks == 1)
3227 export_section = _export_field_values_el(export_section,
3228 helper,
3229 dimension,
3230 interlace,
3231 n_parent_lists,
3232 parent_num_shift,
3233 datatype,
3234 field_values,
3235 f);
3236
3237 } /* End of loop on sections */
3238
3239 } /* End for per element variable */
3240
3241 /* Free helper structures */
3242 /*------------------------*/
3243
3244 fvm_writer_field_helper_destroy(&helper);
3245
3246 BFT_FREE(export_list);
3247
3248 /* Close variable file and update case file */
3249 /*------------------------------------------*/
3250
3251 _free_ensight_file(&f);
3252
3253 fvm_to_ensight_case_write_case(w->case_info, rank);
3254 }
3255
3256 /*----------------------------------------------------------------------------*/
3257
3258 END_C_DECLS
3259