1 /*============================================================================
2 * Manage checkpoint / 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_FCNTL_H)
40 # include <fcntl.h>
41 #endif
42
43 #if defined(HAVE_UNISTD_H)
44 # include <unistd.h>
45 #endif
46
47 #if defined(HAVE_MPI)
48 #include <mpi.h>
49 #endif
50
51 /*----------------------------------------------------------------------------
52 * Local headers
53 *----------------------------------------------------------------------------*/
54
55 #include "bft_mem.h"
56 #include "bft_error.h"
57 #include "bft_printf.h"
58
59 #include "fvm_io_num.h"
60
61 #include "cs_all_to_all.h"
62 #include "cs_base.h"
63 #include "cs_block_dist.h"
64 #include "cs_block_to_part.h"
65 #include "cs_file.h"
66 #include "cs_io.h"
67 #include "cs_mesh.h"
68 #include "cs_mesh_save.h"
69 #include "cs_mesh_location.h"
70 #include "cs_part_to_block.h"
71 #include "cs_parall.h"
72 #include "cs_timer.h"
73 #include "cs_time_step.h"
74
75 /*----------------------------------------------------------------------------
76 * Header for the current file
77 *----------------------------------------------------------------------------*/
78
79 #include "cs_restart.h"
80
81 /*----------------------------------------------------------------------------*/
82
83 BEGIN_C_DECLS
84
85 /*=============================================================================
86 * Additional doxygen documentation
87 *============================================================================*/
88
89 /*!
90 \file cs_restart.c
91 Manage checkpoint / restart files.
92 */
93
94 /*! \cond DOXYGEN_SHOULD_SKIP_THIS */
95
96 /*============================================================================
97 * Local macro definitions
98 *============================================================================*/
99
100 /*============================================================================
101 * Local type definitions
102 *============================================================================*/
103
104 typedef struct _location_t {
105
106 char *name; /* Location name */
107 size_t id; /* Associated id in file */
108 cs_lnum_t n_ents; /* Local number of entities */
109 cs_gnum_t n_glob_ents_f; /* Global number of entities by file */
110 cs_gnum_t n_glob_ents; /* Global number of entities */
111 const cs_gnum_t *ent_global_num; /* Possibly shared global entity numbers,
112 or NULL */
113 cs_gnum_t *_ent_global_num; /* Private global entity numbers,
114 or NULL */
115
116 } _location_t;
117
118 struct _cs_restart_t {
119
120 char *name; /* Name of restart file */
121
122 cs_io_t *fh; /* Pointer to associated file handle */
123 int rank_step; /* Block rank step for parallel IO */
124 int min_block_size; /* Minimum block size for parallel IO */
125
126 size_t n_locations; /* Number of locations */
127 _location_t *location; /* Location definition array */
128
129 cs_restart_mode_t mode; /* Read or write */
130
131 };
132
133 typedef struct {
134
135 int id; /* Id of the writer */
136
137 char *name; /* Name of the checkpoint file */
138 char *path; /* Full path to the checkpoint file */
139
140 int n_prev_files; /* Number of times this file has already
141 been written */
142 int n_prev_files_tot; /* Total number of times this file has already
143 been written */
144 char **prev_files; /* Names of the previous versions */
145
146 } _restart_multiwriter_t;
147
148 /*============================================================================
149 * Prototypes for private functions
150 *============================================================================*/
151
152 static int
153 _check_section(cs_restart_t *restart,
154 void *context,
155 const char *sec_name,
156 int location_id,
157 int n_location_vals,
158 cs_restart_val_type_t val_type);
159
160 static int
161 _read_section(cs_restart_t *restart,
162 void *context,
163 const char *sec_name,
164 int location_id,
165 int n_location_vals,
166 cs_restart_val_type_t val_type,
167 void *val);
168
169 static void
170 _write_section(cs_restart_t *restart,
171 void *context,
172 const char *sec_name,
173 int location_id,
174 int n_location_vals,
175 cs_restart_val_type_t val_type,
176 const void *val);
177
178 /*============================================================================
179 * Static global variables
180 *============================================================================*/
181
182 #if defined(WIN32) || defined(_WIN32)
183 static const char _dir_separator = '\\';
184 #else
185 static const char _dir_separator = '/';
186 #endif
187
188 /* Monitoring info */
189
190 static int _restart_n_opens[2] = {0, 0};
191 static double _restart_wtime[2] = {0.0, 0.0};
192
193 /* Do we have a restart directory ? */
194
195 static int _restart_present = -1;
196
197 /* Restart time steps and frequency */
198
199 static int _checkpoint_mesh = 1; /* checkpoint mesh if possible */
200 /* time step interval */
201 static int _checkpoint_nt_interval = CS_RESTART_INTERVAL_ONLY_AT_END;
202 static int _checkpoint_nt_next = -1; /* next forced time step */
203 static int _checkpoint_nt_last = -1; /* last checkpoint time step */
204 static double _checkpoint_t_interval = -1.; /* physical time interval */
205 static double _checkpoint_t_next = -1.; /* next forced time value */
206 static double _checkpoint_t_last = 0.; /* last forced time value */
207 static double _checkpoint_wt_interval = -1.; /* wall-clock interval */
208 static double _checkpoint_wt_next = -1.; /* next forced wall-clock value */
209 static double _checkpoint_wt_last = 0.; /* wall-clock time of last
210 checkpointing */
211 /* Are we restarting from a NCFD file ? */
212
213 static int _restart_from_ncfd = 0;
214
215 /* Restart modification */
216
217 static void *_restart_context = NULL;
218 static cs_restart_check_section_t *_check_section_f = _check_section;
219 static cs_restart_read_section_t *_read_section_f = _read_section;
220 static cs_restart_write_section_t *_write_section_f = _write_section;
221
222 static size_t _n_locations_ref; /* Number of locations */
223 static _location_t *_location_ref; /* Location definition array */
224
225
226 /* Restart multi writer */
227
228 static int _n_restart_directories_to_write = 1;
229 static int _n_restart_multiwriters = 0;
230 static _restart_multiwriter_t **_restart_multiwriter = NULL;
231
232 /*============================================================================
233 * Private function definitions
234 *============================================================================*/
235
236 /*----------------------------------------------------------------------------
237 * Compute number of values in a record
238 *
239 * parameters:
240 * r <-- associated restart file pointer
241 * location_id <-- location id
242 * n_location_vals <-- number of values per location
243 *----------------------------------------------------------------------------*/
244
245 static cs_gnum_t
_compute_n_ents(const cs_restart_t * r,size_t location_id,size_t n_location_vals)246 _compute_n_ents(const cs_restart_t *r,
247 size_t location_id,
248 size_t n_location_vals)
249 {
250 cs_gnum_t retval = 0;
251
252 if (location_id == 0)
253 retval = n_location_vals;
254
255 else if (location_id > 0 && location_id <= r->n_locations)
256 retval = r->location[location_id-1].n_glob_ents_f
257 * (cs_gnum_t)n_location_vals;
258
259 else
260 bft_error(__FILE__, __LINE__, 0,
261 _("Location number %d given for restart file\n"
262 "\"%s\" is not valid."),
263 (int)location_id, r->name);
264
265 return retval;
266 }
267
268 /*----------------------------------------------------------------------------
269 * Analyze the content of a restart file to build locations
270 *
271 * parameters:
272 * r <-> associated restart file pointer
273 *----------------------------------------------------------------------------*/
274
275 static void
_locations_from_index(cs_restart_t * r)276 _locations_from_index(cs_restart_t *r)
277 {
278 cs_io_sec_header_t h;
279
280 size_t rec_id = 0;
281 size_t index_size = 0;
282
283 /* Initialization */
284
285 index_size = cs_io_get_index_size(r->fh);
286
287 /* Analyze records to determine locations */
288
289 for (rec_id = 0; rec_id < index_size; rec_id++) {
290
291 h = cs_io_get_indexed_sec_header(r->fh, rec_id);
292
293 if (h.location_id > r->n_locations) {
294
295 _location_t *loc = NULL;
296
297 if (h.location_id != r->n_locations + 1)
298 bft_error(__FILE__, __LINE__, 0,
299 _("Restart file \"%s\" declares a location number %d\n"
300 "but no location %d has been declared."),
301 r->name, (int)(h.location_id),
302 (int)(r->n_locations + 1));
303
304 BFT_REALLOC(r->location, r->n_locations + 1, _location_t);
305
306 loc = r->location + r->n_locations;
307 BFT_MALLOC(loc->name, strlen(h.sec_name) + 1, char);
308 strcpy(loc->name, h.sec_name);
309
310 loc->id = h.location_id;
311 loc->n_ents = 0;
312 loc->n_glob_ents = 0;
313
314 cs_io_set_indexed_position(r->fh, &h, rec_id);
315 cs_io_set_cs_gnum(&h, r->fh);
316 cs_io_read_global(&h, &(loc->n_glob_ents_f), r->fh);
317
318 loc->ent_global_num = NULL;
319 loc->_ent_global_num = NULL;
320
321 r->n_locations += 1;
322 }
323
324 }
325 }
326
327 /*----------------------------------------------------------------------------
328 * Initialize a checkpoint / restart file management structure;
329 *
330 * parameters:
331 * r <-> associated restart file pointer
332 *----------------------------------------------------------------------------*/
333
334 static void
_add_file(cs_restart_t * r)335 _add_file(cs_restart_t *r)
336 {
337 double timing[2];
338 cs_file_access_t method;
339
340 const char magic_string[] = "Checkpoint / restart, R0";
341 const long echo = CS_IO_ECHO_NONE;
342
343 timing[0] = cs_timer_wtime();
344
345 /* In read mode, open file to detect header first */
346
347 #if defined(HAVE_MPI)
348 {
349 MPI_Info hints;
350 MPI_Comm block_comm, comm;
351
352 cs_file_get_default_comm(NULL, &block_comm, &comm);
353
354 r->rank_step = 1;
355 r->min_block_size = cs_parall_get_min_coll_buf_size();
356 assert(comm == cs_glob_mpi_comm || comm == MPI_COMM_NULL);
357
358 if (r->mode == CS_RESTART_MODE_READ) {
359 cs_file_get_default_access(CS_FILE_MODE_READ, &method, &hints);
360 r->fh = cs_io_initialize_with_index(r->name,
361 magic_string,
362 method,
363 echo,
364 hints,
365 block_comm,
366 comm);
367 _locations_from_index(r);
368 }
369 else {
370 cs_file_get_default_access(CS_FILE_MODE_WRITE, &method, &hints);
371 r->fh = cs_io_initialize(r->name,
372 magic_string,
373 CS_IO_MODE_WRITE,
374 method,
375 echo,
376 hints,
377 block_comm,
378 comm);
379 }
380 }
381 #else
382 {
383 if (r->mode == CS_RESTART_MODE_READ) {
384 cs_file_get_default_access(CS_FILE_MODE_READ, &method);
385 r->fh = cs_io_initialize_with_index(r->name,
386 magic_string,
387 method,
388 echo);
389 _locations_from_index(r);
390 }
391 else {
392 cs_file_get_default_access(CS_FILE_MODE_WRITE, &method);
393 r->fh = cs_io_initialize(r->name,
394 magic_string,
395 CS_IO_MODE_WRITE,
396 method,
397 echo);
398 }
399 }
400 #endif
401
402 timing[1] = cs_timer_wtime();
403 _restart_wtime[r->mode] += timing[1] - timing[0];
404
405 _restart_n_opens[r->mode] += 1;
406 }
407
408 #if defined(HAVE_MPI)
409
410 /*----------------------------------------------------------------------------
411 * Read variable values defined on a mesh location.
412 *
413 * parameters:
414 * r <-> associated restart file pointer
415 * header <-- header associated with current position in file
416 * n_glob_ents <-- global number of entities
417 * n_ents <-- local number of entities
418 * ent_global_num <-- global entity numbers (1 to n numbering)
419 * n_location_vals <-- number of values par location
420 * val_type <-- data type
421 * vals --> array of values
422 *----------------------------------------------------------------------------*/
423
424 static void
_read_ent_values(cs_restart_t * r,cs_io_sec_header_t * header,cs_gnum_t n_glob_ents,cs_lnum_t n_ents,const cs_gnum_t ent_global_num[],int n_location_vals,cs_restart_val_type_t val_type,cs_byte_t vals[])425 _read_ent_values(cs_restart_t *r,
426 cs_io_sec_header_t *header,
427 cs_gnum_t n_glob_ents,
428 cs_lnum_t n_ents,
429 const cs_gnum_t ent_global_num[],
430 int n_location_vals,
431 cs_restart_val_type_t val_type,
432 cs_byte_t vals[])
433 {
434 cs_byte_t *buffer = NULL;
435
436 cs_lnum_t block_buf_size = 0;
437
438 size_t nbr_byte_ent;
439
440 /* Initialization */
441
442 switch (val_type) {
443 case CS_TYPE_char:
444 nbr_byte_ent = n_location_vals;
445 break;
446 case CS_TYPE_int:
447 nbr_byte_ent = n_location_vals * sizeof(int);
448 cs_io_set_cs_lnum(header, r->fh);
449 break;
450 case CS_TYPE_cs_gnum_t:
451 nbr_byte_ent = n_location_vals * sizeof(cs_gnum_t);
452 cs_io_set_cs_gnum(header, r->fh);
453 break;
454 case CS_TYPE_cs_real_t:
455 nbr_byte_ent = n_location_vals * sizeof(cs_real_t);
456 break;
457 default:
458 nbr_byte_ent = 0;
459 assert(0);
460 }
461
462 cs_block_dist_info_t bi
463 = cs_block_dist_compute_sizes(cs_glob_rank_id,
464 cs_glob_n_ranks,
465 r->rank_step,
466 r->min_block_size / nbr_byte_ent,
467 n_glob_ents);
468
469 cs_all_to_all_t *d
470 = cs_all_to_all_create_from_block(n_ents,
471 CS_ALL_TO_ALL_USE_DEST_ID,
472 ent_global_num,
473 bi,
474 cs_glob_mpi_comm);
475
476 /* Read blocks */
477
478 block_buf_size = (bi.gnum_range[1] - bi.gnum_range[0]) * nbr_byte_ent;
479
480 if (block_buf_size > 0)
481 BFT_MALLOC(buffer, block_buf_size, cs_byte_t);
482
483 cs_io_read_block(header,
484 bi.gnum_range[0],
485 bi.gnum_range[1],
486 buffer,
487 r->fh);
488
489 /* Distribute blocks on ranks */
490
491 cs_all_to_all_copy_array(d,
492 header->elt_type,
493 n_location_vals,
494 true, /* reverse */
495 buffer,
496 vals);
497
498 /* Free buffer */
499
500 BFT_FREE(buffer);
501
502 cs_all_to_all_destroy(&d);
503 }
504
505 /*----------------------------------------------------------------------------
506 * Write variable values defined on a mesh location.
507 *
508 * parameters:
509 * r <-> associated restart file pointer
510 * sec_name <-- section name
511 * n_glob_ents <-- global number of entities
512 * n_ents <-- local number of entities
513 * ent_global_num <-- global entity numbers (1 to n numbering)
514 * location_id <-- id of corresponding location
515 * n_location_vals <-- number of values par location
516 * val_type <-- data type
517 * vals --> array of values
518 *----------------------------------------------------------------------------*/
519
520 static void
_write_ent_values(const cs_restart_t * r,const char * sec_name,cs_gnum_t n_glob_ents,cs_lnum_t n_ents,const cs_gnum_t * ent_global_num,int location_id,int n_location_vals,cs_restart_val_type_t val_type,const cs_byte_t * vals)521 _write_ent_values(const cs_restart_t *r,
522 const char *sec_name,
523 cs_gnum_t n_glob_ents,
524 cs_lnum_t n_ents,
525 const cs_gnum_t *ent_global_num,
526 int location_id,
527 int n_location_vals,
528 cs_restart_val_type_t val_type,
529 const cs_byte_t *vals)
530 {
531 cs_lnum_t block_buf_size = 0;
532
533 cs_datatype_t elt_type = CS_DATATYPE_NULL;
534 size_t nbr_byte_ent = 0;
535 cs_byte_t *buffer = NULL;
536
537 cs_block_dist_info_t bi;
538
539 cs_part_to_block_t *d = NULL;
540
541 /* Initialization */
542
543 switch (val_type) {
544 case CS_TYPE_char:
545 nbr_byte_ent = n_location_vals;
546 elt_type = CS_CHAR;
547 break;
548 case CS_TYPE_int:
549 nbr_byte_ent = n_location_vals * sizeof(int);
550 elt_type = (sizeof(int) == 8) ? CS_INT64 : CS_INT32;
551 break;
552 case CS_TYPE_cs_gnum_t:
553 nbr_byte_ent = n_location_vals * sizeof(cs_gnum_t);
554 elt_type = (sizeof(cs_gnum_t) == 8) ? CS_UINT64 : CS_UINT32;
555 break;
556 case CS_TYPE_cs_real_t:
557 nbr_byte_ent = n_location_vals * sizeof(cs_real_t);
558 elt_type = (sizeof(cs_real_t) == cs_datatype_size[CS_DOUBLE])
559 ? CS_DOUBLE : CS_FLOAT;
560 break;
561 default:
562 assert(0);
563 }
564
565 bi = cs_block_dist_compute_sizes(cs_glob_rank_id,
566 cs_glob_n_ranks,
567 r->rank_step,
568 r->min_block_size / nbr_byte_ent,
569 n_glob_ents);
570
571 d = cs_part_to_block_create_by_gnum(cs_glob_mpi_comm,
572 bi,
573 n_ents,
574 ent_global_num);
575
576 /* Distribute to blocks */
577
578 block_buf_size = (bi.gnum_range[1] - bi.gnum_range[0]) * nbr_byte_ent;
579
580 if (block_buf_size > 0)
581 BFT_MALLOC(buffer, block_buf_size, cs_byte_t);
582
583 /* Distribute blocks on ranks */
584
585 cs_part_to_block_copy_array(d,
586 elt_type,
587 n_location_vals,
588 vals,
589 buffer);
590
591 /* Write blocks */
592
593 cs_io_write_block_buffer(sec_name,
594 n_glob_ents,
595 bi.gnum_range[0],
596 bi.gnum_range[1],
597 location_id,
598 0,
599 n_location_vals,
600 elt_type,
601 buffer,
602 r->fh);
603
604 /* Free buffer */
605
606 BFT_FREE(buffer);
607
608 cs_part_to_block_destroy(&d);
609 }
610
611 #endif /* #if defined(HAVE_MPI) */
612
613 /*----------------------------------------------------------------------------
614 * Swap values of a renumbered array when reading
615 *
616 * parameters:
617 * n_ents --> number of local entities
618 * ini_ent_num --> initial entity numbers
619 * n_location_vals --> number of values per entity
620 * val_type --> data type
621 * vals --> array of values
622 *----------------------------------------------------------------------------*/
623
624 static void
_restart_permute_read(cs_lnum_t n_ents,const cs_gnum_t * ini_ent_num,int n_location_vals,cs_restart_val_type_t val_type,cs_byte_t * vals)625 _restart_permute_read(cs_lnum_t n_ents,
626 const cs_gnum_t *ini_ent_num,
627 int n_location_vals,
628 cs_restart_val_type_t val_type,
629 cs_byte_t *vals)
630 {
631 cs_lnum_t ent_id, jj;
632
633 cs_lnum_t ii = 0;
634
635 /* Instructions */
636
637 if (ini_ent_num == NULL)
638 return;
639
640 switch (val_type) {
641
642 case CS_TYPE_char:
643 {
644 char *val_ord;
645 char *val_cur = (char *)vals;
646
647 BFT_MALLOC(val_ord, n_ents * n_location_vals, char);
648
649 for (ent_id = 0; ent_id < n_ents; ent_id++) {
650 for (jj = 0; jj < n_location_vals; jj++)
651 val_ord[ii++]
652 = val_cur[(ini_ent_num[ent_id] - 1) * n_location_vals + jj];
653 }
654
655 for (ii = 0; ii < n_ents * n_location_vals; ii++)
656 val_cur[ii] = val_ord[ii];
657
658 BFT_FREE(val_ord);
659 }
660 break;
661
662 case CS_TYPE_int:
663 {
664 int *val_ord;
665 int *val_cur = (int *)vals;
666
667 BFT_MALLOC(val_ord, n_ents * n_location_vals, int);
668
669 for (ent_id = 0; ent_id < n_ents; ent_id++) {
670 for (jj = 0; jj < n_location_vals; jj++)
671 val_ord[ii++]
672 = val_cur[(ini_ent_num[ent_id] - 1) * n_location_vals + jj];
673 }
674
675 for (ii = 0; ii < n_ents * n_location_vals; ii++)
676 val_cur[ii] = val_ord[ii];
677
678 BFT_FREE(val_ord);
679 }
680 break;
681
682 case CS_TYPE_cs_gnum_t:
683 {
684 cs_gnum_t *val_ord;
685 cs_gnum_t *val_cur = (cs_gnum_t *)vals;
686
687 BFT_MALLOC(val_ord, n_ents * n_location_vals, cs_gnum_t);
688
689 for (ent_id = 0; ent_id < n_ents; ent_id++) {
690 for (jj = 0; jj < n_location_vals; jj++)
691 val_ord[ii++]
692 = val_cur[(ini_ent_num[ent_id] - 1) * n_location_vals + jj];
693 }
694
695 for (ii = 0; ii < n_ents * n_location_vals; ii++)
696 val_cur[ii] = val_ord[ii];
697
698 BFT_FREE(val_ord);
699 }
700 break;
701
702 case CS_TYPE_cs_real_t:
703 {
704 cs_real_t *val_ord;
705 cs_real_t *val_cur = (cs_real_t *)vals;
706
707 BFT_MALLOC (val_ord, n_ents * n_location_vals, cs_real_t);
708
709 for (ent_id = 0; ent_id < n_ents; ent_id++) {
710 for (jj = 0; jj < n_location_vals; jj++)
711 val_ord[ii++]
712 = val_cur[(ini_ent_num[ent_id] - 1) * n_location_vals + jj];
713 }
714
715 for (ii = 0; ii < n_ents * n_location_vals; ii++)
716 val_cur[ii] = val_ord[ii];
717
718 BFT_FREE(val_ord);
719 }
720 break;
721
722 default:
723 assert(0);
724
725 }
726 }
727
728 /*----------------------------------------------------------------------------
729 * Swap values of a renumbered array when writing
730 *
731 * parameters:
732 * n_ents --> number of local entities
733 * ini_ent_num --> initial entity numbers
734 * n_location_vals --> number of values per entity
735 * val_type --> data type
736 * vals --> array of values
737 *
738 * returns:
739 * pointer to array of values in initial entity order
740 *----------------------------------------------------------------------------*/
741
742 static cs_byte_t *
_restart_permute_write(cs_lnum_t n_ents,const cs_gnum_t * ini_ent_num,int n_location_vals,cs_restart_val_type_t val_type,const cs_byte_t * vals)743 _restart_permute_write(cs_lnum_t n_ents,
744 const cs_gnum_t *ini_ent_num,
745 int n_location_vals,
746 cs_restart_val_type_t val_type,
747 const cs_byte_t *vals)
748 {
749 cs_lnum_t ent_id, jj;
750
751 cs_lnum_t ii = 0;
752
753 /* Instructions */
754
755 if (ini_ent_num == NULL)
756 return NULL;
757
758 switch (val_type) {
759
760 case CS_TYPE_char:
761 {
762 char *val_ord;
763 const char *val_cur = (const char *)vals;
764
765 BFT_MALLOC(val_ord, n_ents * n_location_vals, char);
766
767 for (ent_id = 0; ent_id < n_ents; ent_id++) {
768 for (jj = 0; jj < n_location_vals; jj++)
769 val_ord[(ini_ent_num[ent_id] - 1) * n_location_vals + jj]
770 = val_cur[ii++];
771 }
772
773 return (cs_byte_t *)val_ord;
774 }
775 break;
776
777 case CS_TYPE_int:
778 {
779 int *val_ord;
780 const int *val_cur = (const int *)vals;
781
782 BFT_MALLOC(val_ord, n_ents * n_location_vals, int);
783
784 for (ent_id = 0; ent_id < n_ents; ent_id++) {
785 for (jj = 0; jj < n_location_vals; jj++)
786 val_ord[(ini_ent_num[ent_id] - 1) * n_location_vals + jj]
787 = val_cur[ii++];
788 }
789
790 return (cs_byte_t *)val_ord;
791 }
792 break;
793
794 case CS_TYPE_cs_gnum_t:
795 {
796 cs_gnum_t *val_ord;
797 const cs_gnum_t *val_cur = (const cs_gnum_t *)vals;
798
799 BFT_MALLOC(val_ord, n_ents * n_location_vals, cs_gnum_t);
800
801 for (ent_id = 0; ent_id < n_ents; ent_id++) {
802 for (jj = 0; jj < n_location_vals; jj++)
803 val_ord[(ini_ent_num[ent_id] - 1) * n_location_vals + jj]
804 = val_cur[ii++];
805 }
806
807 return (cs_byte_t *)val_ord;
808 }
809 break;
810
811 case CS_TYPE_cs_real_t:
812 {
813 cs_real_t *val_ord;
814 const cs_real_t *val_cur = (const cs_real_t *)vals;
815
816 BFT_MALLOC(val_ord, n_ents * n_location_vals, cs_real_t);
817
818 for (ent_id = 0; ent_id < n_ents; ent_id++) {
819 for (jj = 0; jj < n_location_vals; jj++)
820 val_ord[(ini_ent_num[ent_id] - 1) * n_location_vals + jj]
821 = val_cur[ii++];
822 }
823
824 return (cs_byte_t *)val_ord;
825 }
826 break;
827
828 default:
829 assert(0);
830 return NULL;
831
832 }
833 }
834
835 /*----------------------------------------------------------------------------
836 * Find a given record in an indexed restart file.
837 *
838 * parameters:
839 * restart <-- associated restart file pointer
840 * prefix <-- prefix to name of record
841 * name <-- base name of record
842 * postfix <-- postfix to name of record
843 *
844 * returns:
845 * the id assigned to the record, or -1 if not found
846 *----------------------------------------------------------------------------*/
847
848 static int
_restart_section_id(cs_restart_t * restart,const char * prefix,const char * name,const char * postfix)849 _restart_section_id(cs_restart_t *restart,
850 const char *prefix,
851 const char *name,
852 const char *postfix)
853 {
854 size_t index_size = cs_io_get_index_size(restart->fh);
855
856 char *_sec_name = NULL;
857 const char *sec_name = name;
858
859 int rec_id = -1;
860
861 if (prefix != NULL || postfix != NULL) {
862
863
864 size_t sec_name_l = strlen(name);
865
866 if (prefix != NULL)
867 sec_name_l += strlen(prefix);
868 if (postfix != NULL)
869 sec_name_l += strlen(postfix);
870
871 BFT_MALLOC(_sec_name, sec_name_l + 1, char);
872 sec_name = _sec_name;
873
874 if (prefix != NULL) {
875 strcpy(_sec_name, prefix);
876 strcat(_sec_name, name);
877 }
878 else
879 strcpy(_sec_name, name);
880
881 if (postfix != NULL)
882 strcat(_sec_name, postfix);
883
884 }
885
886 /* Search for the record in the index */
887
888 for (rec_id = 0; rec_id < (int)index_size; rec_id++) {
889 const char * cmp_name = cs_io_get_indexed_sec_name(restart->fh, rec_id);
890 if (strcmp(cmp_name, sec_name) == 0)
891 break;
892 }
893
894 if (rec_id >= (int)index_size) {
895 bft_printf(_(" %s: section \"%s\" not present.\n"),
896 restart->name, sec_name);
897 rec_id = -1;
898 }
899
900 BFT_FREE(_sec_name);
901
902 return rec_id;
903 }
904
905 #if defined(HAVE_MPI)
906
907 /*----------------------------------------------------------------------------
908 * Compute default particle destination rank array in case of untracked
909 * particles.
910 *
911 * For simplicity, those particles are simply distributed among ranks
912 * (as it is possible to define a global numbering based on a space-filling
913 * curve when generating the restart file, this may be made "geometrically"
914 * balanced also).
915 *
916 * parameters:
917 * p_bi <-- pointer to particles bock distribution info
918 * p_cell_num <-- global cell number associated with each particle
919 * (0 for untracked particles)
920 * comm <-- associated MPI communicator
921 *
922 * returns:
923 * default rank array for particles (>= 0 for untracked particles)
924 *----------------------------------------------------------------------------*/
925
926 static int *
_default_p_rank(cs_block_dist_info_t * p_bi,const cs_gnum_t * p_cell_num,MPI_Comm comm)927 _default_p_rank(cs_block_dist_info_t *p_bi,
928 const cs_gnum_t *p_cell_num,
929 MPI_Comm comm)
930 {
931 cs_lnum_t i;
932 cs_block_dist_info_t free_particle_bi;
933
934 int n_ranks = 0, rank_id = -1;
935
936 cs_lnum_t _n_particles = 0, n_free_particles = 0;
937 cs_gnum_t _n_g_free_particles = 0, n_g_free_particles = 0;
938
939 cs_lnum_t *free_particle_ids = NULL;
940
941 fvm_io_num_t *free_particle_io_num = NULL;
942 const cs_gnum_t *free_particle_num = NULL;
943
944 int *default_rank = NULL;
945
946 /* Initialization */
947
948 assert((sizeof(cs_lnum_t) == 4) || (sizeof(cs_lnum_t) == 8));
949
950 /* Count number of untracked particles */
951
952 _n_particles = p_bi->gnum_range[1] - p_bi->gnum_range[0];
953 n_free_particles = 0;
954
955 for (i = 0; i < _n_particles; i++) {
956 if (p_cell_num[i] == 0)
957 n_free_particles += 1;
958 }
959
960 _n_g_free_particles = n_free_particles;
961 MPI_Allreduce(&_n_g_free_particles, &n_g_free_particles, 1,
962 CS_MPI_GNUM, MPI_SUM, comm);
963
964 /* Return if we do not have untracked particles */
965
966 if (n_g_free_particles == 0)
967 return NULL;
968
969 /* Initialize rank info */
970
971 MPI_Comm_size(comm, &n_ranks);
972 MPI_Comm_size(comm, &rank_id);
973 free_particle_bi = cs_block_dist_compute_sizes(rank_id,
974 n_ranks,
975 0,
976 0,
977 n_g_free_particles);
978
979 /* Define distribution of untracked particles based on scan;
980 *
981 * The main objective of this function
982 * is to ensure some measure of load balancing. */
983
984 BFT_MALLOC(default_rank, _n_particles, int);
985 for (i = 0; i < _n_particles; i++)
986 default_rank[i] = -1;
987
988 BFT_MALLOC(free_particle_ids, n_free_particles, cs_lnum_t);
989
990 n_free_particles = 0;
991 for (i = 0; i < _n_particles; i++) {
992 if (p_cell_num[i] == 0)
993 free_particle_ids[n_free_particles++] = i;
994 }
995
996 free_particle_io_num = fvm_io_num_create_from_scan(n_free_particles);
997 free_particle_num = fvm_io_num_get_global_num(free_particle_io_num);
998
999 /* Determine rank based on global numbering */
1000 for (i = 0; i < n_free_particles; i++) {
1001 default_rank[free_particle_ids[i]]
1002 = ((free_particle_num[i] - 1) / free_particle_bi.block_size)
1003 * free_particle_bi.rank_step;
1004 }
1005
1006 free_particle_io_num = fvm_io_num_destroy(free_particle_io_num);
1007 BFT_FREE(free_particle_ids);
1008
1009 return default_rank;
1010 }
1011
1012 #endif /* defined(HAVE_MPI) */
1013
1014 /*----------------------------------------------------------------------------*/
1015 /*!
1016 * \brief Check the presence of a given section in a restart file.
1017 *
1018 * \param[in] restart associated restart file pointer
1019 * \param[in, out] context associated context
1020 * \param[in] sec_name section name
1021 * \param[in] location_id id of corresponding location
1022 * \param[in] n_location_vals number of values per location (interlaced)
1023 * \param[in] val_type value type
1024 *
1025 * \return 0 (CS_RESTART_SUCCESS) in case of success,
1026 * or error code (CS_RESTART_ERR_xxx) in case of error
1027 */
1028 /*----------------------------------------------------------------------------*/
1029
1030 static int
_check_section(cs_restart_t * restart,void * context,const char * sec_name,int location_id,int n_location_vals,cs_restart_val_type_t val_type)1031 _check_section(cs_restart_t *restart,
1032 void *context,
1033 const char *sec_name,
1034 int location_id,
1035 int n_location_vals,
1036 cs_restart_val_type_t val_type)
1037 {
1038 CS_UNUSED(context);
1039
1040 cs_lnum_t n_ents;
1041 cs_gnum_t n_glob_ents;
1042
1043 size_t rec_id;
1044 cs_io_sec_header_t header;
1045
1046 size_t index_size = 0;
1047
1048 index_size = cs_io_get_index_size(restart->fh);
1049
1050 assert(restart != NULL);
1051
1052 /* Check associated location */
1053
1054 if (location_id == 0) {
1055 n_glob_ents = n_location_vals;
1056 n_ents = n_location_vals;
1057 }
1058
1059 else {
1060 if (location_id < 0 || location_id > (int)(restart->n_locations))
1061 return CS_RESTART_ERR_LOCATION;
1062 n_glob_ents = (restart->location[location_id-1]).n_glob_ents;
1063 if ((restart->location[location_id-1]).n_glob_ents_f != n_glob_ents)
1064 return CS_RESTART_ERR_LOCATION;
1065 n_ents = (restart->location[location_id-1]).n_ents;
1066 }
1067
1068 /* Search for the corresponding record in the index */
1069
1070 for (rec_id = 0; rec_id < index_size; rec_id++) {
1071 const char * cmp_name = cs_io_get_indexed_sec_name(restart->fh, rec_id);
1072 if (strcmp(cmp_name, sec_name) == 0)
1073 break;
1074 }
1075
1076 /* If the record was not found */
1077
1078 if (rec_id >= index_size)
1079 return CS_RESTART_ERR_EXISTS;
1080
1081 /*
1082 If the location does not fit: we search for a location of same
1083 name with the correct location.
1084 */
1085
1086 header = cs_io_get_indexed_sec_header(restart->fh, rec_id);
1087
1088 if (header.location_id != (size_t)location_id) {
1089
1090 rec_id++;
1091
1092 while (rec_id < index_size) {
1093 header = cs_io_get_indexed_sec_header(restart->fh, rec_id);
1094 if ( (strcmp(header.sec_name, sec_name) == 0)
1095 && (header.location_id == (size_t)location_id))
1096 break;
1097 rec_id++;
1098 }
1099
1100 if (rec_id >= index_size)
1101 return CS_RESTART_ERR_LOCATION;
1102 }
1103
1104 /* If the number of values per location does not match */
1105
1106 if ( header.location_id > 0
1107 && header.n_location_vals != (size_t)n_location_vals)
1108 return CS_RESTART_ERR_N_VALS;
1109 else if (header.location_id == 0 && header.n_vals != n_ents)
1110 return CS_RESTART_ERR_N_VALS;
1111
1112 /* If the type of value does not match */
1113
1114 if (header.elt_type == CS_CHAR) {
1115 if (val_type != CS_TYPE_char)
1116 return CS_RESTART_ERR_VAL_TYPE;
1117 }
1118 else if (header.elt_type == CS_INT32 || header.elt_type == CS_INT64) {
1119 cs_io_set_cs_lnum(&header, restart->fh);
1120 if (val_type != CS_TYPE_int)
1121 return CS_RESTART_ERR_VAL_TYPE;
1122 }
1123 else if (header.elt_type == CS_UINT32 || header.elt_type == CS_UINT64) {
1124 if (val_type != CS_TYPE_cs_gnum_t && val_type != CS_TYPE_int)
1125 return CS_RESTART_ERR_VAL_TYPE;
1126 }
1127 else if (header.elt_type == CS_FLOAT || header.elt_type == CS_DOUBLE) {
1128 if (val_type != CS_TYPE_cs_real_t)
1129 return CS_RESTART_ERR_VAL_TYPE;
1130 }
1131
1132 /* Return */
1133
1134 return CS_RESTART_SUCCESS;
1135 }
1136
1137 /*----------------------------------------------------------------------------*/
1138 /*!
1139 * \brief Read a section from a restart file.
1140 *
1141 * \param[in] restart associated restart file pointer
1142 * \param[in] sec_name section name
1143 * \param[in] location_id id of corresponding location
1144 * \param[in] n_location_vals number of values per location (interlaced)
1145 * \param[in] val_type value type
1146 * \param[out] val array of values
1147 *
1148 * \return 0 (CS_RESTART_SUCCESS) in case of success,
1149 * or error code (CS_RESTART_ERR_xxx) in case of error
1150 */
1151 /*----------------------------------------------------------------------------*/
1152
1153 static int
_read_section(cs_restart_t * restart,void * context,const char * sec_name,int location_id,int n_location_vals,cs_restart_val_type_t val_type,void * val)1154 _read_section(cs_restart_t *restart,
1155 void *context,
1156 const char *sec_name,
1157 int location_id,
1158 int n_location_vals,
1159 cs_restart_val_type_t val_type,
1160 void *val)
1161 {
1162 CS_UNUSED(context);
1163
1164 cs_lnum_t n_ents;
1165 cs_gnum_t n_glob_ents;
1166
1167 const cs_gnum_t *ent_global_num;
1168
1169 size_t rec_id, rec_id_tmp;
1170 cs_io_sec_header_t header;
1171
1172 cs_lnum_t _n_location_vals = n_location_vals;
1173 size_t index_size = 0;
1174
1175 index_size = cs_io_get_index_size(restart->fh);
1176
1177 assert(restart != NULL);
1178
1179 /* Check associated location */
1180
1181 if (location_id == 0) {
1182 n_glob_ents = n_location_vals;
1183 n_ents = n_location_vals;
1184 _n_location_vals = 1;
1185 ent_global_num = NULL;
1186 }
1187
1188 else {
1189 if (location_id < 0 || location_id > (int)(restart->n_locations)) {
1190 bft_printf(_(" %s: location id %d for \"%s\" does not exist.\n"),
1191 restart->name, location_id, sec_name);
1192 return CS_RESTART_ERR_LOCATION;
1193 }
1194 n_glob_ents = (restart->location[location_id-1]).n_glob_ents;
1195 if ((restart->location[location_id-1]).n_glob_ents_f != n_glob_ents) {
1196 bft_printf
1197 (_(" %s: location id %d for \"%s\" has "
1198 "size %llu, but %llu is expected.\n"),
1199 restart->name, location_id, sec_name,
1200 (unsigned long long)(restart->location[location_id-1]).n_glob_ents_f,
1201 (unsigned long long)n_glob_ents);
1202 return CS_RESTART_ERR_LOCATION;
1203 }
1204 n_ents = (restart->location[location_id-1]).n_ents;
1205 ent_global_num = (restart->location[location_id-1]).ent_global_num;
1206 }
1207
1208 /* Search for the corresponding record in the index */
1209
1210 for (rec_id = 0; rec_id < index_size; rec_id++) {
1211 const char * cmp_name = cs_io_get_indexed_sec_name(restart->fh, rec_id);
1212 if (strcmp(cmp_name, sec_name) == 0)
1213 break;
1214 }
1215
1216 /* If the record was not found */
1217
1218 if (rec_id >= index_size) {
1219 bft_printf(_(" %s: section \"%s\" not present.\n"),
1220 restart->name, sec_name);
1221 return CS_RESTART_ERR_EXISTS;
1222 }
1223
1224 /*
1225 If the location does not fit: we search for a location of same
1226 name with the correct location.
1227 */
1228
1229 header = cs_io_get_indexed_sec_header(restart->fh, rec_id);
1230
1231 if (header.location_id != (size_t)location_id) {
1232
1233 rec_id_tmp = rec_id;
1234 rec_id++;
1235
1236 while (rec_id < index_size) {
1237 header = cs_io_get_indexed_sec_header(restart->fh, rec_id);
1238 if ( (strcmp(header.sec_name, sec_name) == 0)
1239 && (header.location_id == (size_t)location_id))
1240 break;
1241 rec_id++;
1242 }
1243
1244 if (rec_id >= index_size) {
1245 header = cs_io_get_indexed_sec_header(restart->fh, rec_id_tmp);
1246 bft_printf(_(" %s: section \"%s\" at location id %d but not at %d.\n"),
1247 restart->name, sec_name,
1248 (int)(header.location_id), (int)location_id);
1249 return CS_RESTART_ERR_LOCATION;
1250 }
1251 }
1252
1253 /* If the number of values per location does not match */
1254
1255 if ( header.location_id > 0
1256 && header.n_location_vals != (size_t)n_location_vals) {
1257 bft_printf(_(" %s: section \"%s\" has %d values per location and "
1258 " not %d.\n"),
1259 restart->name, sec_name,
1260 (int)header.n_location_vals, (int)n_location_vals);
1261 return CS_RESTART_ERR_N_VALS;
1262 }
1263 else if (header.location_id == 0 && header.n_vals != n_ents) {
1264 bft_printf(_(" %s: section \"%s\" has %d values and not %d.\n"),
1265 restart->name, sec_name, (int)header.n_vals, (int)n_ents);
1266 return CS_RESTART_ERR_N_VALS;
1267 }
1268
1269 /* If the type of value does not match */
1270
1271 if (header.elt_type == CS_CHAR) {
1272 if (val_type != CS_TYPE_char) {
1273 bft_printf(_(" %s: section \"%s\" is not of character type.\n"),
1274 restart->name, sec_name);
1275 return CS_RESTART_ERR_VAL_TYPE;
1276 }
1277 }
1278 else if (header.elt_type == CS_INT32 || header.elt_type == CS_INT64) {
1279 cs_io_set_cs_lnum(&header, restart->fh);
1280 if (val_type != CS_TYPE_int) {
1281 bft_printf(_(" %s: section \"%s\" is not of integer type.\n"),
1282 restart->name, sec_name);
1283 return CS_RESTART_ERR_VAL_TYPE;
1284 }
1285 }
1286 else if (header.elt_type == CS_UINT32 || header.elt_type == CS_UINT64) {
1287 if (val_type != CS_TYPE_cs_gnum_t && val_type != CS_TYPE_int) {
1288 bft_printf(_(" %s: section \"%s\" is not of global number type.\n"),
1289 restart->name, sec_name);
1290 return CS_RESTART_ERR_VAL_TYPE;
1291 }
1292 }
1293 else if (header.elt_type == CS_FLOAT || header.elt_type == CS_DOUBLE) {
1294 if (val_type != CS_TYPE_cs_real_t) {
1295 bft_printf(_(" %s: section \"%s\" is not of floating-point type.\n"),
1296 restart->name, sec_name);
1297 return CS_RESTART_ERR_VAL_TYPE;
1298 }
1299 }
1300
1301 /* Now set position in file to read data */
1302
1303 cs_io_set_indexed_position(restart->fh, &header, rec_id);
1304
1305 /* Now define conversion info */
1306
1307 if (header.elt_type == CS_UINT32 || header.elt_type == CS_UINT64) {
1308 if (val_type == CS_TYPE_cs_gnum_t)
1309 cs_io_set_cs_gnum(&header, restart->fh);
1310 else if (val_type == CS_TYPE_int)
1311 cs_io_set_cs_lnum(&header, restart->fh);
1312 }
1313 else if (header.elt_type == CS_FLOAT || header.elt_type == CS_DOUBLE) {
1314 if (sizeof(cs_real_t) != cs_datatype_size[header.elt_type]) {
1315 if (sizeof(cs_real_t) == cs_datatype_size[CS_FLOAT])
1316 header.elt_type = CS_FLOAT;
1317 else
1318 header.elt_type = CS_DOUBLE;
1319 }
1320 }
1321
1322 /* Section contents */
1323 /*------------------*/
1324
1325 /* In single processor mode or for global values */
1326
1327 if (cs_glob_n_ranks == 1 || location_id == 0) {
1328
1329 cs_io_read_global(&header, val, restart->fh);
1330
1331 if (ent_global_num != NULL)
1332 _restart_permute_read(n_ents,
1333 ent_global_num,
1334 _n_location_vals,
1335 val_type,
1336 val);
1337 }
1338
1339 #if defined(HAVE_MPI)
1340
1341 /* In parallel mode for a distributed mesh location */
1342
1343 else if (n_glob_ents > 0)
1344 _read_ent_values(restart,
1345 &header,
1346 n_glob_ents,
1347 n_ents,
1348 ent_global_num,
1349 _n_location_vals,
1350 val_type,
1351 (cs_byte_t *)val);
1352
1353 #endif /* #if defined(HAVE_MPI) */
1354
1355 /* Return */
1356
1357 return CS_RESTART_SUCCESS;
1358 }
1359
1360 /*----------------------------------------------------------------------------*/
1361 /*!
1362 * \brief Write a section to a restart file.
1363 *
1364 * \param[in] restart associated restart file pointer
1365 * \param[in] sec_name section name
1366 * \param[in] location_id id of corresponding location
1367 * \param[in] n_location_vals number of values per location (interlaced)
1368 * \param[in] val_type value type
1369 * \param[in] val array of values
1370 */
1371 /*----------------------------------------------------------------------------*/
1372
1373 static void
_write_section(cs_restart_t * restart,void * context,const char * sec_name,int location_id,int n_location_vals,cs_restart_val_type_t val_type,const void * val)1374 _write_section(cs_restart_t *restart,
1375 void *context,
1376 const char *sec_name,
1377 int location_id,
1378 int n_location_vals,
1379 cs_restart_val_type_t val_type,
1380 const void *val)
1381 {
1382 CS_UNUSED(context);
1383
1384 cs_lnum_t n_ents;
1385 cs_gnum_t n_tot_vals, n_glob_ents;
1386 cs_datatype_t elt_type = CS_DATATYPE_NULL;
1387
1388 const cs_gnum_t *ent_global_num;
1389
1390 cs_lnum_t _n_location_vals = n_location_vals;
1391
1392 assert(restart != NULL);
1393
1394 n_tot_vals = _compute_n_ents(restart, location_id, n_location_vals);
1395
1396 /* Check associated location */
1397
1398 if (location_id == 0) {
1399 n_glob_ents = n_location_vals;
1400 n_ents = n_location_vals;
1401 _n_location_vals = 1;
1402 ent_global_num = NULL;
1403 }
1404
1405 else {
1406 assert(location_id >= 0 && location_id <= (int)(restart->n_locations));
1407 n_glob_ents = (restart->location[location_id-1]).n_glob_ents;
1408 n_ents = (restart->location[location_id-1]).n_ents;
1409 ent_global_num = (restart->location[location_id-1]).ent_global_num;
1410 }
1411
1412 /* Set val_type */
1413
1414 switch (val_type) {
1415 case CS_TYPE_char:
1416 elt_type = CS_CHAR;
1417 break;
1418 case CS_TYPE_int:
1419 elt_type = (sizeof(int) == 8) ? CS_INT64 : CS_INT32;
1420 break;
1421 case CS_TYPE_cs_gnum_t:
1422 elt_type = (sizeof(cs_gnum_t) == 8) ? CS_UINT64 : CS_UINT32;
1423 break;
1424 case CS_TYPE_cs_real_t:
1425 elt_type = (sizeof(cs_real_t) == cs_datatype_size[CS_DOUBLE])
1426 ? CS_DOUBLE : CS_FLOAT;
1427 break;
1428 default:
1429 assert(0);
1430 }
1431
1432 /* Section contents */
1433 /*------------------*/
1434
1435 /* In single processor mode of for global values */
1436
1437 if (location_id == 0)
1438 cs_io_write_global(sec_name,
1439 n_tot_vals,
1440 location_id,
1441 0,
1442 1,
1443 elt_type,
1444 val,
1445 restart->fh);
1446
1447
1448 else if (cs_glob_n_ranks == 1 || n_glob_ents == 0) {
1449
1450 cs_byte_t *val_tmp = NULL;
1451
1452 if (ent_global_num != NULL)
1453 val_tmp = _restart_permute_write(n_ents,
1454 ent_global_num,
1455 _n_location_vals,
1456 val_type,
1457 val);
1458 cs_io_write_global(sec_name,
1459 n_tot_vals,
1460 location_id,
1461 0,
1462 _n_location_vals,
1463 elt_type,
1464 (val_tmp != NULL) ? val_tmp : val,
1465 restart->fh);
1466
1467 if (val_tmp != NULL)
1468 BFT_FREE (val_tmp);
1469 }
1470
1471 #if defined(HAVE_MPI)
1472
1473 /* In parallel mode for a distributed mesh location */
1474
1475 else
1476 _write_ent_values(restart,
1477 sec_name,
1478 n_glob_ents,
1479 n_ents,
1480 ent_global_num,
1481 location_id,
1482 _n_location_vals,
1483 val_type,
1484 (const cs_byte_t *)val);
1485
1486 #endif /* #if defined(HAVE_MPI) */
1487 }
1488
1489 /*----------------------------------------------------------------------------*/
1490 /*!
1491 * \brief Add a location definition with a possible reference definition.
1492 *
1493 * \param[in] restart associated restart file pointer
1494 * \param[in] location_name name associated with the location
1495 * \param[in] n_glob_ents global number of entities
1496 * \param[in] n_ents local number of entities
1497 * \param[in] ent_global_num global entity numbers, or NULL
1498 *
1499 * \return the location id assigned, or -1 in case of error
1500 */
1501 /*----------------------------------------------------------------------------*/
1502
1503 static int
_add_location_check_ref(cs_restart_t * restart,const char * location_name,cs_gnum_t n_glob_ents,cs_lnum_t n_ents,const cs_gnum_t * ent_global_num)1504 _add_location_check_ref(cs_restart_t *restart,
1505 const char *location_name,
1506 cs_gnum_t n_glob_ents,
1507 cs_lnum_t n_ents,
1508 const cs_gnum_t *ent_global_num)
1509 {
1510 int loc_id;
1511
1512 if (restart->mode == CS_RESTART_MODE_READ) {
1513
1514 /* Search for a reference location with the same name */
1515
1516 for (loc_id = 0; loc_id < (int)(_n_locations_ref); loc_id++) {
1517
1518 if ((strcmp((restart->location[loc_id]).name, location_name) == 0)) {
1519
1520 const _location_t *_loc_ref = _location_ref + loc_id;
1521
1522 n_glob_ents = _loc_ref->n_glob_ents;
1523 n_ents = _loc_ref->n_ents;
1524 ent_global_num =_loc_ref->ent_global_num;
1525
1526 }
1527
1528 }
1529
1530 }
1531
1532 return cs_restart_add_location(restart, location_name,
1533 n_glob_ents, n_ents,
1534 ent_global_num);
1535 }
1536
1537 /*----------------------------------------------------------------------------*/
1538 /*!
1539 * \brief Update checkpoint directory with mesh.
1540 *
1541 * The mesh_output file is moved to restart/mesh_input if present.
1542 * Otherwise, if mesh_input is present and a file (not a directory),
1543 * is linked to restart/mesh_input (using a hard link if possible)
1544 */
1545 /*----------------------------------------------------------------------------*/
1546
1547 static void
_update_mesh_checkpoint(void)1548 _update_mesh_checkpoint(void)
1549 {
1550 if (cs_glob_rank_id < 1 && _checkpoint_mesh > 0) {
1551
1552 if (cs_glob_rank_id < 1) {
1553 const char _checkpoint[] = "checkpoint";
1554 if (cs_file_mkdir_default(_checkpoint) != 0)
1555 bft_error(__FILE__, __LINE__, 0,
1556 _("The %s directory cannot be created"), _checkpoint);
1557 }
1558
1559 /* Move mesh_output if present */
1560
1561 const char opath_i[] = "mesh_input.csm";
1562 const char opath_o[] = "mesh_output.csm";
1563 const char npath[] = "checkpoint/mesh_input.csm";
1564
1565 if (cs_file_isreg(opath_o)) {
1566 int retval = rename(opath_o, npath);
1567 if (retval != 0) {
1568 cs_base_warn(__FILE__, __LINE__);
1569 bft_printf(_("Failure moving %s to %s:\n"
1570 "%s\n"),
1571 opath_o, npath, strerror(errno));
1572 }
1573 }
1574
1575 /* Otherwise link mesh_input if it is a regular file
1576 (in case of a mesh_input directory, do nothing, since a
1577 directory should be created only in case of multiple meshes,
1578 and concatenation leads to a mesh_output being created,
1579 unless the (advanced) user has explicitely deactivated that
1580 output, in which case we abide by that choice) */
1581
1582 else if (cs_glob_mesh->modified < 1 && cs_file_isreg(opath_i)) {
1583
1584 #if defined(HAVE_LINKAT) && defined(HAVE_FCNTL_H)
1585
1586 int retval = linkat(AT_FDCWD, opath_i,
1587 AT_FDCWD, npath, AT_SYMLINK_FOLLOW);
1588
1589 if (retval != 0) {
1590 cs_base_warn(__FILE__, __LINE__);
1591 bft_printf(_("Failure hard-linking %s to %s:\n"
1592 "%s\n"),
1593 opath_i, npath, strerror(errno));
1594
1595 }
1596
1597 #endif
1598
1599 }
1600 }
1601 }
1602
1603 /*----------------------------------------------------------------------------*/
1604 /*!
1605 * \brief Create and allocate a new multiwriter structure.
1606 *
1607 * \return pointer to the allocated object.
1608 */
1609 /*----------------------------------------------------------------------------*/
1610
1611 static _restart_multiwriter_t *
_restart_multiwriter_create(void)1612 _restart_multiwriter_create(void)
1613 {
1614 _restart_multiwriter_t *new_writer = NULL;
1615 BFT_MALLOC(new_writer, 1, _restart_multiwriter_t);
1616
1617 new_writer->id = -1;
1618 new_writer->name = NULL;
1619 new_writer->path = NULL;
1620 new_writer->n_prev_files = -1; /* set at 0 after first (single) output */
1621 new_writer->n_prev_files_tot = 0;
1622 new_writer->prev_files = NULL;
1623
1624 return new_writer;
1625 }
1626
1627 /*----------------------------------------------------------------------------*/
1628 /*!
1629 * \brief Get a multiwriter using its id.
1630 *
1631 * Returns a pointer to the writer object. If no writer has that id,
1632 * a NULL pointer is returned.
1633 *
1634 * \param[in] id id of the multiwriter (int)
1635 *
1636 * \return pointer to the multiwriter object.
1637 */
1638 /*----------------------------------------------------------------------------*/
1639
1640 static _restart_multiwriter_t *
_restart_multiwriter_by_id(const int id)1641 _restart_multiwriter_by_id(const int id)
1642 {
1643 _restart_multiwriter_t *writer = NULL;
1644 if (id >= 0 && id < _n_restart_multiwriters)
1645 writer = _restart_multiwriter[id];
1646
1647 return writer;
1648 }
1649
1650 /*----------------------------------------------------------------------------*/
1651 /*!
1652 * \brief Add a multiwriter based on a file name and path.
1653 *
1654 * Returns the id of the multiwriter.
1655 *
1656 * \param[in] name name of the restart file (char *)
1657 * \param[in] path path to the restart file (char *)
1658 *
1659 * \return id of the newly added multiwriter
1660 */
1661 /*----------------------------------------------------------------------------*/
1662
1663 static int
_add_restart_multiwriter(const char name[],const char path[])1664 _add_restart_multiwriter(const char name[],
1665 const char path[])
1666 {
1667 /* Check if the writer already exists */
1668 for (int i = 0; i < _n_restart_multiwriters; i++) {
1669 if (strcmp(_restart_multiwriter[i]->name, name) == 0)
1670 return i;
1671 }
1672
1673 /* Check that the file name is neither NULL or empty ("") */
1674 if (name == NULL)
1675 bft_error(__FILE__, __LINE__, 0,
1676 _("Null pointer provided as file name.\n"));
1677 else if (strlen(name) == 0)
1678 bft_error(__FILE__, __LINE__, 0,
1679 _("Empty file name was provided.\n"));
1680
1681 /* Allocate or reallocate the array */
1682 if (_n_restart_multiwriters == 0)
1683 BFT_MALLOC(_restart_multiwriter,
1684 _n_restart_multiwriters + 1,
1685 _restart_multiwriter_t *);
1686 else
1687 BFT_REALLOC(_restart_multiwriter,
1688 _n_restart_multiwriters + 1,
1689 _restart_multiwriter_t *);
1690
1691 /* Create the empty structure */
1692 _restart_multiwriter_t *new_writer = _restart_multiwriter_create();
1693
1694 /* Set id */
1695 new_writer->id = _n_restart_multiwriters;
1696
1697 /* Set name */
1698 size_t lname = strlen(name) + 1;
1699 BFT_MALLOC(new_writer->name, lname, char);
1700 strcpy(new_writer->name, name);
1701
1702 /* Set path to subdir, which can be NULL */
1703 const char *_path = path;
1704 if (_path != NULL) {
1705 if (strlen(_path) == 0)
1706 _path = NULL;
1707 }
1708
1709 if (_path != NULL) {
1710 size_t lpath = strlen(_path) + 1;
1711 BFT_MALLOC(new_writer->path, lpath, char);
1712 strcpy(new_writer->path, _path);
1713 }
1714
1715 _restart_multiwriter[_n_restart_multiwriters] = new_writer;
1716
1717 /* Increment the number of writers, and return the newly created
1718 * writer index */
1719
1720 _n_restart_multiwriters++;
1721
1722 return new_writer->id;
1723 }
1724
1725 /*----------------------------------------------------------------------------*/
1726 /*!
1727 * \brief Increment the multiwriter previous files list.
1728 *
1729 * \param[in] mw pointer to the multiwriter object.
1730 * \param[in] fname name of the file to add to the old files list
1731 *
1732 */
1733 /*----------------------------------------------------------------------------*/
1734
1735 static void
_restart_multiwriter_increment(_restart_multiwriter_t * mw,const char fname[])1736 _restart_multiwriter_increment(_restart_multiwriter_t *mw,
1737 const char fname[])
1738 {
1739 mw->n_prev_files++;
1740 mw->n_prev_files_tot++;
1741
1742 if (mw->prev_files == NULL)
1743 BFT_MALLOC(mw->prev_files, mw->n_prev_files, char *);
1744 else
1745 BFT_REALLOC(mw->prev_files, mw->n_prev_files, char *);
1746
1747 mw->prev_files[mw->n_prev_files - 1] = NULL;
1748 size_t lenf = strlen(fname) + 1;
1749 BFT_MALLOC(mw->prev_files[mw->n_prev_files - 1], lenf, char);
1750 strcpy(mw->prev_files[mw->n_prev_files - 1], fname);
1751 }
1752
1753 /*! (DOXYGEN_SHOULD_SKIP_THIS) \endcond */
1754
1755 /*============================================================================
1756 * Public Fortran function definitions
1757 *============================================================================*/
1758
1759 /*----------------------------------------------------------------------------
1760 * Indicate if a restart directory is present.
1761 *
1762 * Fortran interface
1763 *
1764 * subroutine dflsui (ntsuit, ttsuit, wtsuit)
1765 * *****************
1766 *
1767 * integer ntsuit : <-- : > 0: checkpoint time step interval
1768 * : : 0: default interval
1769 * : : -1: checkpoint at end
1770 * : : -2: no checkpoint
1771 * double precision ttsuit : <-- : if> 0, checkpoint time interval
1772 * double precision wtsuit : <-- : if> 0, checkpoint wall time interval
1773 *----------------------------------------------------------------------------*/
1774
CS_PROCF(dflsui,DFLSUI)1775 void CS_PROCF (dflsui, DFLSUI)
1776 (
1777 int *ntsuit,
1778 cs_real_t *ttsuit,
1779 cs_real_t *wtsuit
1780 )
1781 {
1782 cs_restart_checkpoint_set_defaults(*ntsuit, *ttsuit, *wtsuit);
1783 }
1784
1785 /*----------------------------------------------------------------------------
1786 * Check if checkpointing is recommended at a given time.
1787 *
1788 * Fortran interface
1789 *
1790 * subroutine reqsui (iisuit)
1791 * *****************
1792 *
1793 * integer iisuit : --> : 0 if no restart required, 1 otherwise
1794 *----------------------------------------------------------------------------*/
1795
CS_PROCF(reqsui,REQSUI)1796 void CS_PROCF (reqsui, REQSUI)
1797 (
1798 int *iisuit
1799 )
1800 {
1801 if (cs_restart_checkpoint_required(cs_glob_time_step))
1802 *iisuit = 1;
1803 else
1804 *iisuit = 0;
1805 }
1806
1807 /*----------------------------------------------------------------------------
1808 * Indicate checkpointing has been done at a given time.
1809 *
1810 * This updates the status for future checks to determine
1811 * if checkpointing is recommended at a given time.
1812 *
1813 * Fortran interface
1814 *
1815 * subroutine stusui
1816 * *****************
1817 *----------------------------------------------------------------------------*/
1818
CS_PROCF(stusui,STUSUI)1819 void CS_PROCF (stusui, STUSUI)
1820 (
1821 void
1822 )
1823 {
1824 cs_restart_checkpoint_done(cs_glob_time_step);
1825 }
1826
1827 /*----------------------------------------------------------------------------
1828 * Save output mesh for turbomachinery if needed
1829 *
1830 * Fortran interface
1831 *
1832 * subroutine trbsui
1833 * *****************
1834 *----------------------------------------------------------------------------*/
1835
CS_PROCF(trbsui,TRBSUI)1836 void CS_PROCF (trbsui, TRBSUI)
1837 (
1838 void
1839 )
1840 {
1841 cs_mesh_save(cs_glob_mesh, NULL, "checkpoint", "mesh");
1842 }
1843
1844 /*----------------------------------------------------------------------------
1845 * Indicate if a restart directory is present.
1846 *
1847 * Fortran interface
1848 *
1849 * subroutine indsui (isuite)
1850 * *****************
1851 *
1852 * integer isuite : --> : 1 for restart, 0 otherwise
1853 *----------------------------------------------------------------------------*/
1854
CS_PROCF(indsui,INDSUI)1855 void CS_PROCF (indsui, INDSUI)
1856 (
1857 int *isuite
1858 )
1859 {
1860 *isuite = cs_restart_present();
1861 }
1862
1863 /*============================================================================
1864 * Public function definitions
1865 *============================================================================*/
1866
1867 /*----------------------------------------------------------------------------*/
1868 /*!
1869 * \brief Define default checkpoint interval.
1870 *
1871 * \param[in] nt_interval if > 0 time step interval for checkpoint
1872 * if 0, default of 4 checkpoints per run
1873 * if -1, checkpoint at end
1874 * if -2, no checkpointing
1875 * \param[in] t_interval if > 0, time value interval for checkpoint
1876 * \param[in] wt_interval if > 0, wall-clock interval for checkpoints
1877 */
1878 /*----------------------------------------------------------------------------*/
1879
1880 void
cs_restart_checkpoint_set_defaults(int nt_interval,double t_interval,double wt_interval)1881 cs_restart_checkpoint_set_defaults(int nt_interval,
1882 double t_interval,
1883 double wt_interval)
1884 {
1885 _checkpoint_nt_interval = nt_interval;
1886 _checkpoint_t_interval = t_interval;
1887 _checkpoint_wt_interval = wt_interval;
1888 }
1889
1890 /*----------------------------------------------------------------------------*/
1891 /*!
1892 * \brief Define checkpoint behavior for mesh.
1893 *
1894 * If mesh checkpointing is active (default), upon writing the first
1895 * checkpoint file of a computation, a mesh_output file is moved to
1896 * checkpoint/mesh_input if present. If not present but a mesh_input
1897 * file (or link to file) is present, a hard link to that file is
1898 * added as checkpoint/mesh_input.
1899 *
1900 * A mesh_input directory is ignored, as it is normally only created when
1901 * multiple input files are appended, which leads to the output and thus
1902 * presence of a mesh_output file, unless explicitely deactivated by the
1903 * user.
1904 *
1905 * \param[in] mode if 0, do not checkpoint mesh
1906 * if 1, checkpoint mesh_output or mesh_input file
1907 */
1908 /*----------------------------------------------------------------------------*/
1909
1910 void
cs_restart_checkpoint_set_mesh_mode(int mode)1911 cs_restart_checkpoint_set_mesh_mode(int mode)
1912 {
1913 _checkpoint_mesh = mode;
1914 }
1915
1916 /*----------------------------------------------------------------------------*/
1917 /*!
1918 * \brief Define last forced checkpoint time step.
1919 *
1920 * \param[in] nt_last last time step for forced checkpoint
1921 */
1922 /*----------------------------------------------------------------------------*/
1923
1924 void
cs_restart_checkpoint_set_last_ts(int nt_last)1925 cs_restart_checkpoint_set_last_ts(int nt_last)
1926 {
1927 _checkpoint_nt_last = nt_last;
1928 }
1929
1930 /*----------------------------------------------------------------------------*/
1931 /*!
1932 * \brief Define next forced checkpoint time step.
1933 *
1934 * \param[in] nt_next next time step for forced checkpoint
1935 */
1936 /*----------------------------------------------------------------------------*/
1937
1938 void
cs_restart_checkpoint_set_next_ts(int nt_next)1939 cs_restart_checkpoint_set_next_ts(int nt_next)
1940 {
1941 _checkpoint_nt_next = nt_next;
1942 }
1943
1944 /*----------------------------------------------------------------------------*/
1945 /*!
1946 * \brief Define next forced checkpoint time value.
1947 *
1948 * \param[in] t_next next time value for forced checkpoint
1949 */
1950 /*----------------------------------------------------------------------------*/
1951
1952 void
cs_restart_checkpoint_set_next_tv(double t_next)1953 cs_restart_checkpoint_set_next_tv(double t_next)
1954 {
1955 _checkpoint_t_next = t_next;
1956 }
1957
1958 /*----------------------------------------------------------------------------*/
1959 /*!
1960 * \brief Define next forced checkpoint wall-clock time value.
1961 *
1962 * \param[in] wt_next next wall-clock time value for forced checkpoint
1963 */
1964 /*----------------------------------------------------------------------------*/
1965
1966 void
cs_restart_checkpoint_set_next_wt(double wt_next)1967 cs_restart_checkpoint_set_next_wt(double wt_next)
1968 {
1969 _checkpoint_wt_next = wt_next;
1970 }
1971
1972 /*----------------------------------------------------------------------------*/
1973 /*!
1974 * \brief Check if checkpointing is recommended at a given time.
1975 *
1976 * \param[in] ts time step status structure
1977 *
1978 * \return true if checkpointing is recommended, false otherwise
1979 */
1980 /*----------------------------------------------------------------------------*/
1981
1982 bool
cs_restart_checkpoint_required(const cs_time_step_t * ts)1983 cs_restart_checkpoint_required(const cs_time_step_t *ts)
1984 {
1985 assert(ts != NULL);
1986
1987 int nt = ts->nt_cur - ts->nt_prev;
1988 double t = ts->t_cur - ts->t_prev;
1989
1990 bool retval = false;
1991
1992 if (_checkpoint_nt_interval > CS_RESTART_INTERVAL_NONE) {
1993
1994 if (ts->nt_cur == ts->nt_max) /* Output at the last time step */
1995 retval = true;
1996
1997 else if (_checkpoint_nt_interval == CS_RESTART_INTERVAL_DEFAULT) {
1998 /* default interval: current number of expected time_steps for this run,
1999 with a minimum of 10. */
2000 int nt_def = (ts->nt_max - ts->nt_prev)/4;
2001 if (nt_def < 10)
2002 nt_def = 10;
2003 if (nt % nt_def == 0)
2004 retval = true;
2005 }
2006
2007 else if (_checkpoint_nt_interval > CS_RESTART_INTERVAL_DEFAULT &&
2008 nt % _checkpoint_nt_interval == 0)
2009 retval = true;
2010
2011 else if (_checkpoint_nt_interval > CS_RESTART_INTERVAL_DEFAULT &&
2012 _checkpoint_nt_last > -1) {
2013 if (ts->nt_cur >= _checkpoint_nt_interval + _checkpoint_nt_last)
2014 retval = true;
2015 }
2016
2017 }
2018
2019 if (_checkpoint_t_interval > 0
2020 && _checkpoint_t_last + _checkpoint_t_interval <= t)
2021 retval = true;
2022
2023 else if (_checkpoint_wt_next >= 0) {
2024 double wt = cs_timer_wtime();
2025 if (wt >= _checkpoint_wt_next)
2026 retval = true;
2027 }
2028
2029 else if ( (_checkpoint_nt_next >= 0 && _checkpoint_nt_next <= ts->nt_cur)
2030 || (_checkpoint_t_next >= 0 && _checkpoint_t_next <= ts->t_cur))
2031 retval = true;
2032
2033 else if (_checkpoint_wt_interval >= 0) {
2034 double wt = cs_timer_wtime();
2035 if (wt - _checkpoint_wt_last >= _checkpoint_wt_interval)
2036 retval = true;
2037 }
2038
2039 return retval;
2040 }
2041
2042 /*----------------------------------------------------------------------------*/
2043 /*!
2044 * \brief Indicate checkpointing has been done at a given time.
2045 *
2046 * This updates the status for future checks to determine
2047 * if checkpointing is recommended at a given time.
2048 *
2049 * \param[in] ts time step status structure
2050 */
2051 /*----------------------------------------------------------------------------*/
2052
2053 void
cs_restart_checkpoint_done(const cs_time_step_t * ts)2054 cs_restart_checkpoint_done(const cs_time_step_t *ts)
2055 {
2056 assert(ts != NULL);
2057
2058 double t = ts->t_cur - ts->t_prev;
2059
2060 if (_checkpoint_nt_next >= 0 && _checkpoint_nt_next <= ts->nt_cur)
2061 _checkpoint_nt_next = -1;
2062
2063 if (_checkpoint_t_next >= 0 && _checkpoint_t_next <= ts->t_cur)
2064 _checkpoint_t_next = -1.;
2065
2066 if (_checkpoint_wt_next >= 0) {
2067 double wt = cs_timer_wtime();
2068 if (wt >= _checkpoint_wt_next)
2069 _checkpoint_wt_next = -1.;
2070 }
2071
2072 if (_checkpoint_t_interval > 0
2073 && _checkpoint_t_last + _checkpoint_t_interval <= t)
2074 _checkpoint_t_last = ts->t_cur;
2075
2076 if (_checkpoint_wt_interval >= 0) {
2077 double wt = cs_timer_wtime();
2078 if (wt - _checkpoint_wt_last >= _checkpoint_wt_interval)
2079 _checkpoint_wt_last = cs_timer_wtime();
2080 }
2081 }
2082
2083 /*----------------------------------------------------------------------------*/
2084 /*!
2085 * \brief Check if we have a restart directory.
2086 *
2087 * \return 1 if a restart directory is present, 0 otherwise
2088 */
2089 /*----------------------------------------------------------------------------*/
2090
2091 int
cs_restart_present(void)2092 cs_restart_present(void)
2093 {
2094 if (_restart_present < 0) {
2095 if (cs_glob_rank_id < 1) {
2096 if (cs_file_isdir("restart"))
2097 _restart_present = 1;
2098 else
2099 _restart_present = 0;
2100 }
2101 cs_parall_bcast(0, 1, CS_INT_TYPE, &_restart_present);
2102 }
2103
2104 return _restart_present;
2105 }
2106
2107 /*----------------------------------------------------------------------------*/
2108 /*!
2109 * \brief Initialize a restart file.
2110 *
2111 * \param[in] name file name
2112 * \param[in] path optional directory name for output, or NULL for default
2113 * (directory automatically created if necessary)
2114 * \param[in] mode read or write
2115 *
2116 * \returns pointer to initialized restart file structure
2117 */
2118 /*----------------------------------------------------------------------------*/
2119
2120 cs_restart_t *
cs_restart_create(const char * name,const char * path,cs_restart_mode_t mode)2121 cs_restart_create(const char *name,
2122 const char *path,
2123 cs_restart_mode_t mode)
2124 {
2125 cs_restart_t * restart;
2126
2127 double timing[2];
2128
2129 char *_name = NULL;
2130 size_t ldir, lname, lext;
2131
2132 const char *_path = path;
2133 const char _restart[] = "restart";
2134 const char _checkpoint[] = "checkpoint";
2135 const char _extension[] = ".csc";
2136
2137 const cs_mesh_t *mesh = cs_glob_mesh;
2138
2139 /* Ensure mesh checkpoint is updated on first call */
2140
2141 if ( mode == CS_RESTART_MODE_WRITE
2142 && _restart_n_opens[mode] == 0) {
2143 _update_mesh_checkpoint();
2144 }
2145
2146 /* Initializations */
2147
2148 timing[0] = cs_timer_wtime();
2149
2150 if (_path != NULL) {
2151 if (strlen(_path) == 0)
2152 _path = NULL;
2153 }
2154
2155 if (_path == NULL) {
2156 if (mode == CS_RESTART_MODE_WRITE)
2157 _path = _checkpoint;
2158 else
2159 _path = _restart;
2160 }
2161
2162 /* Create 'checkpoint' directory or read from 'restart' directory */
2163
2164 if (cs_glob_rank_id < 1) {
2165 if (mode == CS_RESTART_MODE_WRITE) {
2166 if (cs_file_mkdir_default(_path) != 0)
2167 bft_error(__FILE__, __LINE__, 0,
2168 _("The %s directory cannot be created"), _path);
2169 }
2170 else if (mode == CS_RESTART_MODE_READ) {
2171 if (cs_file_isdir(_path) == 0) {
2172 bft_error(__FILE__, __LINE__, 0,
2173 _("The %s directory cannot be found"), _path);
2174 }
2175 }
2176 }
2177 #if defined(HAVE_MPI)
2178 if (cs_glob_n_ranks > 1)
2179 MPI_Barrier(cs_glob_mpi_comm);
2180 #endif
2181
2182 ldir = strlen(_path);
2183 lname = strlen(name);
2184
2185 BFT_MALLOC(_name, ldir + lname + 2, char);
2186
2187 strcpy(_name, _path);
2188 _name[ldir] = _dir_separator;
2189 _name[ldir+1] = '\0';
2190 strcat(_name, name);
2191 _name[ldir+lname+1] = '\0';
2192
2193 /* Following the addition of an extension, we check for READ mode
2194 * if a file exists without the extension */
2195
2196 if (mode == CS_RESTART_MODE_READ) {
2197
2198 if (cs_file_isreg(_name) == 0 && cs_file_endswith(name, _extension)) {
2199 BFT_FREE(_name);
2200
2201 lext = strlen(_extension);
2202 BFT_MALLOC(_name, ldir + lname + 2 - lext, char);
2203 strcpy(_name, _path);
2204 _name[ldir] = _dir_separator;
2205 _name[ldir+1] = '\0';
2206 strncat(_name, name, lname-lext);
2207 _name[ldir+lname-lext+1] = '\0';
2208 }
2209
2210 } else if (mode == CS_RESTART_MODE_WRITE) {
2211
2212 /* Check if file already exists, and if so rename and delete if needed */
2213 int writer_id = _add_restart_multiwriter(name, _name);
2214 _restart_multiwriter_t *mw = _restart_multiwriter_by_id(writer_id);
2215
2216 /* Rename an already existing file */
2217 if (cs_file_isreg(_name) && mw->n_prev_files > -1) {
2218
2219 char _subdir[19];
2220 sprintf(_subdir, "previous_dump_%04d", mw->n_prev_files_tot);
2221 size_t lsdir = strlen(_subdir);
2222
2223 char *_re_name = NULL;
2224 BFT_MALLOC(_re_name, ldir + lsdir + lname + 3, char);
2225
2226 strcpy(_re_name, _path);
2227 _re_name[ldir] = _dir_separator;
2228 _re_name[ldir+1] = '\0';
2229
2230 strcat(_re_name, _subdir);
2231
2232 /* Check that the sub-directory exists or can be created */
2233 if (cs_glob_rank_id < 1) {
2234 if (cs_file_mkdir_default(_re_name) != 0)
2235 bft_error(__FILE__, __LINE__, 0,
2236 _("The %s directory cannot be created"), _re_name);
2237 }
2238 #if defined(HAVE_MPI)
2239 if (cs_glob_n_ranks > 1)
2240 MPI_Barrier(cs_glob_mpi_comm);
2241 #endif
2242
2243 _re_name[ldir+lsdir+1] = _dir_separator;
2244 _re_name[ldir+lsdir+2] = '\0';
2245 strcat(_re_name, name);
2246 _re_name[ldir+lsdir+lname+2] = '\0';
2247
2248 rename(_name, _re_name);
2249
2250 _restart_multiwriter_increment(mw, _re_name);
2251
2252 BFT_FREE(_re_name);
2253 }
2254 else
2255 mw->n_prev_files = 0;
2256 }
2257
2258 /* Allocate and initialize base structure */
2259
2260 BFT_MALLOC(restart, 1, cs_restart_t);
2261
2262 BFT_MALLOC(restart->name, strlen(_name) + 1, char);
2263
2264 strcpy(restart->name, _name);
2265
2266 BFT_FREE(_name);
2267
2268 /* Initialize other fields */
2269
2270 restart->mode = mode;
2271
2272 restart->fh = NULL;
2273
2274 restart->rank_step = 1;
2275 restart->min_block_size = 0;
2276
2277 /* Initialize location data */
2278
2279 restart->n_locations = 0;
2280 restart->location = NULL;
2281
2282 /* Open associated file, and build an index of sections in read mode */
2283
2284 _add_file(restart);
2285
2286 /* Add basic location definitions */
2287
2288 _add_location_check_ref(restart, "cells",
2289 mesh->n_g_cells, mesh->n_cells,
2290 mesh->global_cell_num);
2291 _add_location_check_ref(restart, "interior_faces",
2292 mesh->n_g_i_faces, mesh->n_i_faces,
2293 mesh->global_i_face_num);
2294 _add_location_check_ref(restart, "boundary_faces",
2295 mesh->n_g_b_faces, mesh->n_b_faces,
2296 mesh->global_b_face_num);
2297 _add_location_check_ref(restart, "vertices",
2298 mesh->n_g_vertices, mesh->n_vertices,
2299 mesh->global_vtx_num);
2300
2301 timing[1] = cs_timer_wtime();
2302 _restart_wtime[mode] += timing[1] - timing[0];
2303
2304 return restart;
2305 }
2306
2307 /*----------------------------------------------------------------------------*/
2308 /*!
2309 * \brief Destroy structure associated with a restart file (and close the file).
2310 *
2311 * \param[in, out] restart pointer to restart file structure pointer
2312 */
2313 /*----------------------------------------------------------------------------*/
2314
2315 void
cs_restart_destroy(cs_restart_t ** restart)2316 cs_restart_destroy(cs_restart_t **restart)
2317 {
2318 cs_restart_mode_t mode;
2319
2320 cs_restart_t *r = *restart;
2321
2322 double timing[2];
2323
2324 timing[0] = cs_timer_wtime();
2325
2326 assert(restart != NULL);
2327
2328 mode = r->mode;
2329
2330 if (r->fh != NULL)
2331 cs_io_finalize(&(r->fh));
2332
2333 /* Free locations array */
2334
2335 if (r->n_locations > 0) {
2336 size_t loc_id;
2337 for (loc_id = 0; loc_id < r->n_locations; loc_id++) {
2338 BFT_FREE((r->location[loc_id]).name);
2339 BFT_FREE((r->location[loc_id])._ent_global_num);
2340 }
2341 }
2342 if (r->location != NULL)
2343 BFT_FREE(r->location);
2344
2345 /* Free remaining memory */
2346
2347 BFT_FREE(r->name);
2348
2349 BFT_FREE(*restart);
2350
2351 timing[1] = cs_timer_wtime();
2352 _restart_wtime[mode] += timing[1] - timing[0];
2353 }
2354
2355 /*----------------------------------------------------------------------------*/
2356 /*!
2357 * \brief Check the locations associated with a restart file.
2358 *
2359 * For each type of entity, the corresponding flag is set to true if the
2360 * associated number of entities matches the current value (and so that we
2361 * consider the mesh locations are the same), false otherwise.
2362 *
2363 * \param[out] restart associated restart file pointer
2364 * \param[out] match_cell matching cells flag
2365 * \param[out] match_i_face matching interior faces flag
2366 * \param[out] match_b_face matching boundary faces flag
2367 * \param[out] match_vertex matching vertices flag
2368 */
2369 /*----------------------------------------------------------------------------*/
2370
2371 void
cs_restart_check_base_location(const cs_restart_t * restart,bool * match_cell,bool * match_i_face,bool * match_b_face,bool * match_vertex)2372 cs_restart_check_base_location(const cs_restart_t *restart,
2373 bool *match_cell,
2374 bool *match_i_face,
2375 bool *match_b_face,
2376 bool *match_vertex)
2377 {
2378 size_t location_id;
2379
2380 *match_cell = false;
2381 *match_i_face = false;
2382 *match_b_face = false;
2383 *match_vertex = false;
2384
2385 assert(restart != NULL);
2386
2387 for (location_id = 0; location_id < 4; location_id++) {
2388
2389 const _location_t *loc = restart->location + location_id;
2390
2391 if (loc->n_glob_ents_f == loc->n_glob_ents) {
2392 if (location_id == 0)
2393 *match_cell = true;
2394 else if (location_id == 1)
2395 *match_i_face = true;
2396 else if (location_id == 2)
2397 *match_b_face = true;
2398 else if (location_id == 3)
2399 *match_vertex = true;
2400 }
2401
2402 else if (cs_glob_rank_id <= 0) {
2403 cs_base_warn(__FILE__, __LINE__);
2404 bft_printf(_("The size of location \"%s\" associated with\n"
2405 "the restart file \"%s\" is %llu and does not\n"
2406 "correspond to that of the current mesh (%llu).\n"),
2407 loc->name, restart->name,
2408 (unsigned long long)loc->n_glob_ents_f,
2409 (unsigned long long)loc->n_glob_ents);
2410 }
2411
2412 }
2413 }
2414
2415 /*----------------------------------------------------------------------------*/
2416 /*!
2417 * \brief Add a location definition.
2418 *
2419 * \param[in] restart associated restart file pointer
2420 * \param[in] location_name name associated with the location
2421 * \param[in] n_glob_ents global number of entities
2422 * \param[in] n_ents local number of entities
2423 * \param[in] ent_global_num global entity numbers, or NULL
2424 *
2425 * \return the location id assigned, or -1 in case of error
2426 */
2427 /*----------------------------------------------------------------------------*/
2428
2429 int
cs_restart_add_location(cs_restart_t * restart,const char * location_name,cs_gnum_t n_glob_ents,cs_lnum_t n_ents,const cs_gnum_t * ent_global_num)2430 cs_restart_add_location(cs_restart_t *restart,
2431 const char *location_name,
2432 cs_gnum_t n_glob_ents,
2433 cs_lnum_t n_ents,
2434 const cs_gnum_t *ent_global_num)
2435 {
2436 double timing[2];
2437
2438 int loc_id;
2439
2440 timing[0] = cs_timer_wtime();
2441
2442 if (restart->mode == CS_RESTART_MODE_READ) {
2443
2444 /* Search for a location with the same name */
2445
2446 for (loc_id = 0; loc_id < (int)(restart->n_locations); loc_id++) {
2447
2448 if ((strcmp((restart->location[loc_id]).name, location_name) == 0)) {
2449
2450 (restart->location[loc_id]).n_glob_ents = n_glob_ents;
2451
2452 (restart->location[loc_id]).n_ents = n_ents;
2453 (restart->location[loc_id]).ent_global_num = ent_global_num;
2454 (restart->location[loc_id])._ent_global_num = NULL;
2455
2456 timing[1] = cs_timer_wtime();
2457 _restart_wtime[restart->mode] += timing[1] - timing[0];
2458
2459 return loc_id + 1;
2460
2461 }
2462 }
2463
2464 if (loc_id >= ((int)(restart->n_locations)))
2465 bft_error(__FILE__, __LINE__, 0,
2466 _("The restart file \"%s\" references no\n"
2467 "location named \"%s\"."),
2468 restart->name, location_name);
2469
2470 }
2471
2472 else {
2473
2474 cs_datatype_t gnum_type
2475 = (sizeof(cs_gnum_t) == 8) ? CS_UINT64 : CS_UINT32;
2476
2477 /* Create a new location */
2478
2479 restart->n_locations += 1;
2480
2481 BFT_REALLOC(restart->location, restart->n_locations, _location_t);
2482 BFT_MALLOC((restart->location[restart->n_locations-1]).name,
2483 strlen(location_name)+1,
2484 char);
2485
2486 strcpy((restart->location[restart->n_locations-1]).name, location_name);
2487
2488 (restart->location[restart->n_locations-1]).id = restart->n_locations;
2489 (restart->location[restart->n_locations-1]).n_glob_ents = n_glob_ents;
2490 (restart->location[restart->n_locations-1]).n_glob_ents_f = n_glob_ents;
2491 (restart->location[restart->n_locations-1]).n_ents = n_ents;
2492 (restart->location[restart->n_locations-1]).ent_global_num = ent_global_num;
2493 (restart->location[restart->n_locations-1])._ent_global_num = NULL;
2494
2495 cs_io_write_global(location_name, 1, restart->n_locations, 0, 0,
2496 gnum_type, &n_glob_ents,
2497 restart->fh);
2498
2499 timing[1] = cs_timer_wtime();
2500 _restart_wtime[restart->mode] += timing[1] - timing[0];
2501
2502 return restart->n_locations;
2503 }
2504
2505 timing[1] = cs_timer_wtime();
2506 _restart_wtime[restart->mode] += timing[1] - timing[0];
2507
2508 return -1;
2509 }
2510
2511 /*----------------------------------------------------------------------------*/
2512 /*!
2513 * \brief Add a reference location definition with a private copy.
2514 *
2515 * \param[in] location_name name associated with the location
2516 * \param[in] n_glob_ents global number of entities
2517 * \param[in] n_ents local number of entities
2518 * \param[in] ent_global_num global entity numbers, or NULL
2519 *
2520 * \return the location id assigned, or -1 in case of error
2521 */
2522 /*----------------------------------------------------------------------------*/
2523
2524 void
cs_restart_add_location_ref(const char * location_name,cs_gnum_t n_glob_ents,cs_lnum_t n_ents,const cs_gnum_t * ent_global_num)2525 cs_restart_add_location_ref(const char *location_name,
2526 cs_gnum_t n_glob_ents,
2527 cs_lnum_t n_ents,
2528 const cs_gnum_t *ent_global_num)
2529 {
2530 /* Create a new location */
2531
2532 _n_locations_ref += 1;
2533
2534 BFT_REALLOC(_location_ref, _n_locations_ref, _location_t);
2535 BFT_MALLOC((_location_ref[_n_locations_ref-1]).name,
2536 strlen(location_name)+1,
2537 char);
2538
2539 strcpy((_location_ref[_n_locations_ref-1]).name, location_name);
2540
2541 if (ent_global_num != NULL) {
2542 BFT_MALLOC((_location_ref[_n_locations_ref-1])._ent_global_num,
2543 n_ents, cs_gnum_t);
2544 for (cs_lnum_t i = 0; i < n_ents; i++) {
2545 (_location_ref[_n_locations_ref-1])._ent_global_num[i]
2546 = ent_global_num[i];
2547 }
2548 }
2549 else
2550 (_location_ref[_n_locations_ref-1])._ent_global_num = NULL;
2551
2552 (_location_ref[_n_locations_ref-1]).id = _n_locations_ref;
2553 (_location_ref[_n_locations_ref-1]).n_glob_ents = n_glob_ents;
2554 (_location_ref[_n_locations_ref-1]).n_glob_ents_f = n_glob_ents;
2555 (_location_ref[_n_locations_ref-1]).n_ents = n_ents;
2556 (_location_ref[_n_locations_ref-1]).ent_global_num
2557 = (_location_ref[_n_locations_ref-1])._ent_global_num;
2558 }
2559
2560 /*----------------------------------------------------------------------------*/
2561 /*!
2562 * \brief Clear reference location definitions with a private copy.
2563 */
2564 /*----------------------------------------------------------------------------*/
2565
2566 void
cs_restart_clear_locations_ref(void)2567 cs_restart_clear_locations_ref(void)
2568 {
2569 for (size_t loc_id = 0; loc_id < _n_locations_ref; loc_id++) {
2570 BFT_FREE((_location_ref[loc_id]).name);
2571 BFT_FREE((_location_ref[loc_id])._ent_global_num);
2572 }
2573 BFT_FREE(_location_ref);
2574 _n_locations_ref = 0;
2575 }
2576
2577 /*----------------------------------------------------------------------------*/
2578 /*!
2579 * \brief Associate a context to restart section check operations.
2580 *
2581 * This context may be used by the \ref cs_restart_check_section_t,
2582 * \ref cs_restart_read_section_t, and \ref cs_restart_write_section_t
2583 * type functions.
2584 *
2585 * Note that the lifecycle of the data pointed to must be handled separately
2586 * (and the pointer must remain valid until the matching restart structure
2587 * is destroyed.
2588 *
2589 * \param[in] context pointer to associated data, or NULL
2590 */
2591 /*----------------------------------------------------------------------------*/
2592
2593 void
cs_restart_set_context(void * context)2594 cs_restart_set_context(void *context)
2595 {
2596 _restart_context = context;
2597 }
2598
2599 /*----------------------------------------------------------------------------*/
2600 /*!
2601 * \brief Associate a function to restart section check operations.
2602 *
2603 * This allows defining alternate operations when checking restart sections.
2604 *
2605 * \param[in] func associated function
2606 *
2607 * \return pointer to previous function
2608 */
2609 /*----------------------------------------------------------------------------*/
2610
2611 cs_restart_check_section_t *
cs_restart_set_check_section_func(cs_restart_check_section_t * func)2612 cs_restart_set_check_section_func(cs_restart_check_section_t *func)
2613 {
2614 cs_restart_check_section_t *p = _check_section_f;
2615 _check_section_f = func;
2616
2617 return p;
2618 }
2619
2620 /*----------------------------------------------------------------------------*/
2621 /*!
2622 * \brief Associate a function and its input to all restart section
2623 * read operations.
2624 *
2625 * This allows defining alternate operations when reading restart sections.
2626 *
2627 * \param[in] func associated function
2628 *
2629 * \return pointer to previous function
2630 */
2631 /*----------------------------------------------------------------------------*/
2632
2633 cs_restart_read_section_t *
cs_restart_set_read_section_func(cs_restart_read_section_t * func)2634 cs_restart_set_read_section_func(cs_restart_read_section_t *func)
2635 {
2636 cs_restart_read_section_t *p = _read_section_f;
2637 _read_section_f = func;
2638
2639 return p;
2640 }
2641
2642 /*----------------------------------------------------------------------------*/
2643 /*!
2644 * \brief Associate a function and its input to all restart section
2645 * write operations.
2646 *
2647 * This allows defining alternate operations when writing restart sections.
2648 *
2649 * \param[in] func associated hook function
2650 *
2651 * \return pointer to previous function
2652 */
2653 /*----------------------------------------------------------------------------*/
2654
2655 cs_restart_write_section_t *
cs_restart_set_write_section_func(cs_restart_write_section_t * func)2656 cs_restart_set_write_section_func(cs_restart_write_section_t *func)
2657 {
2658 cs_restart_write_section_t *p = _write_section_f;
2659 _write_section_f = func;
2660
2661 return p;
2662 }
2663
2664 /*----------------------------------------------------------------------------*/
2665 /*!
2666 * \brief Return name of restart file
2667 *
2668 * \param[in] restart associated restart file pointer
2669 *
2670 * \return base name of restart file
2671 */
2672 /*----------------------------------------------------------------------------*/
2673
2674 const char *
cs_restart_get_name(const cs_restart_t * restart)2675 cs_restart_get_name(const cs_restart_t *restart)
2676 {
2677 assert(restart != NULL);
2678
2679 return restart->name;
2680 }
2681
2682 /*----------------------------------------------------------------------------*/
2683 /*!
2684 * \brief Return local number of elements associated with a
2685 * given restart location.
2686 *
2687 * \param[in] restart associated restart file pointer
2688 * \param[in] location_id id of corresponding location
2689 *
2690 * \return number of elements associated with location.
2691 */
2692 /*----------------------------------------------------------------------------*/
2693
2694 cs_lnum_t
cs_restart_get_n_location_elts(const cs_restart_t * restart,int location_id)2695 cs_restart_get_n_location_elts(const cs_restart_t *restart,
2696 int location_id)
2697 {
2698 cs_lnum_t retval = 0;
2699
2700 assert(restart != NULL);
2701
2702 if (location_id > 0 && location_id <= (int)(restart->n_locations))
2703 retval = restart->location[location_id-1].n_ents;
2704
2705 return retval;
2706 }
2707
2708 /*----------------------------------------------------------------------------*/
2709 /*!
2710 * \brief Print the index associated with a restart file in read mode
2711 *
2712 * \param[in] restart associated restart file pointer
2713 */
2714 /*----------------------------------------------------------------------------*/
2715
2716 void
cs_restart_dump_index(const cs_restart_t * restart)2717 cs_restart_dump_index(const cs_restart_t *restart)
2718 {
2719 size_t loc_id;
2720
2721 assert(restart != NULL);
2722
2723 for (loc_id = 0; loc_id < restart->n_locations; loc_id++) {
2724 const _location_t *loc = &(restart->location[loc_id]);
2725 bft_printf(_(" Location: %s\n"
2726 " (number: %03d, n_glob_ents: %llu)\n"),
2727 loc->name, (int)(loc->id),
2728 (unsigned long long)(loc->n_glob_ents));
2729 }
2730 if (restart->n_locations > 0)
2731 bft_printf("\n");
2732
2733 /* Dump general file info, including index */
2734
2735 bft_printf(_(" General information associated with the restart file:\n"));
2736
2737 cs_io_dump(restart->fh);
2738 }
2739
2740 /*----------------------------------------------------------------------------*/
2741 /*!
2742 * \brief Check the presence of a given section in a restart file.
2743 *
2744 * \param[in] restart associated restart file pointer
2745 * \param[in] sec_name section name
2746 * \param[in] location_id id of corresponding location
2747 * \param[in] n_location_vals number of values per location (interlaced)
2748 * \param[in] val_type value type
2749 *
2750 * \return 0 (CS_RESTART_SUCCESS) in case of success,
2751 * or error code (CS_RESTART_ERR_xxx) in case of error
2752 */
2753 /*----------------------------------------------------------------------------*/
2754
2755 int
cs_restart_check_section(cs_restart_t * restart,const char * sec_name,int location_id,int n_location_vals,cs_restart_val_type_t val_type)2756 cs_restart_check_section(cs_restart_t *restart,
2757 const char *sec_name,
2758 int location_id,
2759 int n_location_vals,
2760 cs_restart_val_type_t val_type)
2761 {
2762 assert(restart != NULL);
2763
2764 return _check_section_f(restart,
2765 _restart_context,
2766 sec_name,
2767 location_id,
2768 n_location_vals,
2769 val_type);
2770 }
2771
2772 /*----------------------------------------------------------------------------*/
2773 /*!
2774 * \brief Read a section from a restart file.
2775 *
2776 * \param[in] restart associated restart file pointer
2777 * \param[in] sec_name section name
2778 * \param[in] location_id id of corresponding location
2779 * \param[in] n_location_vals number of values per location (interlaced)
2780 * \param[in] val_type value type
2781 * \param[out] val array of values
2782 *
2783 * \return 0 (CS_RESTART_SUCCESS) in case of success,
2784 * or error code (CS_RESTART_ERR_xxx) in case of error
2785 */
2786 /*----------------------------------------------------------------------------*/
2787
2788 int
cs_restart_read_section(cs_restart_t * restart,const char * sec_name,int location_id,int n_location_vals,cs_restart_val_type_t val_type,void * val)2789 cs_restart_read_section(cs_restart_t *restart,
2790 const char *sec_name,
2791 int location_id,
2792 int n_location_vals,
2793 cs_restart_val_type_t val_type,
2794 void *val)
2795 {
2796 double timing[2];
2797
2798 timing[0] = cs_timer_wtime();
2799
2800 assert(restart != NULL);
2801
2802 int retval = _read_section_f(restart,
2803 _restart_context,
2804 sec_name,
2805 location_id,
2806 n_location_vals,
2807 val_type,
2808 val);
2809
2810 timing[1] = cs_timer_wtime();
2811 _restart_wtime[restart->mode] += timing[1] - timing[0];
2812
2813 /* Return */
2814
2815 return retval;
2816 }
2817
2818 /*----------------------------------------------------------------------------*/
2819 /*!
2820 * \brief Write a section to a restart file.
2821 *
2822 * \param[in] restart associated restart file pointer
2823 * \param[in] sec_name section name
2824 * \param[in] location_id id of corresponding location
2825 * \param[in] n_location_vals number of values per location (interlaced)
2826 * \param[in] val_type value type
2827 * \param[in] val array of values
2828 */
2829 /*----------------------------------------------------------------------------*/
2830
2831 void
cs_restart_write_section(cs_restart_t * restart,const char * sec_name,int location_id,int n_location_vals,cs_restart_val_type_t val_type,const void * val)2832 cs_restart_write_section(cs_restart_t *restart,
2833 const char *sec_name,
2834 int location_id,
2835 int n_location_vals,
2836 cs_restart_val_type_t val_type,
2837 const void *val)
2838 {
2839 double timing[2];
2840
2841 timing[0] = cs_timer_wtime();
2842
2843 assert(restart != NULL);
2844
2845 _write_section_f(restart,
2846 _restart_context,
2847 sec_name,
2848 location_id,
2849 n_location_vals,
2850 val_type,
2851 val);
2852
2853 timing[1] = cs_timer_wtime();
2854 _restart_wtime[restart->mode] += timing[1] - timing[0];
2855 }
2856
2857 /*----------------------------------------------------------------------------*/
2858 /*!
2859 * \brief Read basic particles information from a restart file.
2860 *
2861 * This includes building a matching location and associated global numbering.
2862 *
2863 * \param[in] restart associated restart file pointer
2864 * \param[in] name name of particles set
2865 * \param[out] n_particles number of particles, or NULL
2866 *
2867 * \return the location id assigned to the particles, or -1 in case of error
2868 */
2869 /*----------------------------------------------------------------------------*/
2870
2871 int
cs_restart_read_particles_info(cs_restart_t * restart,const char * name,cs_lnum_t * n_particles)2872 cs_restart_read_particles_info(cs_restart_t *restart,
2873 const char *name,
2874 cs_lnum_t *n_particles)
2875 {
2876 double timing[2];
2877
2878 cs_lnum_t c_id;
2879 int rec_id;
2880 cs_io_sec_header_t header;
2881
2882 cs_lnum_t block_buf_size = 0;
2883 size_t nbr_byte_ent = sizeof(cs_gnum_t);
2884 cs_lnum_t n_cells = (restart->location[CS_MESH_LOCATION_CELLS-1]).n_ents;
2885 cs_gnum_t n_glob_cells
2886 = (restart->location[CS_MESH_LOCATION_CELLS-1]).n_glob_ents;
2887 cs_gnum_t n_glob_particles = 0;
2888
2889 int *default_p_rank = NULL;
2890 const cs_gnum_t *g_cell_num
2891 = restart->location[CS_MESH_LOCATION_CELLS-1].ent_global_num;
2892 const cs_datatype_t int_type
2893 = (sizeof(int) == 8) ? CS_INT64 : CS_INT32;
2894
2895 int loc_id = -1;
2896
2897 timing[0] = cs_timer_wtime();
2898
2899 if (n_particles != NULL)
2900 *n_particles = 0;
2901
2902 /* Search for location with the same name */
2903
2904 for (loc_id = 0; loc_id < (int)(restart->n_locations); loc_id++) {
2905 if ((strcmp((restart->location[loc_id]).name, name) == 0))
2906 break;
2907 }
2908
2909 if (loc_id >= (int)(restart->n_locations))
2910 return -1;
2911
2912 n_glob_particles = (restart->location[loc_id]).n_glob_ents_f;
2913
2914 /* Search for the corresponding cell_num record in the index */
2915
2916 rec_id = _restart_section_id(restart, NULL, name, "_cell_num");
2917
2918 if (rec_id < 0)
2919 return -1;
2920
2921 #if defined(HAVE_MPI)
2922
2923 if (cs_glob_n_ranks > 1) {
2924
2925 int *b_cell_rank, *p_cell_rank;
2926 cs_gnum_t *part_cell_num = NULL;
2927 cs_part_to_block_t *pbd = NULL;
2928 cs_all_to_all_t *d = NULL;
2929
2930 /* Now read matching cell_num to an arbitrary block distribution */
2931
2932 cs_block_dist_info_t cell_bi
2933 = cs_block_dist_compute_sizes(cs_glob_rank_id,
2934 cs_glob_n_ranks,
2935 restart->rank_step,
2936 restart->min_block_size / nbr_byte_ent,
2937 n_glob_cells);
2938
2939 cs_block_dist_info_t part_bi
2940 = cs_block_dist_compute_sizes(cs_glob_rank_id,
2941 cs_glob_n_ranks,
2942 restart->rank_step,
2943 restart->min_block_size / nbr_byte_ent,
2944 n_glob_particles);
2945
2946 /* Read data to blocks */
2947
2948 block_buf_size = (part_bi.gnum_range[1] - part_bi.gnum_range[0]);
2949
2950 if (block_buf_size > 0)
2951 BFT_MALLOC(part_cell_num, block_buf_size, cs_gnum_t);
2952
2953 header = cs_io_get_indexed_sec_header(restart->fh, rec_id);
2954
2955 cs_io_set_indexed_position(restart->fh, &header, rec_id);
2956
2957 cs_io_read_block(&header,
2958 part_bi.gnum_range[0],
2959 part_bi.gnum_range[1],
2960 part_cell_num,
2961 restart->fh);
2962
2963 /* Build block distribution cell rank info */
2964
2965 BFT_MALLOC(b_cell_rank,
2966 (cell_bi.gnum_range[1] - cell_bi.gnum_range[0]),
2967 int);
2968
2969 BFT_MALLOC(p_cell_rank, n_cells, int);
2970
2971 pbd = cs_part_to_block_create_by_gnum(cs_glob_mpi_comm,
2972 cell_bi,
2973 n_cells,
2974 g_cell_num);
2975
2976 for (c_id = 0; c_id < n_cells; c_id++)
2977 p_cell_rank[c_id] = cs_glob_rank_id;
2978
2979 cs_part_to_block_copy_array(pbd,
2980 int_type,
2981 1,
2982 p_cell_rank,
2983 b_cell_rank);
2984
2985 cs_part_to_block_destroy(&pbd);
2986
2987 BFT_FREE(p_cell_rank);
2988
2989 /* Now build distribution structure */
2990
2991 default_p_rank = _default_p_rank(&part_bi,
2992 part_cell_num,
2993 cs_glob_mpi_comm);
2994
2995 cs_lnum_t n_part_ents = 0;
2996 cs_gnum_t *ent_global_num = NULL;
2997
2998 d = cs_block_to_part_create_by_adj_s(cs_glob_mpi_comm,
2999 part_bi,
3000 cell_bi,
3001 1,
3002 part_cell_num,
3003 b_cell_rank,
3004 default_p_rank,
3005 &n_part_ents,
3006 &ent_global_num);
3007
3008 if (default_p_rank != NULL)
3009 BFT_FREE(default_p_rank);
3010
3011 BFT_FREE(b_cell_rank);
3012
3013 (restart->location[loc_id])._ent_global_num = ent_global_num;
3014 (restart->location[loc_id]).ent_global_num
3015 = (restart->location[loc_id])._ent_global_num;
3016
3017 (restart->location[loc_id]).n_glob_ents = n_glob_particles;
3018 (restart->location[loc_id]).n_ents = n_part_ents;
3019
3020 cs_all_to_all_destroy(&d);
3021
3022 BFT_FREE(part_cell_num);
3023
3024 }
3025
3026 #endif /* #if defined(HAVE_MPI) */
3027
3028 if (cs_glob_n_ranks == 1) {
3029
3030 (restart->location[loc_id]).n_glob_ents = n_glob_particles;
3031 (restart->location[loc_id]).n_ents = n_glob_particles;
3032
3033 }
3034
3035 if (n_particles != NULL)
3036 *n_particles = (restart->location[loc_id]).n_ents;
3037
3038 timing[1] = cs_timer_wtime();
3039 _restart_wtime[restart->mode] += timing[1] - timing[0];
3040
3041 return loc_id + 1;
3042 }
3043
3044 /*----------------------------------------------------------------------------*/
3045 /*!
3046 * \brief Read basic particles information from a restart file.
3047 *
3048 * \param[in] restart associated restart file pointer
3049 * \param[in] particles_location_id location id of particles set
3050 * \param[out] particle_cell_id local cell id to which particles belong
3051 * \param[out] particle_coords local particle coordinates (interleaved)
3052 *
3053 * \return 0 (CS_RESTART_SUCCESS) in case of success,
3054 * or error code (CS_RESTART_ERR_xxx) in case of error
3055 */
3056 /*----------------------------------------------------------------------------*/
3057
3058 int
cs_restart_read_particles(cs_restart_t * restart,int particles_location_id,cs_lnum_t * particle_cell_id,cs_real_t * particle_coords)3059 cs_restart_read_particles(cs_restart_t *restart,
3060 int particles_location_id,
3061 cs_lnum_t *particle_cell_id,
3062 cs_real_t *particle_coords)
3063 {
3064 double timing[2];
3065 char *sec_name = NULL;
3066
3067 cs_lnum_t n_cells = (restart->location[CS_MESH_LOCATION_CELLS-1]).n_ents;
3068 const cs_gnum_t *g_cells_num
3069 = (restart->location[CS_MESH_LOCATION_CELLS-1]).ent_global_num;
3070
3071 const char *name = (restart->location[particles_location_id - 1]).name;
3072 const char *cell_num_postfix = "_cell_num";
3073 const char *coords_postfix = "_coords";
3074
3075 int retcode = CS_RESTART_SUCCESS;
3076
3077 cs_lnum_t n_particles = (restart->location[particles_location_id - 1]).n_ents;
3078
3079 /* Read particle coordinates */
3080
3081 BFT_MALLOC(sec_name, strlen(name) + strlen(coords_postfix) + 1, char);
3082 strcpy(sec_name, name);
3083 strcat(sec_name, coords_postfix);
3084
3085 retcode = cs_restart_read_section(restart,
3086 sec_name,
3087 particles_location_id,
3088 3,
3089 CS_TYPE_cs_real_t,
3090 particle_coords);
3091
3092 BFT_FREE(sec_name);
3093
3094 if (retcode != CS_RESTART_SUCCESS)
3095 return retcode;
3096
3097 /* Read particle cell id */
3098
3099 BFT_MALLOC(sec_name, strlen(name) + strlen(cell_num_postfix) + 1, char);
3100 strcpy(sec_name, name);
3101 strcat(sec_name, cell_num_postfix);
3102
3103 #if defined(HAVE_MPI)
3104
3105 if (cs_glob_n_ranks > 1) {
3106
3107 cs_gnum_t *g_part_cell_num;
3108
3109 BFT_MALLOC(g_part_cell_num, n_particles, cs_gnum_t);
3110
3111 retcode = cs_restart_read_section(restart,
3112 sec_name,
3113 particles_location_id,
3114 1,
3115 CS_TYPE_cs_gnum_t,
3116 g_part_cell_num);
3117
3118 timing[0] = cs_timer_wtime();
3119
3120 cs_block_to_part_global_to_local(n_particles,
3121 0,
3122 n_cells,
3123 false,
3124 g_cells_num,
3125 g_part_cell_num,
3126 particle_cell_id);
3127
3128 BFT_FREE(g_part_cell_num);
3129
3130 timing[1] = cs_timer_wtime();
3131 _restart_wtime[restart->mode] += timing[1] - timing[0];
3132
3133 }
3134
3135 #endif /* #if defined(HAVE_MPI) */
3136
3137 if (cs_glob_n_ranks == 1) {
3138 retcode = cs_restart_read_section(restart,
3139 sec_name,
3140 particles_location_id,
3141 1,
3142 CS_TYPE_int,
3143 particle_cell_id);
3144 for (cs_lnum_t i = 0; i < n_particles; i++)
3145 particle_cell_id[i] -= 1;
3146 }
3147
3148 BFT_FREE(sec_name);
3149
3150 return retcode;
3151 }
3152
3153 /*----------------------------------------------------------------------------*/
3154 /*!
3155 * \brief Write basic particles information to a restart file.
3156 *
3157 * This includes defining a matching location and associated global numbering,
3158 * then writing particle coordinates and cell ids.
3159 *
3160 * \param[in] restart associated restart file pointer
3161 * \param[in] name name of particles set
3162 * \param[in] number_by_coords if true, numbering is based on current
3163 * coordinates; otherwise, it is simply based
3164 * on local numbers, plus the sum of particles
3165 * on lower MPI ranks
3166 * \param[in] n_particles local number of particles
3167 * \param[in] particle_cell_id local cell id (0 to n-1) to which particles
3168 * belong
3169 * \param[in] particle_coords local particle coordinates (interleaved)
3170 *
3171 * \return the location id assigned to the particles
3172 */
3173 /*----------------------------------------------------------------------------*/
3174
3175 int
cs_restart_write_particles(cs_restart_t * restart,const char * name,bool number_by_coords,cs_lnum_t n_particles,const cs_lnum_t * particle_cell_id,const cs_real_t * particle_coords)3176 cs_restart_write_particles(cs_restart_t *restart,
3177 const char *name,
3178 bool number_by_coords,
3179 cs_lnum_t n_particles,
3180 const cs_lnum_t *particle_cell_id,
3181 const cs_real_t *particle_coords)
3182 {
3183 double timing[2];
3184 cs_lnum_t i;
3185
3186 cs_gnum_t n_glob_particles = n_particles;
3187 cs_gnum_t *global_particle_num = NULL;
3188 cs_gnum_t *global_part_cell_num = NULL;
3189 fvm_io_num_t *io_num = NULL;
3190 char *sec_name = NULL;
3191
3192 const char *cell_num_postfix = "_cell_num";
3193 const char *coords_postfix = "_coords";
3194
3195 int loc_id = -1;
3196
3197 timing[0] = cs_timer_wtime();
3198
3199 /* Build global numbering */
3200
3201 cs_parall_counter(&n_glob_particles, 1);
3202
3203 if (number_by_coords)
3204 io_num = fvm_io_num_create_from_sfc(particle_coords,
3205 3,
3206 n_particles,
3207 FVM_IO_NUM_SFC_MORTON_BOX);
3208
3209 else
3210 io_num = fvm_io_num_create_from_scan(n_particles);
3211
3212 global_particle_num = fvm_io_num_transfer_global_num(io_num);
3213 fvm_io_num_destroy(io_num);
3214
3215 /* Create a new location, with ownership of global numbers */
3216
3217 loc_id = cs_restart_add_location(restart,
3218 name,
3219 n_glob_particles,
3220 n_particles,
3221 global_particle_num);
3222
3223 (restart->location[loc_id-1])._ent_global_num = global_particle_num;
3224 assert((restart->location[loc_id-1]).ent_global_num == global_particle_num);
3225
3226 /* Write particle coordinates */
3227
3228 BFT_MALLOC(sec_name, strlen(name) + strlen(coords_postfix) + 1, char);
3229 strcpy(sec_name, name);
3230 strcat(sec_name, coords_postfix);
3231
3232 timing[1] = cs_timer_wtime();
3233 _restart_wtime[restart->mode] += timing[1] - timing[0];
3234
3235 cs_restart_write_section(restart,
3236 sec_name,
3237 loc_id,
3238 3,
3239 CS_TYPE_cs_real_t,
3240 particle_coords);
3241
3242 timing[0] = cs_timer_wtime();
3243
3244 BFT_FREE(sec_name);
3245
3246 /* Write particle cell location information */
3247
3248 BFT_MALLOC(global_part_cell_num, n_particles, cs_gnum_t);
3249
3250 if (restart->location[CS_MESH_LOCATION_CELLS-1].ent_global_num != NULL) {
3251 const cs_gnum_t *g_cell_num
3252 = restart->location[CS_MESH_LOCATION_CELLS-1].ent_global_num;
3253 for (i = 0; i < n_particles; i++) {
3254 if (particle_cell_id[i] > -1)
3255 global_part_cell_num[i] = g_cell_num[particle_cell_id[i]];
3256 else
3257 global_part_cell_num[i] = 0;
3258 }
3259 }
3260 else {
3261 for (i = 0; i < n_particles; i++)
3262 global_part_cell_num[i] = particle_cell_id[i] + 1;
3263 }
3264
3265 BFT_MALLOC(sec_name, strlen(name) + strlen(cell_num_postfix) + 1, char);
3266 strcpy(sec_name, name);
3267 strcat(sec_name, cell_num_postfix);
3268
3269 timing[1] = cs_timer_wtime();
3270 _restart_wtime[restart->mode] += timing[1] - timing[0];
3271
3272 cs_restart_write_section(restart,
3273 sec_name,
3274 loc_id,
3275 1,
3276 CS_TYPE_cs_gnum_t,
3277 global_part_cell_num);
3278
3279 BFT_FREE(sec_name);
3280
3281 BFT_FREE(global_part_cell_num);
3282
3283 return loc_id;
3284 }
3285
3286 /*----------------------------------------------------------------------------*/
3287 /*!
3288 * \brief Read a referenced location id section from a restart file.
3289 *
3290 * The section read from file contains the global ids matching the local
3291 * element ids of a given location. Global id's are transformed to local
3292 * ids by this function.
3293 *
3294 * In case global referenced ids read do not match those of local elements,
3295 * id_base - 1 is assigned to the corresponding local ids.
3296 *
3297 * \param[in] restart associated restart file pointer
3298 * \param[in] sec_name section name
3299 * \param[in] location_id id of location on which id_ref is defined
3300 * \param[in] ref_location_id id of referenced location
3301 * \param[in] ref_id_base base of location entity id numbers
3302 * (usually 0 or 1)
3303 * \param[out] ref_id array of location entity ids
3304 *
3305 * \return 0 (CS_RESTART_SUCCESS) in case of success,
3306 * or error code (CS_RESTART_ERR_xxx) in case of error
3307 */
3308 /*----------------------------------------------------------------------------*/
3309
3310 int
cs_restart_read_ids(cs_restart_t * restart,const char * sec_name,int location_id,int ref_location_id,cs_lnum_t ref_id_base,cs_lnum_t * ref_id)3311 cs_restart_read_ids(cs_restart_t *restart,
3312 const char *sec_name,
3313 int location_id,
3314 int ref_location_id,
3315 cs_lnum_t ref_id_base,
3316 cs_lnum_t *ref_id)
3317 {
3318 cs_lnum_t n_ents = 0;
3319 cs_gnum_t *g_num;
3320
3321 _location_t *ref_location = NULL;
3322
3323 int retcode = CS_RESTART_SUCCESS;
3324
3325 assert(restart != NULL);
3326
3327 /* Local number of elements for location id */
3328
3329 if (location_id == 0)
3330 n_ents = 1;
3331
3332 else if (location_id > 0 && location_id <= (int)(restart->n_locations))
3333 n_ents = restart->location[location_id-1].n_ents;
3334
3335 else
3336 bft_error(__FILE__, __LINE__, 0,
3337 _("Location number %d given for restart file\n"
3338 "\"%s\" is not valid."),
3339 (int)location_id, restart->name);
3340
3341 if (ref_location_id > 0 && ref_location_id <= (int)(restart->n_locations))
3342 ref_location = restart->location + ref_location_id-1;
3343
3344 else if (ref_location_id != 0)
3345 bft_error(__FILE__, __LINE__, 0,
3346 _("Location number %d given for restart file\n"
3347 "\"%s\" is not valid."),
3348 (int)location_id, restart->name);
3349
3350 /* Transform local to global ids */
3351
3352 BFT_MALLOC(g_num, n_ents, cs_gnum_t);
3353
3354 /* Read associated data */
3355
3356 retcode = cs_restart_read_section(restart,
3357 sec_name,
3358 location_id,
3359 1,
3360 CS_TYPE_cs_gnum_t,
3361 g_num);
3362
3363 if (retcode == CS_RESTART_SUCCESS) {
3364
3365 double timing[2];
3366
3367 timing[0] = cs_timer_wtime();
3368
3369 if (ref_location_id == 0 || ref_location->ent_global_num == NULL) {
3370 cs_lnum_t i;
3371 for (i = 0; i < n_ents; i++)
3372 ref_id[i] = g_num[i] + ref_id_base - 1;
3373 }
3374
3375 else /* if location_id > 0 */
3376 cs_block_to_part_global_to_local(n_ents,
3377 ref_id_base,
3378 ref_location->n_ents,
3379 false,
3380 ref_location->ent_global_num,
3381 g_num,
3382 ref_id);
3383
3384 timing[1] = cs_timer_wtime();
3385 _restart_wtime[restart->mode] += timing[1] - timing[0];
3386
3387 }
3388
3389 BFT_FREE(g_num);
3390
3391 /* Return */
3392
3393 return retcode;
3394 }
3395
3396 /*----------------------------------------------------------------------------*/
3397 /*!
3398 * \brief Write a referenced location id section to a restart file.
3399 *
3400 * The section written to file contains the global ids matching the
3401 * local element ids of a given location.
3402 *
3403 * \param[in] restart associated restart file pointer
3404 * \param[in] sec_name section name
3405 * \param[in] location_id id of location on which id_ref is defined
3406 * \param[in] ref_location_id id of referenced location
3407 * \param[in] ref_id_base base of location entity id numbers
3408 * (usually 0 or 1)
3409 * \param[in] ref_id array of location entity ids
3410 */
3411 /*----------------------------------------------------------------------------*/
3412
3413 void
cs_restart_write_ids(cs_restart_t * restart,const char * sec_name,int location_id,int ref_location_id,cs_lnum_t ref_id_base,const cs_lnum_t * ref_id)3414 cs_restart_write_ids(cs_restart_t *restart,
3415 const char *sec_name,
3416 int location_id,
3417 int ref_location_id,
3418 cs_lnum_t ref_id_base,
3419 const cs_lnum_t *ref_id)
3420 {
3421 double timing[2];
3422 cs_lnum_t i, n_ents = 0;
3423 cs_gnum_t *g_num;
3424
3425 _location_t *ref_location = NULL;
3426
3427 assert(restart != NULL);
3428
3429 /* Local number of elements for location id */
3430
3431 if (location_id == 0)
3432 n_ents = 1;
3433
3434 else if (location_id > 0 && location_id <= (int)(restart->n_locations))
3435 n_ents = restart->location[location_id-1].n_ents;
3436
3437 else
3438 bft_error(__FILE__, __LINE__, 0,
3439 _("Location number %d given for restart file\n"
3440 "\"%s\" is not valid."),
3441 (int)location_id, restart->name);
3442
3443 if (ref_location_id > 0 && ref_location_id <= (int)(restart->n_locations))
3444 ref_location = restart->location + ref_location_id-1;
3445
3446 else if (ref_location_id != 0)
3447 bft_error(__FILE__, __LINE__, 0,
3448 _("Location number %d given for restart file\n"
3449 "\"%s\" is not valid."),
3450 (int)location_id, restart->name);
3451
3452 /* Transform local to global ids */
3453
3454 timing[0] = cs_timer_wtime();
3455
3456 BFT_MALLOC(g_num, n_ents, cs_gnum_t);
3457
3458 if (ref_location_id == 0) {
3459 for (i = 0; i < n_ents; i++)
3460 g_num[0] = ref_id[0] - ref_id_base + 1;
3461 }
3462
3463 else { /* if location_id > 0 */
3464 if (ref_location->ent_global_num != NULL) {
3465 for (i = 0; i < n_ents; i++) {
3466 if (ref_id[i] >= ref_id_base)
3467 g_num[i] = ref_location->ent_global_num[ref_id[i] - ref_id_base];
3468 else
3469 g_num[i] = 0;
3470 }
3471 }
3472 else {
3473 for (i = 0; i < n_ents; i++) {
3474 if (ref_id[i] >= ref_id_base)
3475 g_num[i] = ref_id[i] - ref_id_base + 1;
3476 else
3477 g_num[i] = 0;
3478 }
3479 }
3480
3481 }
3482
3483 timing[1] = cs_timer_wtime();
3484 _restart_wtime[restart->mode] += timing[1] - timing[0];
3485
3486 /* Write associated data */
3487
3488 cs_restart_write_section(restart,
3489 sec_name,
3490 location_id,
3491 1,
3492 CS_TYPE_cs_gnum_t,
3493 g_num);
3494
3495 BFT_FREE(g_num);
3496 }
3497
3498 /*----------------------------------------------------------------------------*/
3499 /*!
3500 * \brief Read a section from a restart file, when that section may have used
3501 * a different name in a previous version.
3502 *
3503 * \param[in] restart associated restart file pointer
3504 * \param[in] sec_name section name
3505 * \param[in] old_name old name
3506 * \param[in] location_id id of corresponding location
3507 * \param[in] n_location_vals number of values per location (interlaced)
3508 * \param[in] val_type value type
3509 * \param[out] val array of values
3510 *
3511 * \return 0 (CS_RESTART_SUCCESS) in case of success,
3512 * or error code (CS_RESTART_ERR_xxx) in case of error
3513 */
3514 /*----------------------------------------------------------------------------*/
3515
3516 int
cs_restart_read_section_compat(cs_restart_t * restart,const char * sec_name,const char * old_name,int location_id,int n_location_vals,cs_restart_val_type_t val_type,void * val)3517 cs_restart_read_section_compat(cs_restart_t *restart,
3518 const char *sec_name,
3519 const char *old_name,
3520 int location_id,
3521 int n_location_vals,
3522 cs_restart_val_type_t val_type,
3523 void *val)
3524 {
3525 int retval = CS_RESTART_SUCCESS;
3526
3527 assert(location_id >= 0);
3528
3529 /* Check for section with current name */
3530
3531 retval = cs_restart_check_section(restart,
3532 sec_name,
3533 location_id,
3534 n_location_vals,
3535 val_type);
3536
3537 /* Check for older name, read and return if present */
3538
3539 if (retval == CS_RESTART_ERR_N_VALS || retval == CS_RESTART_ERR_EXISTS) {
3540
3541 retval =cs_restart_check_section(restart,
3542 old_name,
3543 location_id,
3544 n_location_vals,
3545 val_type);
3546
3547 if (retval == CS_RESTART_SUCCESS) {
3548 retval = cs_restart_read_section(restart,
3549 old_name,
3550 location_id,
3551 n_location_vals,
3552 val_type,
3553 val);
3554 return retval;
3555 }
3556
3557 }
3558
3559 /* Read with current name (if the section is not found,
3560 logging will refer to the current name) */
3561
3562 retval = cs_restart_read_section(restart,
3563 sec_name,
3564 location_id,
3565 n_location_vals,
3566 val_type,
3567 val);
3568
3569 return retval;
3570 }
3571
3572 /*----------------------------------------------------------------------------*/
3573 /*!
3574 * \brief Read a cs_real_3_t vector section from a restart file, when that
3575 * section may have used a different name and been non-interleaved
3576 * in a previous version.
3577 *
3578 * This function assumes a mesh-base location (i.e. location_id > 0)
3579 *
3580 * \param[in] restart associated restart file pointer
3581 * \param[in] sec_name section name
3582 * \param[in] old_name_x old name, x component
3583 * \param[in] old_name_y old name, y component
3584 * \param[in] old_name_z old name, z component
3585 * \param[in] location_id id of corresponding location
3586 * \param[out] val array of values
3587 *
3588 * \return 0 (CS_RESTART_SUCCESS) in case of success,
3589 * or error code (CS_RESTART_ERR_xxx) in case of error
3590 */
3591 /*----------------------------------------------------------------------------*/
3592
3593 int
cs_restart_read_real_3_t_compat(cs_restart_t * restart,const char * sec_name,const char * old_name_x,const char * old_name_y,const char * old_name_z,int location_id,cs_real_3_t * val)3594 cs_restart_read_real_3_t_compat(cs_restart_t *restart,
3595 const char *sec_name,
3596 const char *old_name_x,
3597 const char *old_name_y,
3598 const char *old_name_z,
3599 int location_id,
3600 cs_real_3_t *val)
3601 {
3602 int retval = CS_RESTART_SUCCESS;
3603
3604 assert(location_id > 0);
3605
3606 /* Check for section with current name */
3607
3608 retval = cs_restart_check_section(restart,
3609 sec_name,
3610 location_id,
3611 3,
3612 CS_TYPE_cs_real_t);
3613
3614 /* Check for older name series, read and return if present */
3615
3616 if (retval == CS_RESTART_ERR_N_VALS || retval == CS_RESTART_ERR_EXISTS) {
3617
3618 retval = cs_restart_check_section(restart,
3619 old_name_x,
3620 location_id,
3621 1,
3622 CS_TYPE_cs_real_t);
3623
3624 if (retval == CS_RESTART_SUCCESS) {
3625
3626 cs_real_t *buffer = NULL;
3627 cs_lnum_t i;
3628 cs_lnum_t n_ents = (restart->location[location_id-1]).n_ents;
3629
3630 BFT_MALLOC(buffer, n_ents*3, cs_real_t);
3631
3632 retval = cs_restart_read_section(restart,
3633 old_name_x,
3634 location_id,
3635 1,
3636 CS_TYPE_cs_real_t,
3637 buffer);
3638
3639 if (retval == CS_RESTART_SUCCESS)
3640 retval = cs_restart_read_section(restart,
3641 old_name_y,
3642 location_id,
3643 1,
3644 CS_TYPE_cs_real_t,
3645 buffer + n_ents);
3646
3647 if (retval == CS_RESTART_SUCCESS)
3648 retval = cs_restart_read_section(restart,
3649 old_name_z,
3650 location_id,
3651 1,
3652 CS_TYPE_cs_real_t,
3653 buffer + n_ents*2);
3654
3655 if (retval == CS_RESTART_SUCCESS) {
3656 for (i = 0; i < n_ents; i++) {
3657 val[i][0] = buffer[i];
3658 val[i][1] = buffer[i + n_ents];
3659 val[i][2] = buffer[i + n_ents*2];
3660 }
3661 }
3662
3663 BFT_FREE(buffer);
3664
3665 return retval;
3666
3667 }
3668 }
3669
3670 /* Read with current name (if the section is not found,
3671 logging will refer to the current name) */
3672
3673 retval = cs_restart_read_section(restart,
3674 sec_name,
3675 location_id,
3676 3,
3677 CS_TYPE_cs_real_t,
3678 val);
3679
3680 return retval;
3681 }
3682
3683 /*----------------------------------------------------------------------------*/
3684 /*!
3685 * \brief Read a cs_real_6_t tensor section from a restart file, when that
3686 * section may have used a different name and been non-interleaved
3687 * in a previous version.
3688 *
3689 * This function assumes a mesh-base location (i.e. location_id > 0)
3690 *
3691 * \param[in] restart associated restart file pointer
3692 * \param[in] sec_name section name
3693 * \param[in] old_name_xx old name, xx component
3694 * \param[in] old_name_yy old name, yy component
3695 * \param[in] old_name_zz old name, zz component
3696 * \param[in] old_name_xy old name, xy component
3697 * \param[in] old_name_yz old name, yz component
3698 * \param[in] old_name_xz old name, xz component
3699 * \param[in] location_id id of corresponding location
3700 * \param[out] val array of values
3701 *
3702 * \return 0 (CS_RESTART_SUCCESS) in case of success,
3703 * or error code (CS_RESTART_ERR_xxx) in case of error
3704 */
3705 /*----------------------------------------------------------------------------*/
3706
3707 int
cs_restart_read_real_6_t_compat(cs_restart_t * restart,const char * sec_name,const char * old_name_xx,const char * old_name_yy,const char * old_name_zz,const char * old_name_xy,const char * old_name_yz,const char * old_name_xz,int location_id,cs_real_6_t * val)3708 cs_restart_read_real_6_t_compat(cs_restart_t *restart,
3709 const char *sec_name,
3710 const char *old_name_xx,
3711 const char *old_name_yy,
3712 const char *old_name_zz,
3713 const char *old_name_xy,
3714 const char *old_name_yz,
3715 const char *old_name_xz,
3716 int location_id,
3717 cs_real_6_t *val)
3718 {
3719 int retval = CS_RESTART_SUCCESS;
3720
3721 assert(location_id > 0);
3722
3723 /* Check for section with current name */
3724
3725 retval = cs_restart_check_section(restart,
3726 sec_name,
3727 location_id,
3728 6,
3729 CS_TYPE_cs_real_t);
3730
3731 /* Check for older name series, read and return if present */
3732
3733 if (retval == CS_RESTART_ERR_N_VALS || retval == CS_RESTART_ERR_EXISTS) {
3734
3735 retval = cs_restart_check_section(restart,
3736 old_name_xx,
3737 location_id,
3738 1,
3739 CS_TYPE_cs_real_t);
3740
3741 if (retval == CS_RESTART_SUCCESS) {
3742
3743 cs_real_t *buffer = NULL;
3744 cs_lnum_t i;
3745 cs_lnum_t n_ents = (restart->location[location_id-1]).n_ents;
3746
3747 BFT_MALLOC(buffer, n_ents*6, cs_real_t);
3748
3749 retval = cs_restart_read_section(restart,
3750 old_name_xx,
3751 location_id,
3752 1,
3753 CS_TYPE_cs_real_t,
3754 buffer);
3755
3756 if (retval == CS_RESTART_SUCCESS)
3757 retval = cs_restart_read_section(restart,
3758 old_name_yy,
3759 location_id,
3760 1,
3761 CS_TYPE_cs_real_t,
3762 buffer + n_ents);
3763
3764 if (retval == CS_RESTART_SUCCESS)
3765 retval = cs_restart_read_section(restart,
3766 old_name_zz,
3767 location_id,
3768 1,
3769 CS_TYPE_cs_real_t,
3770 buffer + n_ents*2);
3771 if (retval == CS_RESTART_SUCCESS)
3772 retval = cs_restart_read_section(restart,
3773 old_name_xy,
3774 location_id,
3775 1,
3776 CS_TYPE_cs_real_t,
3777 buffer + n_ents*3);
3778
3779 if (retval == CS_RESTART_SUCCESS)
3780 retval = cs_restart_read_section(restart,
3781 old_name_yz,
3782 location_id,
3783 1,
3784 CS_TYPE_cs_real_t,
3785 buffer + n_ents*4);
3786 if (retval == CS_RESTART_SUCCESS)
3787 retval = cs_restart_read_section(restart,
3788 old_name_xz,
3789 location_id,
3790 1,
3791 CS_TYPE_cs_real_t,
3792 buffer + n_ents*5);
3793
3794 if (retval == CS_RESTART_SUCCESS) {
3795 for (i = 0; i < n_ents; i++) {
3796 val[i][0] = buffer[i];
3797 val[i][1] = buffer[i + n_ents];
3798 val[i][2] = buffer[i + n_ents*2];
3799 val[i][3] = buffer[i + n_ents*3];
3800 val[i][4] = buffer[i + n_ents*4];
3801 val[i][5] = buffer[i + n_ents*5];
3802 }
3803 }
3804
3805 BFT_FREE(buffer);
3806
3807 return retval;
3808
3809 }
3810 }
3811
3812 /* Read with current name (if the section is not found,
3813 logging will refer to the current name) */
3814
3815 retval = cs_restart_read_section(restart,
3816 sec_name,
3817 location_id,
3818 6,
3819 CS_TYPE_cs_real_t,
3820 val);
3821
3822 return retval;
3823 }
3824
3825 /*----------------------------------------------------------------------------*/
3826 /*!
3827 * \brief Read a cs_real_66_t tensor section from a restart file, when that
3828 * section may have used a different name and been non-interleaved
3829 * in a previous version.
3830 *
3831 * This function assumes a mesh-base location (i.e. location_id > 0)
3832 *
3833 * \param[in] restart associated restart file pointer
3834 * \param[in] sec_name section name
3835 * \param[in] old_name_xx old name, xx component
3836 * \param[in] old_name_yy old name, yy component
3837 * \param[in] old_name_zz old name, zz component
3838 * \param[in] old_name_xy old name, xy component
3839 * \param[in] old_name_yz old name, yz component
3840 * \param[in] old_name_xz old name, xz component
3841 * \param[in] location_id id of corresponding location
3842 * \param[out] val array of values
3843 *
3844 * \return 0 (CS_RESTART_SUCCESS) in case of success,
3845 * or error code (CS_RESTART_ERR_xxx) in case of error
3846 */
3847 /*----------------------------------------------------------------------------*/
3848
3849 int
cs_restart_read_real_66_t_compat(cs_restart_t * restart,const char * sec_name,const char * old_name_xx,const char * old_name_yy,const char * old_name_zz,const char * old_name_xy,const char * old_name_yz,const char * old_name_xz,int location_id,cs_real_66_t * val)3850 cs_restart_read_real_66_t_compat(cs_restart_t *restart,
3851 const char *sec_name,
3852 const char *old_name_xx,
3853 const char *old_name_yy,
3854 const char *old_name_zz,
3855 const char *old_name_xy,
3856 const char *old_name_yz,
3857 const char *old_name_xz,
3858 int location_id,
3859 cs_real_66_t *val)
3860 {
3861 int retval = CS_RESTART_SUCCESS;
3862
3863 assert(location_id > 0);
3864
3865 /* Check for section with current name */
3866
3867 retval = cs_restart_check_section(restart,
3868 sec_name,
3869 location_id,
3870 6,
3871 CS_TYPE_cs_real_t);
3872
3873 /* Check for older name series, read and return if present */
3874
3875 if (retval == CS_RESTART_ERR_N_VALS || retval == CS_RESTART_ERR_EXISTS) {
3876
3877 retval = cs_restart_check_section(restart,
3878 old_name_xx,
3879 location_id,
3880 1,
3881 CS_TYPE_cs_real_t);
3882
3883 if (retval == CS_RESTART_SUCCESS) {
3884
3885 cs_real_t *buffer = NULL;
3886 cs_lnum_t i;
3887 cs_lnum_t n_ents = (restart->location[location_id-1]).n_ents;
3888
3889 BFT_MALLOC(buffer, n_ents*6, cs_real_t);
3890
3891 retval = cs_restart_read_section(restart,
3892 old_name_xx,
3893 location_id,
3894 1,
3895 CS_TYPE_cs_real_t,
3896 buffer);
3897
3898 if (retval == CS_RESTART_SUCCESS)
3899 retval = cs_restart_read_section(restart,
3900 old_name_yy,
3901 location_id,
3902 1,
3903 CS_TYPE_cs_real_t,
3904 buffer + n_ents);
3905
3906 if (retval == CS_RESTART_SUCCESS)
3907 retval = cs_restart_read_section(restart,
3908 old_name_zz,
3909 location_id,
3910 1,
3911 CS_TYPE_cs_real_t,
3912 buffer + n_ents*2);
3913 if (retval == CS_RESTART_SUCCESS)
3914 retval = cs_restart_read_section(restart,
3915 old_name_xy,
3916 location_id,
3917 1,
3918 CS_TYPE_cs_real_t,
3919 buffer + n_ents*3);
3920
3921 if (retval == CS_RESTART_SUCCESS)
3922 retval = cs_restart_read_section(restart,
3923 old_name_yz,
3924 location_id,
3925 1,
3926 CS_TYPE_cs_real_t,
3927 buffer + n_ents*4);
3928 if (retval == CS_RESTART_SUCCESS)
3929 retval = cs_restart_read_section(restart,
3930 old_name_xz,
3931 location_id,
3932 1,
3933 CS_TYPE_cs_real_t,
3934 buffer + n_ents*5);
3935
3936 if (retval == CS_RESTART_SUCCESS) {
3937 for (i = 0; i < n_ents; i++) {
3938 val[i][0][0] = buffer[i];
3939 val[i][1][1] = buffer[i + n_ents*7];
3940 val[i][2][2] = buffer[i + n_ents*14];
3941 val[i][3][3] = buffer[i + n_ents*21];
3942 val[i][4][4] = buffer[i + n_ents*28];
3943 val[i][5][5] = buffer[i + n_ents*35];
3944 }
3945 }
3946
3947 BFT_FREE(buffer);
3948
3949 return retval;
3950
3951 }
3952 }
3953
3954 /* Read with current name (if the section is not found,
3955 logging will refer to the current name) */
3956 retval = cs_restart_read_section(restart,
3957 sec_name,
3958 location_id,
3959 3,
3960 CS_TYPE_cs_real_t,
3961 val);
3962
3963 return retval;
3964 }
3965
3966 /*----------------------------------------------------------------------------*/
3967 /*!
3968 * \brief Print statistics associated with restart files
3969 */
3970 /*----------------------------------------------------------------------------*/
3971
3972 void
cs_restart_print_stats(void)3973 cs_restart_print_stats(void)
3974 {
3975 bft_printf(_("\n"
3976 "Checkpoint / restart files summary:\n"
3977 "\n"
3978 " Number of files read: %3d\n"
3979 " Number of files written: %3d\n"
3980 "\n"
3981 " Elapsed time for reading: %12.3f\n"
3982 " Elapsed time for writing: %12.3f\n"),
3983 _restart_n_opens[0], _restart_n_opens[1],
3984 _restart_wtime[0], _restart_wtime[1]);
3985 }
3986
3987 /*----------------------------------------------------------------------------*/
3988 /*!
3989 * \brief Checks if restart is done from a NCFD checkpoint file
3990 *
3991 * \return 0 if no, 1 if yes
3992 */
3993 /*----------------------------------------------------------------------------*/
3994
3995 int
cs_restart_check_if_restart_from_ncfd(cs_restart_t * r)3996 cs_restart_check_if_restart_from_ncfd(cs_restart_t *r)
3997 {
3998 int inttmp[1000];
3999 int ierror
4000 = cs_restart_read_section_compat(r,
4001 "neptune_cfd:checkpoint:main:version",
4002 "version_fichier_suite_principal",
4003 CS_MESH_LOCATION_NONE,
4004 1,
4005 CS_TYPE_int,
4006 inttmp);
4007
4008 if (ierror == 0) {
4009 bft_printf(_("Remark: restarting based on a NEPTUNE_CFD computation.\n"));
4010 _restart_from_ncfd = 1;
4011 }
4012
4013 return _restart_from_ncfd;
4014 }
4015
4016 /*----------------------------------------------------------------------------*/
4017 /*!
4018 * \brief Returns if restart is done from a NCFD checkpoint file
4019 *
4020 * \return 0 if no, 1 if yes
4021 */
4022 /*----------------------------------------------------------------------------*/
4023
4024 int
cs_restart_is_from_ncfd(void)4025 cs_restart_is_from_ncfd(void)
4026 {
4027 return _restart_from_ncfd;
4028 }
4029
4030 /*----------------------------------------------------------------------------*/
4031 /*!
4032 * \brief Set the number of checkpoints to keep.
4033 *
4034 * This function sets the maximum number of checkpoint files to keep.
4035 * Beyond this maximum, the oldest checkpoints are deleted to limit the
4036 * number of saved checkpoints. The default value is 1.
4037 *
4038 * If more than one file is kept, last one is always named "<prefix>.csc",
4039 * while others are names "<prefix_%04d.csc".
4040 *
4041 * %04 provides an id using 4 digits (0 padding), and value is the order
4042 * of writing, starting with 0. Hence: prefix_0000.csc, prefix_0001.csc, ...
4043 *
4044 * \param[in] n_checkpoints number of checkpoints to save.
4045 */
4046 /*----------------------------------------------------------------------------*/
4047
4048 void
cs_restart_set_n_max_checkpoints(int n_checkpoints)4049 cs_restart_set_n_max_checkpoints(int n_checkpoints)
4050 {
4051 if (n_checkpoints <= 0) {
4052
4053 /* deactivate checkpointing */
4054 _checkpoint_nt_interval = CS_RESTART_INTERVAL_NONE;
4055 _n_restart_directories_to_write = 0;
4056
4057 }
4058 else
4059 _n_restart_directories_to_write = n_checkpoints;
4060
4061 return;
4062 }
4063
4064 /*----------------------------------------------------------------------------*/
4065 /*!
4066 * \brief Remove all previous checkpoints which are not to be retained.
4067 */
4068 /*----------------------------------------------------------------------------*/
4069
4070 void
cs_restart_clean_multiwriters_history(void)4071 cs_restart_clean_multiwriters_history(void)
4072 {
4073 /* Check that the structure is allocated */
4074 if ( _restart_multiwriter == NULL
4075 || _n_restart_directories_to_write < 0)
4076 return;
4077
4078 for (int i = 0; i < _n_restart_multiwriters; i++) {
4079 _restart_multiwriter_t *mw = _restart_multiwriter_by_id(i);
4080
4081 int n_files_to_remove
4082 = mw->n_prev_files - _n_restart_directories_to_write + 1;
4083
4084 if (n_files_to_remove > 0) {
4085 for (int ii = 0; ii < n_files_to_remove; ii++) {
4086
4087 if (cs_glob_rank_id <= 0) {
4088 char *path = mw->prev_files[ii];
4089 if (cs_glob_rank_id <= 0)
4090 cs_file_remove(path);
4091
4092 /* Try to remove directory (if it is empty) */
4093 for (int j = strlen(path)-1; j > -1; j--) {
4094 if (path[j] == _dir_separator) {
4095 if (j > 0) {
4096 path[j] = '\0';
4097 cs_file_remove(path);
4098 }
4099 break;
4100 }
4101 }
4102 }
4103
4104 BFT_FREE(mw->prev_files[ii]);
4105
4106 }
4107
4108 /* Rotate available paths */
4109
4110 int ii = 0;
4111 for (int jj = n_files_to_remove; jj < mw->n_prev_files; jj++) {
4112 mw->prev_files[ii] = mw->prev_files[jj];
4113 mw->prev_files[jj] = NULL;
4114 }
4115
4116 mw->n_prev_files -= n_files_to_remove;
4117 /* No need for extra reallocation of mw->prev_files */
4118 }
4119
4120 }
4121 }
4122
4123 /*----------------------------------------------------------------------------*/
4124 /*!
4125 * \brief Destroy the multiwriter structure at the end of the computation.
4126 */
4127 /*----------------------------------------------------------------------------*/
4128
4129 void
cs_restart_multiwriters_destroy_all(void)4130 cs_restart_multiwriters_destroy_all(void)
4131 {
4132 if (_restart_multiwriter != NULL) {
4133 for (int i = 0; i < _n_restart_multiwriters; i++) {
4134 _restart_multiwriter_t *w = _restart_multiwriter[i];
4135
4136 BFT_FREE(w->name);
4137 BFT_FREE(w->path);
4138
4139 for (int j = 0; j < w->n_prev_files; j++)
4140 BFT_FREE(w->prev_files[j]);
4141 BFT_FREE(w->prev_files);
4142
4143 BFT_FREE(w);
4144
4145 }
4146 BFT_FREE(_restart_multiwriter);
4147 }
4148 }
4149
4150 /*----------------------------------------------------------------------------*/
4151
4152 END_C_DECLS
4153