1 /*============================================================================
2 * Write a nodal representation associated with a mesh and associated
3 * variables to time plot files
4 *============================================================================*/
5
6 /*
7 This file is part of Code_Saturne, a general-purpose CFD tool.
8
9 Copyright (C) 1998-2021 EDF S.A.
10
11 This program is free software; you can redistribute it and/or modify it under
12 the terms of the GNU General Public License as published by the Free Software
13 Foundation; either version 2 of the License, or (at your option) any later
14 version.
15
16 This program is distributed in the hope that it will be useful, but WITHOUT
17 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
18 FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
19 details.
20
21 You should have received a copy of the GNU General Public License along with
22 this program; if not, write to the Free Software Foundation, Inc., 51 Franklin
23 Street, Fifth Floor, Boston, MA 02110-1301, USA.
24 */
25
26 /*----------------------------------------------------------------------------*/
27
28 #include "cs_defs.h"
29
30 /*----------------------------------------------------------------------------
31 * Standard C library headers
32 *----------------------------------------------------------------------------*/
33
34 #include <assert.h>
35 #include <errno.h>
36 #include <stdio.h>
37 #include <stdlib.h>
38 #include <string.h>
39
40 /*----------------------------------------------------------------------------
41 * Local headers
42 *----------------------------------------------------------------------------*/
43
44 #include "bft_error.h"
45 #include "bft_mem.h"
46
47 #include "fvm_defs.h"
48 #include "fvm_convert_array.h"
49 #include "fvm_io_num.h"
50 #include "fvm_nodal.h"
51 #include "fvm_nodal_priv.h"
52 #include "fvm_writer_helper.h"
53 #include "fvm_writer_priv.h"
54
55 #include "cs_block_dist.h"
56 #include "cs_file.h"
57 #include "cs_map.h"
58 #include "cs_parall.h"
59 #include "cs_part_to_block.h"
60 #include "cs_time_plot.h"
61
62 /*----------------------------------------------------------------------------
63 * Header for the current file
64 *----------------------------------------------------------------------------*/
65
66 #include "fvm_to_time_plot.h"
67
68 /*----------------------------------------------------------------------------*/
69
70 BEGIN_C_DECLS
71
72 /*! \cond DOXYGEN_SHOULD_SKIP_THIS */
73
74 /*============================================================================
75 * Local Type Definitions
76 *============================================================================*/
77
78 /*----------------------------------------------------------------------------
79 * time plot writer structure
80 *----------------------------------------------------------------------------*/
81
82 typedef struct {
83
84 char *name; /* Writer name */
85 char *prefix; /* Plot file prefix */
86
87 int rank; /* Rank of current process in communicator */
88 int n_ranks; /* Number of processes in communicator */
89
90 cs_time_plot_format_t format; /* Time plot format */
91
92 float flush_wtime; /* Elapsed time interval between file
93 flushes (if < 0, no forced flush) */
94 int n_buf_steps; /* Number of time steps in output
95 buffer if file is not to be kept open */
96
97 bool use_iteration; /* Should we use the iteration number
98 instead of the physical time ? */
99
100 int nt; /* Time step */
101 double t; /* Time value */
102
103 int n_plots; /* Number of associated plots */
104
105 cs_map_name_to_id_t *f_map; /* field names to plots mapping */
106 cs_time_plot_t **tp; /* Associated plots */
107
108 #if defined(HAVE_MPI)
109 MPI_Comm comm; /* Associated MPI communicator */
110 #endif
111
112 } fvm_to_time_plot_writer_t;
113
114 /*----------------------------------------------------------------------------
115 * Context structure for fvm_writer_field_helper_output_* functions.
116 *----------------------------------------------------------------------------*/
117
118 typedef struct {
119
120 fvm_to_time_plot_writer_t *writer; /* Pointer to writer structure */
121
122 const fvm_nodal_t *mesh; /* Pointer to mesh */
123 const char *name; /* Field name */
124
125 } _time_plot_context_t;
126
127 /*============================================================================
128 * Static global variables
129 *============================================================================*/
130
131 /*============================================================================
132 * Private function definitions
133 *============================================================================*/
134
135 /*----------------------------------------------------------------------------
136 * Output function for coordinates files
137 *
138 * This function is passed to fvm_writer_field_helper_output_* functions.
139 *
140 * parameters:
141 * context <-> pointer to writer and field context
142 * datatype <-- output datatype
143 * dimension <-- output field dimension
144 * component_id <-- output component id (if non-interleaved)
145 * block_start <-- start global number of element for current block
146 * block_end <-- past-the-end global number of element for current block
147 * buffer <-> associated output buffer
148 *----------------------------------------------------------------------------*/
149
150 static void
_coords_output(void * context,cs_datatype_t datatype,int dimension,int component_id,cs_gnum_t block_start,cs_gnum_t block_end,void * buffer)151 _coords_output(void *context,
152 cs_datatype_t datatype,
153 int dimension,
154 int component_id,
155 cs_gnum_t block_start,
156 cs_gnum_t block_end,
157 void *buffer)
158 {
159 if (dimension > 3 || buffer == NULL)
160 return;
161
162 _time_plot_context_t *c = context;
163
164 fvm_to_time_plot_writer_t *w = c->writer;
165
166 assert(datatype == CS_REAL_TYPE);
167 assert(component_id == 0);
168
169 char *file_name;
170 FILE *_f;
171
172 const cs_real_t *coords = buffer;
173 const int n_coords = block_end - block_start;
174
175 char t_stamp[32];
176 if (w->nt >= 0)
177 sprintf(t_stamp, "_%.4i", w->nt);
178 else
179 t_stamp[0] = '\0';
180
181 size_t l = strlen(w->prefix) + strlen("coords") + strlen(t_stamp) + 4 + 1;
182
183 BFT_MALLOC(file_name, l, char);
184
185 if (w->format == CS_TIME_PLOT_DAT)
186 sprintf(file_name, "%scoords%s.dat", w->prefix, t_stamp);
187 else if (w->format == CS_TIME_PLOT_CSV)
188 sprintf(file_name, "%scoords%s.csv", w->prefix, t_stamp);
189
190 _f = fopen(file_name, "w");
191 if (_f == NULL) {
192 bft_error(__FILE__, __LINE__, errno,
193 _("Error opening file: \"%s\""), file_name);
194 return;
195 }
196
197 /* whitespace-separated format */
198
199 if (w->format == CS_TIME_PLOT_DAT) {
200
201 const char **probe_names = fvm_nodal_get_global_vertex_labels(c->mesh);
202
203 if (probe_names != NULL) {
204 fprintf(_f, "# Monitoring point names:\n");
205 for (int i = 0; i < n_coords; i++)
206 fprintf(_f, "# %6i %16s\n",
207 i+1, probe_names[i]);
208 fprintf(_f, "#\n");
209 }
210
211 fprintf(_f, _("# Monitoring point coordinates:\n"));
212
213 switch(dimension) {
214 case 3:
215 for (int i = 0; i < n_coords; i++)
216 fprintf(_f, "# %6i %14.7e %14.7e %14.7e\n",
217 i + 1,
218 coords[i*3],
219 coords[i*3 + 1],
220 coords[i*3 + 2]);
221 break;
222 case 2:
223 for (int i = 0; i < n_coords; i++)
224 fprintf(_f, "# %6i %14.7e %14.7e\n",
225 i + 1,
226 coords[i*2],
227 coords[i*2 + 1]);
228 break;
229 case 1:
230 for (int i = 0; i < n_coords; i++)
231 fprintf(_f, "# %6i %14.7e\n",
232 i + 1,
233 coords[i]);
234 break;
235 }
236 fprintf(_f, "#\n");
237
238 }
239
240 /* CSV format */
241
242 else if (w->format == CS_TIME_PLOT_CSV) {
243
244 switch(dimension) {
245 case 3:
246 fprintf(_f, "x, y, z\n");
247 for (int i = 0; i < n_coords; i++) {
248 fprintf(_f, "%14.7e, %14.7e, %14.7e\n",
249 coords[i*3], coords[i*3 + 1], coords[i*3 + 2]);
250 }
251 break;
252 case 2:
253 fprintf(_f, "x, y\n");
254 for (int i = 0; i < n_coords; i++) {
255 fprintf(_f, "%14.7e, %14.7e\n",
256 coords[i*2], coords[i*2 + 1]);
257 }
258 break;
259 case 1:
260 fprintf(_f, "x\n");
261 for (int i = 0; i < n_coords; i++) {
262 fprintf(_f, "%14.7e\n", coords[i]);
263 }
264 break;
265 }
266
267 }
268
269 /* Close file */
270
271 if (fclose(_f) != 0)
272 bft_error(__FILE__, __LINE__, errno,
273 _("Error closing file: \"%s\""), file_name);
274
275 BFT_FREE(file_name);
276 }
277
278 /*----------------------------------------------------------------------------
279 * Output function for field values.
280 *
281 * This function is passed to fvm_writer_field_helper_output_* functions.
282 *
283 * parameters:
284 * context <-> pointer to writer and field context
285 * datatype <-- output datatype
286 * dimension <-- output field dimension
287 * component_id <-- output component id (if non-interleaved)
288 * block_start <-- start global number of element for current block
289 * block_end <-- past-the-end global number of element for current block
290 * buffer <-> associated output buffer
291 *----------------------------------------------------------------------------*/
292
293 static void
_field_output(void * context,cs_datatype_t datatype,int dimension,int component_id,cs_gnum_t block_start,cs_gnum_t block_end,void * buffer)294 _field_output(void *context,
295 cs_datatype_t datatype,
296 int dimension,
297 int component_id,
298 cs_gnum_t block_start,
299 cs_gnum_t block_end,
300 void *buffer)
301 {
302 CS_UNUSED(datatype);
303
304 cs_time_plot_t *p;
305
306 _time_plot_context_t *c = context;
307 fvm_to_time_plot_writer_t *w = c->writer;
308
309 if (buffer == NULL) return;
310
311 assert(component_id == 0);
312
313 int n_vals = block_end - block_start;
314
315 cs_real_t *_vals = NULL;
316
317 if (dimension > 1)
318 BFT_MALLOC(_vals, n_vals, cs_real_t);
319
320 for (int _component_id = 0; _component_id < dimension; _component_id++) {
321
322 /* Build plot name */
323
324 char tmpn[128], tmpe[6];
325
326 char *plot_name = tmpn;
327
328 fvm_writer_field_component_name(tmpe, 6, false, dimension, _component_id);
329
330 size_t lce = strlen(tmpe);
331 size_t l = strlen(c->name) + 1;
332
333 if (lce > 0)
334 l += 2 + lce;
335
336 if (l > 128)
337 BFT_MALLOC(plot_name, l, char);
338
339 if (lce > 0)
340 sprintf(plot_name, "%s[%s]", c->name, tmpe);
341 else
342 strcpy(plot_name, c->name);
343
344 int p_id = cs_map_name_to_id(w->f_map, plot_name);
345
346 if (p_id >= w->n_plots) {
347
348 w->n_plots += 1;
349 BFT_REALLOC(w->tp, w->n_plots, cs_time_plot_t *);
350
351 int n_probes = (block_end > block_start) ? block_end - block_start : 0;
352
353 const char **probe_names = fvm_nodal_get_global_vertex_labels(c->mesh);
354
355 w->tp[p_id] = cs_time_plot_init_probe(plot_name,
356 w->prefix,
357 w->format,
358 w->use_iteration,
359 w->flush_wtime,
360 w->n_buf_steps,
361 n_probes,
362 NULL,
363 NULL, /* probe_coords */
364 probe_names);
365
366 }
367
368 if (plot_name != tmpn)
369 BFT_FREE(plot_name);
370
371 p = w->tp[p_id];
372
373 if (p != NULL) {
374 const cs_real_t *vals = (const cs_real_t *)buffer;
375 if (dimension > 1) {
376 for (int i = 0; i < n_vals; i++)
377 _vals[i] = vals[i*dimension + _component_id];
378 vals = _vals;
379 }
380 cs_time_plot_vals_write(p,
381 w->nt,
382 w->t,
383 n_vals,
384 vals);
385 }
386 }
387
388 BFT_FREE(_vals);
389 }
390
391 /*! (DOXYGEN_SHOULD_SKIP_THIS) \endcond */
392
393 /*============================================================================
394 * Public function definitions
395 *============================================================================*/
396
397 /*----------------------------------------------------------------------------
398 * Initialize FVM to time plot file writer.
399 *
400 * Options are:
401 * csv output CSV (comma-separated-values) files
402 * dat output dat (space-separated) files
403 * use_iteration use time step id instead of time value for
404 * first column
405 * flush_wtime=<wt> flush output file every 'wt' seconds
406 * n_buf_steps=<n> write output to file every 'n' output steps
407 *
408 * parameters:
409 * name <-- base output case name.
410 * options <-- whitespace separated, lowercase options list
411 * time_dependecy <-- indicates if and how meshes will change with time
412 * comm <-- associated MPI communicator.
413 *
414 * returns:
415 * pointer to opaque time plot writer structure.
416 *----------------------------------------------------------------------------*/
417
418 #if defined(HAVE_MPI)
419 void *
fvm_to_time_plot_init_writer(const char * name,const char * path,const char * options,fvm_writer_time_dep_t time_dependency,MPI_Comm comm)420 fvm_to_time_plot_init_writer(const char *name,
421 const char *path,
422 const char *options,
423 fvm_writer_time_dep_t time_dependency,
424 MPI_Comm comm)
425 #else
426 void *
427 fvm_to_time_plot_init_writer(const char *name,
428 const char *path,
429 const char *options,
430 fvm_writer_time_dep_t time_dependency)
431 #endif
432 {
433 CS_NO_WARN_IF_UNUSED(time_dependency);
434
435 fvm_to_time_plot_writer_t *w = NULL;
436
437 /* Initialize writer */
438
439 BFT_MALLOC(w, 1, fvm_to_time_plot_writer_t);
440
441 BFT_MALLOC(w->name, strlen(name) + 1, char);
442 strcpy(w->name, name);
443
444 if (strlen(name) > 0) {
445 BFT_MALLOC(w->prefix, strlen(path) + 1 + strlen(name) + 1, char);
446 sprintf(w->prefix, "%s%s_", path, name);
447 }
448 else {
449 BFT_MALLOC(w->prefix, strlen(path) + 1, char);
450 strcpy(w->prefix, path);
451 }
452
453 w->rank = 0;
454 w->n_ranks = 1;
455
456 #if defined(HAVE_MPI)
457 {
458 int mpi_flag, rank, n_ranks;
459 w->comm = MPI_COMM_NULL;
460 MPI_Initialized(&mpi_flag);
461 if (mpi_flag && comm != MPI_COMM_NULL) {
462 w->comm = comm;
463 MPI_Comm_rank(w->comm, &rank);
464 MPI_Comm_size(w->comm, &n_ranks);
465 w->rank = rank;
466 w->n_ranks = n_ranks;
467 }
468 }
469 #endif /* defined(HAVE_MPI) */
470
471 /* Defaults */
472
473 w->format = CS_TIME_PLOT_CSV;
474
475 cs_time_plot_get_flush_default(&(w->flush_wtime),
476 &(w->n_buf_steps));
477
478 w->use_iteration = false;
479
480 w->nt = -1;
481 w->t = -1;
482
483 w->n_plots = 0;
484 w->f_map = (w->rank > 0) ? NULL : cs_map_name_to_id_create();
485 w->tp = NULL;
486
487 /* Parse options */
488
489 if (options != NULL) {
490
491 int i1, i2, l_opt;
492 int l_tot = strlen(options);
493
494 i1 = 0; i2 = 0;
495 while (i1 < l_tot) {
496
497 for (i2 = i1; i2 < l_tot && options[i2] != ' '; i2++);
498 l_opt = i2 - i1;
499
500 if ((l_opt == 3) && (strncmp(options + i1, "csv", l_opt) == 0))
501 w->format = CS_TIME_PLOT_CSV;
502 else if ((l_opt == 3) && (strncmp(options + i1, "dat", l_opt) == 0))
503 w->format = CS_TIME_PLOT_DAT;
504 else if ((l_opt == 13) && (strcmp(options + i1, "use_iteration") == 0))
505 w->use_iteration = true;
506 else if (strncmp(options + i1, "n_buf_steps=", 12) == 0) {
507 const char *s = options + i1 + 12;
508 int nb, nr;
509 nr = sscanf(s, "%d", &nb);
510 if (nr == 1)
511 w->n_buf_steps = nb;
512 }
513 else if (strncmp(options + i1, "flush_wtime=", 12) == 0) {
514 const char *s = options + i1 + 12;
515 int nr;
516 float wt;
517 nr = sscanf(s, "%g", &wt);
518 if (nr == 1)
519 w->flush_wtime = wt;
520 }
521
522 for (i1 = i2 + 1 ; i1 < l_tot && options[i1] == ' ' ; i1++);
523
524 }
525
526 }
527
528 /* Return writer */
529
530 return w;
531 }
532
533 /*----------------------------------------------------------------------------
534 * Finalize FVM to time plot file writer.
535 *
536 * parameters:
537 * writer <-- pointer to opaque time plot writer structure.
538 *
539 * returns:
540 * NULL pointer
541 *----------------------------------------------------------------------------*/
542
543 void *
fvm_to_time_plot_finalize_writer(void * writer)544 fvm_to_time_plot_finalize_writer(void *writer)
545 {
546 fvm_to_time_plot_writer_t *w = (fvm_to_time_plot_writer_t *)writer;
547
548 BFT_FREE(w->name);
549 BFT_FREE(w->prefix);
550
551 if (w->rank <= 0) {
552 for (int i = 0; i < w->n_plots; i++)
553 cs_time_plot_finalize(&(w->tp[i]));
554 BFT_FREE(w->tp);
555 cs_map_name_to_id_destroy(&(w->f_map));
556 }
557
558 BFT_FREE(w);
559
560 return NULL;
561 }
562
563 /*----------------------------------------------------------------------------
564 * Associate new time step with a time plot geometry.
565 *
566 * parameters:
567 * writer <-- pointer to associated writer
568 * time_step <-- time step number
569 * time_value <-- time_value number
570 *----------------------------------------------------------------------------*/
571
572 void
fvm_to_time_plot_set_mesh_time(void * writer,const int time_step,const double time_value)573 fvm_to_time_plot_set_mesh_time(void *writer,
574 const int time_step,
575 const double time_value)
576 {
577 fvm_to_time_plot_writer_t *w = (fvm_to_time_plot_writer_t *)writer;
578
579 w->nt = time_step;
580 w->t = time_value;
581 }
582
583 /*----------------------------------------------------------------------------
584 * Write nodal mesh to a an time plot file
585 *
586 * parameters:
587 * writer <-- pointer to associated writer
588 * mesh <-- pointer to nodal mesh structure that should be written
589 *----------------------------------------------------------------------------*/
590
591 void
fvm_to_time_plot_export_nodal(void * writer,const fvm_nodal_t * mesh)592 fvm_to_time_plot_export_nodal(void *writer,
593 const fvm_nodal_t *mesh)
594 {
595 if (fvm_nodal_get_max_entity_dim(mesh) == 0) {
596
597 fvm_to_time_plot_writer_t *w
598 = (fvm_to_time_plot_writer_t *)writer;
599
600 fvm_writer_field_helper_t *helper
601 = fvm_writer_field_helper_create(mesh,
602 NULL, /* section list */
603 mesh->dim,
604 CS_INTERLACE,
605 CS_REAL_TYPE,
606 FVM_WRITER_PER_NODE);
607
608 #if defined(HAVE_MPI)
609 if (w->n_ranks > 1)
610 fvm_writer_field_helper_init_g(helper,
611 w->n_ranks,
612 0,
613 w->comm);
614 #endif
615
616 _time_plot_context_t c = {.writer = w,
617 .mesh = mesh,
618 .name = NULL};
619
620 int n_parent_lists = (mesh->parent_vertex_num != NULL) ? 1 : 0;
621 cs_lnum_t parent_num_shift[1] = {0};
622 const cs_real_t *coo_ptr[1] = {mesh->vertex_coords};
623
624 fvm_writer_field_helper_output_n(helper,
625 &c,
626 mesh,
627 mesh->dim,
628 CS_INTERLACE,
629 NULL,
630 n_parent_lists,
631 parent_num_shift,
632 CS_COORD_TYPE,
633 (const void **)coo_ptr,
634 _coords_output);
635
636 fvm_writer_field_helper_destroy(&helper);
637 }
638 }
639
640 /*----------------------------------------------------------------------------
641 * Write field associated with a nodal mesh to an time plot file.
642 *
643 * Assigning a negative value to the time step indicates a time-independent
644 * field (in which case the time_value argument is unused).
645 *
646 * parameters:
647 * writer <-- pointer to associated writer
648 * mesh <-- pointer to associated nodal mesh structure
649 * name <-- variable name
650 * location <-- variable definition location (nodes or elements)
651 * dimension <-- variable dimension (0: constant, 1: scalar,
652 * 3: vector, 6: sym. tensor, 9: asym. tensor)
653 * interlace <-- indicates if variable in memory is interlaced
654 * n_parent_lists <-- indicates if variable values are to be obtained
655 * directly through the local entity index (when 0) or
656 * through the parent entity numbers (when 1 or more)
657 * parent_num_shift <-- parent number to value array index shifts;
658 * size: n_parent_lists
659 * datatype <-- indicates the data type of (source) field values
660 * time_step <-- number of the current time step
661 * time_value <-- associated time value
662 * field_values <-- array of associated field value arrays
663 *----------------------------------------------------------------------------*/
664
665 void
fvm_to_time_plot_export_field(void * writer,const fvm_nodal_t * mesh,const char * name,fvm_writer_var_loc_t location,int dimension,cs_interlace_t interlace,int n_parent_lists,const cs_lnum_t parent_num_shift[],cs_datatype_t datatype,int time_step,double time_value,const void * const field_values[])666 fvm_to_time_plot_export_field(void *writer,
667 const fvm_nodal_t *mesh,
668 const char *name,
669 fvm_writer_var_loc_t location,
670 int dimension,
671 cs_interlace_t interlace,
672 int n_parent_lists,
673 const cs_lnum_t parent_num_shift[],
674 cs_datatype_t datatype,
675 int time_step,
676 double time_value,
677 const void *const field_values[])
678 {
679 fvm_to_time_plot_writer_t *w = (fvm_to_time_plot_writer_t *)writer;
680
681 /* If time step changes, update it */
682
683 if (time_step != w->nt)
684 fvm_to_time_plot_set_mesh_time(writer,
685 time_step,
686 time_value);
687
688 /* Initialize writer helper */
689
690 cs_datatype_t dest_datatype = CS_REAL_TYPE;
691
692 if (datatype >= CS_INT32 && datatype <= CS_UINT64)
693 dest_datatype = CS_INT64;
694
695 fvm_writer_field_helper_t *helper
696 = fvm_writer_field_helper_create(mesh,
697 NULL, /* section list */
698 dimension,
699 CS_INTERLACE,
700 dest_datatype,
701 location);
702
703 #if defined(HAVE_MPI)
704
705 if (w->n_ranks > 1)
706 fvm_writer_field_helper_init_g(helper,
707 w->n_ranks,
708 0,
709 w->comm);
710
711 #endif
712
713 /* Per node variable only */
714
715 if (location == FVM_WRITER_PER_NODE) {
716
717 _time_plot_context_t c = {.writer = w,
718 .mesh = mesh,
719 .name = name};
720
721 fvm_writer_field_helper_output_n(helper,
722 &c,
723 mesh,
724 dimension,
725 interlace,
726 NULL,
727 n_parent_lists,
728 parent_num_shift,
729 datatype,
730 field_values,
731 _field_output);
732
733 }
734
735 /* Free helper structures */
736
737 fvm_writer_field_helper_destroy(&helper);
738 }
739
740 /*----------------------------------------------------------------------------*/
741
742 END_C_DECLS
743