1 /*============================================================================
2 * Low level file I/O utility functions for Preprocessor and restart files
3 *============================================================================*/
4
5 /*
6 This file is part of Code_Saturne, a general-purpose CFD tool.
7
8 Copyright (C) 1998-2021 EDF S.A.
9
10 This program is free software; you can redistribute it and/or modify it under
11 the terms of the GNU General Public License as published by the Free Software
12 Foundation; either version 2 of the License, or (at your option) any later
13 version.
14
15 This program is distributed in the hope that it will be useful, but WITHOUT
16 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
17 FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
18 details.
19
20 You should have received a copy of the GNU General Public License along with
21 this program; if not, write to the Free Software Foundation, Inc., 51 Franklin
22 Street, Fifth Floor, Boston, MA 02110-1301, USA.
23 */
24
25 /*----------------------------------------------------------------------------*/
26
27 #include "cs_defs.h"
28
29 /*----------------------------------------------------------------------------
30 * Standard C library headers
31 *----------------------------------------------------------------------------*/
32
33 #include <assert.h>
34 #include <errno.h>
35 #include <stdio.h>
36 #include <stdlib.h>
37 #include <string.h>
38
39 #if defined(HAVE_MPI)
40 #include <mpi.h>
41 #endif
42
43 #undef HAVE_STDINT_H
44 #if defined(__STDC_VERSION__)
45 # if (__STDC_VERSION__ >= 199901L)
46 # define HAVE_STDINT_H
47 # include <stdint.h>
48 # endif
49 #endif
50
51 /*----------------------------------------------------------------------------
52 * Local headers
53 *----------------------------------------------------------------------------*/
54
55 #include "bft_error.h"
56 #include "bft_mem.h"
57 #include "bft_printf.h"
58
59 #include "cs_base.h"
60 #include "cs_log.h"
61 #include "cs_map.h"
62 #include "cs_file.h"
63 #include "cs_timer.h"
64
65 /*----------------------------------------------------------------------------
66 * Header for the current file
67 *----------------------------------------------------------------------------*/
68
69 #include "cs_io.h"
70
71 /*----------------------------------------------------------------------------*/
72
73 BEGIN_C_DECLS
74
75 /*! \cond DOXYGEN_SHOULD_SKIP_THIS */
76
77 /*============================================================================
78 * Local types and structures
79 *============================================================================*/
80
81 /* Basic per linear system options and logging */
82 /*---------------------------------------------*/
83
84 typedef struct {
85
86 unsigned n_opens; /* Number of times file opened */
87
88 double wtimes[3]; /* Wall-clock time for headers,
89 data, and open time */
90
91 unsigned long long data_size[2]; /* Cumulative header and data
92 size */
93
94 } cs_io_log_t;
95
96 /* Structure used to index cs_io_file contents when reading */
97 /*----------------------------------------------------------*/
98
99 typedef struct {
100
101 size_t size; /* Current number of entries */
102 size_t max_size; /* Maximum number of entries */
103
104 /* For each entry, we need 8 values, which we store in h_vals :
105 * 0: number of values in section
106 * 1: location_id
107 * 2: index id
108 * 3: number of values per location
109 * 4: index of section name in names array
110 * 5: index of embedded data in data array + 1 if data is
111 * embedded, 0 otherwise
112 * 6: datatype id in file
113 * 7: associated file id (in case of multiple files)
114 */
115
116 cs_file_off_t *h_vals; /* Base values associated
117 with each header */
118
119 cs_file_off_t *offset; /* Position of associated data
120 in file (-1 if embedded) */
121
122 size_t max_names_size; /* Maximum size of names array */
123 size_t names_size; /* Current size of names array */
124 char *names; /* Array containing section names */
125
126 size_t max_data_size; /* Maximum size of embedded data array */
127 size_t data_size; /* Current size of data array */
128 unsigned char *data; /* Array containing embedded data */
129
130 } cs_io_sec_index_t;
131
132 /* Main kernel IO state structure */
133 /*--------------------------------*/
134
135 struct _cs_io_t {
136
137 /* File information */
138
139 cs_file_t *f; /* Pointer to associated file */
140
141 char contents[64]; /* String describing file contents */
142
143 cs_io_mode_t mode; /* File access mode */
144
145 size_t header_size; /* Header default size */
146 size_t header_align; /* Header alignment */
147 size_t body_align; /* Body alignment */
148
149 cs_io_sec_index_t *index; /* Optional section index (on read) */
150
151 /* Current section buffer state */
152
153 size_t buffer_size; /* Current size of header buffer */
154 unsigned char *buffer; /* Header buffer */
155
156 cs_file_off_t n_vals; /* Number of values in section header */
157 size_t location_id; /* Id of location, or 0 */
158 size_t index_id; /* Id of index, or 0 */
159 size_t n_loc_vals; /* Number of values per location */
160 size_t type_size; /* Size of current type */
161 char *sec_name; /* Pointer to name in section header */
162 char *type_name; /* Pointer to type in section header */
163 void *data; /* Pointer to data in section header
164 (if embedded; NULL otherwise) */
165
166 /* Other flags */
167
168 long echo; /* Data echo level (verbosity) */
169 int log_id; /* Id of log entry, or -1 */
170 double start_time; /* Wall-clock time at open */
171
172 #if defined(HAVE_MPI)
173 MPI_Comm comm; /* Assigned communicator */
174 #endif
175 };
176
177 /*============================================================================
178 * Constants and Macros
179 *============================================================================*/
180
181 #define CS_IO_MPI_TAG 'C'+'S'+'_'+'I'+'O'
182
183 /*============================================================================
184 * Static global variables
185 *============================================================================*/
186
187 static char _type_name_none[] = " ";
188 static char _type_name_char[] = "c "; /* Character string */
189 static char _type_name_i4[] = "i4"; /* Signed 32 bit integer */
190 static char _type_name_i8[] = "i8"; /* Signed 64 bit integer */
191 static char _type_name_u4[] = "u4"; /* Unsigned 32 bit integer */
192 static char _type_name_u8[] = "u8"; /* Unsigned 64 bit integer */
193 static char _type_name_r4[] = "r4"; /* Single precision real */
194 static char _type_name_r8[] = "r8"; /* Double precsision real */
195
196 /* Logging */
197
198 static int _cs_io_map_size[2] = {0, 0};
199 static int _cs_io_map_size_max[2] = {0, 0};
200 static cs_map_name_to_id_t *_cs_io_map[2] = {NULL, NULL};
201 static cs_io_log_t *_cs_io_log[2] = {NULL, NULL};
202
203 /*============================================================================
204 * Private function definitions
205 *============================================================================*/
206
207 /*----------------------------------------------------------------------------
208 * Convert data from "little-endian" to "big-endian" or the reverse.
209 *
210 * parameters:
211 * buf <-- pointer to converted data location.
212 * size <-- size of each item of data in bytes.
213 * ni <-- number of data items.
214 */
215 /*----------------------------------------------------------------------------*/
216
217 static void
_swap_endian(void * buf,size_t size,size_t ni)218 _swap_endian(void *buf,
219 size_t size,
220 size_t ni)
221 {
222 size_t i, ib, shift;
223 unsigned char tmpswap;
224 unsigned char *pbuf = (unsigned char *)buf;
225
226 for (i = 0; i < ni; i++) {
227
228 shift = i * size;
229
230 for (ib = 0; ib < (size / 2); ib++) {
231
232 tmpswap = *(pbuf + shift + ib);
233 *(pbuf + shift + ib) = *(pbuf + shift + (size - 1) - ib);
234 *(pbuf + shift + (size - 1) - ib) = tmpswap;
235
236 }
237
238 }
239 }
240
241 /*----------------------------------------------------------------------------
242 * Default conversion rule from type in file to type in memory.
243 *
244 * parameters:
245 * type_read <-- type in file
246 *
247 * returns:
248 * default corresponding type in memory (may need conversion)
249 *----------------------------------------------------------------------------*/
250
251 static cs_datatype_t
_type_read_to_elt_type(cs_datatype_t type_read)252 _type_read_to_elt_type(cs_datatype_t type_read)
253 {
254 cs_datatype_t elt_type = CS_DATATYPE_NULL;
255
256 if (type_read == CS_INT32 || type_read == CS_INT64) {
257 assert(sizeof(cs_lnum_t) == 4 || sizeof(cs_lnum_t) == 8);
258 if (sizeof(cs_lnum_t) == 4)
259 elt_type = CS_INT32;
260 else
261 elt_type = CS_INT64;
262 }
263
264 else if (type_read == CS_UINT32 || type_read == CS_UINT64) {
265 assert(sizeof(cs_gnum_t) == 4 || sizeof(cs_gnum_t) == 8);
266 if (sizeof(cs_gnum_t) == 4)
267 elt_type = CS_UINT32;
268 else
269 elt_type = CS_UINT64;
270 }
271
272 else if (type_read == CS_FLOAT || type_read == CS_DOUBLE) {
273 if (sizeof(cs_real_t) == 4)
274 elt_type = CS_FLOAT;
275 else
276 elt_type = CS_DOUBLE;
277 }
278
279 else if (type_read == CS_CHAR)
280 elt_type = CS_CHAR;
281
282 return elt_type;
283 }
284
285 /*----------------------------------------------------------------------------
286 * Convert a buffer of type uint64_t to cs_file_off_t
287 *
288 * parameters:
289 * buf <-- buffer
290 * val --> array to which values are converted
291 * n <-- number of values to convert
292 *----------------------------------------------------------------------------*/
293
294 static void
_convert_to_offset(const unsigned char buf[],cs_file_off_t val[],size_t n)295 _convert_to_offset(const unsigned char buf[],
296 cs_file_off_t val[],
297 size_t n)
298 {
299 size_t i;
300
301 #if defined(HAVE_STDINT_H)
302
303 #pragma _NEC novector
304 for (i = 0; i < n; i++)
305 val[i] = ((const uint64_t *)buf)[i];
306
307 #else
308
309 if (sizeof(size_t) == 8) {
310 for (i = 0; i < n; i++)
311 val[i] = ((const size_t *)buf)[i];
312 }
313 else if (sizeof(unsigned long long) == 8) {
314 for (i = 0; i < n; i++)
315 val[i] = ((const unsigned long long *)buf)[i];
316 }
317 else
318 bft_error(__FILE__, __LINE__, 0,
319 "Compilation configuration / porting error:\n"
320 "Unable to determine a 64-bit unsigned int type.\n"
321 "size_t is %d bits, unsigned long long %d bits",
322 sizeof(size_t)*8, sizeof(unsigned long long)*8);
323
324 #endif
325 }
326
327 /*----------------------------------------------------------------------------
328 * Convert a buffer of type cs_file_off_t to uint64_t
329 *
330 * parameters:
331 * buf --> buffer
332 * val <-- array from which values are converted
333 * n <-- number of values to convert
334 *----------------------------------------------------------------------------*/
335
336 static void
_convert_from_offset(unsigned char buf[],const cs_file_off_t val[],size_t n)337 _convert_from_offset(unsigned char buf[],
338 const cs_file_off_t val[],
339 size_t n)
340 {
341 size_t i;
342
343 #if defined(HAVE_STDINT_H)
344
345 for (i = 0; i < n; i++)
346 ((uint64_t *)buf)[i]= val[i];
347
348 #else
349
350 if (sizeof(size_t) == 8) {
351 for (i = 0; i < n; i++)
352 ((size_t *)buf)[i] = val[i];
353 }
354 else if (sizeof(unsigned long long) == 8) {
355 for (i = 0; i < n; i++)
356 ((unsigned long long *)buf)[i] = val[i];
357 }
358 else
359 bft_error(__FILE__, __LINE__, 0,
360 "Compilation configuration / porting error:\n"
361 "Unable to determine a 64-bit unsigned int type.\n"
362 "size_t is %d bits, unsigned long long %d bits",
363 sizeof(size_t)*8, sizeof(unsigned long long)*8);
364
365 #endif
366 }
367
368 /*----------------------------------------------------------------------------
369 * Return an empty kernel IO file structure.
370 *
371 * parameters:
372 * mode --> read or write
373 * echo --> echo on main output (< 0 if none, header if 0,
374 * n first and last elements if n > 0)
375 *
376 * returns:
377 * pointer to kernel IO structure
378 *----------------------------------------------------------------------------*/
379
380 static cs_io_t *
_cs_io_create(cs_io_mode_t mode,size_t echo)381 _cs_io_create(cs_io_mode_t mode,
382 size_t echo)
383 {
384 cs_io_t *cs_io = NULL;
385
386 BFT_MALLOC(cs_io, 1, cs_io_t);
387
388 /* Set structure fields */
389
390 cs_io->mode = mode;
391
392 cs_io->f = NULL;
393
394 memset(cs_io->contents, 0, 64);
395
396 cs_io->header_size = 0;
397 cs_io->header_align = 0;
398 cs_io->body_align = 0;
399
400 cs_io->index = NULL;
401
402 /* Current section buffer state */
403
404 cs_io->buffer_size = 0;
405 cs_io->buffer = NULL;
406
407 cs_io->n_vals = 0;
408 cs_io->type_size = 0;
409 cs_io->sec_name = NULL;
410 cs_io->type_name = NULL;
411 cs_io->data = NULL;
412
413 /* Verbosity and logging */
414
415 cs_io->echo = echo;
416 cs_io->log_id = -1;
417 cs_io->start_time = 0;
418
419 #if defined(HAVE_MPI)
420 cs_io->comm = MPI_COMM_NULL;
421 #endif
422
423 return cs_io;
424 }
425
426 /*----------------------------------------------------------------------------
427 * Add an empty index structure to a cs_io_t structure.
428 *
429 * parameters:
430 * inp <-> pointer to cs_io_t structure
431 *----------------------------------------------------------------------------*/
432
433 static void
_create_index(cs_io_t * inp)434 _create_index(cs_io_t *inp)
435 {
436 cs_io_sec_index_t *idx = NULL;
437
438 BFT_MALLOC(idx, 1, cs_io_sec_index_t);
439
440 /* Set structure fields */
441
442 idx->size = 0;
443 idx->max_size = 32;
444
445 BFT_MALLOC(idx->h_vals, idx->max_size*7, cs_file_off_t);
446 BFT_MALLOC(idx->offset, idx->max_size, cs_file_off_t);
447
448 idx->max_names_size = 256;
449 idx->names_size = 0;
450
451 BFT_MALLOC(idx->names, idx->max_names_size, char);
452
453 idx->max_data_size = 256;
454 idx->data_size = 0;
455
456 BFT_MALLOC(idx->data, idx->max_data_size, unsigned char);
457
458 /* Add structure */
459
460 inp->index = idx;
461 }
462
463 /*----------------------------------------------------------------------------
464 * Destroy a cs_io_t structure's optional index structure.
465 *
466 * parameters:
467 * inp <-> pointer to cs_io_t structure
468 *----------------------------------------------------------------------------*/
469
470 static void
_destroy_index(cs_io_t * inp)471 _destroy_index(cs_io_t *inp)
472 {
473 cs_io_sec_index_t *idx = inp->index;
474
475 if (idx == NULL)
476 return;
477
478 BFT_FREE(idx->h_vals);
479 BFT_FREE(idx->offset);
480 BFT_FREE(idx->names);
481 BFT_FREE(idx->data);
482
483 BFT_FREE(inp->index);
484 }
485
486 /*----------------------------------------------------------------------------
487 * Update an index structure with info from the last header read
488 *
489 * Also sets the file position for the next read
490 *
491 * parameters:
492 * inp <-> input kernel IO structure
493 * header <-- header structure
494 *----------------------------------------------------------------------------*/
495
496 static void
_update_index_and_shift(cs_io_t * inp,cs_io_sec_header_t * header)497 _update_index_and_shift(cs_io_t *inp,
498 cs_io_sec_header_t *header)
499 {
500 size_t id = 0;
501 size_t new_names_size = 0;
502 size_t new_data_size = 0;
503
504 cs_io_sec_index_t *idx = inp->index;
505
506 if (idx == NULL)
507 return;
508
509 /* Reallocate if necessary */
510
511 if (idx->size + 1 == idx->max_size) {
512 if (idx->max_size == 0)
513 idx->max_size = 32;
514 else
515 idx->max_size *= 2;
516 BFT_REALLOC(idx->h_vals, idx->max_size*7, cs_file_off_t);
517 BFT_REALLOC(idx->offset, idx->max_size, cs_file_off_t);
518 };
519
520 new_names_size = idx->names_size + strlen(inp->sec_name) + 1;
521
522 if (inp->data != NULL)
523 new_data_size
524 = idx->data_size + ( inp->n_vals
525 * cs_datatype_size[header->type_read]);
526
527 if (new_names_size > idx->max_names_size) {
528 if (idx->max_names_size == 0)
529 idx->max_names_size = 128;
530 while (new_names_size > idx->max_names_size)
531 idx->max_names_size *= 2;
532 BFT_REALLOC(idx->names, idx->max_names_size, char);
533 }
534
535 if (new_data_size > idx->max_data_size) {
536 if (idx->max_data_size == 0)
537 idx->max_data_size = 128;
538 while (new_data_size > idx->max_data_size)
539 idx->max_data_size *= 2;
540 BFT_REALLOC(idx->data, idx->max_data_size, unsigned char);
541 }
542
543 /* Set values */
544
545 id = idx->size;
546
547 idx->h_vals[id*7] = inp->n_vals;
548 idx->h_vals[id*7 + 1] = inp->location_id;
549 idx->h_vals[id*7 + 2] = inp->index_id;
550 idx->h_vals[id*7 + 3] = inp->n_loc_vals;
551 idx->h_vals[id*7 + 4] = idx->names_size;
552 idx->h_vals[id*7 + 5] = 0;
553 idx->h_vals[id*7 + 6] = header->type_read;
554
555 strcpy(idx->names + idx->names_size, inp->sec_name);
556 idx->names[new_names_size - 1] = '\0';
557 idx->names_size = new_names_size;
558
559 if (inp->data == NULL) {
560 cs_file_off_t offset = cs_file_tell(inp->f);
561 cs_file_off_t data_shift = inp->n_vals * inp->type_size;
562 if (inp->body_align > 0) {
563 size_t ba = inp->body_align;
564 idx->offset[id] = offset + (ba - (offset % ba)) % ba;
565 }
566 else
567 idx->offset[id] = offset;
568 cs_file_seek(inp->f, idx->offset[id] + data_shift, CS_FILE_SEEK_SET);
569 }
570 else {
571 idx->h_vals[id*7 + 5] = idx->data_size + 1;
572 memcpy(idx->data + idx->data_size,
573 inp->data,
574 new_data_size - idx->data_size);
575 idx->data_size = new_data_size;
576 idx->offset[id] = -1;
577 }
578
579 idx->size += 1;
580 }
581
582 /*----------------------------------------------------------------------------
583 * Open the interface file descriptor and initialize the file by writing
584 * or reading a "magic string" used to check the file content type.
585 *
586 * parameters:
587 * cs_io <-> kernel IO structure
588 * name <-- file name
589 * magic_string <-- magic string associated with file content type
590 * method <-- file access method options
591 * hints <-- MPI-IO hints
592 * block_comm <-- associated MPI communicator for block IO
593 * comm <-- associated MPI communicator
594 *----------------------------------------------------------------------------*/
595
596 #if defined(HAVE_MPI)
597 static void
_file_open(cs_io_t * cs_io,const char * name,const char * magic_string,cs_file_access_t method,MPI_Info hints,MPI_Comm block_comm,MPI_Comm comm)598 _file_open(cs_io_t *cs_io,
599 const char *name,
600 const char *magic_string,
601 cs_file_access_t method,
602 MPI_Info hints,
603 MPI_Comm block_comm,
604 MPI_Comm comm)
605 #else
606 static void
607 _file_open(cs_io_t *cs_io,
608 const char *name,
609 const char *magic_string,
610 cs_file_access_t method)
611 #endif
612 {
613 cs_file_mode_t f_mode;
614 char header_data[128 + 24];
615 cs_file_off_t header_vals[3];
616
617 char base_header[] = "Code_Saturne I/O, BE, R0";
618
619 /* Prepare file open */
620
621 switch(cs_io->mode) {
622
623 case CS_IO_MODE_READ:
624 f_mode = CS_FILE_MODE_READ;
625 break;
626
627 case CS_IO_MODE_WRITE:
628 f_mode = CS_FILE_MODE_WRITE;
629 break;
630
631 default:
632 assert( cs_io->mode == CS_IO_MODE_READ
633 || cs_io->mode == CS_IO_MODE_WRITE);
634 return;
635 }
636
637 /* Start logging if active */
638
639 if (_cs_io_map[cs_io->mode] != NULL) {
640
641 int mode = cs_io->mode;
642 int i = cs_map_name_to_id(_cs_io_map[mode], name);
643 cs_io_log_t *l = NULL;
644
645 if (i >= _cs_io_map_size_max[mode]) {
646 _cs_io_map_size_max[mode] *= 2;
647 BFT_REALLOC(_cs_io_log[mode], _cs_io_map_size_max[mode], cs_io_log_t);
648 }
649 cs_io->log_id = i;
650
651 l = _cs_io_log[mode] + i;
652
653 if (i >= _cs_io_map_size[mode]) {
654
655 /* Add new log entry */
656
657 int j;
658
659 l->n_opens = 1;
660 for (j = 0; j < 3; j++)
661 l->wtimes[j] = 0.0;
662 for (j = 0; j < 2; j++)
663 l->data_size[j] = 0;
664
665 _cs_io_map_size[mode] += 1;
666 }
667 else
668 l->n_opens += 1;
669 }
670
671 cs_io->start_time = cs_timer_wtime();
672
673 /* Create interface file descriptor */
674
675 #if defined(HAVE_MPI)
676 cs_io->f = cs_file_open(name, f_mode, method, hints, block_comm, comm);
677 #else
678 cs_io->f = cs_file_open(name, f_mode, method);
679 #endif
680
681 cs_file_set_big_endian(cs_io->f);
682
683 #if defined(HAVE_MPI)
684 cs_io->comm = comm;
685 #endif
686
687 /* Write or read a magic string */
688 /*------------------------------*/
689
690 if (cs_io->mode == CS_IO_MODE_READ) {
691
692 cs_file_read_global(cs_io->f, header_data, 1, 128 + 24);
693
694 header_data[63] = '\0';
695 header_data[127] = '\0';
696
697 /* If the format does not correspond, we have an error */
698
699 if (strncmp(header_data, base_header, 64) != 0) {
700
701 bft_error(__FILE__, __LINE__, 0,
702 _("Error reading file: \"%s\".\n"
703 "File format is not the correct version.\n"
704 "The first 64 bytes expected contain:\n"
705 "\"%s\"\n"
706 "The first 64 bytes read contain:\n"
707 "\"%s\"\n"),
708 cs_file_get_name(cs_io->f), base_header, header_data);
709
710 }
711
712 /* Copy magic string */
713
714 strncpy(cs_io->contents, header_data + 64, 64);
715 cs_io->contents[63] = '\0';
716
717 /* If the magic string does not correspond, we have an error */
718
719 if (magic_string != NULL) {
720 if (strncmp(cs_io->contents, magic_string, 64) != 0)
721 bft_error(__FILE__, __LINE__, 0,
722 _("Error reading file: \"%s\".\n"
723 "The file contents are not of the expected type.\n"
724 "\"%s\" was expected,\n"
725 "\"%s\" was read."),
726 cs_file_get_name(cs_io->f), magic_string, cs_io->contents);
727 }
728
729 /* Now decode the sizes */
730
731 if (cs_file_get_swap_endian(cs_io->f) == 1)
732 _swap_endian(header_data + 128, 8, 3);
733
734 _convert_to_offset((unsigned char *)(header_data + 128), header_vals, 3);
735
736 cs_io->header_size = header_vals[0];
737 cs_io->header_align = header_vals[1];
738 cs_io->body_align = header_vals[2];
739
740 }
741 else if (cs_io->mode == CS_IO_MODE_WRITE) {
742
743 size_t n_written = 0;
744
745 memset(header_data, 0, sizeof(header_data));
746 strcpy(header_data, base_header);
747 strncpy(header_data + 64, magic_string, 64);
748 header_data[127] = '\0';
749
750 /* Set default header size and alignments */
751
752 cs_io->header_size = 128;
753 cs_io->header_align = 64;
754 cs_io->body_align = 64;
755
756 header_vals[0] = cs_io->header_size;
757 header_vals[1] = cs_io->header_align;
758 header_vals[2] = cs_io->body_align;
759
760 _convert_from_offset((unsigned char *)(header_data + 128), header_vals, 3);
761
762 if (cs_file_get_swap_endian(cs_io->f) == 1)
763 _swap_endian(header_data + 128, 8, 3);
764
765 n_written = cs_file_write_global(cs_io->f, header_data, 1, 128 + 24);
766
767 if (n_written < 128 + 24)
768 bft_error(__FILE__, __LINE__, 0,
769 _("Error writing the header of file: \"%s\".\n"),
770 cs_file_get_name(cs_io->f));
771 }
772
773 cs_io->buffer_size = cs_io->header_size;
774 BFT_MALLOC(cs_io->buffer, cs_io->buffer_size, unsigned char);
775 }
776
777 /*----------------------------------------------------------------------------
778 * Re-open the interface file descriptor when building an index.
779 *
780 * parameters:
781 * inp <-> kernel IO structure
782 * method <-- file access method options
783 * hints <-- MPI-IO hints
784 * block_comm <-- associated MPI communicator for block IO
785 * comm <-- associated MPI communicator
786 *----------------------------------------------------------------------------*/
787
788 #if defined(HAVE_MPI)
789 static void
_file_reopen_read(cs_io_t * inp,cs_file_access_t method,MPI_Info hints,MPI_Comm block_comm,MPI_Comm comm)790 _file_reopen_read(cs_io_t *inp,
791 cs_file_access_t method,
792 MPI_Info hints,
793 MPI_Comm block_comm,
794 MPI_Comm comm)
795 #else
796 static void
797 _file_reopen_read(cs_io_t *inp,
798 cs_file_access_t method)
799 #endif /* HAVE_MPI */
800 {
801 char _tmpname[128];
802 char *tmpname = _tmpname;
803
804 assert(inp->index != NULL);
805
806 if (inp->f == NULL)
807 return;
808
809 const char *filename = cs_file_get_name(inp->f);
810
811 if (strlen(filename) >= 128)
812 BFT_MALLOC(tmpname, strlen(filename) + 1, char);
813 strcpy(tmpname, filename);
814
815 inp->f = cs_file_free(inp->f);
816
817 #if defined(HAVE_MPI)
818 inp->f = cs_file_open(tmpname,
819 CS_FILE_MODE_READ,
820 method,
821 hints,
822 block_comm,
823 comm);
824 #else
825 inp->f = cs_file_open(tmpname, CS_FILE_MODE_READ, method);
826 #endif
827
828 cs_file_set_big_endian(inp->f);
829
830 if (tmpname != _tmpname)
831 BFT_FREE(tmpname);
832 }
833
834 /*----------------------------------------------------------------------------
835 * Close the interface file.
836 *
837 * parameters:
838 * cs_io <-> kernel IO structure
839 *----------------------------------------------------------------------------*/
840
841 static void
_file_close(cs_io_t * cs_io)842 _file_close(cs_io_t *cs_io)
843 {
844 if (cs_io->f != NULL)
845 cs_io->f = cs_file_free(cs_io->f);
846
847 if (cs_io->log_id > -1) {
848 double t_end = cs_timer_wtime();
849 cs_io_log_t *log = _cs_io_log[cs_io->mode] + cs_io->log_id;
850 log->wtimes[2] += t_end - cs_io->start_time;
851 }
852 }
853
854 /*----------------------------------------------------------------------------
855 * Echo pending section read or write
856 *
857 * parameters:
858 * cs_io --> kernel IO structure
859 *----------------------------------------------------------------------------*/
860
861 static void
_echo_pre(const cs_io_t * cs_io)862 _echo_pre(const cs_io_t *cs_io)
863 {
864 assert(cs_io != NULL);
865
866 switch(cs_io->mode) {
867
868 case CS_IO_MODE_READ:
869 bft_printf(_("\n Section read on \"%s\":\n"),
870 cs_file_get_name(cs_io->f));
871 break;
872
873 case CS_IO_MODE_WRITE:
874 bft_printf(_("\n Section written on \"%s\":\n"),
875 cs_file_get_name(cs_io->f));
876 break;
877
878 default:
879 assert( cs_io->mode == CS_IO_MODE_READ
880 || cs_io->mode == CS_IO_MODE_WRITE);
881 }
882
883 bft_printf_flush();
884 }
885
886 /*----------------------------------------------------------------------------
887 * Echo a section header
888 *
889 * parameters:
890 * sec_name --> section name
891 * n_elts --> number of elements
892 * elt_type --> expected (final) element type
893 * type --> element type in file
894 *----------------------------------------------------------------------------*/
895
896 static void
_echo_header(const char * sec_name,cs_gnum_t n_elts,cs_datatype_t type)897 _echo_header(const char *sec_name,
898 cs_gnum_t n_elts,
899 cs_datatype_t type)
900 {
901 /* Instructions */
902
903 bft_printf(_(" section name: \"%s\"\n"
904 " number of elements: %llu\n"),
905 sec_name, (unsigned long long)n_elts);
906
907 if (n_elts > 0) {
908
909 char *type_name;
910
911 switch(type) {
912 case CS_CHAR:
913 type_name = _type_name_char;
914 break;
915 case CS_INT32:
916 type_name = _type_name_i4;
917 break;
918 case CS_INT64:
919 type_name = _type_name_i8;
920 break;
921 case CS_UINT32:
922 type_name = _type_name_u4;
923 break;
924 case CS_UINT64:
925 type_name = _type_name_u8;
926 break;
927 case CS_FLOAT:
928 type_name = _type_name_r4;
929 break;
930 case CS_DOUBLE:
931 type_name = _type_name_r8;
932 break;
933 default:
934 assert(0);
935 type_name = _type_name_none;
936 }
937
938 bft_printf(_(" element type name: \"%s\"\n"), type_name);
939
940 }
941
942 bft_printf_flush();
943 }
944
945 /*----------------------------------------------------------------------------
946 * Partial echo of a section's contents.
947 *
948 * FVM datatypes must have been converted to the corresponding
949 * Code_Saturne compatible datatype before calling this function:
950 * CS_CHAR -> char
951 * CS_INT32 / CS_INT64 -> cs_lnum_t
952 * CS_UINT32 / CS_UINT64 -> cs_gnum_t
953 * CS_REAL / CS_FLOAT -> double / cs_real_t
954 *
955 * If global_num_start and global_num_end are > 0, the echo shows that
956 * a different block is assigned to each processor. Otherwise, the full
957 * data is replicated for each processor, and this is also indicated.
958 *
959 * parameters:
960 * echo --> echo level
961 * n_elts --> number of elements
962 * global_num_start <-- global number of first block item (1 to n numbering)
963 * global_num_end <-- global number of past-the end block item
964 * (1 to n numbering)
965 * elt_type --> element type
966 * elts --> element buffer
967 *----------------------------------------------------------------------------*/
968
969 static void
_echo_data(size_t echo,cs_file_off_t n_elts,cs_gnum_t global_num_start,cs_gnum_t global_num_end,cs_datatype_t elt_type,const void * elts)970 _echo_data(size_t echo,
971 cs_file_off_t n_elts,
972 cs_gnum_t global_num_start,
973 cs_gnum_t global_num_end,
974 cs_datatype_t elt_type,
975 const void *elts)
976 {
977 cs_file_off_t i;
978 cs_gnum_t num_shift = 1;
979 size_t _n_elts = n_elts;
980 cs_file_off_t echo_start = 0;
981 cs_file_off_t echo_end = 0;
982 const char *_loc_glob[] = {N_(" (local)"), ""};
983 const char *loc_glob = _loc_glob[1];
984
985 /* Instructions */
986
987 if (n_elts == 0) return;
988
989 if (cs_glob_n_ranks == 1)
990 loc_glob = _loc_glob[0];
991 else if (global_num_start > 0) {
992 loc_glob = _loc_glob[0];
993 num_shift = global_num_start;
994 }
995
996 if (global_num_start > 0 && global_num_end > 0) {
997 assert(global_num_end >= global_num_start);
998 _n_elts = global_num_end - global_num_start;
999 }
1000
1001 if (echo * 2 < _n_elts) {
1002 echo_end = echo;
1003 bft_printf(_(" %d first and last elements%s:\n"),
1004 (int)echo, loc_glob);
1005 }
1006 else {
1007 echo_end = _n_elts;
1008 bft_printf(_(" elements%s:\n"), _(loc_glob));
1009 }
1010
1011 /* Note that FVM datatypes will have been converted to
1012 the corresponding datatype if necessary before
1013 calling this function, hence the identical treatment
1014 for different cases. */
1015
1016 do {
1017
1018 switch (elt_type) {
1019
1020 case CS_INT32:
1021 case CS_INT64:
1022 {
1023 const cs_lnum_t *_elts = elts;
1024
1025 for (i = echo_start ; i < echo_end ; i++)
1026 bft_printf(" %10llu : %12ld\n",
1027 (unsigned long long)(i + num_shift),
1028 (long)*(_elts + i));
1029 }
1030 break;
1031
1032 case CS_UINT32:
1033 case CS_UINT64:
1034 {
1035 const cs_gnum_t *_elts = elts;
1036
1037 for (i = echo_start ; i < echo_end ; i++)
1038 bft_printf(" %10llu : %12llu\n",
1039 (unsigned long long)(i + num_shift),
1040 (unsigned long long)*(_elts + i));
1041 }
1042 break;
1043
1044 case CS_FLOAT:
1045 case CS_DOUBLE:
1046 {
1047 const cs_real_t *_elts = elts;
1048
1049 for (i = echo_start ; i < echo_end ; i++)
1050 bft_printf(" %10llu : %12.5e\n",
1051 (unsigned long long)(i + num_shift), *(_elts + i));
1052 }
1053 break;
1054
1055 case CS_CHAR:
1056 {
1057 const char *_elts = elts;
1058
1059 for (i = echo_start ; i < echo_end ; i++) {
1060 if (*(_elts + i) != '\0')
1061 bft_printf(" %10llu : '%c'\n",
1062 (unsigned long long)(i + num_shift), *(_elts + i));
1063 else
1064 bft_printf(" %10llu : '\\0'\n",
1065 (unsigned long long)(i + num_shift));
1066 }
1067 }
1068 break;
1069
1070 default:
1071 assert(0);
1072
1073 }
1074
1075 if (echo_end < (cs_file_off_t)_n_elts) {
1076 bft_printf(" .......... ............\n");
1077 echo_start = _n_elts - echo;
1078 echo_end = _n_elts;
1079 }
1080 else {
1081 assert(echo_end == (cs_file_off_t)_n_elts);
1082 echo_end = _n_elts + 1;
1083 }
1084
1085 } while (echo_end <= (cs_file_off_t)_n_elts);
1086
1087 bft_printf_flush();
1088 }
1089
1090 /*----------------------------------------------------------------------------
1091 * Convert read data.
1092 *
1093 * dest_type must have been set to the corresponding
1094 * Code_Saturne compatible datatype before calling this function and
1095 * conversion will be done accordingly:
1096 * CS_INT32 / CS_INT64 -> cs_lnum_t
1097 * CS_UINT32 / CS_UINT64 -> cs_gnum_t
1098 * CS_REAL / CS_FLOAT -> double / cs_real_t
1099 *
1100 * parameters:
1101 * buffer --> input buffer
1102 * dest <-- output buffer
1103 * n_elts --> number of elements
1104 * buffer_type --> input buffer element type
1105 * dest_type --> output buffer element type
1106 *----------------------------------------------------------------------------*/
1107
1108 static void
_cs_io_convert_read(void * buffer,void * dest,cs_file_off_t n_elts,cs_datatype_t buffer_type,cs_datatype_t dest_type)1109 _cs_io_convert_read(void *buffer,
1110 void *dest,
1111 cs_file_off_t n_elts,
1112 cs_datatype_t buffer_type,
1113 cs_datatype_t dest_type)
1114 {
1115 cs_file_off_t ii;
1116
1117 assert(dest_type != buffer_type);
1118
1119 /* Note that dest will have been set to the corresponding datatype
1120 before calling this function, hence the identical treatment
1121 for different cases. */
1122
1123 switch(dest_type) {
1124
1125 case CS_INT32:
1126 {
1127 int32_t *_dest = dest;
1128
1129 if (buffer_type == CS_INT32) {
1130 int32_t * _buffer = buffer;
1131 for (ii = 0; ii < n_elts; ii++)
1132 _dest[ii] = _buffer[ii];
1133 }
1134 else if (buffer_type == CS_INT64) {
1135 int64_t * _buffer = buffer;
1136 for (ii = 0; ii < n_elts; ii++)
1137 _dest[ii] = _buffer[ii];
1138 }
1139 if (buffer_type == CS_UINT32) {
1140 uint32_t * _buffer = buffer;
1141 for (ii = 0; ii < n_elts; ii++)
1142 _dest[ii] = _buffer[ii];
1143 }
1144 else if (buffer_type == CS_UINT64) {
1145 uint64_t * _buffer = buffer;
1146 for (ii = 0; ii < n_elts; ii++)
1147 _dest[ii] = _buffer[ii];
1148 }
1149 }
1150 break;
1151
1152 case CS_INT64:
1153 {
1154 int64_t *_dest = dest;
1155
1156 if (buffer_type == CS_INT32) {
1157 int32_t * _buffer = buffer;
1158 for (ii = 0; ii < n_elts; ii++)
1159 _dest[ii] = _buffer[ii];
1160 }
1161 else if (buffer_type == CS_INT64) {
1162 int64_t * _buffer = buffer;
1163 for (ii = 0; ii < n_elts; ii++)
1164 _dest[ii] = _buffer[ii];
1165 }
1166 if (buffer_type == CS_UINT32) {
1167 uint32_t * _buffer = buffer;
1168 for (ii = 0; ii < n_elts; ii++)
1169 _dest[ii] = _buffer[ii];
1170 }
1171 else if (buffer_type == CS_UINT64) {
1172 uint64_t * _buffer = buffer;
1173 for (ii = 0; ii < n_elts; ii++)
1174 _dest[ii] = _buffer[ii];
1175 }
1176 }
1177 break;
1178
1179 case CS_UINT32:
1180 {
1181 uint32_t *_dest = dest;
1182
1183 if (buffer_type == CS_INT32) {
1184 int32_t * _buffer = buffer;
1185 for (ii = 0; ii < n_elts; ii++)
1186 _dest[ii] = _buffer[ii];
1187 }
1188 else if (buffer_type == CS_INT64) {
1189 int64_t * _buffer = buffer;
1190 for (ii = 0; ii < n_elts; ii++)
1191 _dest[ii] = _buffer[ii];
1192 }
1193 if (buffer_type == CS_UINT32) {
1194 uint32_t * _buffer = buffer;
1195 for (ii = 0; ii < n_elts; ii++)
1196 _dest[ii] = _buffer[ii];
1197 }
1198 else if (buffer_type == CS_UINT64) {
1199 uint64_t * _buffer = buffer;
1200 for (ii = 0; ii < n_elts; ii++)
1201 _dest[ii] = _buffer[ii];
1202 }
1203 }
1204 break;
1205
1206 case CS_UINT64:
1207 {
1208 uint64_t *_dest = dest;
1209
1210 if (buffer_type == CS_INT32) {
1211 int32_t * _buffer = buffer;
1212 for (ii = 0; ii < n_elts; ii++)
1213 _dest[ii] = _buffer[ii];
1214 }
1215 else if (buffer_type == CS_INT64) {
1216 int64_t * _buffer = buffer;
1217 for (ii = 0; ii < n_elts; ii++)
1218 _dest[ii] = _buffer[ii];
1219 }
1220 if (buffer_type == CS_UINT32) {
1221 uint32_t * _buffer = buffer;
1222 for (ii = 0; ii < n_elts; ii++)
1223 _dest[ii] = _buffer[ii];
1224 }
1225 else if (buffer_type == CS_UINT64) {
1226 uint64_t * _buffer = buffer;
1227 for (ii = 0; ii < n_elts; ii++)
1228 _dest[ii] = _buffer[ii];
1229 }
1230 }
1231 break;
1232
1233 case CS_FLOAT:
1234 {
1235 cs_real_t *_dest = dest;
1236 double * _buffer = buffer;
1237
1238 assert(buffer_type == CS_DOUBLE);
1239
1240 for (ii = 0; ii < n_elts; ii++)
1241 _dest[ii] = _buffer[ii];
1242 }
1243 break;
1244
1245 case CS_DOUBLE:
1246 {
1247 cs_real_t *_dest = dest;
1248 float * _buffer = buffer;
1249
1250 assert(buffer_type == CS_FLOAT);
1251
1252 for (ii = 0; ii < n_elts; ii++)
1253 _dest[ii] = _buffer[ii];
1254 }
1255 break;
1256
1257 default:
1258 assert(0);
1259 }
1260 }
1261
1262 /*----------------------------------------------------------------------------
1263 * Read a section body.
1264 *
1265 * If global_num_start and global_num_end are > 0, a different block is
1266 * assigned to each processor. Otherwise, the full data is replicated
1267 * for each processor.
1268 *
1269 * If location_id > 0 and header->n_location_vals > 1, then
1270 * global_num_start and global_num_end will be based on location element
1271 * numbers, so the total number of values read equals
1272 * (global_num_end - global_num_start) * header->n_location_vals.
1273 *
1274 * If the array intended to receive the data already exists, we pass an
1275 * "elt" pointer to this array; this same pointer is then returned.
1276 * Otherwise, if this pointer is passed as NULL, memory is allocated
1277 * by this function, and the corresponding pointer is returned. It is
1278 * the caller's responsibility to free this array.
1279 *
1280 * parameters:
1281 * header <-- header structure
1282 * global_num_start <-- global number of first block item (1 to n numbering)
1283 * global_num_end <-- global number of past-the end block item
1284 * (1 to n numbering)
1285 * elts <-> pointer to data array, or NULL
1286 * inp --> input kernel IO structure
1287 *
1288 * returns:
1289 * elts if non NULL, or pointer to allocated array otherwise
1290 *----------------------------------------------------------------------------*/
1291
1292 static void *
_cs_io_read_body(const cs_io_sec_header_t * header,cs_gnum_t global_num_start,cs_gnum_t global_num_end,void * elts,cs_io_t * inp)1293 _cs_io_read_body(const cs_io_sec_header_t *header,
1294 cs_gnum_t global_num_start,
1295 cs_gnum_t global_num_end,
1296 void *elts,
1297 cs_io_t *inp)
1298 {
1299 double t_start = 0.;
1300 size_t type_size = 0;
1301 cs_file_off_t n_vals = inp->n_vals;
1302 cs_io_log_t *log = NULL;
1303 bool convert_type = false;
1304 void *_elts = NULL;
1305 void *_buf = NULL;
1306 size_t stride = 1;
1307
1308 assert(inp != NULL);
1309
1310 assert(header->n_vals == inp->n_vals);
1311
1312 assert(global_num_end <= (cs_gnum_t)header->n_vals + 1);
1313
1314 if (inp->log_id > -1) {
1315 log = _cs_io_log[inp->mode] + inp->log_id;
1316 t_start = cs_timer_wtime();
1317 }
1318
1319 /* Choose global or block mode */
1320
1321 if (header->n_location_vals > 1)
1322 stride = header->n_location_vals;
1323
1324 if (global_num_start > 0 && global_num_end > 0) {
1325 assert(global_num_end >= global_num_start);
1326 n_vals = (global_num_end - global_num_start)*stride;
1327 }
1328
1329 /* Datatype size given by FVM datatype */
1330
1331 type_size = cs_datatype_size[header->type_read];
1332
1333 /* Assign or allocate */
1334
1335 _elts = elts;
1336
1337 if (_elts == NULL && n_vals != 0) {
1338 if (header->elt_type == CS_CHAR && header->location_id == 0)
1339 BFT_MALLOC(_elts, n_vals + 1, char);
1340 else
1341 BFT_MALLOC(_elts, n_vals*type_size, char);
1342 }
1343
1344 /* Element values */
1345
1346 if (n_vals != 0 && header->elt_type != header->type_read)
1347 convert_type = true;
1348
1349 if (inp->data != NULL)
1350 _buf = NULL;
1351 else if (convert_type == true
1352 && ( cs_datatype_size[header->type_read]
1353 != cs_datatype_size[header->elt_type]))
1354 BFT_MALLOC(_buf, n_vals*type_size, char);
1355 else
1356 _buf = _elts;
1357
1358 /* Read data from file */
1359
1360 if (inp->data == NULL) {
1361
1362 /* Position read pointer if necessary */
1363
1364 if (inp->body_align > 0) {
1365 cs_file_off_t offset = cs_file_tell(inp->f);
1366 size_t ba = inp->body_align;
1367 offset += (ba - (offset % ba)) % ba;
1368 cs_file_seek(inp->f, offset, CS_FILE_SEEK_SET);
1369 }
1370
1371 /* Read local or global values */
1372
1373 if (global_num_start > 0 && global_num_end > 0) {
1374 cs_file_read_block(inp->f,
1375 _buf,
1376 type_size,
1377 stride,
1378 global_num_start,
1379 global_num_end);
1380 if (log != NULL)
1381 log->data_size[1] += (global_num_end - global_num_start)*type_size;
1382 }
1383
1384 else if (n_vals > 0) {
1385 cs_file_read_global(inp->f,
1386 _buf,
1387 type_size,
1388 n_vals);
1389 if (log != NULL)
1390 log->data_size[0] += n_vals*type_size;
1391 }
1392
1393 }
1394
1395 /* If data is embedded in header, simply point to it */
1396
1397 else {
1398
1399 if (global_num_start > 0 && global_num_end > 0)
1400 _buf = ((unsigned char *)inp->data)
1401 + ( (global_num_start - 1) * stride
1402 * cs_datatype_size[header->type_read]);
1403 else
1404 _buf = inp->data;
1405
1406 }
1407
1408 /* Convert data if necessary */
1409
1410 if (convert_type == true) {
1411 _cs_io_convert_read(_buf,
1412 _elts,
1413 n_vals,
1414 header->type_read,
1415 header->elt_type);
1416 if ( inp->data == NULL
1417 && _buf != _elts)
1418 BFT_FREE(_buf);
1419 }
1420 else if (inp->data != NULL) {
1421 memcpy(_elts, _buf, n_vals*cs_datatype_size[header->type_read]);
1422 _buf = NULL;
1423 }
1424
1425 if (inp->data != NULL) /* Reset for next read */
1426 inp->data = NULL;
1427
1428 /* Add null character at end of string to ensure C-type string */
1429
1430 if (n_vals != 0 && header->elt_type == CS_CHAR && header->location_id == 0)
1431 ((char *)_elts)[header->n_vals] = '\0';
1432
1433 if (log != NULL) {
1434 double t_end = cs_timer_wtime();
1435 int t_id = (global_num_start > 0 && global_num_end > 0) ? 1 : 0;
1436 log->wtimes[t_id] += t_end - t_start;
1437 }
1438
1439 /* Optional echo */
1440
1441 if (header->n_vals != 0 && inp->echo > CS_IO_ECHO_HEADERS)
1442 _echo_data(inp->echo,
1443 n_vals,
1444 (global_num_start-1)*stride + 1,
1445 (global_num_end-1)*stride + 1,
1446 header->elt_type,
1447 _elts);
1448
1449 /* Return read values */
1450
1451 return _elts;
1452 }
1453
1454 /*----------------------------------------------------------------------------
1455 * Build an index for a kernel IO file structure in read mode.
1456 *
1457 * The magic string may be NULL, if we choose to ignore it.
1458 *
1459 * parameters:
1460 * inp <-> empty input kernel IO file structure
1461 * name <-- file name
1462 * magic_string <-- magic string associated with file type
1463 * method <-- file access method options
1464 * hints <-- MPI-IO hints
1465 * block_comm <-- associated MPI communicator for block IO
1466 * comm <-- associated MPI communicator
1467 *----------------------------------------------------------------------------*/
1468
1469 #if defined(HAVE_MPI)
1470 static void
_cs_io_initialize_with_index(cs_io_t * inp,const char * file_name,const char * magic_string,cs_file_access_t method,MPI_Info hints,MPI_Comm block_comm,MPI_Comm comm)1471 _cs_io_initialize_with_index(cs_io_t *inp,
1472 const char *file_name,
1473 const char *magic_string,
1474 cs_file_access_t method,
1475 MPI_Info hints,
1476 MPI_Comm block_comm,
1477 MPI_Comm comm)
1478 #else
1479 static void
1480 _cs_io_initialize_with_index(cs_io_t *inp,
1481 const char *file_name,
1482 const char *magic_string,
1483 cs_file_access_t method)
1484 #endif /* HAVE_MPI */
1485 {
1486 cs_io_sec_header_t h;
1487 int end_reached = 0;
1488
1489 #if defined(HAVE_MPI)
1490 _file_open(inp,
1491 file_name,
1492 magic_string,
1493 method,
1494 hints,
1495 block_comm,
1496 comm);
1497 #else
1498 _file_open(inp, file_name, magic_string, method);
1499 #endif
1500
1501 /* Read headers to build index index */
1502
1503 while (end_reached == 0) {
1504
1505 end_reached = cs_io_read_header(inp, &h);
1506
1507 if (end_reached == 0)
1508 _update_index_and_shift(inp, &h);
1509
1510 }
1511 }
1512
1513 /*----------------------------------------------------------------------------
1514 * Write padding with zero bytes if necessary to ensure alignment.
1515 *
1516 * Under MPI, data is only written by the associated communicator's root
1517 * rank. The section data on other ranks is ignored, though the file offset
1518 * is updated (i.e. the call to this function is collective).
1519 *
1520 * parameters:
1521 * align <-- required alignment
1522 * outp --> output kernel IO structure
1523 *----------------------------------------------------------------------------*/
1524
1525 static void
_write_padding(size_t align,cs_io_t * outp)1526 _write_padding(size_t align,
1527 cs_io_t *outp)
1528 {
1529 if (align > 0) {
1530
1531 cs_file_off_t offset = cs_file_tell(outp->f);
1532 cs_file_off_t add_offset = (align - (offset % align)) % align;
1533
1534 if (add_offset > 0) {
1535
1536 size_t pad_size = add_offset;
1537 size_t n_written = 0;
1538
1539 if (pad_size > outp->buffer_size) {
1540 while (pad_size > outp->buffer_size)
1541 outp->buffer_size *=2;
1542 BFT_REALLOC(outp->buffer, outp->buffer_size, unsigned char);
1543 }
1544
1545 memset(outp->buffer, 0, pad_size);
1546
1547 n_written = cs_file_write_global(outp->f,
1548 outp->buffer,
1549 1,
1550 pad_size);
1551
1552 if (pad_size != n_written)
1553 bft_error(__FILE__, __LINE__, 0,
1554 _("Error writing %llu bytes to file \"%s\"."),
1555 (unsigned long long)pad_size, cs_file_get_name(outp->f));
1556 }
1557 }
1558 }
1559
1560 /*----------------------------------------------------------------------------
1561 * Write a section header, with possibly embedded data.
1562 *
1563 * Under MPI, data is only written by the associated communicator's root
1564 * rank. The section data on other ranks is ignored, though the file offset
1565 * is updated (i.e. the call to this function is collective).
1566 *
1567 * parameters:
1568 * section_name <-- section name
1569 * n_vals <-- total number of values
1570 * location_id <-- id of associated location, or 0
1571 * index_id <-- id of associated index, or 0
1572 * n_location_vals <-- number of values per location
1573 * elt_type <-- element type
1574 * elts <-- pointer to element data, if it may be embedded
1575 * outp --> output kernel IO structure
1576 *
1577 * returns:
1578 * true if element data is embedded, false otherwise
1579 *----------------------------------------------------------------------------*/
1580
1581 static bool
_write_header(const char * sec_name,cs_gnum_t n_vals,size_t location_id,size_t index_id,size_t n_location_vals,cs_datatype_t elt_type,const void * elts,cs_io_t * outp)1582 _write_header(const char *sec_name,
1583 cs_gnum_t n_vals,
1584 size_t location_id,
1585 size_t index_id,
1586 size_t n_location_vals,
1587 cs_datatype_t elt_type,
1588 const void *elts,
1589 cs_io_t *outp)
1590 {
1591 cs_file_off_t header_vals[6];
1592
1593 double t_start = 0.;
1594 cs_file_off_t write_size = 0;
1595 cs_file_off_t data_size = n_vals * cs_datatype_size[elt_type];
1596 size_t name_size = 0, name_pad_size = 0;
1597 size_t n_written = 0;
1598 cs_io_log_t *log = NULL;
1599
1600 bool embed = false;
1601
1602 assert(outp != NULL);
1603
1604 if (outp->echo >= CS_IO_ECHO_HEADERS)
1605 _echo_pre(outp);
1606
1607 if (outp->log_id > -1) {
1608 log = _cs_io_log[outp->mode] + outp->log_id;
1609 t_start = cs_timer_wtime();
1610 }
1611
1612 _write_padding(outp->header_align, outp); /* Pad if necessary */
1613
1614 /* Prepare header data */
1615 /*---------------------*/
1616
1617 header_vals[0] = 56;
1618 header_vals[1] = n_vals;
1619 header_vals[2] = location_id;
1620 header_vals[3] = index_id;
1621
1622 header_vals[4] = n_location_vals;
1623
1624 name_size = strlen(sec_name);
1625 name_pad_size = 8-(name_size%8); /* At least 1 NULL
1626 character with this rule */
1627
1628 header_vals[5] = name_size + name_pad_size;
1629 header_vals[0] += (name_size + name_pad_size);
1630
1631 /* Decide if data is to be embedded */
1632
1633 if ( n_vals > 0
1634 && elts != NULL
1635 && (header_vals[0] + data_size <= (cs_file_off_t)(outp->header_size))) {
1636 header_vals[0] += data_size;
1637 embed = true;
1638 }
1639
1640 /* Ensure buffer is big enough for data */
1641
1642 if (header_vals[0] > (cs_file_off_t)(outp->buffer_size)) {
1643 while (header_vals[0] > (cs_file_off_t)(outp->buffer_size))
1644 outp->buffer_size *=2;
1645 BFT_REALLOC(outp->buffer, outp->buffer_size, unsigned char);
1646 }
1647
1648 memset(outp->buffer, 0, outp->buffer_size);
1649
1650 /* Now build buffer */
1651
1652 _convert_from_offset(outp->buffer, header_vals, 6);
1653
1654 if (cs_file_get_swap_endian(outp->f) == 1)
1655 _swap_endian(outp->buffer, 8, 6);
1656
1657 /* Element type name */
1658
1659 outp->type_name = (char *)outp->buffer + 48;
1660
1661 switch(elt_type) {
1662 case CS_FLOAT:
1663 outp->type_name[0] = 'r';
1664 if (sizeof(float) == 4)
1665 outp->type_name[1] = '4';
1666 else
1667 outp->type_name[1] = '8';
1668 assert(sizeof(float) == 4 || sizeof(float) == 8);
1669 break;
1670 case CS_DOUBLE:
1671 outp->type_name[0] = 'r';
1672 if (sizeof(double) == 8)
1673 outp->type_name[1] = '8';
1674 else {
1675 outp->type_name[1] = '1';
1676 outp->type_name[2] = '6';
1677 }
1678 assert(sizeof(double) == 8 || sizeof(float) == 16);
1679 break;
1680 case CS_INT32:
1681 outp->type_name[0] = 'i';
1682 outp->type_name[1] = '4';
1683 break;
1684 case CS_INT64:
1685 outp->type_name[0] = 'i';
1686 outp->type_name[1] = '8';
1687 break;
1688 case CS_UINT32:
1689 outp->type_name[0] = 'u';
1690 outp->type_name[1] = '4';
1691 break;
1692 case CS_UINT64:
1693 outp->type_name[0] = 'u';
1694 outp->type_name[1] = '8';
1695 break;
1696 case CS_CHAR:
1697 outp->type_name[0] = 'c';
1698 outp->type_name[1] = ' ';
1699 break;
1700 default:
1701 break;
1702 }
1703
1704 if (embed == true)
1705 outp->type_name[7] = 'e';
1706
1707 /* Section name */
1708
1709 strcpy((char *)(outp->buffer) + 56, sec_name);
1710
1711 if (embed == true) {
1712
1713 unsigned char *data = (unsigned char *)(outp->buffer)
1714 + (56 + name_size + name_pad_size);
1715
1716 memcpy(data, elts, data_size);
1717
1718 if ( cs_file_get_swap_endian(outp->f) == 1
1719 && cs_datatype_size[elt_type] > 1)
1720 _swap_endian(data, cs_datatype_size[elt_type], n_vals);
1721 }
1722
1723 /* Now write header data */
1724
1725 write_size = CS_MAX((cs_file_off_t)(outp->header_size), header_vals[0]);
1726
1727 n_written = cs_file_write_global(outp->f,
1728 outp->buffer,
1729 1,
1730 write_size);
1731
1732 if (write_size != (cs_file_off_t)n_written)
1733 bft_error(__FILE__, __LINE__, 0,
1734 _("Error writing %llu bytes to file \"%s\"."),
1735 (unsigned long long)write_size, cs_file_get_name(outp->f));
1736
1737 if (log != NULL) {
1738 double t_end = cs_timer_wtime();
1739 log->wtimes[0] += t_end - t_start;
1740 log->data_size[0] += write_size;
1741 }
1742
1743 if (outp->echo >= CS_IO_ECHO_HEADERS)
1744 _echo_header(sec_name, n_vals, elt_type);
1745
1746 return embed;
1747 }
1748
1749 /*----------------------------------------------------------------------------
1750 * Dump a kernel IO file handle's metadata.
1751 *
1752 * parameters:
1753 * idx <-- kernel IO index
1754 *----------------------------------------------------------------------------*/
1755
1756 static void
_dump_index(const cs_io_sec_index_t * idx)1757 _dump_index(const cs_io_sec_index_t *idx)
1758 {
1759 size_t ii;
1760
1761 assert(idx != NULL);
1762
1763 bft_printf(_(" %llu indexed records:\n"
1764 " (name, n_vals, location_id, index_id, n_loc_vals, type, "
1765 "embed, file_id, offset)\n\n"),
1766 (unsigned long long)(idx->size));
1767
1768 for (ii = 0; ii < idx->size; ii++) {
1769
1770 char embed = 'n';
1771 cs_file_off_t *h_vals = idx->h_vals + ii*7;
1772 const char *name = idx->names + h_vals[4];
1773
1774 if (h_vals[5] > 0)
1775 embed = 'y';
1776
1777 bft_printf(_(" %40s %10llu %2u %2u %2u %6s %c %2u %ld\n"),
1778 name, (unsigned long long)(h_vals[0]),
1779 (unsigned)(h_vals[1]), (unsigned)(h_vals[2]),
1780 (unsigned)(h_vals[3]), cs_datatype_name[h_vals[6]],
1781 embed, (unsigned)(h_vals[7]),
1782 (long)(idx->offset[ii]));
1783
1784 }
1785
1786 bft_printf("\n");
1787 }
1788
1789 /*! (DOXYGEN_SHOULD_SKIP_THIS) \endcond */
1790
1791 /*============================================================================
1792 * Public function definitions
1793 *============================================================================*/
1794
1795 /*----------------------------------------------------------------------------
1796 * Initialize a kernel IO file structure.
1797 *
1798 * The magic string may be NULL only in read mode, if we choose to ignore it.
1799 *
1800 * parameters:
1801 * name <-- file name
1802 * magic_string <-- magic string associated with file type
1803 * mode <-- read or write
1804 * method <-- file access method
1805 * echo <-- echo on main output (< 0 if none, header if 0,
1806 * n first and last elements if n > 0)
1807 * hints <-- associated hints for MPI-IO, or MPI_INFO_NULL
1808 * block_comm <-- handle to MPI communicator used for distributed file
1809 * block access (may be a subset of comm if some ranks do
1810 * not directly access distributed data blocks)
1811 * comm <-- handle to main MPI communicator
1812 *
1813 * returns:
1814 * pointer to kernel IO structure
1815 *----------------------------------------------------------------------------*/
1816
1817 #if defined(HAVE_MPI)
1818 cs_io_t *
cs_io_initialize(const char * file_name,const char * magic_string,cs_io_mode_t mode,cs_file_access_t method,long echo,MPI_Info hints,MPI_Comm block_comm,MPI_Comm comm)1819 cs_io_initialize(const char *file_name,
1820 const char *magic_string,
1821 cs_io_mode_t mode,
1822 cs_file_access_t method,
1823 long echo,
1824 MPI_Info hints,
1825 MPI_Comm block_comm,
1826 MPI_Comm comm)
1827 #else
1828 cs_io_t *
1829 cs_io_initialize(const char *file_name,
1830 const char *magic_string,
1831 cs_io_mode_t mode,
1832 cs_file_access_t method,
1833 long echo)
1834 #endif /* HAVE_MPI */
1835 {
1836 cs_io_t *cs_io =_cs_io_create(mode, echo);
1837
1838 /* Info on interface creation */
1839
1840 if (echo >= CS_IO_ECHO_OPEN_CLOSE) {
1841 if (mode == CS_IO_MODE_READ)
1842 bft_printf(_("\n Reading file: %s\n"), file_name);
1843 else
1844 bft_printf(_("\n Writing file: %s\n"), file_name);
1845 bft_printf_flush();
1846 }
1847
1848 /* Create interface file descriptor */
1849
1850 #if defined(HAVE_MPI)
1851 _file_open(cs_io, file_name, magic_string, method, hints, block_comm, comm);
1852 #else
1853 _file_open(cs_io, file_name, magic_string, method);
1854 #endif
1855
1856 return cs_io;
1857 }
1858
1859 /*----------------------------------------------------------------------------
1860 * Initialize a kernel IO file structure in read mode, building an index.
1861 *
1862 * The magic string may be NULL, if we choose to ignore it.
1863 *
1864 * parameters:
1865 * name <-- file name
1866 * magic_string <-- magic string associated with file type
1867 * method <-- file access method
1868 * echo <-- echo on main output (< 0 if none, header if 0,
1869 * n first and last elements if n > 0)
1870 * hints <-- associated hints for MPI-IO, or MPI_INFO_NULL
1871 * block_comm <-- handle to MPI communicator used for distributed file
1872 * block access (may be a subset of comm if some ranks do
1873 * not directly access distributed data blocks)
1874 * comm <-- handle to main MPI communicator
1875 *
1876 * returns:
1877 * pointer to kernel IO structure
1878 *----------------------------------------------------------------------------*/
1879
1880 #if defined(HAVE_MPI)
1881 cs_io_t *
cs_io_initialize_with_index(const char * file_name,const char * magic_string,cs_file_access_t method,long echo,MPI_Info hints,MPI_Comm block_comm,MPI_Comm comm)1882 cs_io_initialize_with_index(const char *file_name,
1883 const char *magic_string,
1884 cs_file_access_t method,
1885 long echo,
1886 MPI_Info hints,
1887 MPI_Comm block_comm,
1888 MPI_Comm comm)
1889 #else
1890 cs_io_t *
1891 cs_io_initialize_with_index(const char *file_name,
1892 const char *magic_string,
1893 cs_file_access_t method,
1894 long echo)
1895 #endif /* HAVE_MPI */
1896 {
1897 cs_io_t *inp =_cs_io_create(CS_IO_MODE_READ, echo);
1898
1899 /* Info on interface creation */
1900
1901 if (echo >= CS_IO_ECHO_OPEN_CLOSE) {
1902 bft_printf(_("\n Reading file: %s\n"), file_name);
1903 bft_printf_flush();
1904 }
1905
1906 /* Initialize index */
1907
1908 _create_index(inp);
1909
1910 /* Create interface file descriptor; do not use MPI-IO at this
1911 stage, as we only read global headers of limited size, and
1912 a "lighter" method than MPI-IO should be well adapted. */
1913
1914 cs_file_access_t _method = CS_FILE_STDIO_SERIAL;
1915
1916 #if defined(HAVE_MPI)
1917
1918 MPI_Info _hints = MPI_INFO_NULL;
1919
1920 _cs_io_initialize_with_index(inp,
1921 file_name,
1922 magic_string,
1923 _method,
1924 _hints,
1925 block_comm,
1926 comm);
1927
1928 #else
1929
1930 _cs_io_initialize_with_index(inp, file_name, magic_string, _method);
1931
1932 #endif
1933
1934 /* Now reopen all indexed files using hints */
1935
1936 #if defined(HAVE_MPI)
1937 _file_reopen_read(inp, method, hints, block_comm, comm);
1938 #else
1939 _file_reopen_read(inp, method);
1940 #endif
1941
1942 return inp;
1943 }
1944
1945 /*----------------------------------------------------------------------------
1946 * Free a kernel IO file structure, closing the associated file.
1947 *
1948 * parameters:
1949 * cs_io <-> kernel IO structure
1950 *----------------------------------------------------------------------------*/
1951
1952 void
cs_io_finalize(cs_io_t ** cs_io)1953 cs_io_finalize(cs_io_t **cs_io)
1954 {
1955 cs_io_t *_cs_io = *cs_io;
1956
1957 if(_cs_io->mode == CS_IO_MODE_WRITE)
1958 cs_io_write_global("EOF", 0, 0, 0, 0, CS_DATATYPE_NULL, NULL, _cs_io);
1959
1960 /* Info on closing of interface file */
1961
1962 if (_cs_io->echo >= CS_IO_ECHO_OPEN_CLOSE) {
1963 if (_cs_io->mode == CS_IO_MODE_READ)
1964 bft_printf(_(" Finished reading: %s\n"),
1965 cs_file_get_name(_cs_io->f));
1966 else
1967 bft_printf(_(" Finished writing: %s\n"),
1968 cs_file_get_name(_cs_io->f));
1969 bft_printf_flush();
1970 }
1971
1972 if (_cs_io->index != NULL)
1973 _destroy_index(_cs_io);
1974
1975 _file_close(_cs_io);
1976
1977 _cs_io->buffer_size = 0;
1978 BFT_FREE(_cs_io->buffer);
1979
1980 BFT_FREE(*cs_io);
1981 }
1982
1983 /*----------------------------------------------------------------------------
1984 * Return a pointer to a kernel IO structure's name.
1985 *
1986 * parameters:
1987 * cs_io <-- kernel IO structure
1988 *----------------------------------------------------------------------------*/
1989
1990 const char *
cs_io_get_name(const cs_io_t * cs_io)1991 cs_io_get_name(const cs_io_t *cs_io)
1992 {
1993 assert(cs_io != NULL);
1994
1995 return(cs_file_get_name(cs_io->f));
1996 }
1997
1998 /*----------------------------------------------------------------------------
1999 * Return the number of indexed entries in a kernel IO structure.
2000 *
2001 * parameters:
2002 * inp <-- input kernel IO structure
2003 *
2004 * returns:
2005 * size of index if present, 0 otherwise,
2006 *----------------------------------------------------------------------------*/
2007
2008 size_t
cs_io_get_index_size(const cs_io_t * inp)2009 cs_io_get_index_size(const cs_io_t *inp)
2010 {
2011 size_t retval = 0;
2012
2013 if (inp != NULL && inp->index != NULL)
2014 retval = inp->index->size;
2015
2016 return retval;
2017 }
2018
2019 /*----------------------------------------------------------------------------
2020 * Return the name of an indexed section in a kernel IO structure.
2021 *
2022 * parameters:
2023 * inp <-- input kernel IO structure
2024 * id <-- id of section in index (0 to n-1 numbering)
2025 *
2026 * returns:
2027 * pointer to section name if id in index range, NULL otherwise
2028 *----------------------------------------------------------------------------*/
2029
2030 const char *
cs_io_get_indexed_sec_name(const cs_io_t * inp,size_t id)2031 cs_io_get_indexed_sec_name(const cs_io_t *inp,
2032 size_t id)
2033 {
2034 const char *retval = NULL;
2035
2036 if (inp != NULL && inp->index != NULL) {
2037 if (id < inp->index->size) {
2038 size_t name_id = inp->index->h_vals[7*id + 4];
2039 retval = inp->index->names + name_id;
2040 }
2041 }
2042
2043 return retval;
2044 }
2045
2046 /*----------------------------------------------------------------------------
2047 * Return header data for an indexed section in a kernel IO structure.
2048 *
2049 * parameters:
2050 * inp <-- input kernel IO structure
2051 * id <-- id of section in index (0 to n-1 numbering)
2052 *
2053 * returns:
2054 * section header data (if id not in index range, fields set to zero)
2055 *----------------------------------------------------------------------------*/
2056
2057 cs_io_sec_header_t
cs_io_get_indexed_sec_header(const cs_io_t * inp,size_t id)2058 cs_io_get_indexed_sec_header(const cs_io_t *inp,
2059 size_t id)
2060 {
2061 cs_io_sec_header_t h;
2062
2063 h.sec_name = NULL;
2064
2065 if (inp != NULL && inp->index != NULL) {
2066 if (id < inp->index->size) {
2067
2068 size_t name_id = inp->index->h_vals[7*id + 4];
2069
2070 h.sec_name = inp->index->names + name_id;
2071
2072 h.n_vals = inp->index->h_vals[7*id];
2073 h.location_id = inp->index->h_vals[7*id + 1];
2074 h.index_id = inp->index->h_vals[7*id + 2];
2075 h.n_location_vals = inp->index->h_vals[7*id + 3];
2076 h.type_read = (cs_datatype_t)(inp->index->h_vals[7*id + 6]);
2077 h.elt_type = _type_read_to_elt_type(h.type_read);
2078 }
2079 }
2080
2081 if (h.sec_name == NULL) {
2082 h.n_vals = 0;
2083 h.location_id = 0;
2084 h.index_id = 0;
2085 h.n_location_vals = 0;
2086 h.type_read = CS_DATATYPE_NULL;
2087 h.elt_type = h.type_read;
2088 }
2089
2090 return h;
2091 }
2092
2093 /*----------------------------------------------------------------------------
2094 * Return a kernel IO structure's echo (verbosity) level.
2095 *
2096 * parameters:
2097 * cs_io --> kernel IO structure
2098 *----------------------------------------------------------------------------*/
2099
2100 size_t
cs_io_get_echo(const cs_io_t * cs_io)2101 cs_io_get_echo(const cs_io_t *cs_io)
2102 {
2103 assert(cs_io != NULL);
2104
2105 return (size_t)(cs_io->echo);
2106 }
2107
2108 /*----------------------------------------------------------------------------
2109 * Read a section header.
2110 *
2111 * The header values remain valid until the next section header read
2112 * or the file is closed.
2113 *
2114 * parameters:
2115 * inp --> input kernel IO structure
2116 * header <-- header structure
2117 *
2118 * returns:
2119 * 0 if a header was read, 1 in case of error or end-of-file
2120 *----------------------------------------------------------------------------*/
2121
2122 int
cs_io_read_header(cs_io_t * inp,cs_io_sec_header_t * header)2123 cs_io_read_header(cs_io_t *inp,
2124 cs_io_sec_header_t *header)
2125 {
2126 cs_file_off_t header_vals[6];
2127
2128 double t_start = 0.;
2129 int type_name_error = 0;
2130 cs_io_log_t *log = NULL;
2131 size_t n_read = 0, n_add = 0;
2132
2133 assert(inp != NULL);
2134
2135 if (inp->echo >= CS_IO_ECHO_HEADERS)
2136 _echo_pre(inp);
2137
2138 if (inp->log_id > -1) {
2139 log = _cs_io_log[inp->mode] + inp->log_id;
2140 t_start = cs_timer_wtime();
2141 }
2142
2143 /* Position read pointer if necessary */
2144 /*------------------------------------*/
2145
2146 if (inp->header_align > 0) {
2147 size_t ha = inp->header_align;
2148 cs_file_off_t offset = cs_file_tell(inp->f);
2149 cs_file_off_t add_offset = (ha - (offset % ha)) % ha;
2150 if (add_offset > 0) {
2151 int errcode = 0;
2152 errcode = cs_file_seek(inp->f, add_offset, CS_FILE_SEEK_CUR);
2153 if (errcode != 0)
2154 return 1;
2155 }
2156 }
2157
2158 inp->n_vals = 0;
2159
2160 /* Read header */
2161 /*-------------*/
2162
2163 n_read = cs_file_read_global(inp->f, inp->buffer, 1, inp->header_size);
2164
2165 if (n_read < inp->header_size)
2166 return 1;
2167
2168 if (cs_file_get_swap_endian(inp->f) == 1)
2169 _swap_endian(inp->buffer, 8, 6);
2170
2171 _convert_to_offset(inp->buffer, header_vals, 6);
2172
2173 if (header_vals[0] > (cs_file_off_t)(inp->header_size)) {
2174
2175 n_add = header_vals[0] - inp->header_size;
2176
2177 if (header_vals[0] > (cs_file_off_t)(inp->buffer_size)) {
2178 while (header_vals[0] > (cs_file_off_t)(inp->buffer_size))
2179 inp->buffer_size *=2;
2180 BFT_REALLOC(inp->buffer, inp->buffer_size, unsigned char);
2181 }
2182
2183 n_read = cs_file_read_global(inp->f,
2184 inp->buffer + inp->header_size,
2185 1,
2186 n_add);
2187
2188 if (n_read < n_add)
2189 return 1;
2190 }
2191
2192 /* Set pointers to data fields */
2193
2194 inp->n_vals = header_vals[1];
2195 inp->location_id = header_vals[2];
2196 inp->index_id = header_vals[3];
2197 inp->n_loc_vals = header_vals[4];
2198 inp->type_size = 0;
2199 inp->data = NULL;
2200 inp->type_name = (char *)(inp->buffer + 48);
2201 inp->sec_name = (char *)(inp->buffer + 56);
2202
2203 if (header_vals[1] > 0 && inp->type_name[7] == 'e')
2204 inp->data = inp->buffer + 56 + header_vals[5];
2205
2206 inp->type_size = 0;
2207
2208 /* Return immediately if we have an end-of file marker */
2209
2210 if ((inp->n_vals == 0) && (strcmp(inp->sec_name, "EOF") == 0))
2211 return 1;
2212
2213 if (inp->n_vals > 0) {
2214
2215 /* Check type name and compute size of data */
2216
2217 if (inp->type_name[0] == 'c') {
2218 if (inp->type_name[1] != ' ')
2219 type_name_error = 1;
2220 else
2221 inp->type_size = 1;
2222 }
2223 else if ( inp->type_name[0] == 'i'
2224 || inp->type_name[0] == 'u'
2225 || inp->type_name[0] == 'r') {
2226
2227 if (inp->type_name[1] == '4')
2228 inp->type_size = 4;
2229 else if (inp->type_name[1] == '8')
2230 inp->type_size = 8;
2231 else
2232 type_name_error = 1;
2233
2234 }
2235 else
2236 type_name_error = 1;
2237
2238 if (type_name_error)
2239 bft_error(__FILE__, __LINE__, 0,
2240 _("Type \"%s\" is not known\n"
2241 "Known types: \"c \", \"i4\", \"i8\", \"u4\", \"u8\", "
2242 "\"r4\", \"r8\"."), inp->type_name);
2243
2244 else if ( inp->data != NULL
2245 && cs_file_get_swap_endian(inp->f) == 1 && inp->type_size > 1)
2246 _swap_endian(inp->data, inp->type_size, inp->n_vals);
2247 }
2248
2249 /* Set externally visible header values */
2250 /*--------------------------------------*/
2251
2252 header->sec_name = inp->sec_name;
2253 header->n_vals = inp->n_vals;
2254 header->location_id = inp->location_id;
2255 header->index_id = inp->index_id;
2256 header->n_location_vals = inp->n_loc_vals;
2257
2258 /* Initialize data type */
2259 /*----------------------*/
2260
2261 assert(sizeof(unsigned long) == 8 || sizeof(unsigned long long) == 8);
2262
2263 if (header->n_vals != 0) {
2264
2265 const char *elt_type_name = inp->type_name;
2266
2267 if ( strcmp(elt_type_name, _type_name_i4) == 0
2268 || strcmp(elt_type_name, "i ") == 0)
2269 header->type_read = CS_INT32;
2270
2271 else if (strcmp(elt_type_name, _type_name_i8) == 0)
2272 header->type_read = CS_INT64;
2273
2274 else if (strcmp(elt_type_name, _type_name_u4) == 0)
2275 header->type_read = CS_UINT32;
2276
2277 else if (strcmp(elt_type_name, _type_name_u8) == 0)
2278 header->type_read = CS_UINT64;
2279
2280 else if (strcmp(elt_type_name, _type_name_r4) == 0)
2281 header->type_read = CS_FLOAT;
2282
2283 else if (strcmp(elt_type_name, _type_name_r8) == 0)
2284 header->type_read = CS_DOUBLE;
2285
2286 else if (strcmp(elt_type_name, _type_name_char) == 0)
2287 header->type_read = CS_CHAR;
2288
2289 else
2290 bft_error(__FILE__, __LINE__, 0,
2291 _("Error reading file: \"%s\".\n"
2292 "Data type \"%s\" is not recognized."),
2293 cs_file_get_name(inp->f), elt_type_name);
2294
2295 header->elt_type = _type_read_to_elt_type(header->type_read);
2296 }
2297 else {
2298 header->type_read = CS_DATATYPE_NULL;
2299 header->elt_type = CS_DATATYPE_NULL;
2300 }
2301
2302 /* Possible echo */
2303
2304 if (log != NULL) {
2305 double t_end = cs_timer_wtime();
2306 log->wtimes[0] += t_end - t_start;
2307 log->data_size[0] += (inp->header_size + n_add);
2308 }
2309
2310 if (inp->echo >= CS_IO_ECHO_HEADERS)
2311 _echo_header(header->sec_name,
2312 header->n_vals,
2313 header->type_read);
2314
2315 return 0;
2316 }
2317
2318 /*----------------------------------------------------------------------------
2319 * Set a kernel IO's state so as to be ready to read an indexed section.
2320 *
2321 * The header values and position in the file are set so as to be equivalent
2322 * to those they would have if the corresponding header had just been read.
2323 *
2324 * parameters:
2325 * inp <-> input kernel IO structure
2326 * header --> associated header
2327 * id <-- id of section in index (0 to n-1 numbering)
2328 *
2329 * returns:
2330 * 0 in case of success, 1 in case of error
2331 *----------------------------------------------------------------------------*/
2332
2333 int
cs_io_set_indexed_position(cs_io_t * inp,cs_io_sec_header_t * header,size_t id)2334 cs_io_set_indexed_position(cs_io_t *inp,
2335 cs_io_sec_header_t *header,
2336 size_t id)
2337 {
2338 int retval = 0;
2339
2340 /* Return immediately if out of range */
2341
2342 if (inp == NULL || inp->index == NULL)
2343 return 1;
2344 if (id >= inp->index->size)
2345 return 1;
2346
2347 header->sec_name = inp->index->names + inp->index->h_vals[7*id + 4];
2348
2349 header->n_vals = inp->index->h_vals[7*id];
2350 header->location_id = inp->index->h_vals[7*id + 1];
2351 header->index_id = inp->index->h_vals[7*id + 2];
2352 header->n_location_vals = inp->index->h_vals[7*id + 3];
2353 header->type_read = (cs_datatype_t)(inp->index->h_vals[7*id + 6]);
2354 header->elt_type = _type_read_to_elt_type(header->type_read);
2355
2356 inp->n_vals = header->n_vals;
2357 inp->location_id = header->location_id;
2358 inp->index_id = header->index_id;
2359 inp->n_loc_vals = header->n_location_vals;
2360 inp->type_size = cs_datatype_size[header->type_read];
2361
2362 /* The following values are not taken from the header buffer as
2363 usual, but are base on the index */
2364
2365 strcpy((char *)(inp->buffer + 56), header->sec_name);
2366 inp->sec_name = (char *)(inp->buffer + 56);
2367 inp->type_name = NULL; /* should not be needed now that datatype is known */
2368
2369 /* Non-embedded values */
2370
2371 if (inp->index->h_vals[7*id + 5] == 0) {
2372 cs_file_off_t offset = inp->index->offset[id];
2373 retval = cs_file_seek(inp->f, offset, CS_FILE_SEEK_SET);
2374 }
2375
2376 /* Embedded values */
2377
2378 else {
2379 size_t data_id = inp->index->h_vals[7*id + 5] - 1;
2380 unsigned char *_data = inp->index->data + data_id;
2381 inp->data = _data;
2382 }
2383
2384 return retval;
2385 }
2386
2387 /*----------------------------------------------------------------------------
2388 * Set a section's final data type to int.
2389 *
2390 * It the datatype is not compatible, throw an error.
2391 *
2392 * parameters:
2393 * header <-- header structure
2394 * cs_io --> kernel IO structure
2395 *----------------------------------------------------------------------------*/
2396
2397 void
cs_io_set_int(cs_io_sec_header_t * header,const cs_io_t * cs_io)2398 cs_io_set_int(cs_io_sec_header_t *header,
2399 const cs_io_t *cs_io)
2400 {
2401 assert(header != NULL);
2402
2403 if ( header->type_read != CS_INT32
2404 && header->type_read != CS_INT64
2405 && header->type_read != CS_UINT32
2406 && header->type_read != CS_UINT64)
2407 bft_error(__FILE__, __LINE__, 0,
2408 _("Error reading file: \"%s\".\n"
2409 "Type expected for section: "
2410 "\"%s\" is a signed integer\n"
2411 "and is not convertible from type read: \"%s\"."),
2412 cs_file_get_name(cs_io->f),
2413 header->sec_name,
2414 cs_io->type_name);
2415
2416 assert(sizeof(int) == 4 || sizeof(int) == 8);
2417
2418 if (sizeof(int) == 4)
2419 header->elt_type = CS_INT32;
2420 else
2421 header->elt_type = CS_INT64;
2422 }
2423
2424 /*----------------------------------------------------------------------------
2425 * Set a section's final data type to cs_lnum_t.
2426 *
2427 * It the datatype is not compatible, throw an error.
2428 *
2429 * parameters:
2430 * header <-- header structure
2431 * cs_io --> kernel IO structure
2432 *----------------------------------------------------------------------------*/
2433
2434 void
cs_io_set_cs_lnum(cs_io_sec_header_t * header,const cs_io_t * cs_io)2435 cs_io_set_cs_lnum(cs_io_sec_header_t *header,
2436 const cs_io_t *cs_io)
2437 {
2438 assert(header != NULL);
2439
2440 if ( header->type_read != CS_INT32
2441 && header->type_read != CS_INT64
2442 && header->type_read != CS_UINT32
2443 && header->type_read != CS_UINT64)
2444 bft_error(__FILE__, __LINE__, 0,
2445 _("Error reading file: \"%s\".\n"
2446 "Type expected for section: "
2447 "\"%s\" is a signed integer\n"
2448 "and is not convertible from type read: \"%s\"."),
2449 cs_file_get_name(cs_io->f),
2450 header->sec_name,
2451 cs_io->type_name);
2452
2453 assert(sizeof(cs_lnum_t) == 4 || sizeof(cs_lnum_t) == 8);
2454
2455 if (sizeof(cs_lnum_t) == 4)
2456 header->elt_type = CS_INT32;
2457 else
2458 header->elt_type = CS_INT64;
2459 }
2460
2461 /*----------------------------------------------------------------------------
2462 * Set a section's final data type to cs_gnum_t.
2463 *
2464 * It the datatype is not compatible, throw an error.
2465 *
2466 * parameters:
2467 * header <-- header structure
2468 * cs_io --> kernel IO structure
2469 *----------------------------------------------------------------------------*/
2470
2471 void
cs_io_set_cs_gnum(cs_io_sec_header_t * header,const cs_io_t * cs_io)2472 cs_io_set_cs_gnum(cs_io_sec_header_t *header,
2473 const cs_io_t *cs_io)
2474 {
2475 assert(header != NULL);
2476
2477 if ( header->type_read != CS_INT32
2478 && header->type_read != CS_INT64
2479 && header->type_read != CS_UINT32
2480 && header->type_read != CS_UINT64)
2481 bft_error(__FILE__, __LINE__, 0,
2482 _("Error reading file: \"%s\".\n"
2483 "Type expected for section: "
2484 "\"%s\" is an unsigned integer\n"
2485 "and is not convertible from type read: \"%s\"."),
2486 cs_file_get_name(cs_io->f), header->sec_name,
2487 cs_io->type_name);
2488
2489 assert(sizeof(cs_gnum_t) == 4 || sizeof(cs_gnum_t) == 8);
2490
2491 if (sizeof(cs_gnum_t) == 4)
2492 header->elt_type = CS_UINT32;
2493 else
2494 header->elt_type = CS_UINT64;
2495 }
2496
2497 /*----------------------------------------------------------------------------
2498 * Check that a section's final data type corresponds to cs_real_t.
2499 *
2500 * parameters:
2501 * header <-- header structure
2502 * cs_io --> kernel IO structure
2503 *----------------------------------------------------------------------------*/
2504
2505 void
cs_io_assert_cs_real(const cs_io_sec_header_t * header,const cs_io_t * cs_io)2506 cs_io_assert_cs_real(const cs_io_sec_header_t *header,
2507 const cs_io_t *cs_io)
2508 {
2509 assert(header != NULL);
2510
2511 if ( header->elt_type != CS_FLOAT
2512 && header->elt_type != CS_DOUBLE)
2513 bft_error(__FILE__, __LINE__, 0,
2514 _("Error reading file: \"%s\".\n"
2515 "Type expected for section: \"%s\"\n"
2516 "is \"r4\" or \"r8\" (real), and not \"%s\"."),
2517 cs_file_get_name(cs_io->f),
2518 header->sec_name,
2519 cs_io->type_name);
2520 }
2521
2522 /*----------------------------------------------------------------------------
2523 * Read a section body and replicate it to all processors.
2524 *
2525 * If the array intended to receive the data already exists, we pass an
2526 * "elt" pointer to this array; this same pointer is then returned.
2527 * Otherwise, if this pointer is passed as NULL, memory is allocated
2528 * by this function, and the corresponding pointer is returned. It is
2529 * the caller's responsibility to free this array.
2530 *
2531 * parameters:
2532 * header <-- header structure
2533 * elts <-> pointer to data array, or NULL
2534 * cs_io --> kernel IO structure
2535 *
2536 * returns:
2537 * elts if non NULL, or pointer to allocated array otherwise
2538 *----------------------------------------------------------------------------*/
2539
2540 void *
cs_io_read_global(const cs_io_sec_header_t * header,void * elts,cs_io_t * cs_io)2541 cs_io_read_global(const cs_io_sec_header_t *header,
2542 void *elts,
2543 cs_io_t *cs_io)
2544 {
2545 return _cs_io_read_body(header, 0, 0, elts, cs_io);
2546 }
2547
2548 /*----------------------------------------------------------------------------
2549 * Read a section body, assigning a different block to each processor.
2550 *
2551 * If location_id > 0 and header->n_location_vals > 1, then
2552 * global_num_start and global_num_end will be based on location element
2553 * numbers, so the total number of values read equals
2554 * (global_num_end - global_num_start) * header->n_location_vals.
2555 *
2556 * If the array intended to receive the data already exists, we pass an
2557 * "elt" pointer to this array; this same pointer is then returned.
2558 * Otherwise, if this pointer is passed as NULL, memory is allocated
2559 * by this function, and the corresponding pointer is returned. It is
2560 * the caller's responsibility to free this array.
2561 *
2562 * parameters:
2563 * header <-- header structure
2564 * global_num_start <-- global number of first block item (1 to n numbering)
2565 * global_num_end <-- global number of past-the end block item
2566 * (1 to n numbering)
2567 * elts <-> pointer to data array, or NULL
2568 * cs_io --> kernel IO structure
2569 *
2570 * returns:
2571 * elts if non NULL, or pointer to allocated array otherwise
2572 *----------------------------------------------------------------------------*/
2573
2574 void *
cs_io_read_block(const cs_io_sec_header_t * header,cs_gnum_t global_num_start,cs_gnum_t global_num_end,void * elts,cs_io_t * cs_io)2575 cs_io_read_block(const cs_io_sec_header_t *header,
2576 cs_gnum_t global_num_start,
2577 cs_gnum_t global_num_end,
2578 void *elts,
2579 cs_io_t *cs_io)
2580 {
2581 assert(global_num_start > 0);
2582 assert(global_num_end >= global_num_start);
2583
2584 return _cs_io_read_body(header,
2585 global_num_start,
2586 global_num_end,
2587 elts,
2588 cs_io);
2589 }
2590
2591 /*----------------------------------------------------------------------------
2592 * Read a section body, assigning a different block to each processor,
2593 * when the body corresponds to an index.
2594 *
2595 * In serial mode, this function behaves just like cs_io_read_block(),
2596 * except that it allows only unsigned integer values (cs_gnum_t).
2597 *
2598 * In parallel mode, global_num_end should be set to the past-the-end value
2599 * of the base data block, the same as for regular data (and not increased
2600 * by 1 for the last rank, as this will be handled internally).
2601 * On each rank, the buffer size should be:
2602 * global_num_end - global_num_start + 1, as the past-the end index
2603 * for the local block is added automatically.
2604 *
2605 * If the array intended to receive the data already exists, we pass an
2606 * "elt" pointer to this array; this same pointer is then returned.
2607 * Otherwise, if this pointer is passed as NULL, memory is allocated
2608 * by this function, and the corresponding pointer is returned. It is
2609 * the caller's responsibility to free this array.
2610 *
2611 * parameters:
2612 * header <-- header structure
2613 * global_num_start <-- global number of first block item (1 to n numbering)
2614 * global_num_end <-- global number of past-the end block item
2615 * (1 to n numbering)
2616 * elts <-> pointer to data array, or NULL
2617 * cs_io --> kernel IO structure
2618 *
2619 * returns:
2620 * elts if non NULL, or pointer to allocated array otherwise
2621 *----------------------------------------------------------------------------*/
2622
2623 void *
cs_io_read_index_block(cs_io_sec_header_t * header,cs_gnum_t global_num_start,cs_gnum_t global_num_end,cs_gnum_t * elts,cs_io_t * cs_io)2624 cs_io_read_index_block(cs_io_sec_header_t *header,
2625 cs_gnum_t global_num_start,
2626 cs_gnum_t global_num_end,
2627 cs_gnum_t *elts,
2628 cs_io_t *cs_io)
2629 {
2630 cs_gnum_t _global_num_start = global_num_start;
2631 cs_gnum_t _global_num_end = global_num_end;
2632
2633 cs_gnum_t *retval = NULL;
2634
2635 #if defined(HAVE_MPI)
2636 int rank_id = 0;
2637 int n_ranks = 1;
2638 MPI_Comm comm = cs_io->comm;
2639
2640 if (comm != MPI_COMM_NULL) {
2641 MPI_Comm_rank(comm, &rank_id);
2642 MPI_Comm_size(comm, &n_ranks);
2643 }
2644 #endif
2645
2646 assert(global_num_start > 0);
2647 assert(global_num_end >= global_num_start);
2648
2649 /* Check type */
2650
2651 cs_io_set_cs_gnum(header, cs_io);
2652
2653 /* Increase _global_num_end by 1 for the last rank containing data */
2654
2655 if (global_num_end == (cs_gnum_t)header->n_vals) {
2656
2657 if (global_num_start < global_num_end)
2658 _global_num_end += 1;
2659
2660 /* Also shift start values for possibly empty
2661 blocks past the last rank reading data */
2662
2663 else {
2664 _global_num_start += 1;
2665 _global_num_end += 1;
2666 }
2667
2668 }
2669
2670 retval = _cs_io_read_body(header,
2671 _global_num_start,
2672 _global_num_end,
2673 elts,
2674 cs_io);
2675
2676 /* Ensure element value initialized in case of empty block */
2677
2678 if (retval == NULL)
2679 BFT_MALLOC(retval, 1, cs_gnum_t);
2680
2681 if (_global_num_start == _global_num_end)
2682 retval[0] = 0;
2683
2684 /* Exchange past-the-end values */
2685
2686 #if defined(HAVE_MPI)
2687
2688 if (n_ranks > 1) {
2689
2690 cs_gnum_t past_last_max = 0;
2691 cs_gnum_t past_last_max_0 = 0;
2692 cs_gnum_t past_last = 0;
2693 cs_gnum_t *past_last_0 = NULL;
2694
2695 if ( _global_num_end > global_num_end
2696 && _global_num_end > _global_num_start)
2697 past_last_max = retval[_global_num_end - _global_num_start - 1];
2698
2699 MPI_Reduce(&past_last_max, &past_last_max_0, 1, CS_MPI_GNUM, MPI_MAX,
2700 0, comm);
2701
2702 /* Initially, past_last values contain the first index value
2703 for a given rank; they will later be shifted to define the
2704 past-the-last value for the previous rank) */
2705
2706 if (retval != NULL)
2707 past_last = retval[0];
2708
2709 if (rank_id == 0)
2710 BFT_MALLOC(past_last_0, n_ranks, cs_gnum_t);
2711
2712 MPI_Gather(&past_last, 1, CS_MPI_GNUM,
2713 past_last_0, 1, CS_MPI_GNUM,
2714 0, comm);
2715
2716 if (rank_id == 0) {
2717
2718 int i;
2719 int last_data_rank = n_ranks - 1;
2720
2721 while (last_data_rank > 0 && past_last_0[last_data_rank] == 0)
2722 last_data_rank -= 1;
2723
2724 /* fill empty values before last rank with data */
2725 for (i = last_data_rank; i > 0; i--) {
2726 if (past_last_0[i-1] == 0)
2727 past_last_0[i-1] = past_last_0[i];
2728 }
2729
2730 /* Shift values one rank (to really define past-the-last values) */
2731
2732 for (i = 0; i < last_data_rank; i++)
2733 past_last_0[i] = past_last_0[i+1];
2734
2735 /* Define for all other ranks from last rank with data onwards */
2736
2737 for (i = last_data_rank; i < n_ranks; i++)
2738 past_last_0[i] = past_last_max_0;
2739 }
2740
2741 MPI_Scatter(past_last_0, 1, CS_MPI_GNUM,
2742 &past_last, 1, CS_MPI_GNUM,
2743 0, comm);
2744
2745 if (rank_id == 0)
2746 BFT_FREE(past_last_0);
2747
2748 if (retval != NULL)
2749 retval[global_num_end - global_num_start] = past_last;
2750 }
2751
2752 if ( retval != NULL
2753 && header->n_vals != 0 && (cs_gnum_t)header->n_vals != global_num_end
2754 && cs_io->echo > CS_IO_ECHO_HEADERS)
2755 bft_printf(_(" first element for next rank:\n"
2756 " %10llu : %12llu\n"),
2757 (unsigned long long)(global_num_end),
2758 (unsigned long long)retval[global_num_end - global_num_start]);
2759
2760 #endif /* defined(HAVE_MPI) */
2761
2762 return retval;
2763 }
2764
2765 /*----------------------------------------------------------------------------
2766 * Write a global section.
2767 *
2768 * Under MPI, data is only written by the associated communicator's root
2769 * rank. The section data on other ranks is ignored, though the file offset
2770 * is updated (i.e. the call to this function is collective).
2771 *
2772 * parameters:
2773 * section_name <-- section name
2774 * n_vals <-- total number of values
2775 * location_id <-- id of associated location, or 0
2776 * index_id <-- id of associated index, or 0
2777 * n_location_vals <-- number of values per location
2778 * elt_type <-- element type
2779 * elts <-- pointer to element data
2780 * outp <-> output kernel IO structure
2781 *----------------------------------------------------------------------------*/
2782
2783 void
cs_io_write_global(const char * sec_name,cs_gnum_t n_vals,size_t location_id,size_t index_id,size_t n_location_vals,cs_datatype_t elt_type,const void * elts,cs_io_t * outp)2784 cs_io_write_global(const char *sec_name,
2785 cs_gnum_t n_vals,
2786 size_t location_id,
2787 size_t index_id,
2788 size_t n_location_vals,
2789 cs_datatype_t elt_type,
2790 const void *elts,
2791 cs_io_t *outp)
2792 {
2793 bool embed = false;
2794
2795 if (outp->echo >= CS_IO_ECHO_HEADERS)
2796 _echo_header(sec_name, n_vals, elt_type);
2797
2798 embed = _write_header(sec_name,
2799 n_vals,
2800 location_id,
2801 index_id,
2802 n_location_vals,
2803 elt_type,
2804 elts,
2805 outp);
2806
2807 if (n_vals > 0 && embed == false) {
2808
2809 double t_start = 0.;
2810 cs_io_log_t *log = NULL;
2811 size_t n_written = 0;
2812
2813 if (outp->log_id > -1) {
2814 log = _cs_io_log[outp->mode] + outp->log_id;
2815 t_start = cs_timer_wtime();
2816 }
2817
2818 _write_padding(outp->body_align, outp);
2819
2820 n_written = cs_file_write_global(outp->f,
2821 elts,
2822 cs_datatype_size[elt_type],
2823 n_vals);
2824
2825 if (n_vals != (cs_gnum_t)n_written)
2826 bft_error(__FILE__, __LINE__, 0,
2827 _("Error writing %llu bytes to file \"%s\"."),
2828 (unsigned long long)n_vals, cs_file_get_name(outp->f));
2829
2830 if (log != NULL) {
2831 double t_end = cs_timer_wtime();
2832 log->wtimes[0] += t_end - t_start;
2833 log->data_size[0] += n_written*cs_datatype_size[elt_type];
2834 }
2835 }
2836
2837 if (n_vals != 0 && outp->echo > CS_IO_ECHO_HEADERS)
2838 _echo_data(outp->echo, n_vals, 1, n_vals + 1, elt_type, elts);
2839 }
2840
2841 /*----------------------------------------------------------------------------
2842 * Write a section to file, each associated process providing a contiguous
2843 * of the section's body.
2844 *
2845 * Each process should provide a (possibly empty) block of the body,
2846 * and we should have:
2847 * global_num_start at rank 0 = 1
2848 * global_num_start at rank i+1 = global_num_end at rank i.
2849 * Otherwise, behavior (especially positioning for future reads) is undefined.
2850 *
2851 * If location_id > 0 and n_location_vals > 1, then global_num_start
2852 * and global_num_end will be based on location element numbers, so the
2853 * total number of values read equals
2854 * (global_num_end - global_num_start) * header->n_location_vals.
2855 *
2856 * This function does not modify the values in its input buffer (notably,
2857 * a copy is used to convert from little-endian to big-endian or vice-versa
2858 * if necessary).
2859 *
2860 * parameters:
2861 * section_name <-- section name
2862 * n_g_elts <-- number of global elements (locations)
2863 * global_num_start <-- global number of first block item (1 to n numbering)
2864 * global_num_end <-- global number of past-the end block item
2865 * location_id <-- id of associated location, or 0
2866 * index_id <-- id of associated index, or 0
2867 * n_location_vals <-- number of values per location
2868 * elt_type <-- element type
2869 * (1 to n numbering)
2870 * elts <-- pointer to element data
2871 * outp <-> output kernel IO structure
2872 *----------------------------------------------------------------------------*/
2873
2874 void
cs_io_write_block(const char * sec_name,cs_gnum_t n_g_elts,cs_gnum_t global_num_start,cs_gnum_t global_num_end,size_t location_id,size_t index_id,size_t n_location_vals,cs_datatype_t elt_type,const void * elts,cs_io_t * outp)2875 cs_io_write_block(const char *sec_name,
2876 cs_gnum_t n_g_elts,
2877 cs_gnum_t global_num_start,
2878 cs_gnum_t global_num_end,
2879 size_t location_id,
2880 size_t index_id,
2881 size_t n_location_vals,
2882 cs_datatype_t elt_type,
2883 const void *elts,
2884 cs_io_t *outp)
2885 {
2886 double t_start = 0.;
2887 size_t n_written = 0;
2888 size_t n_g_vals = n_g_elts;
2889 size_t n_vals = global_num_end - global_num_start;
2890 size_t stride = 1;
2891 cs_io_log_t *log = NULL;
2892
2893 if (n_location_vals > 1) {
2894 stride = n_location_vals;
2895 n_g_vals *= n_location_vals;
2896 n_vals *= n_location_vals;
2897 }
2898
2899 _write_header(sec_name,
2900 n_g_vals,
2901 location_id,
2902 index_id,
2903 n_location_vals,
2904 elt_type,
2905 NULL,
2906 outp);
2907
2908 if (outp->log_id > -1) {
2909 log = _cs_io_log[outp->mode] + outp->log_id;
2910 t_start = cs_timer_wtime();
2911 }
2912
2913 _write_padding(outp->body_align, outp);
2914
2915 n_written = cs_file_write_block(outp->f,
2916 elts,
2917 cs_datatype_size[elt_type],
2918 stride,
2919 global_num_start,
2920 global_num_end);
2921
2922 if (n_vals != (cs_gnum_t)n_written)
2923 bft_error(__FILE__, __LINE__, 0,
2924 _("Error writing %llu bytes to file \"%s\"."),
2925 (unsigned long long)n_vals, cs_file_get_name(outp->f));
2926
2927 if (log != NULL) {
2928 double t_end = cs_timer_wtime();
2929 log->wtimes[1] += t_end - t_start;
2930 log->data_size[1] += n_written*cs_datatype_size[elt_type];
2931 }
2932
2933 if (n_vals != 0 && outp->echo > CS_IO_ECHO_HEADERS)
2934 _echo_data(outp->echo, n_g_vals,
2935 (global_num_start-1)*stride + 1,
2936 (global_num_end -1)*stride + 1,
2937 elt_type, elts);
2938 }
2939
2940 /*----------------------------------------------------------------------------
2941 * Write a section to file, each associated process providing a contiguous
2942 * of the section's body.
2943 *
2944 * Each process should provide a (possibly empty) block of the body,
2945 * and we should have:
2946 * global_num_start at rank 0 = 1
2947 * global_num_start at rank i+1 = global_num_end at rank i.
2948 * Otherwise, behavior (especially positioning for future reads) is undefined.
2949 *
2950 * If location_id > 0 and n_location_vals > 1, then global_num_start
2951 * and global_num_end will be based on location element numbers, so the
2952 * total number of values read equals
2953 * (global_num_end - global_num_start) * header->n_location_vals.
2954 *
2955 * This function is intended to be used mainly on data that is already of
2956 * copy of original data (such as data that has been redistributed across
2957 * processors just for the sake of output), or that is to be deleted after
2958 * writing, so it may modify the values in its input buffer (notably to
2959 * convert from little-endian to big-endian or vice-versa if necessary).
2960 *
2961 * parameters:
2962 * section_name <-- section name
2963 * n_g_elts <-- number of global elements (locations)
2964 * global_num_start <-- global number of first block item (1 to n numbering)
2965 * global_num_end <-- global number of past-the end block item
2966 * location_id <-- id of associated location, or 0
2967 * index_id <-- id of associated index, or 0
2968 * n_location_vals <-- number of values per location
2969 * elt_type <-- element type
2970 * (1 to n numbering)
2971 * elts <-- pointer to element data
2972 * outp <-> output kernel IO structure
2973 *----------------------------------------------------------------------------*/
2974
2975 void
cs_io_write_block_buffer(const char * sec_name,cs_gnum_t n_g_elts,cs_gnum_t global_num_start,cs_gnum_t global_num_end,size_t location_id,size_t index_id,size_t n_location_vals,cs_datatype_t elt_type,void * elts,cs_io_t * outp)2976 cs_io_write_block_buffer(const char *sec_name,
2977 cs_gnum_t n_g_elts,
2978 cs_gnum_t global_num_start,
2979 cs_gnum_t global_num_end,
2980 size_t location_id,
2981 size_t index_id,
2982 size_t n_location_vals,
2983 cs_datatype_t elt_type,
2984 void *elts,
2985 cs_io_t *outp)
2986 {
2987 double t_start = 0.;
2988 size_t n_written = 0;
2989 size_t n_g_vals = n_g_elts;
2990 size_t n_vals = global_num_end - global_num_start;
2991 size_t stride = 1;
2992 cs_io_log_t *log = NULL;
2993
2994 if (n_location_vals > 1) {
2995 stride = n_location_vals;
2996 n_g_vals *= n_location_vals;
2997 n_vals *= n_location_vals;
2998 }
2999
3000 _write_header(sec_name,
3001 n_g_vals,
3002 location_id,
3003 index_id,
3004 n_location_vals,
3005 elt_type,
3006 NULL,
3007 outp);
3008
3009 if (outp->log_id > -1) {
3010 log = _cs_io_log[outp->mode] + outp->log_id;
3011 t_start = cs_timer_wtime();
3012 }
3013
3014 _write_padding(outp->body_align, outp);
3015
3016 n_written = cs_file_write_block_buffer(outp->f,
3017 elts,
3018 cs_datatype_size[elt_type],
3019 stride,
3020 global_num_start,
3021 global_num_end);
3022
3023 if (n_vals != (cs_gnum_t)n_written)
3024 bft_error(__FILE__, __LINE__, 0,
3025 _("Error writing %llu bytes to file \"%s\"."),
3026 (unsigned long long)n_vals, cs_file_get_name(outp->f));
3027
3028 if (log != NULL) {
3029 double t_end = cs_timer_wtime();
3030 log->wtimes[1] += t_end - t_start;
3031 log->data_size[1] += n_written*cs_datatype_size[elt_type];
3032 }
3033
3034 if (n_vals != 0 && outp->echo > CS_IO_ECHO_HEADERS)
3035 _echo_data(outp->echo, n_g_vals,
3036 (global_num_start-1)*stride + 1,
3037 (global_num_end -1)*stride + 1,
3038 elt_type, elts);
3039 }
3040
3041 /*----------------------------------------------------------------------------
3042 * Skip a message.
3043 *
3044 * parameters:
3045 * header <-- header structure
3046 * pp_io --> kernel IO structure
3047 *----------------------------------------------------------------------------*/
3048
3049 void
cs_io_skip(const cs_io_sec_header_t * header,cs_io_t * pp_io)3050 cs_io_skip(const cs_io_sec_header_t *header,
3051 cs_io_t *pp_io)
3052 {
3053 double t_start = 0.;
3054 size_t type_size = 0;
3055 cs_file_off_t n_vals = pp_io->n_vals;
3056 cs_io_log_t *log = NULL;
3057
3058 assert(pp_io != NULL);
3059 assert(header->n_vals == pp_io->n_vals);
3060
3061 if (pp_io->log_id > -1) {
3062 log = _cs_io_log[pp_io->mode] + pp_io->log_id;
3063 t_start = cs_timer_wtime();
3064 }
3065
3066 /* Datatype size given by FVM datatype */
3067
3068 type_size = cs_datatype_size[header->type_read];
3069
3070 /* If data is present in file, skip it */
3071
3072 if (pp_io->data == NULL) {
3073
3074 /* Position read pointer if necessary */
3075
3076 if (pp_io->body_align > 0) {
3077 cs_file_off_t offset = cs_file_tell(pp_io->f);
3078 size_t ba = pp_io->body_align;
3079 offset += (ba - (offset % ba)) % ba;
3080 offset += n_vals*type_size;
3081 cs_file_seek(pp_io->f, offset, CS_FILE_SEEK_SET);
3082 }
3083
3084 pp_io->data = NULL; /* Reset for next read */
3085 }
3086
3087 if (log != NULL) {
3088 double t_end = cs_timer_wtime();
3089 log->wtimes[0] += t_end - t_start;
3090 }
3091 }
3092
3093 /*----------------------------------------------------------------------------
3094 * Return the position of the file pointer for an open kernel IO file.
3095 *
3096 * parameters:
3097 * inp <-- input kernel IO structure
3098 *
3099 * returns:
3100 * offset in file
3101 *----------------------------------------------------------------------------*/
3102
3103 cs_file_off_t
cs_io_get_offset(cs_io_t * inp)3104 cs_io_get_offset(cs_io_t *inp)
3105 {
3106 assert(inp != NULL);
3107 assert(inp->f != NULL);
3108
3109 return cs_file_tell(inp->f);
3110 }
3111
3112 /*----------------------------------------------------------------------------
3113 * Set the position of the file pointer for an open kernel IO file.
3114 *
3115 * parameters:
3116 * inp <-- input kernel IO structure
3117 * offset <-- offset in file
3118 *----------------------------------------------------------------------------*/
3119
3120 void
cs_io_set_offset(cs_io_t * inp,cs_file_off_t offset)3121 cs_io_set_offset(cs_io_t *inp,
3122 cs_file_off_t offset)
3123 {
3124 assert(inp != NULL);
3125 assert(inp->f != NULL);
3126
3127 cs_file_seek(inp->f,
3128 offset,
3129 CS_FILE_SEEK_SET);
3130 }
3131
3132 /*----------------------------------------------------------------------------
3133 * Initialize performance logging for cs_io_t structures.
3134 *----------------------------------------------------------------------------*/
3135
3136 void
cs_io_log_initialize(void)3137 cs_io_log_initialize(void)
3138 {
3139 int i = 0;
3140
3141 for (i = 0; i < 2; i++) {
3142 _cs_io_map_size[i] = 0;
3143 _cs_io_map_size_max[i] = 1;
3144 _cs_io_map[i] = cs_map_name_to_id_create();
3145 BFT_MALLOC(_cs_io_log[i], _cs_io_map_size_max[i], cs_io_log_t);
3146 }
3147 }
3148
3149 /*----------------------------------------------------------------------------
3150 * Finalize performance logging for cs_io_t structures.
3151 *----------------------------------------------------------------------------*/
3152
3153 void
cs_io_log_finalize(void)3154 cs_io_log_finalize(void)
3155 {
3156 int i;
3157
3158 const char unit[8] = {'K', 'M', 'G', 'T', 'P', 'E', 'Z', 'Y'};
3159
3160 for (i = 0; i < 2; i++) {
3161
3162 size_t j;
3163 size_t map_size = cs_map_name_to_id_size(_cs_io_map[i]);
3164
3165 /* Print info */
3166
3167 if (map_size > 0) {
3168 if (i == 0)
3169 cs_log_printf(CS_LOG_PERFORMANCE,
3170 _("\nCode_Saturne IO files read:\n\n"));
3171 else
3172 cs_log_printf(CS_LOG_PERFORMANCE,
3173 _("\nCode_Saturne IO files written:\n\n"));
3174 }
3175
3176 for (j = 0; j < map_size; j++) {
3177
3178 const char *key = cs_map_name_to_id_key(_cs_io_map[i], j);
3179 cs_io_log_t *log = ( _cs_io_log[i]
3180 + cs_map_name_to_id(_cs_io_map[i], key));
3181
3182 #if defined(HAVE_MPI)
3183
3184 if (cs_glob_n_ranks > 1) {
3185
3186 int k, l;
3187 double _wtimes[3];
3188 double _data_size[2];
3189 memcpy(_wtimes, log->wtimes, 3*sizeof(double));
3190 int _data_mult[2] = {0, 0};
3191 unsigned long long data_size_loc = log->data_size[1];
3192
3193 MPI_Allreduce(_wtimes, log->wtimes, 3, MPI_DOUBLE, MPI_MAX,
3194 cs_glob_mpi_comm);
3195 #if defined(MPI_UNSIGNED_LONG_LONG) && !defined(MSMPI_VER) /* By-pass MS-MPI */
3196 MPI_Allreduce(&data_size_loc, log->data_size + 1, 1,
3197 MPI_UNSIGNED_LONG_LONG, MPI_SUM, cs_glob_mpi_comm);
3198 #elif defined(MPI_LONG_LONG)
3199 MPI_Allreduce(&data_size_loc, log->data_size + 1, 1,
3200 MPI_LONG_LONG, MPI_SUM, cs_glob_mpi_comm);
3201 #endif
3202 for (k = 0; k < 2; k++) {
3203 _data_size[k] = (double)(log->data_size[k]) / 1024.;
3204 for (l = 0; _data_size[k] > 1024. && l < 8; l++)
3205 _data_size[k] /= 1024.;
3206 _data_mult[k] = l;
3207 }
3208
3209 cs_log_printf(CS_LOG_PERFORMANCE,
3210 _(" %s\n"
3211 " global: %12.5f s, %12.3f %ciB\n"
3212 " local: %12.5f s, %12.3f %ciB\n"
3213 " open: %12.5f s, %u open(s)\n"),
3214 key,
3215 log->wtimes[0], _data_size[0], unit[_data_mult[0]],
3216 log->wtimes[1], _data_size[1], unit[_data_mult[1]],
3217 log->wtimes[2], log->n_opens);
3218 }
3219 #endif
3220
3221 if (cs_glob_n_ranks == 1) {
3222
3223 int k = 0;
3224 double _data_size = (double)(log->data_size[0] + log->data_size[1])
3225 / 1024.;
3226
3227 for (k = 0; _data_size > 1024. && k < 8; k++)
3228 _data_size /= 1024.;
3229
3230 cs_log_printf(CS_LOG_PERFORMANCE,
3231 _(" %s\n"
3232 " data: %12.5f s, %12.3f %ciB\n"
3233 " open: %12.5f s, %u open(s)\n"),
3234 key,
3235 log->wtimes[0] + log->wtimes[1], _data_size, unit[k],
3236 log->wtimes[2], log->n_opens);
3237 }
3238 }
3239
3240 /* Free structures */
3241
3242 _cs_io_map_size[i] = 0;
3243 _cs_io_map_size_max[i] = 0;
3244 cs_map_name_to_id_destroy(&(_cs_io_map[i]));
3245 BFT_FREE(_cs_io_log[i]);
3246 }
3247
3248 cs_log_printf(CS_LOG_PERFORMANCE, "\n");
3249 cs_log_separator(CS_LOG_PERFORMANCE);
3250 }
3251
3252 /*----------------------------------------------------------------------------
3253 * Dump a kernel IO file handle's metadata.
3254 *
3255 * parameters:
3256 * cs_io <-- kernel IO structure
3257 *----------------------------------------------------------------------------*/
3258
3259 void
cs_io_dump(const cs_io_t * cs_io)3260 cs_io_dump(const cs_io_t *cs_io)
3261 {
3262 assert(cs_io != NULL);
3263
3264 bft_printf(_("\n\n file contents:\n\n"));
3265
3266 if (cs_io->f != NULL)
3267 bft_printf(_(" file: %s\n"), cs_file_get_name(cs_io->f));
3268
3269 bft_printf(_(" contents: \"%s\"\n"), cs_io->contents);
3270 if (cs_io->mode == CS_IO_MODE_READ)
3271 bft_printf(_(" mode: CS_IO_MODE_READ\n"));
3272 else if (cs_io->mode == CS_IO_MODE_WRITE)
3273 bft_printf(_(" mode: CS_IO_MODE_WRITE\n"));
3274
3275 #if defined(HAVE_MPI)
3276 bft_printf(_(" MPI communicator: %ld\n"), (long)(cs_io->comm));
3277 #endif
3278
3279 bft_printf(_(" default header size: %lu\n"
3280 " header alignment: %lu\n"
3281 " body alignment: %lu\n"
3282 " verbosity level: %ld\n\n"),
3283 cs_io->header_size, cs_io->header_align, cs_io->body_align,
3284 cs_io->echo);
3285
3286 if (cs_io->index != NULL)
3287 _dump_index(cs_io->index);
3288 }
3289
3290 /*----------------------------------------------------------------------------*/
3291
3292 END_C_DECLS
3293