1 /*============================================================================
2 * Write a nodal representation associated with a mesh and associated
3 * variables to 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 "bft_printf.h"
53 #include "fvm_writer_helper.h"
54 #include "fvm_writer_priv.h"
55
56 #include "cs_block_dist.h"
57 #include "cs_file.h"
58 #include "cs_parall.h"
59 #include "cs_part_to_block.h"
60
61 /*----------------------------------------------------------------------------
62 * Header for the current file
63 *----------------------------------------------------------------------------*/
64
65 #include "fvm_to_plot.h"
66
67 /*----------------------------------------------------------------------------*/
68
69 BEGIN_C_DECLS
70
71 /*! \cond DOXYGEN_SHOULD_SKIP_THIS */
72
73 /*============================================================================
74 * Local Type Definitions
75 *============================================================================*/
76
77 /* Type of 1D plot file format */
78
79 typedef enum {
80 CS_PLOT_DAT, /* .dat file (usable by Qtplot or Grace) */
81 CS_PLOT_CSV /* .csv file (readable by ParaView or spreadsheat) */
82 } cs_plot_format_t;
83
84 /*----------------------------------------------------------------------------
85 * Plot writer structure
86 *----------------------------------------------------------------------------*/
87
88 typedef struct {
89
90 char *name; /* Writer name */
91 char *path; /* Path prefix */
92
93 int rank; /* Rank of current process in communicator */
94 int n_ranks; /* Number of processes in communicator */
95
96 cs_plot_format_t format; /* Plot format */
97
98 bool no_time_step; /* Do not append time step to file name */
99
100 int nt; /* Time step */
101 double t; /* Time value */
102
103 int n_cols; /* Number of columns */
104 int n_cols_max; /* Max. number of columns */
105 int n_rows; /* Number of rows */
106
107 cs_real_t *buffer; /* Values buffer */
108
109 char *file_name; /* File name */
110 FILE *f; /* Associated plots */
111
112 #if defined(HAVE_MPI)
113 MPI_Comm comm; /* Associated MPI communicator */
114 #endif
115
116 } fvm_to_plot_writer_t;
117
118 /*----------------------------------------------------------------------------
119 * Context structure for fvm_writer_field_helper_output_* functions.
120 *----------------------------------------------------------------------------*/
121
122 typedef struct {
123
124 fvm_to_plot_writer_t *writer; /* Pointer to writer structure */
125
126 const char *name; /* Field name */
127
128 } _plot_context_t;
129
130 /*============================================================================
131 * Static global variables
132 *============================================================================*/
133
134 /*============================================================================
135 * Private function definitions
136 *============================================================================*/
137
138 /*----------------------------------------------------------------------------
139 * Output function for field values.
140 *
141 * This function is passed to fvm_writer_field_helper_output_* functions.
142 *
143 * parameters:
144 * context <-> pointer to writer and field context
145 * datatype <-- output datatype
146 * dimension <-- output field dimension
147 * component_id <-- output component id (if non-interleaved)
148 * block_start <-- start global number of element for current block
149 * block_end <-- past-the-end global number of element for current block
150 * buffer <-> associated output buffer
151 *----------------------------------------------------------------------------*/
152
153 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)154 _field_output(void *context,
155 cs_datatype_t datatype,
156 int dimension,
157 int component_id,
158 cs_gnum_t block_start,
159 cs_gnum_t block_end,
160 void *buffer)
161 {
162 CS_UNUSED(datatype);
163 CS_UNUSED(component_id);
164
165 _plot_context_t *c = context;
166
167 fvm_to_plot_writer_t *w = c->writer;
168
169 if (w->rank > 0)
170 return;
171
172 const int n_rows = (block_end > block_start) ? block_end - block_start : 0;
173
174 if (w->n_cols == 0)
175 w->n_rows = n_rows;
176 else if (w->n_rows != n_rows) {
177 const char e[] = "";
178 const char *name = (c->name != NULL) ? c->name : e;
179 bft_printf(_("Warning: inconsistent data size for plot \"%s\" between\n"
180 "field \"%s\" and previous outputs; values dropped.\n"),
181 w->name, name);
182 return;
183 }
184
185 /* Create file and write headers on first output */
186
187 if (w->f == NULL) {
188
189 char t_stamp[32];
190 if (w->no_time_step || w->nt < 0)
191 t_stamp[0] = '\0';
192 else
193 sprintf(t_stamp, "_%.4i", w->nt);
194 size_t l = strlen(w->path) + strlen(w->name)
195 + strlen(t_stamp) + 4 + 1;
196 BFT_REALLOC(w->file_name, l, char);
197
198 if (w->format == CS_PLOT_DAT)
199 sprintf(w->file_name, "%s%s%s.dat", w->path, w->name, t_stamp);
200 else
201 sprintf(w->file_name, "%s%s%s.csv", w->path, w->name, t_stamp);
202
203 w->f = fopen(w->file_name, "w");
204
205 if (w->f == NULL) {
206 bft_error(__FILE__, __LINE__, errno,
207 _("Error opening file: \"%s\""), w->file_name);
208 return;
209 }
210
211 if (w->format == CS_PLOT_DAT) {
212 fprintf(w->f, "# Code_Saturne plot output\n#\n");
213 if (w->nt < 0)
214 fprintf(w->f, "# time independent\n");
215 else {
216 fprintf(w->f, "# time step id: %i\n", w->nt);
217 fprintf(w->f, "# time: %12.5e\n#\n", w->t);
218 }
219 fprintf(w->f, "#COLUMN_TITLES: ");
220 }
221
222 }
223
224 assert(component_id == 0);
225
226 if (w->n_cols + dimension > w->n_cols_max) {
227 while (w->n_cols + dimension > w->n_cols_max) {
228 if (w->n_cols_max == 0)
229 w->n_cols_max = 4;
230 else
231 w->n_cols_max *= 2;
232 }
233 BFT_REALLOC(w->buffer, w->n_rows*w->n_cols_max, cs_real_t);
234 }
235
236 for (int i = 0; i < dimension; i++) {
237
238 /* Add header info */
239
240 char name_buf[64];
241
242 if (c->name != NULL)
243 strncpy(name_buf, c->name, 63);
244 else
245 name_buf[0] = '\0';
246 name_buf[63] = '\0';
247
248 if (dimension > 1) {
249 size_t l = strlen(name_buf);
250 if (l > 59)
251 l = 59;
252 if (l > 0)
253 name_buf[l++] = '_';
254 fvm_writer_field_component_name(name_buf+l, 3, true, dimension, i);
255 }
256
257 if (w->format == CS_PLOT_DAT) {
258 if (w->n_cols > 0)
259 fprintf(w->f, " | %s", name_buf);
260 else
261 fprintf(w->f, " %s", name_buf);
262 }
263 else if (w->format == CS_PLOT_CSV) {
264 if (w->n_cols > 0)
265 fprintf(w->f, ", %s", name_buf);
266 else
267 fprintf(w->f, "%s", name_buf);
268 }
269
270 /* Update buffer */
271
272 const cs_real_t *src = buffer;
273 cs_real_t *dest = w->buffer + w->n_rows*w->n_cols;
274
275 for (cs_lnum_t j = 0; j < w->n_rows; j++)
276 dest[j] = src[j*dimension + i];
277 w->n_cols++;
278
279 }
280
281 }
282
283 /*----------------------------------------------------------------------------
284 * Close the file associated with a given writer.
285 *
286 * This assumes w->f != NULL on entry.
287 *
288 * parameters:
289 * w <-- pointer to associated writer
290 *----------------------------------------------------------------------------*/
291
292 static void
_file_close(fvm_to_plot_writer_t * w)293 _file_close(fvm_to_plot_writer_t *w)
294 {
295 assert(w->f != NULL);
296
297 if (fclose(w->f) != 0)
298 bft_error(__FILE__, __LINE__, errno,
299 _("Error closing file: \"%s\""), w->file_name);
300
301 w->f = NULL;
302 }
303
304 /*! (DOXYGEN_SHOULD_SKIP_THIS) \endcond */
305
306 /*============================================================================
307 * Public function definitions
308 *============================================================================*/
309
310 /*----------------------------------------------------------------------------
311 * Initialize FVM to plot file writer.
312 *
313 * Options are:
314 * csv output CSV (comma-separated-values) files
315 * dat output dat (space-separated) files
316 * no_time_step do not add time step id to file name
317 *
318 * parameters:
319 * name <-- base output case name.
320 * options <-- whitespace separated, lowercase options list
321 * time_dependecy <-- indicates if and how meshes will change with time
322 * comm <-- associated MPI communicator.
323 *
324 * returns:
325 * pointer to opaque plot writer structure.
326 *----------------------------------------------------------------------------*/
327
328 #if defined(HAVE_MPI)
329 void *
fvm_to_plot_init_writer(const char * name,const char * path,const char * options,fvm_writer_time_dep_t time_dependency,MPI_Comm comm)330 fvm_to_plot_init_writer(const char *name,
331 const char *path,
332 const char *options,
333 fvm_writer_time_dep_t time_dependency,
334 MPI_Comm comm)
335 #else
336 void *
337 fvm_to_plot_init_writer(const char *name,
338 const char *path,
339 const char *options,
340 fvm_writer_time_dep_t time_dependency)
341 #endif
342 {
343 CS_UNUSED(time_dependency);
344
345 fvm_to_plot_writer_t *w = NULL;
346
347 /* Initialize writer */
348
349 BFT_MALLOC(w, 1, fvm_to_plot_writer_t);
350
351 BFT_MALLOC(w->name, strlen(name) + 1, char);
352 strcpy(w->name, name);
353
354 BFT_MALLOC(w->path, strlen(path) + 1, char);
355 strcpy(w->path, path);
356
357 w->rank = 0;
358 w->n_ranks = 1;
359
360 #if defined(HAVE_MPI)
361 {
362 int mpi_flag, rank, n_ranks;
363 w->comm = MPI_COMM_NULL;
364 MPI_Initialized(&mpi_flag);
365 if (mpi_flag && comm != MPI_COMM_NULL) {
366 w->comm = comm;
367 MPI_Comm_rank(w->comm, &rank);
368 MPI_Comm_size(w->comm, &n_ranks);
369 w->rank = rank;
370 w->n_ranks = n_ranks;
371 }
372 }
373 #endif /* defined(HAVE_MPI) */
374
375 /* Defaults */
376
377 w->format = CS_PLOT_CSV;
378 w->no_time_step = false;
379
380 w->nt = -1;
381 w->t = -1;
382
383 w->n_cols = 0;
384 w->n_cols_max = 0;
385 w->n_rows = 0;
386
387 w->buffer = NULL;
388
389 w->file_name = NULL;
390 w->f = NULL;
391
392 /* Parse options */
393
394 if (options != NULL) {
395
396 int i1, i2, l_opt;
397 int l_tot = strlen(options);
398
399 i1 = 0; i2 = 0;
400 while (i1 < l_tot) {
401
402 for (i2 = i1; i2 < l_tot && options[i2] != ' '; i2++);
403 l_opt = i2 - i1;
404
405 if ((l_opt == 3) && (strncmp(options + i1, "csv", l_opt) == 0))
406 w->format = CS_PLOT_CSV;
407 else if ((l_opt == 3) && (strncmp(options + i1, "dat", l_opt) == 0))
408 w->format = CS_PLOT_DAT;
409 else if ((l_opt == 12) && (strncmp(options + i1, "no_time_step", l_opt) == 0))
410 w->no_time_step = true;
411
412 for (i1 = i2 + 1 ; i1 < l_tot && options[i1] == ' ' ; i1++);
413
414 }
415
416 }
417
418 /* Return writer */
419
420 return w;
421 }
422
423 /*----------------------------------------------------------------------------
424 * Finalize FVM to plot file writer.
425 *
426 * parameters:
427 * writer <-- pointer to opaque plot Gold writer structure.
428 *
429 * returns:
430 * NULL pointer
431 *----------------------------------------------------------------------------*/
432
433 void *
fvm_to_plot_finalize_writer(void * writer)434 fvm_to_plot_finalize_writer(void *writer)
435 {
436 fvm_to_plot_writer_t *w
437 = (fvm_to_plot_writer_t *)writer;
438
439 BFT_FREE(w->name);
440 BFT_FREE(w->path);
441
442 fvm_to_plot_flush(writer);
443
444 if (w->f != NULL)
445 _file_close(w);
446
447 BFT_FREE(w->file_name);
448
449 BFT_FREE(w);
450
451 return NULL;
452 }
453
454 /*----------------------------------------------------------------------------
455 * Associate new time step with an plot geometry.
456 *
457 * parameters:
458 * writer <-- pointer to associated writer
459 * time_step <-- time step number
460 * time_value <-- time_value number
461 *----------------------------------------------------------------------------*/
462
463 void
fvm_to_plot_set_mesh_time(void * writer,int time_step,double time_value)464 fvm_to_plot_set_mesh_time(void *writer,
465 int time_step,
466 double time_value)
467 {
468 fvm_to_plot_writer_t *w = (fvm_to_plot_writer_t *)writer;
469
470 if (w->nt != time_step) {
471 if (w->n_cols > 0)
472 fvm_to_plot_flush(writer);
473 if (w->f != NULL)
474 _file_close(w);
475 }
476 w->nt = time_step;
477 w->t = time_value;
478
479 }
480
481 /*----------------------------------------------------------------------------
482 * Write nodal mesh to a an plot file
483 *
484 * parameters:
485 * writer <-- pointer to associated writer
486 * mesh <-- pointer to nodal mesh structure that should be written
487 *----------------------------------------------------------------------------*/
488
489 void
fvm_to_plot_export_nodal(void * writer,const fvm_nodal_t * mesh)490 fvm_to_plot_export_nodal(void *writer,
491 const fvm_nodal_t *mesh)
492 {
493 CS_UNUSED(writer);
494 CS_UNUSED(mesh);
495
496 #if 0
497 if (fvm_nodal_get_max_entity_dim(mesh) == 0) {
498
499 fvm_to_plot_writer_t *w
500 = (fvm_to_plot_writer_t *)writer;
501
502 fvm_writer_field_helper_t *helper
503 = fvm_writer_field_helper_create(mesh,
504 NULL, /* section list */
505 mesh->dim,
506 CS_INTERLACE,
507 CS_REAL_TYPE,
508 FVM_WRITER_PER_NODE);
509
510 #if defined(HAVE_MPI)
511 if (w->n_ranks > 1)
512 fvm_writer_field_helper_init_g(helper,
513 w->n_ranks,
514 0,
515 w->comm);
516 #endif
517
518 _plot_context_t c = {.writer = w};
519
520 int n_parent_lists = (mesh->parent_vertex_num != NULL) ? 1 : 0;
521 cs_lnum_t parent_num_shift[1] = {0};
522 const cs_real_t *coo_ptr[1] = {mesh->vertex_coords};
523
524 fvm_writer_field_helper_output_n(helper,
525 &c,
526 mesh,
527 mesh->dim,
528 CS_INTERLACE,
529 NULL,
530 n_parent_lists,
531 parent_num_shift,
532 CS_COORD_TYPE,
533 (const void **)coo_ptr,
534 _field_output);
535
536 }
537
538 /* Free helper structures */
539
540 fvm_writer_field_helper_destroy(&helper);
541 #endif
542 }
543
544 /*----------------------------------------------------------------------------
545 * Write field associated with a nodal mesh to an plot file.
546 *
547 * Assigning a negative value to the time step indicates a time-independent
548 * field (in which case the time_value argument is unused).
549 *
550 * parameters:
551 * writer <-- pointer to associated writer
552 * mesh <-- pointer to associated nodal mesh structure
553 * name <-- variable name
554 * location <-- variable definition location (nodes or elements)
555 * dimension <-- variable dimension (0: constant, 1: scalar,
556 * 3: vector, 6: sym. tensor, 9: asym. tensor)
557 * interlace <-- indicates if variable in memory is interlaced
558 * n_parent_lists <-- indicates if variable values are to be obtained
559 * directly through the local entity index (when 0) or
560 * through the parent entity numbers (when 1 or more)
561 * parent_num_shift <-- parent number to value array index shifts;
562 * size: n_parent_lists
563 * datatype <-- indicates the data type of (source) field values
564 * time_step <-- number of the current time step
565 * time_value <-- associated time value
566 * field_values <-- array of associated field value arrays
567 *----------------------------------------------------------------------------*/
568
569 void
fvm_to_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[])570 fvm_to_plot_export_field(void *writer,
571 const fvm_nodal_t *mesh,
572 const char *name,
573 fvm_writer_var_loc_t location,
574 int dimension,
575 cs_interlace_t interlace,
576 int n_parent_lists,
577 const cs_lnum_t parent_num_shift[],
578 cs_datatype_t datatype,
579 int time_step,
580 double time_value,
581 const void *const field_values[])
582 {
583 fvm_to_plot_writer_t *w = (fvm_to_plot_writer_t *)writer;
584
585 /* If time step changes, update it */
586
587 if (time_step != w->nt)
588 fvm_to_plot_set_mesh_time(writer,
589 time_step,
590 time_value);
591
592 /* Initialize writer helper */
593
594 cs_datatype_t dest_datatype = CS_REAL_TYPE;
595
596 if (datatype >= CS_INT32 && datatype <= CS_UINT64)
597 dest_datatype = CS_INT64;
598
599 fvm_writer_field_helper_t *helper
600 = fvm_writer_field_helper_create(mesh,
601 NULL, /* section list */
602 dimension,
603 CS_INTERLACE,
604 dest_datatype,
605 location);
606
607 #if defined(HAVE_MPI)
608
609 if (w->n_ranks > 1)
610 fvm_writer_field_helper_init_g(helper,
611 w->n_ranks,
612 0,
613 w->comm);
614
615 #endif
616
617 /* Per node variable only */
618
619 if (location == FVM_WRITER_PER_NODE) {
620
621 _plot_context_t c = {.writer = w, .name = name};
622
623 fvm_writer_field_helper_output_n(helper,
624 &c,
625 mesh,
626 dimension,
627 interlace,
628 NULL,
629 n_parent_lists,
630 parent_num_shift,
631 datatype,
632 field_values,
633 _field_output);
634
635 }
636
637 /* Free helper structures */
638
639 fvm_writer_field_helper_destroy(&helper);
640 }
641
642 /*----------------------------------------------------------------------------
643 * Flush files associated with a given writer.
644 *
645 * In this case, the effective writing to file is done.
646 *
647 * parameters:
648 * writer <-- pointer to associated writer
649 *----------------------------------------------------------------------------*/
650
651 void
fvm_to_plot_flush(void * writer)652 fvm_to_plot_flush(void *writer)
653 {
654 fvm_to_plot_writer_t *w = (fvm_to_plot_writer_t *)writer;
655
656 if (w->f != NULL && w->buffer != NULL) {
657
658 /* Transpose output on write */
659
660 int n_cols = w->n_cols;
661 int n_rows = w->n_rows;
662
663 if (w->format == CS_PLOT_DAT) {
664
665 fprintf(w->f, "\n");
666
667 for (int i = 0; i < n_rows; i++) {
668 for (int j = 0; j < n_cols - 1; j++)
669 fprintf(w->f, "%12.5e ", w->buffer[w->n_rows*j + i]);
670 if (n_cols > 0) {
671 int j = n_cols -1;
672 fprintf(w->f, "%12.5e\n", w->buffer[w->n_rows*j + i]);
673 }
674 }
675
676 }
677 else if (w->format == CS_PLOT_CSV) {
678
679 fprintf(w->f, "\n");
680
681 for (int i = 0; i < n_rows; i++) {
682 for (int j = 0; j < n_cols-1; j++)
683 fprintf(w->f, "%12.5e, ", w->buffer[w->n_rows*j + i]);
684 if (n_cols > 0) {
685 int j = n_cols -1;
686 fprintf(w->f, "%12.5e\n", w->buffer[w->n_rows*j + i]);
687 }
688 }
689
690 }
691
692 w->n_rows = 0;
693 w->n_cols = 0;
694 w->n_cols_max = 0;
695
696 _file_close(w);
697
698 }
699
700 BFT_FREE(w->buffer);
701 }
702
703 /*----------------------------------------------------------------------------*/
704
705 END_C_DECLS
706