1
/*============================================================================
2
* Write a nodal representation associated with a mesh and associated
3
* variables to EnSight Gold files
4
*============================================================================*/
7
This file is part of Code_Saturne, a general-purpose CFD tool.
9
Copyright (C) 1998-2011 EDF S.A.
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
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
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.
26
/*----------------------------------------------------------------------------*/
28
#if defined(HAVE_CONFIG_H)
29
#include "cs_config.h"
32
/*----------------------------------------------------------------------------
33
* Standard C library headers
34
*----------------------------------------------------------------------------*/
42
/*----------------------------------------------------------------------------
44
*----------------------------------------------------------------------------*/
46
#include <bft_error.h>
49
#include <bft_printf.h>
51
/*----------------------------------------------------------------------------
53
*----------------------------------------------------------------------------*/
55
#include "fvm_config_defs.h"
57
#include "fvm_convert_array.h"
58
#include "fvm_part_to_block.h"
59
#include "fvm_block_to_part.h"
60
#include "fvm_gather.h"
61
#include "fvm_io_num.h"
62
#include "fvm_nodal.h"
63
#include "fvm_nodal_priv.h"
64
#include "fvm_parall.h"
65
#include "fvm_to_ensight_case.h"
66
#include "fvm_writer_helper.h"
67
#include "fvm_writer_priv.h"
70
/*----------------------------------------------------------------------------
71
* Header for the current file
72
*----------------------------------------------------------------------------*/
74
#include "fvm_to_ensight.h"
76
/*----------------------------------------------------------------------------*/
81
} /* Fake brace to force back Emacs auto-indentation back to column 0 */
83
#endif /* __cplusplus */
85
/*============================================================================
86
* Local Type Definitions
87
*============================================================================*/
89
/*----------------------------------------------------------------------------
90
* EnSight Gold writer structure
91
*----------------------------------------------------------------------------*/
95
char *name; /* Writer name */
97
int rank; /* Rank of current process in communicator */
98
int n_ranks; /* Number of processes in communicator */
100
_Bool text_mode; /* true if using text output */
101
_Bool swap_endian; /* true if binary file endianness must
104
_Bool discard_polygons; /* Option to discard polygonal elements */
105
_Bool discard_polyhedra; /* Option to discard polyhedral elements */
107
_Bool divide_polygons; /* Option to tesselate polygonal elements */
108
_Bool divide_polyhedra; /* Option to tesselate polyhedral elements */
110
fvm_to_ensight_case_t *case_info; /* Associated case structure */
112
#if defined(HAVE_MPI)
113
MPI_Comm comm; /* Associated MPI communicator */
116
} fvm_to_ensight_writer_t;
118
/*----------------------------------------------------------------------------
119
* Indirect file structure to handle both text and binary files
120
*----------------------------------------------------------------------------*/
124
bft_file_t *tf; /* Text file handing structure */
125
fvm_file_t *bf; /* Binary file handling structure */
129
/*============================================================================
130
* Static global variables
131
*============================================================================*/
133
static const char *_ensight_type_name[FVM_N_ELEMENT_TYPES] = {"bar2",
143
/*============================================================================
144
* Private function definitions
145
*============================================================================*/
147
/*----------------------------------------------------------------------------
148
* Open an EnSight Gold geometry or variable file
151
* this_writer <-- pointer to Ensight Gold writer structure.
152
* filename <-- name of file to open.
153
* apend <-- if true, append to file instead of overwriting
154
*----------------------------------------------------------------------------*/
156
static _ensight_file_t
157
_open_ensight_file(const fvm_to_ensight_writer_t *this_writer,
158
const char *filename,
161
_ensight_file_t f = {NULL, NULL};
163
if (this_writer->text_mode == true) {
164
bft_file_mode_t mode = append ? BFT_FILE_MODE_APPEND : BFT_FILE_MODE_WRITE;
165
if (this_writer->rank == 0)
166
f.tf = bft_file_open(filename, mode, BFT_FILE_TYPE_TEXT);
169
fvm_file_mode_t mode = append ? FVM_FILE_MODE_APPEND : FVM_FILE_MODE_WRITE;
171
#if defined(HAVE_MPI)
172
f.bf = fvm_file_open(filename, mode, hints, this_writer->comm);
174
f.bf = fvm_file_open(filename, mode, hints);
176
if (this_writer->swap_endian == true)
177
fvm_file_set_swap_endian(f.bf, 1);
183
/*----------------------------------------------------------------------------
184
* close an EnSight Gold geometry or variable file
187
* f <-- pointer to file handler structure.
188
*----------------------------------------------------------------------------*/
191
_free_ensight_file(_ensight_file_t *f)
194
f->tf = bft_file_free(f->tf);
196
else if (f->bf != NULL)
197
f->bf = fvm_file_free(f->bf);
200
/*----------------------------------------------------------------------------
201
* Count number of extra vertices when tesselations are present
204
* mesh <-- pointer to nodal mesh structure
205
* divide_polyhedra <-- true if polyhedra are tesselated
206
* n_extra_vertices_g --> global number of extra vertices (optional)
207
* n_extra_vertices --> local number of extra vertices (optional)
208
*----------------------------------------------------------------------------*/
211
_count_extra_vertices(const fvm_nodal_t *mesh,
212
_Bool divide_polyhedra,
213
fvm_gnum_t *n_extra_vertices_g,
214
fvm_lnum_t *n_extra_vertices)
218
const int export_dim = fvm_nodal_get_max_entity_dim(mesh);
220
/* Initial count and allocation */
222
if (n_extra_vertices_g != NULL)
223
*n_extra_vertices_g = 0;
224
if (n_extra_vertices != NULL)
225
*n_extra_vertices = 0;
227
if (divide_polyhedra) {
229
for (i = 0; i < mesh->n_sections; i++) {
231
const fvm_nodal_section_t *section = mesh->sections[i];
233
/* Output if entity dimension equal to highest in mesh
234
(i.e. no output of faces if cells present, or edges
235
if cells or faces) */
237
if ( section->entity_dim == export_dim
238
&& section->type == FVM_CELL_POLY
239
&& section->tesselation != NULL) {
241
if (n_extra_vertices_g != NULL)
243
+= fvm_tesselation_n_g_vertices_add(section->tesselation);
245
if (n_extra_vertices != NULL)
247
+= fvm_tesselation_n_vertices_add(section->tesselation);
254
/*----------------------------------------------------------------------------
255
* Write string to a text or C binary EnSight Gold file
258
* f <-- file to write to
259
* s <-- string to write
260
*----------------------------------------------------------------------------*/
263
_write_string(_ensight_file_t f,
272
bft_file_printf(f.tf, "%s\n", buf);
275
else if (f.bf != NULL) {
278
for (i = strlen(buf); i < 80; i++)
280
fvm_file_write_global(f.bf, buf, 1, 80);
284
/*----------------------------------------------------------------------------
285
* Write integer to a text or C binary EnSight Gold file
288
* f <-- file to write to
289
* n <-- integer value to write
290
*----------------------------------------------------------------------------*/
293
_write_int(_ensight_file_t f,
297
bft_file_printf(f.tf, "%10d\n", (int)n);
299
else if (f.bf != NULL) {
301
fvm_file_write_global(f.bf, &_n, sizeof(int32_t), 1);
305
/*----------------------------------------------------------------------------
306
* Write headers to an EnSight Gold geometry file.
309
* this_writer <-- pointer to Ensight Gold writer structure.
310
* f <-- pointer to file to initialize.
311
*----------------------------------------------------------------------------*/
314
_write_geom_headers(fvm_to_ensight_writer_t *this_writer,
318
_write_string(f, "C Binary");
320
/* 1st description line */
323
if (this_writer->name != NULL)
324
strncpy(buf, this_writer->name, 80);
326
_write_string(f, buf);
328
/* 2nd description line */
329
_write_string(f, "Output by Code_Saturne version "VERSION);
330
_write_string(f, "node id assign");
331
_write_string(f, "element id assign");
334
#if defined(HAVE_MPI)
336
/*----------------------------------------------------------------------------
337
* Write block of a vector of floats to an EnSight Gold file.
339
* This variant is called in parallel mode, where the values are already
340
* in a temporary block buffer and may be discarded.
343
* num_start <-- global number of first element for this block
344
* num_end <-- global number of past the last element for this block
345
* values <-- pointer to values block array
346
* comm <-- associated MPI communicator
347
* f <-- file to write to
348
*----------------------------------------------------------------------------*/
351
_write_block_floats_g(fvm_gnum_t num_start,
360
/* In Binary mode, all ranks have a file structure,
361
we may use use a collective call */
364
fvm_file_write_block_buffer(f.bf,
371
/* If all ranks do not have a binary file structure pointer, then
372
we are using a text file, open only on rank 0 */
376
float *_values = NULL;
377
fvm_file_serializer_t *s
378
= fvm_file_serializer_create(sizeof(float),
388
fvm_gnum_t range[2] = {num_start, num_end};
390
_values = fvm_file_serializer_advance(s, range);
392
if (_values != NULL) { /* only possible on rank 0 */
393
assert(f.tf != NULL);
394
for (i = 0, j = range[0]; j < range[1]; i++, j++)
395
bft_file_printf(f.tf, "%12.5e\n", _values[i]);
398
} while (_values != NULL);
400
fvm_file_serializer_destroy(&s);
405
#endif /* defined(HAVE_MPI) */
407
/*----------------------------------------------------------------------------
408
* Write block of a vector of floats to an EnSight Gold file
411
* n_values <-- number of values to write
412
* num_end <-- global number of past the last element for this block
413
* values <-- pointer to values block array
414
* f <-- file to write to
415
*----------------------------------------------------------------------------*/
418
_write_block_floats_l(size_t n_values,
419
const float values[],
425
for (i = 0; i < n_values; i++)
426
bft_file_printf(f.tf, "%12.5e\n", values[i]);
428
else if (f.bf != NULL)
429
fvm_file_write_global(f.bf, values, sizeof(float), n_values);
432
/*----------------------------------------------------------------------------
433
* Return extra vertex coordinates when tesselations are present
436
* mesh <-- pointer to nodal mesh structure
437
* n_extra_vertices <-- number of extra vertices
440
* array containing all extra vertex coordinates
441
*----------------------------------------------------------------------------*/
444
_extra_vertex_coords(const fvm_nodal_t *mesh,
445
fvm_lnum_t n_extra_vertices)
448
fvm_lnum_t n_extra_vertices_section;
450
size_t coord_shift = 0;
451
fvm_coord_t *coords = NULL;
453
if (n_extra_vertices > 0) { /* This implies divide_polyhedra = true */
455
BFT_MALLOC(coords, n_extra_vertices * 3, fvm_coord_t);
457
for (i = 0; i < mesh->n_sections; i++) {
459
const fvm_nodal_section_t *const section = mesh->sections[i];
461
if ( section->type == FVM_CELL_POLY
462
&& section->tesselation != NULL) {
464
n_extra_vertices_section
465
= fvm_tesselation_n_vertices_add(section->tesselation);
467
if (n_extra_vertices_section > 0) {
469
fvm_tesselation_vertex_coords(section->tesselation,
470
coords + coord_shift);
472
coord_shift += n_extra_vertices_section * 3;
483
#if defined(HAVE_MPI)
485
/*----------------------------------------------------------------------------
486
* Get extra vertex global numbers when tesselations are present
489
* mesh <-- pointer to nodal mesh structure
490
* n_extra_vertices <-- number of extra vertices
491
* vtx_gnum --> extra vertex global numbers (size: n_extra_vertices)
492
*----------------------------------------------------------------------------*/
495
_extra_vertex_get_gnum(const fvm_nodal_t *mesh,
496
fvm_lnum_t n_extra_vertices,
497
fvm_gnum_t vtx_gnum[])
501
fvm_lnum_t start_id = 0;
502
fvm_gnum_t gnum_shift
503
= fvm_io_num_get_global_count(mesh->global_vertex_num);
505
if (n_extra_vertices > 0) { /* Implies divide_polyhedra */
507
for (i = 0; i < mesh->n_sections; i++) {
509
const fvm_nodal_section_t *const section = mesh->sections[i];
511
if ( section->type == FVM_CELL_POLY
512
&& section->tesselation != NULL) {
514
fvm_lnum_t n_extra_vertices_section
515
= fvm_tesselation_n_vertices_add(section->tesselation);
517
if (n_extra_vertices_section > 0) {
519
const fvm_io_num_t *extra_vertex_num
520
= fvm_tesselation_global_vertex_num(section->tesselation);
521
const fvm_gnum_t *extra_gnum
522
= fvm_io_num_get_global_num(extra_vertex_num);
524
for (j = 0; j < n_extra_vertices_section; j++)
525
vtx_gnum[start_id + j] = extra_gnum[j] + gnum_shift;
527
start_id += n_extra_vertices_section;
532
= fvm_tesselation_n_g_vertices_add(section->tesselation);
539
/*----------------------------------------------------------------------------
540
* Build block info and part to block distribution helper for vertices.
543
* mesh <-- pointer to nodal mesh structure
544
* divide_polyhedra <-- true if polyhedra are tesselated
545
* comm <-- associated MPI communicator
546
* bi --> block information structure
547
* d --> part to bloc distributor
548
*----------------------------------------------------------------------------*/
551
_vertex_part_to_block_create(const fvm_nodal_t *mesh,
552
_Bool divide_polyhedra,
554
fvm_part_to_block_info_t *bi,
555
fvm_part_to_block_t **d)
559
fvm_gnum_t n_g_extra_vertices = 0, n_g_vertices_tot = 0;
560
fvm_lnum_t n_extra_vertices = 0, n_vertices_tot = 0;
562
fvm_gnum_t *_g_num = NULL;
564
fvm_part_to_block_info_t _bi;
565
fvm_part_to_block_t *_d;
567
size_t min_block_size
568
= fvm_parall_get_min_coll_buf_size() / sizeof(float);
570
const fvm_lnum_t n_vertices
571
= fvm_io_num_get_local_count(mesh->global_vertex_num);
572
fvm_gnum_t n_g_vertices
573
= fvm_io_num_get_global_count(mesh->global_vertex_num);
574
const fvm_gnum_t *g_num
575
= fvm_io_num_get_global_num(mesh->global_vertex_num);
577
/* Get info on the current MPI communicator */
579
MPI_Comm_rank(comm, &rank);
580
MPI_Comm_size(comm, &n_ranks);
582
/* Compute extra vertex coordinates if present */
584
_count_extra_vertices(mesh,
589
n_vertices_tot = n_vertices + n_extra_vertices;
590
n_g_vertices_tot = n_g_vertices + n_g_extra_vertices;
592
_bi = fvm_part_to_block_compute_sizes(rank,
598
/* Global vertex numbers */
601
if (n_extra_vertices > 0) {
603
BFT_MALLOC(_g_num, n_vertices_tot, fvm_gnum_t);
605
memcpy(_g_num, g_num, n_vertices*sizeof(fvm_gnum_t));
606
_extra_vertex_get_gnum(mesh, n_extra_vertices, _g_num + n_vertices);
612
/* Build distribution structures */
614
_d = fvm_part_to_block_create_by_gnum(comm, _bi, n_vertices_tot, g_num);
616
if (n_extra_vertices > 0)
617
fvm_part_to_block_transfer_gnum(_d, _g_num);
619
/* Return initialized structures */
628
/*----------------------------------------------------------------------------
629
* Write vertex coordinates to an EnSight Gold file in parallel mode
632
* this_writer <-- pointer to associated writer
633
* mesh <-- pointer to nodal mesh structure
634
* f <-- associated file handle
635
*----------------------------------------------------------------------------*/
638
_export_vertex_coords_g(const fvm_to_ensight_writer_t *this_writer,
639
const fvm_nodal_t *mesh,
644
fvm_part_to_block_info_t bi;
646
fvm_gnum_t n_g_extra_vertices = 0, n_g_vertices_tot = 0;
647
fvm_lnum_t n_extra_vertices = 0, n_vertices_tot = 0;
648
fvm_coord_t *extra_vertex_coords = NULL;
649
float *part_coords = NULL, *block_coords = NULL;
651
fvm_part_to_block_t *d = NULL;
652
size_t block_buf_size = 0;
654
const double *vertex_coords = mesh->vertex_coords;
655
const fvm_lnum_t *parent_vertex_num = mesh->parent_vertex_num;
656
const fvm_lnum_t n_vertices
657
= fvm_io_num_get_local_count(mesh->global_vertex_num);
658
fvm_gnum_t n_g_vertices
659
= fvm_io_num_get_global_count(mesh->global_vertex_num);
661
/* Initialize distribution info */
663
_vertex_part_to_block_create(mesh,
664
this_writer->divide_polyhedra,
669
/* Compute extra vertex coordinates if present */
671
_count_extra_vertices(mesh,
672
this_writer->divide_polyhedra,
676
n_vertices_tot = n_vertices + n_extra_vertices;
677
n_g_vertices_tot = n_g_vertices + n_g_extra_vertices;
679
extra_vertex_coords = _extra_vertex_coords(mesh,
684
block_buf_size = (bi.gnum_range[1] - bi.gnum_range[0]);
685
BFT_MALLOC(block_coords, block_buf_size, float);
686
BFT_MALLOC(part_coords, n_vertices_tot, float);
688
/* Vertex coordinates */
689
/*--------------------*/
691
stride = (size_t)(mesh->dim);
693
_write_string(f, "coordinates");
694
_write_int(f, (int)(n_g_vertices_tot));
696
/* Loop on dimension (de-interlace coordinates, always 3D for EnSight) */
698
for (j = 0; j < 3; j++) {
702
if (parent_vertex_num != NULL) {
703
for (i = 0; i < n_vertices; i++)
704
part_coords[i] = vertex_coords[(parent_vertex_num[i]-1)*stride + j];
707
for (i = 0; i < n_vertices; i++)
708
part_coords[i] = vertex_coords[i*stride + j];
711
for (i = 0; i < n_extra_vertices; i++)
712
part_coords[n_vertices + i] = extra_vertex_coords[(i*stride) + j];
715
for (i = 0; i < n_vertices_tot; i++)
716
part_coords[i] = 0.0;
719
fvm_part_to_block_copy_array(d,
725
_write_block_floats_g(bi.gnum_range[0],
731
} /* end of loop on spatial dimension */
733
fvm_part_to_block_destroy(&d);
735
BFT_FREE(block_coords);
736
BFT_FREE(part_coords);
737
if (extra_vertex_coords != NULL)
738
BFT_FREE(extra_vertex_coords);
741
#endif /* defined(HAVE_MPI) */
743
/*----------------------------------------------------------------------------
744
* Write vertex coordinates to an EnSight Gold file in serial mode
747
* this_writer <-- pointer to associated writer
748
* mesh <-- pointer to nodal mesh structure
749
* f <-- associated file handle
750
*----------------------------------------------------------------------------*/
753
_export_vertex_coords_l(const fvm_to_ensight_writer_t *this_writer,
754
const fvm_nodal_t *mesh,
758
fvm_lnum_t n_extra_vertices = 0;
759
fvm_coord_t *extra_vertex_coords = NULL;
760
float *coords_tmp = NULL;
762
const fvm_lnum_t n_vertices = mesh->n_vertices;
763
const double *vertex_coords = mesh->vertex_coords;
764
const fvm_lnum_t *parent_vertex_num = mesh->parent_vertex_num;
766
const size_t stride = (size_t)(mesh->dim);
768
/* Compute extra vertex coordinates if present */
770
_count_extra_vertices(mesh,
771
this_writer->divide_polyhedra,
775
extra_vertex_coords = _extra_vertex_coords(mesh,
778
/* Vertex coordinates */
779
/*--------------------*/
781
BFT_MALLOC(coords_tmp, FVM_MAX(n_vertices, n_extra_vertices), float);
783
_write_string(f, "coordinates");
784
_write_int(f, n_vertices + n_extra_vertices);
786
/* Loop on dimension (de-interlace coordinates, always 3D for EnSight) */
788
for (j = 0; j < 3; j++) {
790
/* First, handle regular vertices */
793
if (parent_vertex_num != NULL) {
794
for (i = 0; i < n_vertices; i++)
796
= (float)(vertex_coords[(parent_vertex_num[i]-1)*stride + j]);
799
for (i = 0; i < n_vertices; i++)
800
coords_tmp[i] = (float)(vertex_coords[i*stride + j]);
804
for (i = 0; i < (n_vertices); i++)
808
_write_block_floats_l(n_vertices,
812
/* Handle extra vertices (only occur with polyhedra tesselations in 3d) */
814
for (i = 0; i < n_extra_vertices; i++)
815
coords_tmp[i] = (float)(extra_vertex_coords[i*stride + j]);
817
if (n_extra_vertices > 0)
818
_write_block_floats_l(n_extra_vertices,
822
} /* end of loop on mesh dimension */
824
BFT_FREE(coords_tmp);
826
if (extra_vertex_coords != NULL)
827
BFT_FREE(extra_vertex_coords);
830
#if defined(HAVE_MPI)
832
/*----------------------------------------------------------------------------
833
* Write strided connectivity block to an EnSight Gold file in text mode
836
* stride <-- number of vertices per element type
837
* n_elems <-- number of elements
838
* connect <-- connectivity array
839
* tf <-- file to write to
840
*----------------------------------------------------------------------------*/
843
_write_connect_block_gt(int stride,
845
const int32_t connect[],
851
assert(bft_file_get_type(tf) == BFT_FILE_TYPE_TEXT);
856
for (i = 0; i < n_elems; i++)
857
bft_file_printf(tf, "%10d\n",
862
for (i = 0; i < n_elems; i++)
863
bft_file_printf(tf, "%10d%10d\n",
865
(int)connect[i*2+1]);
868
case 3: /* FVM_FACE_TRIA */
869
for (i = 0; i < n_elems; i++)
870
bft_file_printf(tf, "%10d%10d%10d\n",
873
(int)connect[i*3+2]);
876
case 4: /* FVM_FACE_QUAD or FVM_CELL_TETRA */
877
for (i = 0; i < n_elems; i++)
878
bft_file_printf(tf, "%10d%10d%10d%10d\n",
882
(int)connect[i*4+3]);
885
case 5: /* FVM_CELL_PYRAM */
886
for (i = 0; i < n_elems; i++)
887
bft_file_printf(tf, "%10d%10d%10d%10d%10d\n",
892
(int)connect[i*5+4]);
895
case 6: /* FVM_CELL_PRISM */
896
for (i = 0; i < n_elems; i++)
897
bft_file_printf(tf, "%10d%10d%10d%10d%10d%10d\n",
903
(int)connect[i*6+5]);
906
case 8: /* FVM_CELL_HEXA */
907
for (i = 0; i < n_elems; i++)
909
"%10d%10d%10d%10d%10d%10d%10d%10d\n",
917
(int)connect[i*8+7]);
925
/*----------------------------------------------------------------------------
926
* Write strided global connectivity block to an EnSight Gold file
929
* stride <-- number of vertices per element type
930
* num_start <-- global number of first element for this block
931
* num_end <-- global number of past last element for this block
932
* block_connect <-> global connectivity block array
933
* comm <-- associated MPI communicator
934
* f <-- associated file handle
935
*----------------------------------------------------------------------------*/
938
_write_block_connect_g(int stride,
939
fvm_gnum_t num_start,
941
int32_t block_connect[],
945
/* In Binary mode, all ranks have a file structure,
946
we may use use a collective call */
949
fvm_file_write_block(f.bf,
956
/* If all ranks do not have a binary file structure pointer, then
957
we are using a text file, open only on rank 0 */
961
int32_t *_block_connect = NULL;
963
fvm_file_serializer_t *s = fvm_file_serializer_create(sizeof(int32_t),
972
fvm_gnum_t range[2] = {num_start, num_end};
974
_block_connect = fvm_file_serializer_advance(s, range);
976
if (_block_connect != NULL) /* only possible on rank 0 */
977
_write_connect_block_gt(stride,
978
(range[1] - range[0]),
982
} while (_block_connect != NULL);
984
fvm_file_serializer_destroy(&s);
988
#endif /* defined(HAVE_MPI) */
990
/*----------------------------------------------------------------------------
991
* Write strided local connectivity to an EnSight Gold file
994
* stride <-- number of vertices per element type
995
* n_elems <-- number of elements
996
* connect <-- connectivity array
997
* f <-- file to write to
998
*----------------------------------------------------------------------------*/
1001
_write_connect_l(int stride,
1003
const fvm_lnum_t connect[],
1008
if (f.tf != NULL) { /* Text mode */
1013
for (i = 0; i < n_elems; i++)
1014
bft_file_printf(f.tf, "%10d%10d\n",
1016
(int)connect[i*2+1]);
1019
case 3: /* FVM_FACE_TRIA */
1020
for (i = 0; i < n_elems; i++)
1021
bft_file_printf(f.tf, "%10d%10d%10d\n",
1023
(int)connect[i*3+1],
1024
(int)connect[i*3+2]);
1027
case 4: /* FVM_FACE_QUAD or FVM_CELL_TETRA */
1028
for (i = 0; i < n_elems; i++)
1029
bft_file_printf(f.tf, "%10d%10d%10d%10d\n",
1031
(int)connect[i*4+1],
1032
(int)connect[i*4+2],
1033
(int)connect[i*4+3]);
1036
case 5: /* FVM_CELL_PYRAM */
1037
for (i = 0; i < n_elems; i++)
1038
bft_file_printf(f.tf, "%10d%10d%10d%10d%10d\n",
1040
(int)connect[i*5+1],
1041
(int)connect[i*5+2],
1042
(int)connect[i*5+3],
1043
(int)connect[i*5+4]);
1046
case 6: /* FVM_CELL_PRISM */
1047
for (i = 0; i < n_elems; i++)
1048
bft_file_printf(f.tf, "%10d%10d%10d%10d%10d%10d\n",
1050
(int)connect[i*6+1],
1051
(int)connect[i*6+2],
1052
(int)connect[i*6+3],
1053
(int)connect[i*6+4],
1054
(int)connect[i*6+5]);
1057
case 8: /* FVM_CELL_HEXA */
1058
for (i = 0; i < n_elems; i++)
1059
bft_file_printf(f.tf,
1060
"%10d%10d%10d%10d%10d%10d%10d%10d\n",
1062
(int)connect[i*8+1],
1063
(int)connect[i*8+2],
1064
(int)connect[i*8+3],
1065
(int)connect[i*8+4],
1066
(int)connect[i*8+5],
1067
(int)connect[i*8+6],
1068
(int)connect[i*8+7]);
1076
else if (f.bf != NULL) { /* Binary mode */
1080
int32_t *buffer = NULL;
1081
const size_t n_values = n_elems*stride;
1082
const size_t buffer_size = n_values > 64 ? (n_values / 8) : n_values;
1084
BFT_MALLOC(buffer, buffer_size, int32_t);
1086
while (k < n_values) {
1087
for (j = 0; j < buffer_size && k < n_values; j++)
1088
buffer[j] = (int)(connect[k++]);
1089
fvm_file_write_global(f.bf, buffer, sizeof(int32_t), j);
1097
#if defined(HAVE_MPI)
1099
/*----------------------------------------------------------------------------
1100
* Write "trivial" point elements to an EnSight Gold file in parallel mode
1103
* mesh <-- pointer to nodal mesh structure
1104
* comm <-- associated MPI communicator
1105
* f <-- file to write to
1106
*----------------------------------------------------------------------------*/
1109
_export_point_elements_g(const fvm_nodal_t *mesh,
1113
const fvm_gnum_t n_g_vertices
1114
= fvm_io_num_get_global_count(mesh->global_vertex_num);
1116
_write_string(f, "point");
1117
_write_int(f, (int)n_g_vertices);
1119
if (f.tf != NULL) { /* Text mode, rank 0 only */
1124
for (i = 0; i < n_g_vertices; i++)
1125
bft_file_printf(f.tf, "%10d\n", j++);
1128
else if (f.bf != NULL) { /* Binary mode */
1133
fvm_part_to_block_info_t bi;
1135
size_t min_block_size
1136
= fvm_parall_get_min_coll_buf_size() / sizeof(float);
1137
int32_t *connect = NULL;
1139
/* Get info on the current MPI communicator */
1141
MPI_Comm_rank(comm, &rank);
1142
MPI_Comm_size(comm, &n_ranks);
1144
bi = fvm_part_to_block_compute_sizes(rank,
1150
BFT_MALLOC(connect, bi.gnum_range[1] - bi.gnum_range[0], int32_t);
1152
for (i = 0, j = bi.gnum_range[0]; j < bi.gnum_range[1]; i++, j++)
1155
_write_block_connect_g(1,
1166
#endif /* defined(HAVE_MPI) */
1168
/*----------------------------------------------------------------------------
1169
* Write "trivial" point elements to an EnSight Gold file in serial mode
1172
* mesh <-- pointer to nodal mesh structure
1173
* f <-- file to write to
1174
*----------------------------------------------------------------------------*/
1177
_export_point_elements_l(const fvm_nodal_t *mesh,
1180
const fvm_lnum_t n_vertices = mesh->n_vertices;
1182
_write_string(f, "point");
1183
_write_int(f, (int)n_vertices);
1185
if (n_vertices == 0)
1188
if (f.tf != NULL) { /* Text mode */
1190
for (i = 0; i < n_vertices; i++)
1191
bft_file_printf(f.tf, "%10d\n", i+1);
1194
else if (f.bf != NULL) { /* Binary mode */
1198
int32_t *buf = NULL;
1199
const int32_t bufsize = n_vertices > 64 ? (n_vertices / 8) : n_vertices;
1201
BFT_MALLOC(buf, bufsize, int32_t);
1203
j_end = n_vertices + 1;
1205
for (k = 0; j < j_end && k < bufsize; k++)
1207
fvm_file_write_global(f.bf, buf, sizeof(int32_t), k);
1214
#if defined(HAVE_MPI)
1216
/*----------------------------------------------------------------------------
1217
* Write indexed element lengths from a nodal mesh to an EnSight Gold
1218
* file in parallel mode
1221
* global_element_num <-- global element numbering
1222
* vertex_index <-- pointer to element -> vertex index
1223
* comm <-- associated MPI communicator
1224
* n_ranks <-- number of processes in communicator
1225
* f <-- associated file handle
1226
*----------------------------------------------------------------------------*/
1229
_write_lengths_g(const fvm_io_num_t *global_element_num,
1230
const fvm_lnum_t vertex_index[],
1236
fvm_part_to_block_info_t bi;
1238
int32_t *part_lengths = NULL;
1239
int32_t *block_lengths = NULL;
1241
fvm_part_to_block_t *d = NULL;
1243
const size_t min_block_size
1244
= fvm_parall_get_min_coll_buf_size() / sizeof(int32_t);
1245
const fvm_lnum_t n_elements
1246
= fvm_io_num_get_local_count(global_element_num);
1247
const fvm_lnum_t n_g_elements
1248
= fvm_io_num_get_global_count(global_element_num);
1249
const fvm_gnum_t *g_num
1250
= fvm_io_num_get_global_num(global_element_num);
1252
/* Get info on the current MPI communicator */
1254
MPI_Comm_rank(comm, &rank);
1255
MPI_Comm_size(comm, &n_ranks);
1257
/* Allocate block buffer */
1259
bi = fvm_part_to_block_compute_sizes(rank,
1265
/* Build distribution structures */
1267
BFT_MALLOC(block_lengths, bi.gnum_range[1] - bi.gnum_range[0], int);
1268
BFT_MALLOC(part_lengths, n_elements, int32_t);
1270
for (i = 0; i < n_elements; i++)
1271
part_lengths[i] = vertex_index[i+1] - vertex_index[i];
1273
d = fvm_part_to_block_create_by_gnum(comm, bi, n_elements, g_num);
1275
fvm_part_to_block_copy_array(d,
1281
fvm_part_to_block_destroy(&d);
1282
BFT_FREE(part_lengths);
1286
_write_block_connect_g(1,
1293
BFT_FREE(block_lengths);
1296
/*----------------------------------------------------------------------------
1297
* Write block-distributed indexed element (polygons or polyhedra)
1298
* cell -> vertex connectivity to an EnSight Gold file in parallel mode.
1300
* In text mode, zeroes may be used in place of extra vertex numbers
1301
* to indicate extra newlines between face -> vertex definitions.
1304
* num_start <-- global number of first element for this block
1305
* num_end <-- global number of past last element for this block
1306
* block_index <-- global connectivity block array
1307
* block_connect <-> global connectivity block array
1308
* comm <-- associated MPI communicator
1309
* f <-- associated file handle
1310
*----------------------------------------------------------------------------*/
1313
_write_block_indexed(fvm_gnum_t num_start,
1315
const fvm_lnum_t block_index[],
1316
int32_t block_connect[],
1320
fvm_gnum_t block_size = 0, block_start = 0, block_end = 0;
1322
/* Prepare write to file */
1324
block_size = block_index[num_end - num_start];
1326
MPI_Scan(&block_size, &block_end, 1, FVM_MPI_GNUM, MPI_SUM, comm);
1328
block_start = block_end - block_size;
1330
/* In Binary mode, all ranks have a file structure,
1331
we may use use a collective call */
1334
fvm_file_write_block(f.bf,
1341
/* If all ranks do not have a binary file structure pointer, then
1342
we are using a text file, open only on rank 0 */
1346
int32_t *_block_vtx_num = NULL;
1347
fvm_file_serializer_t *s = fvm_file_serializer_create(sizeof(int32_t),
1357
fvm_gnum_t range[2] = {block_start, block_end};
1358
_block_vtx_num = fvm_file_serializer_advance(s, range);
1359
if (_block_vtx_num != NULL) { /* only possible on rank 0 */
1360
assert(f.tf != NULL);
1361
for (i = 0, j = range[0]; j < range[1]; i++, j++) {
1362
if (_block_vtx_num[i] != 0)
1363
bft_file_printf(f.tf, "%10d", _block_vtx_num[i]);
1365
bft_file_printf(f.tf, "\n");
1368
} while (_block_vtx_num != NULL);
1370
fvm_file_serializer_destroy(&s);
1374
/*----------------------------------------------------------------------------
1375
* Write indexed element (polygons or polyhedra) cell -> vertex connectivity
1376
* to an EnSight Gold file in parallel mode.
1378
* In text mode, zeroes may be used in place of extra vertex numbers
1379
* to indicate extra newlines between face -> vertex definitions.
1382
* global_vertex_num <-- vertex global numbering
1383
* global_element_num <-- global element numbering
1384
* vertex_index <-- element -> vertex index
1385
* vertex_num <-- element -> vertex number
1386
* comm <-- associated MPI communicator
1387
* f <-- associated file handle
1388
*----------------------------------------------------------------------------*/
1391
_write_indexed_connect_g(const fvm_io_num_t *global_element_num,
1392
const fvm_lnum_t vertex_index[],
1393
const int32_t vertex_num[],
1398
fvm_part_to_block_info_t bi;
1400
fvm_gnum_t loc_size = 0, tot_size = 0, block_size = 0;
1401
fvm_part_to_block_t *d = NULL;
1402
fvm_lnum_t *block_index = NULL;
1403
int32_t *block_vtx_num = NULL;
1404
size_t min_block_size
1405
= fvm_parall_get_min_coll_buf_size() / sizeof(int32_t);
1407
const fvm_gnum_t n_g_elements
1408
= fvm_io_num_get_global_count(global_element_num);
1409
const fvm_lnum_t n_elements
1410
= fvm_io_num_get_local_count(global_element_num);
1411
const fvm_gnum_t *g_elt_num
1412
= fvm_io_num_get_global_num(global_element_num);
1414
/* Get info on the current MPI communicator */
1416
MPI_Comm_rank(comm, &rank);
1417
MPI_Comm_size(comm, &n_ranks);
1419
/* Adjust min block size based on minimum element size */
1421
loc_size = vertex_index[n_elements];
1422
MPI_Allreduce(&loc_size, &tot_size, 1, FVM_MPI_GNUM, MPI_SUM, comm);
1424
min_block_size /= (tot_size / n_g_elements);
1426
/* Allocate memory for additionnal indexes */
1428
bi = fvm_part_to_block_compute_sizes(rank,
1434
BFT_MALLOC(block_index, bi.gnum_range[1] - bi.gnum_range[0] + 1, fvm_lnum_t);
1436
d = fvm_part_to_block_create_by_gnum(comm, bi, n_elements, g_elt_num);
1438
fvm_part_to_block_copy_index(d,
1442
block_size = block_index[bi.gnum_range[1] - bi.gnum_range[0]];
1444
BFT_MALLOC(block_vtx_num, block_size, int32_t);
1446
fvm_part_to_block_copy_indexed(d,
1455
_write_block_indexed(bi.gnum_range[0],
1464
BFT_FREE(block_vtx_num);
1465
fvm_part_to_block_destroy(&d);
1466
BFT_FREE(block_index);
1469
/*----------------------------------------------------------------------------
1470
* Write polyhedra from a nodal mesh to an EnSight Gold file in parallel mode
1473
* export_section <-- pointer to EnSight section helper structure
1474
* global_vertex_num <-- pointer to vertex global numbering
1475
* comm <-- associated MPI communicator
1476
* f <-- associated file handle
1479
* pointer to next EnSight section helper structure in list
1480
*----------------------------------------------------------------------------*/
1482
static const fvm_writer_section_t *
1483
_export_nodal_polyhedra_g(const fvm_writer_section_t *export_section,
1484
const fvm_io_num_t *global_vertex_num,
1489
fvm_lnum_t i, j, k, l, face_id;
1491
fvm_lnum_t face_length, cell_length;
1492
fvm_part_to_block_info_t bi;
1494
fvm_part_to_block_t *d = NULL;
1495
const fvm_writer_section_t *current_section;
1497
/* Get info on the current MPI communicator */
1499
MPI_Comm_rank(comm, &rank);
1500
MPI_Comm_size(comm, &n_ranks);
1502
/* Export number of faces per polyhedron */
1503
/*---------------------------------------*/
1505
current_section = export_section;
1507
do { /* loop on sections which should be appended */
1509
const fvm_nodal_section_t *section = current_section->section;
1511
_write_lengths_g(section->global_element_num,
1512
section->face_index,
1516
current_section = current_section->next;
1518
} while ( current_section != NULL
1519
&& current_section->continues_previous == true);
1521
/* Export number of vertices per face per polyhedron */
1522
/*---------------------------------------------------*/
1524
current_section = export_section;
1526
do { /* loop on sections which should be appended */
1528
fvm_gnum_t block_size = 0, block_start = 0, block_end = 0;
1529
fvm_lnum_t *block_index = NULL;
1531
size_t min_block_size
1532
= fvm_parall_get_min_coll_buf_size() / sizeof(int32_t);
1533
int32_t *part_face_len = NULL, *block_face_len = NULL;
1535
const fvm_nodal_section_t *section = current_section->section;
1536
const fvm_lnum_t n_elements
1537
= fvm_io_num_get_local_count(section->global_element_num);
1538
const fvm_gnum_t n_g_elements
1539
= fvm_io_num_get_global_count(section->global_element_num);
1540
const fvm_gnum_t *g_elt_num
1541
= fvm_io_num_get_global_num(section->global_element_num);
1543
/* Build local polyhedron face lengths information */
1545
BFT_MALLOC(part_face_len,
1546
section->face_index[section->n_elements],
1550
for (i = 0; i < section->n_elements; i++) {
1551
for (j = section->face_index[i]; j < section->face_index[i+1]; j++) {
1552
face_id = FVM_ABS(section->face_num[j]) - 1;
1553
face_length = ( section->vertex_index[face_id+1]
1554
- section->vertex_index[face_id]);
1555
part_face_len[k++] = face_length;
1558
assert(k == section->face_index[section->n_elements]);
1560
/* Prepare distribution structures */
1562
bi = fvm_part_to_block_compute_sizes(rank,
1568
d = fvm_part_to_block_create_by_gnum(comm,
1573
BFT_MALLOC(block_index, bi.gnum_range[1] - bi.gnum_range[0] + 1, fvm_lnum_t);
1575
fvm_part_to_block_copy_index(d,
1576
section->face_index,
1579
block_size = block_index[bi.gnum_range[1] - bi.gnum_range[0]];
1581
BFT_MALLOC(block_face_len, block_size, int32_t);
1583
fvm_part_to_block_copy_indexed(d,
1585
section->face_index,
1590
MPI_Scan(&block_size, &block_end, 1, FVM_MPI_GNUM, MPI_SUM, comm);
1592
block_start = block_end - block_size;
1594
_write_block_connect_g(1,
1601
BFT_FREE(block_face_len);
1603
fvm_part_to_block_destroy(&d);
1605
BFT_FREE(block_index);
1606
BFT_FREE(part_face_len);
1608
current_section = current_section->next;
1610
} while ( current_section != NULL
1611
&& current_section->continues_previous == true);
1613
/* Export cell->vertex connectivity by blocks */
1614
/*--------------------------------------------*/
1616
current_section = export_section;
1618
do { /* loop on sections which should be appended */
1620
fvm_lnum_t *part_vtx_idx = NULL;
1621
int32_t *part_vtx_num = NULL;
1623
const fvm_nodal_section_t *section = current_section->section;
1624
const fvm_gnum_t *g_vtx_num
1625
= fvm_io_num_get_global_num(global_vertex_num);
1627
BFT_MALLOC(part_vtx_idx, section->n_elements + 1, fvm_lnum_t);
1631
if (f.bf != NULL) { /* In binary mode, build cell -> vertex connectivity */
1633
part_vtx_idx[0] = 0;
1634
for (i = 0; i < section->n_elements; i++) {
1636
for (j = section->face_index[i]; j < section->face_index[i+1]; j++) {
1637
face_id = FVM_ABS(section->face_num[j]) - 1;
1638
face_length = ( section->vertex_index[face_id+1]
1639
- section->vertex_index[face_id]);
1640
cell_length += face_length;
1642
part_vtx_idx[i+1] = part_vtx_idx[i] + cell_length;
1647
/* In text mode, add zeroes to cell vertex connectivity to mark face
1648
bounds (so as to add newlines) */
1650
else { /* we are in text mode if f.bf == NULL */
1652
part_vtx_idx[0] = 0;
1653
for (i = 0; i < section->n_elements; i++) {
1655
for (j = section->face_index[i]; j < section->face_index[i+1]; j++) {
1656
face_id = FVM_ABS(section->face_num[j]) - 1;
1657
face_length = ( section->vertex_index[face_id+1]
1658
- section->vertex_index[face_id]);
1659
cell_length += face_length + 1;
1661
part_vtx_idx[i+1] = part_vtx_idx[i] + cell_length;
1666
BFT_MALLOC(part_vtx_num, part_vtx_idx[section->n_elements], fvm_lnum_t);
1670
for (i = 0; i < section->n_elements; i++) {
1671
for (j = section->face_index[i]; j < section->face_index[i+1]; j++) {
1672
if (section->face_num[j] > 0) {
1673
face_id = section->face_num[j] - 1;
1674
for (k = section->vertex_index[face_id];
1675
k < section->vertex_index[face_id+1];
1677
part_vtx_num[l++] = g_vtx_num[section->vertex_num[k] - 1];
1680
face_id = -section->face_num[j] - 1;
1681
k = section->vertex_index[face_id];
1682
part_vtx_num[l++] = g_vtx_num[section->vertex_num[k] - 1];
1683
for (k = section->vertex_index[face_id+1] - 1;
1684
k > section->vertex_index[face_id];
1686
part_vtx_num[l++] = g_vtx_num[section->vertex_num[k] - 1];
1689
part_vtx_num[l++] = 0; /* mark face limits in text mode */
1694
/* Now distribute and write cells -> vertices connectivity */
1696
_write_indexed_connect_g(section->global_element_num,
1702
BFT_FREE(part_vtx_num);
1703
BFT_FREE(part_vtx_idx);
1705
current_section = current_section->next;
1707
} while ( current_section != NULL
1708
&& current_section->continues_previous == true);
1710
return current_section;
1713
#endif /* defined(HAVE_MPI) */
1715
/*----------------------------------------------------------------------------
1716
* Write polyhedra from a nodal mesh to an EnSight Gold file in serial mode
1719
* export_section <-- pointer to EnSight section helper structure
1720
* f <-- associated file handle
1723
* pointer to next EnSight section helper structure in list
1724
*----------------------------------------------------------------------------*/
1726
static const fvm_writer_section_t *
1727
_export_nodal_polyhedra_l(const fvm_writer_section_t *export_section,
1731
fvm_lnum_t i, j, k, l;
1734
fvm_lnum_t face_length, face_id;
1736
size_t buffer_size = 0;
1737
int32_t *buffer = NULL;
1739
const fvm_writer_section_t *current_section;
1741
/* Write number of faces per cell */
1742
/*--------------------------------*/
1744
current_section = export_section;
1746
do { /* loop on sections which should be appended */
1748
const fvm_nodal_section_t *section = current_section->section;
1750
if (f.tf != NULL) { /* Text mode */
1751
for (i = 0; i < section->n_elements; i++)
1752
bft_file_printf(f.tf, "%10d\n",
1753
(int)( section->face_index[i+1]
1754
- section->face_index[i]));
1756
else { /* binary mode */
1758
/* First, allocate a buffer large enough so that the number of
1759
writes is limited, small enough so that the memory overhead is
1760
minimal; polyhedral connectivity is at least 4 faces x 3 vertices
1761
per cell, usually quite a bit more, so this is 1/3 of the minimum */
1763
if (buffer_size < (size_t)section->n_elements * 4) {
1764
buffer_size = section->n_elements * 4;
1765
BFT_REALLOC(buffer, buffer_size, int32_t);
1768
/* Now fill buffer and write */
1770
for (i = 0, i_buf = 0; i < section->n_elements; i++) {
1771
if (i_buf == buffer_size) {
1772
fvm_file_write_global(f.bf, buffer, sizeof(int32_t), i_buf);
1775
buffer[i_buf++] = (int)( section->face_index[i+1]
1776
- section->face_index[i]);
1779
fvm_file_write_global(f.bf, buffer, sizeof(int32_t), i_buf);
1783
current_section = current_section->next;
1785
} while ( current_section != NULL
1786
&& current_section->continues_previous == true);
1788
/* Write number of vertices/face */
1789
/*-------------------------------*/
1791
current_section = export_section;
1793
do { /* loop on sections which should be appended */
1795
const fvm_nodal_section_t *section = current_section->section;
1797
for (i = 0, i_buf = 0; i < section->n_elements; i++) {
1799
/* Loop on cell faces */
1801
for (j = section->face_index[i];
1802
j < section->face_index[i+1];
1805
if (section->face_num[j] > 0)
1806
face_id = section->face_num[j] - 1;
1808
face_id = -section->face_num[j] - 1;
1810
face_length = ( section->vertex_index[face_id+1]
1811
- section->vertex_index[face_id]);
1814
bft_file_printf(f.tf, "%10d\n",
1817
if (i_buf == buffer_size) {
1818
fvm_file_write_global(f.bf, buffer, sizeof(int32_t), i_buf);
1821
buffer[i_buf++] = (int)face_length;
1828
if (f.bf != NULL && i_buf > 0)
1829
fvm_file_write_global(f.bf, buffer, sizeof(int32_t), i_buf);
1831
current_section = current_section->next;
1833
} while ( current_section != NULL
1834
&& current_section->continues_previous == true);
1836
/* Write cell/vertex connectivity */
1837
/*--------------------------------*/
1839
current_section = export_section;
1841
do { /* loop on sections which should be appended */
1843
const fvm_nodal_section_t *section = current_section->section;
1845
for (i = 0, i_buf = 0; i < section->n_elements; i++) {
1847
/* Loop on cell faces */
1849
for (j = section->face_index[i];
1850
j < section->face_index[i+1];
1853
/* Print face vertex numbers */
1855
if (section->face_num[j] > 0) {
1856
face_id = section->face_num[j] - 1;
1860
face_id = -section->face_num[j] - 1;
1864
face_length = ( section->vertex_index[face_id+1]
1865
- section->vertex_index[face_id]);
1867
if (f.tf != NULL) { /* text mode */
1868
for (k = 0; k < face_length; k++) {
1869
l = section->vertex_index[face_id]
1870
+ (face_length + (k*face_sgn))%face_length;
1871
bft_file_printf(f.tf, "%10d", (int)section->vertex_num[l]);
1873
bft_file_printf(f.tf, "\n");
1875
else { /* binary mode */
1876
for (k = 0; k < face_length; k++) {
1877
l = section->vertex_index[face_id]
1878
+ (face_length + (k*face_sgn))%face_length;
1879
if (i_buf == buffer_size) {
1880
fvm_file_write_global(f.bf, buffer, sizeof(int32_t), i_buf);
1883
buffer[i_buf++] = (int)section->vertex_num[l];
1887
} /* End of loop on cell faces */
1889
} /* End of loop on polyhedral cells */
1891
if (f.bf != NULL && i_buf > 0)
1892
fvm_file_write_global(f.bf, buffer, sizeof(int32_t), i_buf);
1894
current_section = current_section->next;
1896
} while ( current_section != NULL
1897
&& current_section->continues_previous == true);
1902
return current_section;
1905
#if defined(HAVE_MPI)
1907
/*----------------------------------------------------------------------------
1908
* Write polygons from a nodal mesh to an EnSight Gold file in parallel mode
1911
* export_section <-- pointer to EnSight section helper structure
1912
* global_vertex_num <-- pointer to vertex global numbering
1913
* comm <-- associated MPI communicator
1914
* f <-- associated file handle
1917
* pointer to next EnSight section helper structure in list
1918
*----------------------------------------------------------------------------*/
1920
static const fvm_writer_section_t *
1921
_export_nodal_polygons_g(const fvm_writer_section_t *export_section,
1922
const fvm_io_num_t *global_vertex_num,
1926
const fvm_writer_section_t *current_section;
1928
/* Export number of vertices per polygon by blocks */
1929
/*-------------------------------------------------*/
1931
current_section = export_section;
1933
do { /* loop on sections which should be appended */
1935
const fvm_nodal_section_t *section = current_section->section;
1937
_write_lengths_g(section->global_element_num,
1938
section->vertex_index,
1942
current_section = current_section->next;
1944
} while ( current_section != NULL
1945
&& current_section->continues_previous == true);
1947
/* Export face->vertex connectivity by blocks */
1948
/*--------------------------------------------*/
1950
current_section = export_section;
1952
do { /* loop on sections which should be appended */
1955
fvm_lnum_t *_part_vtx_idx = NULL;
1956
const fvm_lnum_t *part_vtx_idx = NULL;
1957
int32_t *part_vtx_num = NULL;
1959
const fvm_nodal_section_t *section = current_section->section;
1960
const fvm_gnum_t *g_vtx_num
1961
= fvm_io_num_get_global_num(global_vertex_num);
1963
if (f.bf != NULL) /* In binary mode, use existing index */
1964
part_vtx_idx = section->vertex_index;
1966
/* In text mode, add zeroes to cell vertex connectivity to mark face
1967
bounds (so as to add newlines) */
1969
else { /* we are in text mode if f.bf == NULL */
1971
BFT_MALLOC(_part_vtx_idx, section->n_elements + 1, fvm_lnum_t);
1973
_part_vtx_idx[0] = 0;
1974
for (i = 0; i < section->n_elements; i++)
1975
_part_vtx_idx[i+1] = _part_vtx_idx[i] + ( section->vertex_index[i+1]
1976
- section->vertex_index[i]) + 1;
1978
part_vtx_idx = _part_vtx_idx;
1981
/* Build connectivity array */
1983
BFT_MALLOC(part_vtx_num, part_vtx_idx[section->n_elements], int32_t);
1985
if (f.bf != NULL) { /* binary mode */
1986
for (i = 0, k = 0; i < section->n_elements; i++) {
1987
for (j = section->vertex_index[i];
1988
j < section->vertex_index[i+1];
1990
part_vtx_num[k++] = g_vtx_num[section->vertex_num[j] - 1];
1994
else { /* text mode */
1995
for (i = 0, k = 0; i < section->n_elements; i++) {
1996
for (j = section->vertex_index[i];
1997
j < section->vertex_index[i+1];
1999
part_vtx_num[k++] = g_vtx_num[section->vertex_num[j] - 1];
2000
part_vtx_num[k++] = 0; /* mark face bounds in text mode */
2004
/* Now distribute and write cell -> vertices connectivity */
2006
_write_indexed_connect_g(section->global_element_num,
2012
BFT_FREE(part_vtx_num);
2013
if (_part_vtx_idx != NULL)
2014
BFT_FREE(_part_vtx_idx);
2016
current_section = current_section->next;
2018
} while ( current_section != NULL
2019
&& current_section->continues_previous == true);
2021
return current_section;
2024
#endif /* defined(HAVE_MPI) */
2026
/*----------------------------------------------------------------------------
2027
* Write polygons from a nodal mesh to a text file in serial mode
2030
* export_section <-- pointer to EnSight section helper structure
2031
* f <-- associated file handle
2034
* pointer to next EnSight section helper structure in list
2035
*----------------------------------------------------------------------------*/
2037
static const fvm_writer_section_t *
2038
_export_nodal_polygons_l(const fvm_writer_section_t *export_section,
2046
size_t buffer_size = 0;
2047
int32_t *buffer = NULL;
2049
const fvm_writer_section_t *current_section = NULL;
2051
/* Print face connectivity directly, without using extra buffers */
2053
/* First loop on all polygonal faces, to write number of vertices */
2054
/*----------------------------------------------------------------*/
2056
current_section = export_section;
2058
do { /* Loop on sections that should be grouped */
2060
const fvm_nodal_section_t *section = current_section->section;
2062
if (f.tf != NULL) { /* Text mode */
2063
for (i = 0; i < section->n_elements; i++)
2064
bft_file_printf(f.tf, "%10d\n", (int)( section->vertex_index[i+1]
2065
- section->vertex_index[i]));
2067
else { /* binary mode */
2069
/* First, allocate a buffer large enough so that the number of
2070
writes is limited, small enough so that the memory overhead is
2071
minimal; polygonal connectivity is at least 3 vertices per face,
2072
usually 5 or more, so this is 1/3 of the minimum */
2074
if (buffer_size < (size_t)section->n_elements) {
2075
buffer_size = section->n_elements;
2076
BFT_REALLOC(buffer, buffer_size, int32_t);
2079
/* Now fill buffer and write */
2081
for (i = 0, i_buf = 0; i < section->n_elements; i++) {
2082
if (i_buf == buffer_size) {
2083
fvm_file_write_global(f.bf, buffer, sizeof(int32_t), i_buf);
2086
buffer[i_buf++] = (int)( section->vertex_index[i+1]
2087
- section->vertex_index[i]);
2090
fvm_file_write_global(f.bf, buffer, sizeof(int32_t), i_buf);
2093
current_section = current_section->next;
2095
} while ( current_section != NULL
2096
&& current_section->continues_previous == true);
2098
/* Loop on all polygonal faces */
2099
/*-----------------------------*/
2101
current_section = export_section;
2103
do { /* Loop on sections that should be grouped */
2105
const fvm_nodal_section_t *section = current_section->section;
2107
for (i = 0, i_buf = 0; i < section->n_elements; i++) {
2109
/* Print face vertex numbers */
2111
if (f.tf != NULL) { /* text mode */
2112
for (j = section->vertex_index[i];
2113
j < section->vertex_index[i+1];
2115
bft_file_printf(f.tf, "%10d", (int)section->vertex_num[j]);
2116
bft_file_printf(f.tf, "\n");
2118
else { /* binary mode */
2119
for (j = section->vertex_index[i];
2120
j < section->vertex_index[i+1];
2122
if (i_buf == buffer_size) {
2123
fvm_file_write_global(f.bf, buffer, sizeof(int32_t), i_buf);
2126
buffer[i_buf++] = (int)section->vertex_num[j];
2130
} /* End of loop on polygonal faces */
2132
if (f.bf != NULL && i_buf > 0)
2133
fvm_file_write_global(f.bf, buffer, sizeof(int32_t), i_buf);
2135
current_section = current_section->next;
2137
} while ( current_section != NULL
2138
&& current_section->continues_previous == true);
2143
return current_section;
2146
#if defined(HAVE_MPI)
2148
/*----------------------------------------------------------------------------
2149
* Write tesselated element cell -> vertex connectivity to an EnSight Gold
2150
* file in parallel mode.
2153
* global_vertex_num <-- vertex global numbering
2154
* global_element_num <-- global element numbering
2155
* tesselation <-- element tesselation description
2156
* type <-- tesselated sub-element type
2157
* extra_vertex_base <-- starting number for added vertices
2158
* comm <-- associated MPI communicator
2159
* f <-- associated file handle
2160
*----------------------------------------------------------------------------*/
2163
_write_tesselated_connect_g(const fvm_io_num_t *global_vertex_num,
2164
const fvm_io_num_t *global_element_num,
2165
const fvm_tesselation_t *tesselation,
2167
const fvm_gnum_t extra_vertex_base,
2173
fvm_part_to_block_info_t bi;
2175
fvm_lnum_t part_size = 0;
2177
fvm_lnum_t end_id = 0;
2178
fvm_gnum_t n_g_sub_elements = 0, global_num_end = 0;
2179
fvm_gnum_t block_size = 0, block_start = 0, block_end = 0;
2181
fvm_part_to_block_t *d = NULL;
2182
fvm_lnum_t *part_index, *block_index = NULL;
2183
int32_t *part_vtx_num = NULL, *block_vtx_num = NULL;
2184
fvm_gnum_t *part_vtx_gnum = NULL;
2186
size_t min_block_size
2187
= fvm_parall_get_min_coll_buf_size() / sizeof(int32_t);
2189
const int stride = fvm_nodal_n_vertices_element[type];
2190
const fvm_lnum_t n_elements = fvm_tesselation_n_elements(tesselation);
2191
const fvm_gnum_t n_g_elements
2192
= fvm_io_num_get_global_count(global_element_num);
2193
const fvm_lnum_t n_sub_elements
2194
= fvm_tesselation_n_sub_elements(tesselation, type);
2195
const fvm_lnum_t *sub_element_idx
2196
= fvm_tesselation_sub_elt_index(tesselation, type);
2197
const fvm_gnum_t *g_elt_num
2198
= fvm_io_num_get_global_num(global_element_num);
2200
/* Get info on the current MPI communicator */
2202
MPI_Comm_rank(comm, &rank);
2203
MPI_Comm_size(comm, &n_ranks);
2205
/* Adjust min block size based on mean number of sub-elements */
2207
fvm_tesselation_get_global_size(tesselation,
2212
min_block_size /= (n_g_sub_elements/n_g_elements) * stride;
2214
/* Decode connectivity */
2216
part_size = n_sub_elements * stride;
2217
assert(sub_element_idx[n_elements]*stride == part_size);
2219
global_num_end = n_g_elements + 1;
2221
if (n_elements > 0) {
2222
BFT_MALLOC(part_vtx_num, part_size, int32_t);
2223
BFT_MALLOC(part_vtx_gnum, part_size, fvm_gnum_t);
2226
end_id = fvm_tesselation_decode_g(tesselation,
2236
assert(end_id == n_elements);
2237
assert(global_num_end == n_g_elements + 1);
2239
/* Convert to write type */
2241
if (n_elements > 0) {
2242
for (i = 0; i < part_size; i++)
2243
part_vtx_num[i] = part_vtx_gnum[i];
2244
BFT_FREE(part_vtx_gnum);
2247
/* Allocate memory for additionnal indexes and decoded connectivity */
2249
bi = fvm_part_to_block_compute_sizes(rank,
2255
BFT_MALLOC(block_index, bi.gnum_range[1] - bi.gnum_range[0] + 1, fvm_lnum_t);
2256
BFT_MALLOC(part_index, n_elements + 1, fvm_lnum_t);
2258
d = fvm_part_to_block_create_by_gnum(comm, bi, n_elements, g_elt_num);
2261
for (i = 0; i < n_elements; i++) {
2262
part_index[i+1] = part_index[i] + ( sub_element_idx[i+1]
2263
- sub_element_idx[i]) * stride;
2268
fvm_part_to_block_copy_index(d,
2272
block_size = (block_index[bi.gnum_range[1] - bi.gnum_range[0]]);
2274
/* Copy connectivity */
2276
BFT_MALLOC(block_vtx_num, block_size, int32_t);
2278
fvm_part_to_block_copy_indexed(d,
2285
fvm_part_to_block_destroy(&d);
2287
BFT_FREE(part_vtx_num);
2288
BFT_FREE(part_index);
2289
BFT_FREE(block_index);
2293
block_size /= stride;
2295
MPI_Scan(&block_size, &block_end, 1, FVM_MPI_GNUM, MPI_SUM, comm);
2297
block_start = block_end - block_size;
2299
_write_block_connect_g(stride,
2306
/* Free remaining memory */
2308
BFT_FREE(block_vtx_num);
2311
/*----------------------------------------------------------------------------
2312
* Write tesselated element connectivity from a nodal mesh to an EnSight Gold
2313
* file in parallel mode
2316
* export_section <-- pointer to EnSight section helper structure
2317
* global_vertex_num <-- pointer to vertex global numbering
2318
* comm <-- associated MPI communicator
2319
* f <-- associated file handle
2322
* pointer to next EnSight section helper structure in list
2323
*----------------------------------------------------------------------------*/
2325
static const fvm_writer_section_t *
2326
_export_nodal_tesselated_g(const fvm_writer_section_t *export_section,
2327
const fvm_io_num_t *global_vertex_num,
2331
const fvm_writer_section_t *current_section;
2333
/* Export face->vertex connectivity by blocks */
2334
/*--------------------------------------------*/
2336
current_section = export_section;
2338
do { /* loop on sections which should be appended */
2340
const fvm_nodal_section_t *section = current_section->section;
2342
_write_tesselated_connect_g(global_vertex_num,
2343
section->global_element_num,
2344
section->tesselation,
2345
current_section->type,
2346
current_section->extra_vertex_base,
2350
current_section = current_section->next;
2352
} while ( current_section != NULL
2353
&& current_section->continues_previous == true);
2355
return current_section;
2358
#endif /* defined(HAVE_MPI) */
2360
/*----------------------------------------------------------------------------
2361
* Write tesselated element connectivity from a nodal mesh to an EnSight Gold
2362
* file in parallel mode
2365
* export_section <-- pointer to EnSight section helper structure
2366
* f <-- associated file handle
2369
* pointer to next EnSight section helper structure in list
2370
*----------------------------------------------------------------------------*/
2372
static const fvm_writer_section_t *
2373
_export_nodal_tesselated_l(const fvm_writer_section_t *export_section,
2376
const fvm_writer_section_t *current_section;
2378
current_section = export_section;
2380
do { /* loop on sections which should be appended */
2382
const fvm_nodal_section_t *section = current_section->section;
2384
fvm_lnum_t start_id, end_id;
2385
fvm_lnum_t n_sub_elements_max;
2386
fvm_lnum_t n_buffer_elements_max = section->n_elements;
2387
fvm_lnum_t *vertex_num = NULL;
2389
const fvm_lnum_t *sub_element_idx
2390
= fvm_tesselation_sub_elt_index(section->tesselation,
2391
export_section->type);
2393
fvm_tesselation_get_global_size(section->tesselation,
2394
export_section->type,
2396
&n_sub_elements_max);
2397
if (n_sub_elements_max > n_buffer_elements_max)
2398
n_buffer_elements_max = n_sub_elements_max;
2400
BFT_MALLOC(vertex_num,
2401
( n_buffer_elements_max
2402
* fvm_nodal_n_vertices_element[export_section->type]),
2406
start_id < section->n_elements;
2407
start_id = end_id) {
2410
= fvm_tesselation_decode(section->tesselation,
2411
current_section->type,
2413
n_buffer_elements_max,
2414
export_section->extra_vertex_base,
2417
_write_connect_l(fvm_nodal_n_vertices_element[export_section->type],
2418
( sub_element_idx[end_id]
2419
- sub_element_idx[start_id]),
2425
BFT_FREE(vertex_num);
2427
current_section = current_section->next;
2429
} while ( current_section != NULL
2430
&& current_section->continues_previous == true);
2432
return current_section;
2435
#if defined(HAVE_MPI)
2437
/*----------------------------------------------------------------------------
2438
* Write strided elements from a nodal mesh to an EnSight Gold file in
2442
* export_section <-- pointer to EnSight section helper structure
2443
* global_vertex_num <-- pointer to vertex global numbering
2444
* comm <-- associated MPI communicator
2445
* f <-- associated file handle
2448
* pointer to next EnSight section helper structure in list
2449
*----------------------------------------------------------------------------*/
2451
static const fvm_writer_section_t *
2452
_export_nodal_strided_g(const fvm_writer_section_t *export_section,
2453
const fvm_io_num_t *global_vertex_num,
2460
const fvm_writer_section_t *current_section;
2462
/* Get info on the current MPI communicator */
2464
MPI_Comm_rank(comm, &rank);
2465
MPI_Comm_size(comm, &n_ranks);
2467
/* Export vertex connectivity */
2469
current_section = export_section;
2471
do { /* loop on sections which should be appended */
2473
fvm_part_to_block_info_t bi;
2475
fvm_lnum_t block_size = 0;
2476
fvm_part_to_block_t *d = NULL;
2477
int32_t *part_vtx_num = NULL, *block_vtx_num = NULL;
2479
const fvm_nodal_section_t *section = current_section->section;
2480
const int stride = fvm_nodal_n_vertices_element[section->type];
2482
const size_t min_block_size
2483
= fvm_parall_get_min_coll_buf_size() / (sizeof(int32_t) * stride);
2485
const fvm_lnum_t n_elements
2486
= fvm_io_num_get_local_count(section->global_element_num);
2487
const fvm_gnum_t n_g_elements
2488
= fvm_io_num_get_global_count(section->global_element_num);
2489
const fvm_gnum_t *g_elt_num
2490
= fvm_io_num_get_global_num(section->global_element_num);
2491
const fvm_gnum_t *g_vtx_num
2492
= fvm_io_num_get_global_num(global_vertex_num);
2494
/* Prepare distribution structures */
2496
bi = fvm_part_to_block_compute_sizes(rank,
2502
d = fvm_part_to_block_create_by_gnum(comm,
2507
/* Build connectivity */
2509
block_size = bi.gnum_range[1] - bi.gnum_range[0];
2511
BFT_MALLOC(block_vtx_num, block_size*stride, int32_t);
2512
BFT_MALLOC(part_vtx_num, n_elements*stride, int32_t);
2514
for (i = 0; i < n_elements; i++) {
2515
for (j = 0; j < stride; j++) {
2516
part_vtx_num[i*stride + j]
2517
= g_vtx_num[section->vertex_num[i*stride + j] - 1];
2521
fvm_part_to_block_copy_array(d,
2527
BFT_FREE(part_vtx_num);
2529
_write_block_connect_g(stride,
2536
BFT_FREE(block_vtx_num);
2538
fvm_part_to_block_destroy(&d);
2540
current_section = current_section->next;
2542
} while ( current_section != NULL
2543
&& current_section->continues_previous == true);
2545
return current_section;
2548
/*----------------------------------------------------------------------------
2549
* Write field values associated with nodal values of a nodal mesh to
2550
* an EnSight Gold file in serial mode.
2552
* Output fields ar either scalar or 3d vectors or scalars, and are
2553
* non interlaced. Input arrays may be less than 2d, in which case the z
2554
* values are set to 0, and may be interlaced or not.
2557
* mesh <-- pointer to nodal mesh structure
2558
* divide_polyhedra <-- true if polyhedra are tesselated
2559
* input_dim <-- input field dimension
2560
* output_dim <-- output field dimension
2561
* interlace <-- indicates if field in memory is interlaced
2562
* n_parent_lists <-- indicates if field values are to be obtained
2563
* directly through the local entity index (when 0) or
2564
* through the parent entity numbers (when 1 or more)
2565
* parent_num_shift <-- parent list to common number index shifts;
2566
* size: n_parent_lists
2567
* datatype <-- input data type (output is real)
2568
* field_values <-- array of associated field value arrays
2569
* comm <-- associated MPI communicator
2570
* f <-- associated file handle
2571
*----------------------------------------------------------------------------*/
2574
_export_field_values_ng(const fvm_nodal_t *mesh,
2575
_Bool divide_polyhedra,
2578
fvm_interlace_t interlace,
2580
const fvm_lnum_t parent_num_shift[],
2581
fvm_datatype_t datatype,
2582
const void *const field_values[],
2587
fvm_part_to_block_info_t bi;
2589
fvm_lnum_t part_size = 0, block_size = 0;
2590
float *part_values = NULL, *block_values = NULL;
2591
fvm_part_to_block_t *d = NULL;
2593
/* Initialize distribution info */
2595
_vertex_part_to_block_create(mesh,
2601
part_size = fvm_part_to_block_get_n_part_ents(d);
2602
block_size = bi.gnum_range[1] - bi.gnum_range[0];
2604
BFT_MALLOC(part_values, part_size, float);
2605
BFT_MALLOC(block_values, block_size, float);
2607
for (i = 0; i < output_dim; i++) {
2609
fvm_lnum_t start_id = 0;
2610
fvm_lnum_t end_id = mesh->n_vertices;
2612
/* Distribute partition to block values */
2614
if (i < input_dim) {
2620
fvm_convert_array(input_dim,
2630
mesh->parent_vertex_num,
2634
/* Additional vertices in case of tesselation
2635
(end_id == part_size with no tesselation or if all tesselated
2636
sections have been accounted for).*/
2638
for (j = 0; end_id < part_size && j < mesh->n_sections; j++) {
2640
const fvm_nodal_section_t *section = mesh->sections[j];
2642
assert(divide_polyhedra == true);
2644
if (section->type == FVM_CELL_POLY && section->tesselation != NULL) {
2646
fvm_lnum_t n_extra_vertices
2647
= fvm_tesselation_n_vertices_add(section->tesselation);
2650
end_id = start_id + n_extra_vertices;
2652
fvm_tesselation_vertex_values(section->tesselation,
2663
mesh->parent_vertex_num,
2665
part_values + start_id);
2669
} /* End of loops on tesselated sections */
2671
assert(end_id == part_size);
2673
fvm_part_to_block_copy_array(d,
2681
/* Zero extra dimensions */
2685
for (j = 0; j < block_size; j++)
2686
block_values[j] = 0.0;
2689
_write_block_floats_g(bi.gnum_range[0],
2696
BFT_FREE(block_values);
2697
BFT_FREE(part_values);
2699
fvm_part_to_block_destroy(&d);
2702
#endif /* defined(HAVE_MPI) */
2704
/*----------------------------------------------------------------------------
2705
* Write field values associated with nodal values of a nodal mesh to
2706
* an EnSight Gold file in serial mode.
2708
* Output fields ar either scalar or 3d vectors or scalars, and are
2709
* non interlaced. Input arrays may be less than 2d, in which case the z
2710
* values are set to 0, and may be interlaced or not.
2713
* n_entities <-- number of entities
2714
* input_dim <-- input field dimension
2715
* output_dim <-- output field dimension
2716
* interlace <-- indicates if field in memory is interlaced
2717
* n_parent_lists <-- indicates if field values are to be obtained
2718
* directly through the local entity index (when 0) or
2719
* through the parent entity numbers (when 1 or more)
2720
* parent_num_shift <-- parent list to common number index shifts;
2721
* size: n_parent_lists
2722
* datatype <-- input data type (output is real)
2723
* field_values <-- array of associated field value arrays
2724
* f <-- associated file handle
2725
*----------------------------------------------------------------------------*/
2728
_export_field_values_nl(const fvm_nodal_t *mesh,
2729
fvm_writer_field_helper_t *helper,
2731
fvm_interlace_t interlace,
2733
const fvm_lnum_t parent_num_shift[],
2734
fvm_datatype_t datatype,
2735
const void *const field_values[],
2740
float *output_buffer;
2742
int output_dim = fvm_writer_field_helper_field_dim(helper);
2744
const size_t output_buffer_size
2745
= mesh->n_vertices > 16 ? (mesh->n_vertices / 4) : mesh->n_vertices;
2747
BFT_MALLOC(output_buffer, output_buffer_size, float);
2749
for (i = 0; i < output_dim; i++) {
2751
while (fvm_writer_field_helper_step_n(helper,
2762
&output_size) == 0) {
2764
_write_block_floats_l(output_size,
2771
BFT_FREE(output_buffer);
2774
#if defined(HAVE_MPI)
2776
/*----------------------------------------------------------------------------
2777
* Write field values associated with element values of a nodal mesh to
2778
* an EnSight Gold file.
2780
* Output fields ar either scalar or 3d vectors or scalars, and are
2781
* non interlaced. Input arrays may be less than 2d, in which case the z
2782
* values are set to 0, and may be interlaced or not.
2785
* export_section <-- pointer to EnSight section helper structure
2786
* helper <-- pointer to general writer helper structure
2787
* input_dim <-- input field dimension
2788
* output_dim <-- output field dimension
2789
* interlace <-- indicates if field in memory is interlaced
2790
* n_parent_lists <-- indicates if field values are to be obtained
2791
* directly through the local entity index (when 0) or
2792
* through the parent entity numbers (when 1 or more)
2793
* parent_num_shift <-- parent list to common number index shifts;
2794
* size: n_parent_lists
2795
* datatype <-- indicates the data type of (source) field values
2796
* field_values <-- array of associated field value arrays
2797
* comm <-- associated MPI communicator
2798
* f <-- associated file handle
2801
* pointer to next EnSight section helper structure in list
2802
*----------------------------------------------------------------------------*/
2804
static const fvm_writer_section_t *
2805
_export_field_values_eg(const fvm_writer_section_t *export_section,
2808
fvm_interlace_t interlace,
2810
const fvm_lnum_t parent_num_shift[],
2811
fvm_datatype_t datatype,
2812
const void *const field_values[],
2820
fvm_part_to_block_info_t bi;
2821
fvm_part_to_block_t *d = NULL;
2824
_Bool have_tesselation = false;
2825
fvm_lnum_t part_size = 0, block_size = 0;
2826
fvm_gnum_t block_sub_size = 0, block_start = 0, block_end = 0;
2827
fvm_gnum_t n_g_elements = 0;
2829
int *part_n_sub = NULL, *block_n_sub = NULL;
2830
float *part_values = NULL, *block_values = NULL, *_block_values = NULL;
2832
fvm_gnum_t *_g_elt_num = NULL;
2833
const fvm_gnum_t *g_elt_num
2834
= fvm_io_num_get_global_num(export_section->section->global_element_num);
2836
const fvm_writer_section_t *current_section = NULL;
2838
size_t min_block_size
2839
= fvm_parall_get_min_coll_buf_size() / sizeof(int32_t);
2841
MPI_Comm_rank(comm, &rank);
2842
MPI_Comm_size(comm, &n_ranks);
2844
/* Loop on sections to count output size */
2846
current_section = export_section;
2849
const fvm_nodal_section_t *section = current_section->section;
2852
n_g_elements += fvm_io_num_get_global_count(section->global_element_num);
2853
part_size += fvm_io_num_get_local_count(section->global_element_num);
2854
if (current_section->type != section->type)
2855
have_tesselation = true;
2857
current_section = current_section->next;
2859
} while ( current_section != NULL
2860
&& current_section->continues_previous == true);
2862
/* Build global numbering if necessary */
2864
if (n_sections > 1) {
2866
fvm_lnum_t start_id = 0;
2868
BFT_MALLOC(_g_elt_num, part_size, fvm_gnum_t);
2869
g_elt_num = _g_elt_num;
2871
/* loop on sections which should be appended */
2873
current_section = export_section;
2876
const fvm_nodal_section_t *section = current_section->section;
2877
const fvm_lnum_t section_size
2878
= fvm_io_num_get_local_count(section->global_element_num);
2880
memcpy(_g_elt_num + start_id,
2881
fvm_io_num_get_global_num(section->global_element_num),
2882
sizeof(fvm_gnum_t)*section_size);
2883
start_id += section_size;
2885
current_section = current_section->next;
2887
} while ( current_section != NULL
2888
&& current_section->continues_previous == true);
2891
/* Build sub-element count if necessary */
2893
if (have_tesselation) {
2895
fvm_lnum_t start_id = 0;
2897
BFT_MALLOC(part_n_sub, part_size, int);
2899
current_section = export_section;
2902
const fvm_nodal_section_t *section = current_section->section;
2903
const fvm_lnum_t section_size
2904
= fvm_io_num_get_local_count(section->global_element_num);
2906
if (current_section->type != section->type) {
2907
const fvm_lnum_t *sub_element_idx
2908
= fvm_tesselation_sub_elt_index(section->tesselation,
2909
current_section->type);
2910
for (j = 0; j < section_size; j++)
2911
part_n_sub[start_id + j] = sub_element_idx[j+1] - sub_element_idx[j];
2914
for (j = 0; j < section_size; j++)
2915
part_n_sub[start_id + j] = 1;
2917
start_id += section_size;
2919
current_section = current_section->next;
2921
} while ( current_section != NULL
2922
&& current_section->continues_previous == true);
2925
/* Build distribution structures */
2927
bi = fvm_part_to_block_compute_sizes(rank,
2933
block_size = bi.gnum_range[1] - bi.gnum_range[0];
2935
d = fvm_part_to_block_create_by_gnum(comm, bi, part_size, g_elt_num);
2937
if (_g_elt_num != NULL)
2938
fvm_part_to_block_transfer_gnum(d, _g_elt_num);
2943
/* Distribute sub-element info in case of tesselation */
2945
if (have_tesselation) {
2947
BFT_MALLOC(block_n_sub, block_size, int);
2948
fvm_part_to_block_copy_array(d,
2953
BFT_FREE(part_n_sub);
2955
for (j = 0; j < block_size; j++)
2956
block_sub_size += block_n_sub[j];
2960
block_sub_size = block_size;
2962
/* To save space, in case of tesselation, part_values and _block_n_sub
2963
point to the same memory space, as they are not needed simultaneously.
2964
Without tesselation, _block_n_sub simply points to block_n_sub */
2966
BFT_MALLOC(part_values,
2967
FVM_MAX(part_size, (fvm_lnum_t)block_sub_size),
2969
BFT_MALLOC(block_values, block_size, float);
2971
if (have_tesselation) {
2972
MPI_Scan(&block_sub_size, &block_end, 1, FVM_MPI_GNUM, MPI_SUM, comm);
2974
block_start = block_end - block_sub_size;
2975
_block_values = part_values;
2978
block_start = bi.gnum_range[0];
2979
block_end = bi.gnum_range[1];
2980
_block_values = block_values;
2983
/* Loop on dimension (de-interlace vectors, always 3D for EnSight) */
2985
for (i = 0; i < output_dim; i++) {
2987
/* Distribute partition to block values */
2989
if (i < input_dim) {
2991
fvm_lnum_t start_id = 0;
2992
fvm_lnum_t src_shift = 0;
2994
/* loop on sections which should be appended */
2996
current_section = export_section;
2999
const fvm_nodal_section_t *section = current_section->section;
3001
if (n_parent_lists == 0)
3002
src_shift = export_section->num_shift;
3004
fvm_convert_array(input_dim,
3008
section->n_elements + src_shift,
3014
section->parent_element_num,
3016
part_values + start_id);
3018
start_id += fvm_io_num_get_local_count(section->global_element_num);
3020
current_section = current_section->next;
3022
} while ( current_section != NULL
3023
&& current_section->continues_previous == true);
3025
/* Distribute part values */
3027
fvm_part_to_block_copy_array(d,
3033
/* Scatter values to sub-elements in case of tesselation */
3035
if (have_tesselation) {
3037
for (j = 0; j < block_size; j++) {
3038
for (k = 0; k < block_n_sub[j]; k++)
3039
_block_values[l++] = block_values[j];
3045
/* Zero extra dimensions */
3048
for (j = 0; j < (fvm_lnum_t)block_sub_size; j++)
3049
block_values[j] = 0.0;
3052
/* Write block values */
3054
_write_block_floats_g(block_start,
3060
} /* end of loop on spatial dimension */
3062
BFT_FREE(block_values);
3063
BFT_FREE(part_values);
3065
fvm_part_to_block_destroy(&d);
3067
if (block_n_sub != NULL)
3068
BFT_FREE(block_n_sub);
3070
return current_section;
3073
#endif /* defined(HAVE_MPI) */
3075
/*----------------------------------------------------------------------------
3076
* Write field values associated with element values of a nodal mesh to
3077
* an EnSight Gold file.
3079
* Output fields ar either scalar or 3d vectors or scalars, and are
3080
* non interlaced. Input arrays may be less than 2d, in which case the z
3081
* values are set to 0, and may be interlaced or not.
3084
* export_section <-- pointer to EnSight section helper structure
3085
* helper <-- pointer to general writer helper structure
3086
* input_dim <-- input field dimension
3087
* interlace <-- indicates if field in memory is interlaced
3088
* n_parent_lists <-- indicates if field values are to be obtained
3089
* directly through the local entity index (when 0) or
3090
* through the parent entity numbers (when 1 or more)
3091
* parent_num_shift <-- parent list to common number index shifts;
3092
* size: n_parent_lists
3093
* datatype <-- indicates the data type of (source) field values
3094
* field_values <-- array of associated field value arrays
3095
* f <-- associated file handle
3098
* pointer to next EnSight section helper structure in list
3099
*----------------------------------------------------------------------------*/
3101
static const fvm_writer_section_t *
3102
_export_field_values_el(const fvm_writer_section_t *export_section,
3103
fvm_writer_field_helper_t *helper,
3105
fvm_interlace_t interlace,
3107
const fvm_lnum_t parent_num_shift[],
3108
fvm_datatype_t datatype,
3109
const void *const field_values[],
3113
size_t input_size = 0, output_size = 0;
3114
size_t min_output_buffer_size = 0, output_buffer_size = 0;
3115
float *output_buffer = NULL;
3117
const fvm_writer_section_t *current_section = NULL;
3119
int output_dim = fvm_writer_field_helper_field_dim(helper);
3121
/* Blocking for arbitrary buffer size, but should be small enough
3122
to add little additional memory requirement (in proportion), large
3123
enough to limit number of write and gather calls. */
3125
fvm_writer_field_helper_get_size(helper,
3129
&min_output_buffer_size);
3131
output_buffer_size = input_size / 4;
3132
output_buffer_size = FVM_MAX(output_buffer_size, min_output_buffer_size);
3133
output_buffer_size = FVM_MAX(output_buffer_size, 128);
3134
output_buffer_size = FVM_MIN(output_buffer_size, output_size);
3136
BFT_MALLOC(output_buffer, output_buffer_size, float);
3138
/* Loop on dimension (de-interlace vectors, always 3D for EnSight) */
3140
for (i = 0; i < output_dim; i++) {
3142
_Bool loop_on_sections = true;
3144
current_section = export_section;
3146
while (loop_on_sections == true) {
3148
while (fvm_writer_field_helper_step_e(helper,
3159
&output_size) == 0) {
3161
_write_block_floats_l(output_size,
3167
current_section = current_section->next;
3169
if ( current_section == NULL
3170
|| current_section->continues_previous == false)
3171
loop_on_sections = false;
3173
} /* while (loop on sections) */
3175
} /* end of loop on spatial dimension */
3177
BFT_FREE(output_buffer);
3179
return current_section;
3182
/*============================================================================
3183
* Public function definitions
3184
*============================================================================*/
3186
/*----------------------------------------------------------------------------
3187
* Initialize FVM to EnSight Gold file writer.
3190
* text output text files
3191
* binary output binary files (default)
3192
* big_endian force binary files to big-endian
3193
* discard_polygons do not output polygons or related values
3194
* discard_polyhedra do not output polyhedra or related values
3195
* divide_polygons tesselate polygons with triangles
3196
* divide_polyhedra tesselate polyhedra with tetrahedra and pyramids
3197
* (adding a vertex near each polyhedron's center)
3200
* name <-- base output case name.
3201
* options <-- whitespace separated, lowercase options list
3202
* time_dependecy <-- indicates if and how meshes will change with time
3203
* comm <-- associated MPI communicator.
3206
* pointer to opaque EnSight Gold writer structure.
3207
*----------------------------------------------------------------------------*/
3209
#if defined(HAVE_MPI)
3211
fvm_to_ensight_init_writer(const char *name,
3213
const char *options,
3214
fvm_writer_time_dep_t time_dependency,
3218
fvm_to_ensight_init_writer(const char *name,
3220
const char *options,
3221
fvm_writer_time_dep_t time_dependency)
3224
fvm_to_ensight_writer_t *this_writer = NULL;
3226
/* Initialize writer */
3228
BFT_MALLOC(this_writer, 1, fvm_to_ensight_writer_t);
3230
BFT_MALLOC(this_writer->name, strlen(name) + 1, char);
3231
strcpy(this_writer->name, name);
3233
this_writer->text_mode = false;
3234
this_writer->swap_endian = false;
3235
this_writer->discard_polygons = false;
3236
this_writer->discard_polyhedra = false;
3237
this_writer->divide_polygons = false;
3238
this_writer->divide_polyhedra = false;
3240
this_writer->rank = 0;
3241
this_writer->n_ranks = 1;
3243
#if defined(HAVE_MPI)
3245
int mpi_flag, rank, n_ranks;
3246
MPI_Initialized(&mpi_flag);
3248
if (mpi_flag && comm != MPI_COMM_NULL) {
3249
this_writer->comm = comm;
3250
MPI_Comm_rank(this_writer->comm, &rank);
3251
MPI_Comm_size(this_writer->comm, &n_ranks);
3252
this_writer->rank = rank;
3253
this_writer->n_ranks = n_ranks;
3256
this_writer->comm = MPI_COMM_NULL;
3258
#endif /* defined(HAVE_MPI) */
3262
if (options != NULL) {
3265
int l_tot = strlen(options);
3268
while (i1 < l_tot) {
3270
for (i2 = i1; i2 < l_tot && options[i2] != ' '; i2++);
3273
if ((l_opt == 4) && (strncmp(options + i1, "text", l_opt) == 0))
3274
this_writer->text_mode = true;
3275
else if ((l_opt == 6) && (strncmp(options + i1, "binary", l_opt) == 0))
3276
this_writer->text_mode = false;
3278
else if ( (l_opt == 10)
3279
&& (strncmp(options + i1, "big_endian", l_opt) == 0)) {
3281
this_writer->text_mode = false;
3282
/* Check if system is "big-endian" or "little-endian" */
3283
*((char *)(&int_endian)) = '\1';
3284
if (int_endian == 1)
3285
this_writer->swap_endian = 1;
3288
else if ( (l_opt == 16)
3289
&& (strncmp(options + i1, "discard_polygons", l_opt) == 0))
3290
this_writer->discard_polygons = true;
3291
else if ( (l_opt == 17)
3292
&& (strncmp(options + i1, "discard_polyhedra", l_opt) == 0))
3293
this_writer->discard_polyhedra = true;
3295
else if ( (l_opt == 15)
3296
&& (strncmp(options + i1, "divide_polygons", l_opt) == 0))
3297
this_writer->divide_polygons = true;
3298
else if ( (l_opt == 16)
3299
&& (strncmp(options + i1, "divide_polyhedra", l_opt) == 0))
3300
this_writer->divide_polyhedra = true;
3302
for (i1 = i2 + 1; i1 < l_tot && options[i1] == ' '; i1++);
3308
this_writer->case_info = fvm_to_ensight_case_create(name,
3317
/*----------------------------------------------------------------------------
3318
* Finalize FVM to EnSight Gold file writer.
3321
* this_writer_p <-- pointer to opaque Ensight Gold writer structure.
3325
*----------------------------------------------------------------------------*/
3328
fvm_to_ensight_finalize_writer(void *this_writer_p)
3330
fvm_to_ensight_writer_t *this_writer
3331
= (fvm_to_ensight_writer_t *)this_writer_p;
3333
BFT_FREE(this_writer->name);
3335
fvm_to_ensight_case_destroy(this_writer->case_info);
3337
BFT_FREE(this_writer);
3342
/*----------------------------------------------------------------------------
3343
* Associate new time step with an EnSight geometry.
3346
* this_writer_p <-- pointer to associated writer
3347
* time_step <-- time step number
3348
* time_value <-- time_value number
3349
*----------------------------------------------------------------------------*/
3352
fvm_to_ensight_set_mesh_time(void *this_writer_p,
3353
const int time_step,
3354
const double time_value)
3356
fvm_to_ensight_writer_t *this_writer
3357
= (fvm_to_ensight_writer_t *)this_writer_p;
3359
fvm_to_ensight_case_set_geom_time(this_writer->case_info,
3364
/*----------------------------------------------------------------------------
3365
* Indicate if a elements of a given type in a mesh associated to a given
3366
* EnSight Gold file writer need to be tesselated.
3369
* this_writer_p <-- pointer to associated writer
3370
* mesh <-- pointer to nodal mesh structure that should be written
3371
* element_type <-- element type we are interested in
3374
* 1 if tesselation of the given element type is needed, 0 otherwise
3375
*----------------------------------------------------------------------------*/
3378
fvm_to_ensight_needs_tesselation(fvm_writer_t *this_writer_p,
3379
const fvm_nodal_t *mesh,
3380
fvm_element_t element_type)
3384
fvm_to_ensight_writer_t *this_writer
3385
= (fvm_to_ensight_writer_t *)this_writer_p;
3387
const int export_dim = fvm_nodal_get_max_entity_dim(mesh);
3389
/* Initial count and allocation */
3391
if ( ( element_type == FVM_FACE_POLY
3392
&& this_writer->divide_polygons == true)
3393
|| ( element_type == FVM_CELL_POLY
3394
&& this_writer->divide_polyhedra == true)) {
3396
for (i = 0; i < mesh->n_sections; i++) {
3398
const fvm_nodal_section_t *const section = mesh->sections[i];
3400
/* Output if entity dimension equal to highest in mesh
3401
(i.e. no output of faces if cells present, or edges
3402
if cells or faces) */
3404
if (section->entity_dim == export_dim) {
3405
if (section->type == element_type)
3416
/*----------------------------------------------------------------------------
3417
* Write nodal mesh to a an EnSight Gold file
3420
* this_writer_p <-- pointer to associated writer
3421
* mesh <-- pointer to nodal mesh structure that should be written
3422
*----------------------------------------------------------------------------*/
3425
fvm_to_ensight_export_nodal(void *this_writer_p,
3426
const fvm_nodal_t *mesh)
3430
const fvm_writer_section_t *export_section = NULL;
3431
fvm_writer_section_t *export_list = NULL;
3432
fvm_to_ensight_writer_t *this_writer
3433
= (fvm_to_ensight_writer_t *)this_writer_p;
3434
_ensight_file_t f = {NULL, NULL};
3436
const int rank = this_writer->rank;
3437
const int n_ranks = this_writer->n_ranks;
3439
/* Initialization */
3440
/*----------------*/
3442
fvm_to_ensight_case_file_info_t file_info;
3444
/* Get part number */
3446
part_num = fvm_to_ensight_case_get_part_num(this_writer->case_info,
3449
part_num = fvm_to_ensight_case_add_part(this_writer->case_info,
3452
/* Open geometry file in append mode */
3454
file_info = fvm_to_ensight_case_get_geom_file(this_writer->case_info);
3456
f = _open_ensight_file(this_writer,
3460
if (file_info.queried == false)
3461
_write_geom_headers(this_writer, f);
3465
_write_string(f, "part");
3466
_write_int(f, part_num);
3467
if (mesh->name != NULL)
3468
_write_string(f, mesh->name);
3470
_write_string(f, "unnamed");
3472
/* Build list of sections that are used here, in order of output */
3474
export_list = fvm_writer_export_list(mesh,
3475
fvm_nodal_get_max_entity_dim(mesh),
3477
this_writer->discard_polygons,
3478
this_writer->discard_polyhedra,
3479
this_writer->divide_polygons,
3480
this_writer->divide_polyhedra);
3482
/* Vertex coordinates */
3483
/*--------------------*/
3485
#if defined(HAVE_MPI)
3487
_export_vertex_coords_g(this_writer, mesh, f);
3491
_export_vertex_coords_l(this_writer, mesh, f);
3493
/* If no sections are present (i.e. we may only have vertices),
3494
add "point" elements */
3496
if (export_list == NULL) {
3498
#if defined(HAVE_MPI)
3500
_export_point_elements_g(mesh, this_writer->comm, f);
3503
_export_point_elements_l(mesh, f);
3507
/* Element connectivity */
3508
/*----------------------*/
3510
export_section = export_list;
3512
while (export_section != NULL) {
3514
const fvm_nodal_section_t *section = export_section->section;
3516
/* Print header if start of corresponding EnSight section */
3518
if (export_section->continues_previous == false) {
3520
fvm_gnum_t n_g_elements = 0;
3521
const fvm_writer_section_t *next_section = export_section;
3525
if (next_section->section->type == export_section->type)
3526
n_g_elements += fvm_nodal_section_n_g_elements(next_section->section);
3529
fvm_gnum_t n_g_sub_elements = 0;
3530
fvm_tesselation_get_global_size(next_section->section->tesselation,
3534
n_g_elements += n_g_sub_elements;
3537
next_section = next_section->next;
3539
} while (next_section != NULL && next_section->continues_previous == true);
3541
_write_string(f, _ensight_type_name[export_section->type]);
3542
_write_int(f, n_g_elements);
3545
/* Output for strided (regular) element types */
3546
/*--------------------------------------------*/
3548
if (section->stride > 0) {
3550
#if defined(HAVE_MPI)
3553
export_section = _export_nodal_strided_g(export_section,
3554
mesh->global_vertex_num,
3558
#endif /* defined(HAVE_MPI) */
3560
if (n_ranks == 1) { /* start of output in serial mode */
3562
_write_connect_l(section->stride,
3563
section->n_elements,
3564
section->vertex_num,
3567
export_section = export_section->next;
3571
} /* end of output for strided element types */
3573
/* Output for tesselated polygons or polyhedra */
3574
/*---------------------------------------------*/
3576
else if (export_section->type != section->type) {
3578
#if defined(HAVE_MPI)
3580
/* output in parallel mode */
3583
export_section = _export_nodal_tesselated_g(export_section,
3584
mesh->global_vertex_num,
3587
#endif /* defined(HAVE_MPI) */
3590
export_section = _export_nodal_tesselated_l(export_section,
3595
/* Output for polygons */
3596
/*---------------------*/
3598
else if (export_section->type == FVM_FACE_POLY) {
3599
#if defined(HAVE_MPI)
3601
/* output in parallel mode */
3604
export_section = _export_nodal_polygons_g(export_section,
3605
mesh->global_vertex_num,
3608
#endif /* defined(HAVE_MPI) */
3611
export_section = _export_nodal_polygons_l(export_section,
3616
/* Output for polyhedra */
3617
/*----------------------*/
3619
else if (export_section->type == FVM_CELL_POLY) {
3621
#if defined(HAVE_MPI)
3623
/* output in parallel mode */
3626
export_section =_export_nodal_polyhedra_g(export_section,
3627
mesh->global_vertex_num,
3631
#endif /* defined(HAVE_MPI) */
3634
export_section = _export_nodal_polyhedra_l(export_section,
3639
} /* End of loop on sections */
3641
BFT_FREE(export_list);
3643
/* Close geometry file and update case file */
3644
/*------------------------------------------*/
3646
_free_ensight_file(&f);
3648
fvm_to_ensight_case_write_case(this_writer->case_info, rank);
3651
/*----------------------------------------------------------------------------
3652
* Write field associated with a nodal mesh to an EnSight Gold file.
3654
* Assigning a negative value to the time step indicates a time-independent
3655
* field (in which case the time_value argument is unused).
3658
* this_writer_p <-- pointer to associated writer
3659
* mesh <-- pointer to associated nodal mesh structure
3660
* name <-- variable name
3661
* location <-- variable definition location (nodes or elements)
3662
* dimension <-- variable dimension (0: constant, 1: scalar,
3663
* 3: vector, 6: sym. tensor, 9: asym. tensor)
3664
* interlace <-- indicates if variable in memory is interlaced
3665
* n_parent_lists <-- indicates if variable values are to be obtained
3666
* directly through the local entity index (when 0) or
3667
* through the parent entity numbers (when 1 or more)
3668
* parent_num_shift <-- parent number to value array index shifts;
3669
* size: n_parent_lists
3670
* datatype <-- indicates the data type of (source) field values
3671
* time_step <-- number of the current time step
3672
* time_value <-- associated time value
3673
* field_values <-- array of associated field value arrays
3674
*----------------------------------------------------------------------------*/
3677
fvm_to_ensight_export_field(void *this_writer_p,
3678
const fvm_nodal_t *mesh,
3680
fvm_writer_var_loc_t location,
3682
fvm_interlace_t interlace,
3684
const fvm_lnum_t parent_num_shift[],
3685
fvm_datatype_t datatype,
3688
const void *const field_values[])
3690
int output_dim, part_num;
3691
fvm_to_ensight_case_file_info_t file_info;
3693
const fvm_writer_section_t *export_section = NULL;
3694
fvm_writer_field_helper_t *helper = NULL;
3695
fvm_writer_section_t *export_list = NULL;
3696
fvm_to_ensight_writer_t *this_writer
3697
= (fvm_to_ensight_writer_t *)this_writer_p;
3698
_ensight_file_t f = {NULL, NULL};
3700
const int rank = this_writer->rank;
3701
const int n_ranks = this_writer->n_ranks;
3703
/* Initialization */
3704
/*----------------*/
3708
output_dim = dimension;
3711
else if (dimension > 3 && dimension != 6 && dimension != 9)
3712
bft_error(__FILE__, __LINE__, 0,
3713
_("Data of dimension %d not handled"), dimension);
3715
/* Get part number */
3717
part_num = fvm_to_ensight_case_get_part_num(this_writer->case_info,
3720
part_num = fvm_to_ensight_case_add_part(this_writer->case_info,
3723
/* Open variable file */
3725
file_info = fvm_to_ensight_case_get_var_file(this_writer->case_info,
3732
f = _open_ensight_file(this_writer, file_info.name, file_info.queried);
3734
if (file_info.queried == false) {
3738
/* New files start with description line */
3741
snprintf(buf, 80, "%s (time values: %d, %g)",
3742
name, time_step, time_value);
3744
strncpy(buf, name, 80);
3746
strncpy(buf, name, 80);
3749
_write_string(f, buf);
3752
/* Initialize writer helper */
3753
/*--------------------------*/
3755
/* Build list of sections that are used here, in order of output */
3757
export_list = fvm_writer_export_list(mesh,
3758
fvm_nodal_get_max_entity_dim(mesh),
3760
this_writer->discard_polygons,
3761
this_writer->discard_polyhedra,
3762
this_writer->divide_polygons,
3763
this_writer->divide_polyhedra);
3766
helper = fvm_writer_field_helper_create(mesh,
3775
_write_string(f, "part");
3776
_write_int(f, part_num);
3778
/* Per node variable */
3779
/*-------------------*/
3781
if (location == FVM_WRITER_PER_NODE) {
3783
_write_string(f, "coordinates");
3785
#if defined(HAVE_MPI)
3788
_export_field_values_ng(mesh,
3789
this_writer->divide_polyhedra,
3800
#endif /* defined(HAVE_MPI) */
3803
_export_field_values_nl(mesh,
3814
/* Per element variable */
3815
/*----------------------*/
3817
else if (location == FVM_WRITER_PER_ELEMENT) {
3819
export_section = export_list;
3821
while (export_section != NULL) {
3823
/* Print header if start of corresponding EnSight section */
3825
if (export_section->continues_previous == false)
3826
_write_string(f, _ensight_type_name[export_section->type]);
3828
/* Output per grouped sections */
3830
#if defined(HAVE_MPI)
3833
export_section = _export_field_values_eg(export_section,
3844
#endif /* defined(HAVE_MPI) */
3847
export_section = _export_field_values_el(export_section,
3857
} /* End of loop on sections */
3859
} /* End for per element variable */
3861
/* Free helper structures */
3862
/*------------------------*/
3865
helper = fvm_writer_field_helper_destroy(helper);
3867
BFT_FREE(export_list);
3869
/* Close variable file and update case file */
3870
/*------------------------------------------*/
3872
_free_ensight_file(&f);
3874
fvm_to_ensight_case_write_case(this_writer->case_info, rank);
3877
/*----------------------------------------------------------------------------*/
3881
#endif /* __cplusplus */