1 /*============================================================================
2 * Write a nodal representation associated with a mesh and associated
3 * variables to histogram 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 <ctype.h>
36 #include <errno.h>
37 #include <float.h>
38 #include <stdio.h>
39 #include <stdlib.h>
40 #include <string.h>
41
42 /*----------------------------------------------------------------------------
43 * Local headers
44 *----------------------------------------------------------------------------*/
45
46 #include "bft_error.h"
47 #include "bft_mem.h"
48 #include "bft_printf.h"
49
50 #include "cs_base.h"
51 #include "fvm_defs.h"
52
53 #include "fvm_convert_array.h"
54 #include "fvm_io_num.h"
55 #include "fvm_nodal.h"
56 #include "fvm_nodal_priv.h"
57 #include "fvm_to_vtk_histogram.h"
58 #include "fvm_writer_helper.h"
59 #include "fvm_writer_priv.h"
60
61 #include "cs_block_dist.h"
62 #include "cs_file.h"
63 #include "cs_parall.h"
64 #include "cs_part_to_block.h"
65
66 /*----------------------------------------------------------------------------
67 * Header for the current file
68 *----------------------------------------------------------------------------*/
69
70 #include "fvm_to_histogram.h"
71
72 /*----------------------------------------------------------------------------*/
73
74 BEGIN_C_DECLS
75
76 /*! \cond DOXYGEN_SHOULD_SKIP_THIS */
77
78 /*=============================================================================
79 * Local Macro Definitions
80 *============================================================================*/
81
82 /*============================================================================
83 * Local Type Definitions
84 *============================================================================*/
85
86 /*----------------------------------------------------------------------------
87 * Context structure for fvm_writer_field_helper_output_* functions.
88 *----------------------------------------------------------------------------*/
89
90 typedef struct {
91
92 fvm_to_histogram_writer_t *writer; /* Pointer to writer structure */
93
94 const char *name; /* Field name */
95
96 } _histogram_context_t;
97
98 /*============================================================================
99 * Static global variables
100 *============================================================================*/
101
102 #if defined(HAVE_PLUGIN_CATALYST)
103 static void * _catalyst_plugin = NULL;
104 #endif
105
106 #if defined(HAVE_CATALYST)
107 static fvm_to_histogram_display_t *_fvm_to_vtk_display_histogram_png = NULL;
108 #endif
109
110 /*============================================================================
111 * Private function definitions
112 *============================================================================*/
113
114 #if defined(HAVE_PLUGIN_CATALYST)
115
116 /*----------------------------------------------------------------------------
117 * Load Catalyst plugin.
118 *----------------------------------------------------------------------------*/
119
120 static void
_load_plugin_catalyst(void)121 _load_plugin_catalyst(void)
122 {
123 /* Open from shared library */
124
125 _catalyst_plugin = cs_base_dlopen_plugin("fvm_catalyst");
126
127 /* Load symbols from shared library */
128
129 /* Function pointers need to be double-casted so as to first convert
130 a (void *) type to a memory address and then convert it back to the
131 original type. Otherwise, the compiler may issue a warning.
132 This is a valid ISO C construction. */
133
134 _fvm_to_vtk_display_histogram_png = (fvm_to_histogram_display_t *) (intptr_t)
135 cs_base_get_dl_function_pointer(_catalyst_plugin,
136 "fvm_to_vtk_display_histogram_png",
137 true);
138 }
139
140 #endif /* defined(HAVE_PLUGIN_CATALYST)*/
141
142 /*----------------------------------------------------------------------------*/
143 /*!
144 * Display histograms in txt format.
145 *
146 * parameters:
147 * var_min <-- minimum variable value
148 * var_max <-- maximum variable value
149 * count <-- count for each histogram slice
150 * w <--> histogram writer
151 * var_name <-- name of the variable
152 */
153 /*----------------------------------------------------------------------------*/
154
155 static void
_display_histogram_txt(cs_real_t var_min,cs_real_t var_max,cs_gnum_t count[],fvm_to_histogram_writer_t * w,char * var_name)156 _display_histogram_txt(cs_real_t var_min,
157 cs_real_t var_max,
158 cs_gnum_t count[],
159 fvm_to_histogram_writer_t *w,
160 char *var_name)
161 {
162 double var_step;
163 int i, j;
164
165 /* Open the txt file */
166 w->f = fopen(w->file_name, "w");
167
168 if (w->f == NULL) {
169 bft_error(__FILE__, __LINE__, errno,
170 _("Error opening file: \"%s\""), w->file_name);
171 return;
172 }
173
174 /* Print header */
175 fprintf(w->f, "# Code_Saturne histograms\n#\n");
176 if (w->nt < 0)
177 fprintf(w->f, "# time independent\n");
178 else {
179 fprintf(w->f, "# time step id: %i\n", w->nt);
180 fprintf(w->f, "# time: %12.5e\n#\n", w->t);
181 }
182
183 /* Print base min, max, and increment */
184 fprintf(w->f, "# Variable : %s\n\n",var_name);
185
186 fprintf(w->f, _(" minimum value = %10.5e\n"), (double)var_min);
187 fprintf(w->f, _(" maximum value = %10.5e\n\n"), (double)var_max);
188
189 var_step = CS_ABS(var_max - var_min) / w->n_sub;
190
191 if (CS_ABS(var_max - var_min) > 0.) {
192
193 /* Number of elements in each subdivision */
194
195 for (i = 0, j = 1; i < w->n_sub - 1; i++, j++)
196 fprintf(w->f, " %3d : [ %10.5e ; %10.5e [ = %10llu\n",
197 i+1, var_min + i*var_step, var_min + j*var_step,
198 (unsigned long long)(count[i]));
199
200 fprintf(w->f, " %3d : [ %10.5e ; %10.5e ] = %10llu\n",
201 w->n_sub,
202 var_min + (w->n_sub - 1)*var_step,
203 var_max,
204 (unsigned long long)(count[w->n_sub - 1]));
205
206 }
207 }
208
209 /*----------------------------------------------------------------------------*/
210 /*!
211 * Display histograms in tex format.
212 *
213 * parameters:
214 * var_min <-- minimum variable value
215 * var_max <-- maximum variable value
216 * count <-- count for each histogram slice
217 * w <--> histogram writer
218 * var_name <-- name of the variable
219 */
220 /*----------------------------------------------------------------------------*/
221
222 static void
_display_histogram_tex(cs_real_t var_min,cs_real_t var_max,cs_gnum_t count[],fvm_to_histogram_writer_t * w,char * var_name)223 _display_histogram_tex(cs_real_t var_min,
224 cs_real_t var_max,
225 cs_gnum_t count[],
226 fvm_to_histogram_writer_t *w,
227 char *var_name)
228 {
229 double var_step = CS_ABS(var_max - var_min) / w->n_sub;
230 int i;
231
232 if (var_step > 0) {
233 /* Open the tex file if non-zero histogram */
234 w->f = fopen(w->file_name, "w");
235
236 if (w->f == NULL) {
237 bft_error(__FILE__, __LINE__, errno,
238 _("Error opening file: \"%s\""), w->file_name);
239 return;
240 }
241
242 fprintf(w->f, "\\begin{center}\n");
243 fprintf(w->f, "\\begin{tikzpicture}[font=\\footnotesize]\n");
244
245 fprintf(w->f, " \\begin{axis}[\n");
246 fprintf(w->f, " ybar,\n");
247 fprintf(w->f, " bar width=18pt,\n");
248 fprintf(w->f, " xlabel={%s},\n",var_name);
249 fprintf(w->f, " ylabel={Number of elements},\n");
250 fprintf(w->f, " ymin=0,\n");
251 fprintf(w->f, " ytick=\\empty,\n");
252 fprintf(w->f, " xtick=data,\n");
253 fprintf(w->f, " axis x line=bottom,\n");
254 fprintf(w->f, " axis y line=left,\n");
255 fprintf(w->f, " enlarge x limits=0.06,\n");
256
257 fprintf(w->f, " symbolic x coords={");
258 for (i = 0 ; i < w->n_sub-1 ; i++)
259 fprintf(w->f, "%.3e,",var_min + (i+0.5)*var_step);
260 fprintf(w->f, "%.3e},\n",var_min + (w->n_sub-0.5)*var_step);
261
262 fprintf(w->f, " xticklabel style={rotate=45,font=\\scriptsize},\n");
263 fprintf(w->f,
264 " nodes near coords={\\pgfmathprintnumber\\pgfplotspointmeta}\n");
265 fprintf(w->f, " ]\n");
266 fprintf(w->f, " \\addplot[fill=blue] coordinates {\n");
267 for (i = 0 ; i < w->n_sub ; i++)
268 fprintf(w->f, " (%.3e,%llu)\n",
269 var_min + (i+0.5)*var_step, (unsigned long long)count[i]);
270 fprintf(w->f, " };\n");
271 fprintf(w->f, " \\end{axis}\n");
272
273 fprintf(w->f, "\\end{tikzpicture}\n");
274 fprintf(w->f, "\\end{center}\n");
275 }
276 }
277
278 /*----------------------------------------------------------------------------
279 * Compute the minimum and the maximum of a vector (locally).
280 *
281 * parameters:
282 * n_vals <-- local number of elements
283 * var <-- pointer to vector
284 * min --> minimum
285 * max --> maximum
286 *----------------------------------------------------------------------------*/
287
288 static void
_compute_local_minmax(cs_lnum_t n_vals,const cs_real_t var[],cs_real_t * min,cs_real_t * max)289 _compute_local_minmax(cs_lnum_t n_vals,
290 const cs_real_t var[],
291 cs_real_t *min,
292 cs_real_t *max)
293 {
294 cs_lnum_t i;
295 cs_real_t _min = DBL_MAX, _max = -DBL_MAX;
296
297 for (i = 0; i < n_vals; i++) {
298 _min = CS_MIN(_min, var[i]);
299 _max = CS_MAX(_max, var[i]);
300 }
301
302 *min = _min;
303 *max = _max;
304 }
305
306 /*----------------------------------------------------------------------------
307 * Display the distribution of values of a real vector.
308 *
309 * parameters:
310 * var_min <-- minimum variable value
311 * var_max <-- maximum variable value
312 * count <-> count for each histogram slice (size: w->n_sub)
313 * local values in, global values out
314 * display_func <-- function pointer to display the histogram
315 * w <--> histogram writer
316 * var_name <-- name of the variable
317 *----------------------------------------------------------------------------*/
318
319 static void
_display_histograms(cs_real_t var_min,cs_real_t var_max,cs_gnum_t count[],fvm_to_histogram_display_t * display_func,fvm_to_histogram_writer_t * w,char * var_name)320 _display_histograms(cs_real_t var_min,
321 cs_real_t var_max,
322 cs_gnum_t count[],
323 fvm_to_histogram_display_t *display_func,
324 fvm_to_histogram_writer_t *w,
325 char *var_name)
326 {
327 int i;
328
329 #if defined(HAVE_MPI)
330
331 if (w->n_ranks > 1) {
332
333 cs_gnum_t *g_count = NULL;
334
335 BFT_MALLOC(g_count, w->n_sub, cs_gnum_t);
336
337 MPI_Allreduce(count, g_count, w->n_sub, CS_MPI_GNUM, MPI_SUM,
338 w->comm);
339
340 for (i = 0; i < w->n_sub; i++)
341 count[i] = g_count[i];
342
343 BFT_FREE(g_count);
344 }
345
346 #endif
347
348 if (w->rank == 0)
349 display_func(var_min, var_max, count, w, var_name);
350 }
351
352 /*----------------------------------------------------------------------------
353 * Display the distribution of values of a real vector on cells or
354 * boundary faces.
355 *
356 * parameters:
357 * n_vals <-- number of values
358 * var <-- variable values
359 * display_func <-- function pointer to display the histogram
360 * w <--> histogram writer
361 * var_name <-- name of the variable
362 *----------------------------------------------------------------------------*/
363
364 static void
_histogram(cs_lnum_t n_vals,const cs_real_t var[],fvm_to_histogram_display_t * display_func,fvm_to_histogram_writer_t * w,char * var_name)365 _histogram(cs_lnum_t n_vals,
366 const cs_real_t var[],
367 fvm_to_histogram_display_t *display_func,
368 fvm_to_histogram_writer_t *w,
369 char *var_name)
370 {
371 cs_lnum_t i;
372 int j, k;
373
374 cs_real_t step;
375 cs_real_t max, min, _max, _min;
376
377 cs_gnum_t *count = NULL;
378
379 BFT_MALLOC(count, w->n_sub, cs_gnum_t);
380
381 assert (sizeof(double) == sizeof(cs_real_t));
382
383 /* Compute global min and max */
384
385 _compute_local_minmax(n_vals, var, &_min, &_max);
386
387 /* Default initialization */
388
389 min = _min;
390 max = _max;
391
392 #if defined(HAVE_MPI)
393
394 if (w->n_ranks > 1) {
395 MPI_Allreduce(&_min, &min, 1, CS_MPI_REAL, MPI_MIN,
396 w->comm);
397
398 MPI_Allreduce(&_max, &max, 1, CS_MPI_REAL, MPI_MAX,
399 w->comm);
400 }
401
402 #endif
403
404 /* Define axis subdivisions */
405
406 for (j = 0; j < w->n_sub; j++)
407 count[j] = 0;
408
409 if (CS_ABS(max - min) > 0.) {
410
411 step = CS_ABS(max - min) / w->n_sub;
412
413 /* Loop on values */
414
415 for (i = 0; i < n_vals; i++) {
416
417 /* Associated subdivision */
418
419 for (j = 0, k = 1; k < w->n_sub; j++, k++) {
420 if (var[i] < min + k*step)
421 break;
422 }
423 count[j] += 1;
424
425 }
426
427 }
428
429 _display_histograms(min, max, count, display_func, w, var_name);
430
431 BFT_FREE(count);
432 }
433
434 /*----------------------------------------------------------------------------
435 * Output function for field values.
436 *
437 * This function is passed to fvm_writer_field_helper_output_* functions.
438 *
439 * parameters:
440 * context <-> pointer to writer and field context
441 * datatype <-- output datatype
442 * dimension <-- output field dimension
443 * component_id <-- output component id (if non-interleaved)
444 * block_start <-- start global number of element for current block
445 * block_end <-- past-the-end global number of element for current block
446 * buffer <-> associated output buffer
447 *----------------------------------------------------------------------------*/
448
449 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)450 _field_output(void *context,
451 cs_datatype_t datatype,
452 int dimension,
453 int component_id,
454 cs_gnum_t block_start,
455 cs_gnum_t block_end,
456 void *buffer)
457 {
458 CS_UNUSED(datatype);
459
460 _histogram_context_t *c = (_histogram_context_t *)context;
461
462 fvm_to_histogram_writer_t *w = c->writer;
463
464 const int n_vals = (block_end > block_start) ? block_end - block_start : 0;
465
466 char tmpn[128], tmpe[6];
467
468 char *var_name = tmpn;
469
470 fvm_writer_field_component_name(tmpe, 6, false, dimension, component_id);
471
472 size_t lce = strlen(tmpe);
473 size_t lv = strlen(c->name) + 1;
474
475 if (lce > 0)
476 lv += 2 + lce;
477
478 if (lv > 128)
479 BFT_MALLOC(var_name, lv, char);
480
481 if (lce > 0)
482 sprintf(var_name, "%s[%s]", c->name, tmpe);
483 else
484 strcpy(var_name, c->name);
485
486 /* File name */
487 char t_stamp[141];
488 sprintf(t_stamp, "_%s_%.4i", var_name, w->nt);
489 size_t l = strlen(w->path) + strlen(w->name)
490 + strlen(t_stamp) + 4 + 1;
491
492 BFT_REALLOC(w->file_name, l, char);
493
494 const cs_real_t *vals = (const cs_real_t *)buffer;
495
496 if (w->format == CS_HISTOGRAM_TXT) {
497
498 sprintf(w->file_name, "%s%s%s.txt", w->path, w->name, t_stamp);
499 _histogram(n_vals, vals, _display_histogram_txt, w, var_name);
500
501 } else if (w->format == CS_HISTOGRAM_TEX) {
502
503 sprintf(w->file_name, "%s%s%s.tex", w->path, w->name, t_stamp);
504 _histogram(n_vals, vals, _display_histogram_tex, w, var_name);
505
506 #if defined(HAVE_CATALYST)
507 } else if (w->format == CS_HISTOGRAM_PNG) {
508
509 sprintf(w->file_name, "%s%s%s.png", w->path, w->name, t_stamp);
510 _histogram(n_vals, vals, _fvm_to_vtk_display_histogram_png, w, var_name);
511
512 #endif
513 }
514 }
515
516 /*! (DOXYGEN_SHOULD_SKIP_THIS) \endcond */
517
518 /*============================================================================
519 * Public function definitions
520 *============================================================================*/
521
522 /*----------------------------------------------------------------------------
523 * Initialize FVM to histogram file writer.
524 *
525 * Options are:
526 * txt output txt (space-separated) files
527 * tex output TeX (TixZ) files
528 * png output PNG files
529 * [n_sub] number of subdivisions
530 *
531 * parameters:
532 * name <-- base output case name.
533 * options <-- whitespace separated, lowercase options list
534 * time_dependecy <-- indicates if and how meshes will change with time
535 * comm <-- associated MPI communicator.
536 *
537 * returns:
538 * pointer to opaque histogram writer structure.
539 *----------------------------------------------------------------------------*/
540
541 #if defined(HAVE_MPI)
542 void *
fvm_to_histogram_init_writer(const char * name,const char * path,const char * options,fvm_writer_time_dep_t time_dependency,MPI_Comm comm)543 fvm_to_histogram_init_writer(const char *name,
544 const char *path,
545 const char *options,
546 fvm_writer_time_dep_t time_dependency,
547 MPI_Comm comm)
548 #else
549 void *
550 fvm_to_histogram_init_writer(const char *name,
551 const char *path,
552 const char *options,
553 fvm_writer_time_dep_t time_dependency)
554 #endif
555 {
556 CS_UNUSED(time_dependency);
557
558 fvm_to_histogram_writer_t *w = NULL;
559
560 /* Initialize writer */
561
562 BFT_MALLOC(w, 1, fvm_to_histogram_writer_t);
563
564 BFT_MALLOC(w->name, strlen(name) + 1, char);
565 strcpy(w->name, name);
566
567 BFT_MALLOC(w->path, strlen(path) + 1, char);
568 strcpy(w->path, path);
569
570 w->rank = 0;
571 w->n_ranks = 1;
572
573 #if defined(HAVE_MPI)
574 {
575 int mpi_flag, rank, n_ranks;
576 w->comm = MPI_COMM_NULL;
577 MPI_Initialized(&mpi_flag);
578 if (mpi_flag && comm != MPI_COMM_NULL) {
579 w->comm = comm;
580 MPI_Comm_rank(w->comm, &rank);
581 MPI_Comm_size(w->comm, &n_ranks);
582 w->rank = rank;
583 w->n_ranks = n_ranks;
584 }
585 }
586 #endif /* defined(HAVE_MPI) */
587
588 /* Defaults */
589
590 w->format = CS_HISTOGRAM_TXT;
591
592 w->nt = -1;
593 w->t = -1;
594 w->buffer = NULL;
595
596 w->file_name = NULL;
597 w->f = NULL;
598
599 w->n_sub = 5; /* default */
600
601 /* Parse options */
602
603 if (options != NULL) {
604
605 int i1, i2, l_opt;
606 int l_tot = strlen(options);
607 bool n_sub_read = false;
608
609 i1 = 0; i2 = 0;
610 while (i1 < l_tot) {
611
612 for (i2 = i1; i2 < l_tot && options[i2] != ' '; i2++);
613 l_opt = i2 - i1;
614
615 if ((l_opt == 3) && (strncmp(options + i1, "txt", l_opt) == 0))
616 w->format = CS_HISTOGRAM_TXT;
617 else if ((l_opt == 3) && (strncmp(options + i1, "tex", l_opt) == 0)) {
618 w->format = CS_HISTOGRAM_TEX;
619 if (!n_sub_read)
620 w->n_sub = 10;
621 }
622 #if defined(HAVE_CATALYST)
623 else if ((l_opt == 3) && (strncmp(options + i1, "png", l_opt) == 0)) {
624 w->format = CS_HISTOGRAM_PNG;
625 if (!n_sub_read)
626 w->n_sub = 10;
627 #if defined(HAVE_PLUGIN_CATALYST)
628 _load_plugin_catalyst();
629 #else
630 _fvm_to_vtk_display_histogram_png = fvm_to_vtk_display_histogram_png;
631 #endif
632 }
633 #endif
634 else {
635 const char *n_sub_opt = options+i1;
636 while (*n_sub_opt != '\0' && !isdigit(*n_sub_opt))
637 n_sub_opt++;
638 if (atoi(n_sub_opt) > 0) {
639 w->n_sub = atoi(n_sub_opt);
640 n_sub_read = true;
641 }
642 }
643
644 for (i1 = i2 + 1 ; i1 < l_tot && options[i1] == ' ' ; i1++);
645
646 }
647
648 }
649
650 /* Return writer */
651
652 return w;
653 }
654
655 /*----------------------------------------------------------------------------
656 * Finalize FVM to histogram file writer.
657 *
658 * parameters:
659 * writer <-- pointer to opaque histogram Gold writer structure.
660 *
661 * returns:
662 * NULL pointer
663 *----------------------------------------------------------------------------*/
664
665 void *
fvm_to_histogram_finalize_writer(void * writer)666 fvm_to_histogram_finalize_writer(void *writer)
667 {
668 fvm_to_histogram_writer_t *w
669 = (fvm_to_histogram_writer_t *)writer;
670
671 BFT_FREE(w->name);
672 BFT_FREE(w->path);
673
674 fvm_to_histogram_flush(writer);
675
676 BFT_FREE(w->file_name);
677
678 #if defined(HAVE_PLUGIN_CATALYST)
679 if (w->format == CS_HISTOGRAM_PNG)
680 cs_base_dlclose("fvm_catalyst",
681 _catalyst_plugin); /* decrease reference count or free */
682 #endif
683
684 BFT_FREE(w);
685
686 return NULL;
687 }
688
689 /*----------------------------------------------------------------------------
690 * Associate new time step with a histogram geometry.
691 *
692 * parameters:
693 * writer <-- pointer to associated writer
694 * time_step <-- time step number
695 * time_value <-- time_value number
696 *----------------------------------------------------------------------------*/
697
698 void
fvm_to_histogram_set_mesh_time(void * writer,int time_step,double time_value)699 fvm_to_histogram_set_mesh_time(void *writer,
700 int time_step,
701 double time_value)
702 {
703 fvm_to_histogram_writer_t *w = (fvm_to_histogram_writer_t *)writer;
704
705 w->nt = time_step;
706 w->t = time_value;
707
708 fvm_to_histogram_flush(writer);
709 }
710
711 /*----------------------------------------------------------------------------
712 * Write field associated with a nodal mesh to a histogram file.
713 *
714 * Assigning a negative value to the time step indicates a time-independent
715 * field (in which case the time_value argument is unused).
716 *
717 * parameters:
718 * writer <-- pointer to associated writer
719 * mesh <-- pointer to associated nodal mesh structure
720 * name <-- variable name
721 * location <-- variable definition location (nodes or elements)
722 * dimension <-- variable dimension (0: constant, 1: scalar,
723 * 3: vector, 6: sym. tensor, 9: asym. tensor)
724 * interlace <-- indicates if variable in memory is interlaced
725 * n_parent_lists <-- indicates if variable values are to be obtained
726 * directly through the local entity index (when 0) or
727 * through the parent entity numbers (when 1 or more)
728 * parent_num_shift <-- parent number to value array index shifts;
729 * size: n_parent_lists
730 * datatype <-- indicates the data type of (source) field values
731 * time_step <-- number of the current time step
732 * time_value <-- associated time value
733 * field_values <-- array of associated field value arrays
734 *----------------------------------------------------------------------------*/
735
736 void
fvm_to_histogram_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[])737 fvm_to_histogram_export_field(void *writer,
738 const fvm_nodal_t *mesh,
739 const char *name,
740 fvm_writer_var_loc_t location,
741 int dimension,
742 cs_interlace_t interlace,
743 int n_parent_lists,
744 const cs_lnum_t parent_num_shift[],
745 cs_datatype_t datatype,
746 int time_step,
747 double time_value,
748 const void *const field_values[])
749 {
750 fvm_to_histogram_writer_t *w = (fvm_to_histogram_writer_t *)writer;
751
752 /* If time step changes, update it */
753
754 if (time_step != w->nt)
755 fvm_to_histogram_set_mesh_time(writer,
756 time_step,
757 time_value);
758
759 /* Initialize writer helper */
760
761 cs_datatype_t dest_datatype = CS_REAL_TYPE;
762
763 if (datatype >= CS_INT32 && datatype <= CS_UINT64)
764 dest_datatype = CS_INT64;
765
766 int export_dim = fvm_nodal_get_max_entity_dim(mesh);
767
768 fvm_writer_section_t *export_list
769 = fvm_writer_export_list(mesh,
770 export_dim,
771 export_dim,
772 -1,
773 true,
774 true,
775 false,
776 false,
777 false,
778 true);
779
780 fvm_writer_field_helper_t *helper
781 = fvm_writer_field_helper_create(mesh,
782 export_list,
783 dimension,
784 CS_NO_INTERLACE,
785 dest_datatype,
786 location);
787
788 #if defined(HAVE_MPI)
789
790 if (w->n_ranks > 1)
791 fvm_writer_field_helper_init_g(helper,
792 w->n_ranks,
793 0,
794 w->comm);
795
796 #endif
797
798 _histogram_context_t c = {.writer = w, .name = name};
799
800 fvm_writer_field_helper_output_e(helper,
801 &c,
802 export_list,
803 dimension,
804 interlace,
805 NULL,
806 n_parent_lists,
807 parent_num_shift,
808 datatype,
809 field_values,
810 _field_output);
811
812 BFT_FREE(export_list);
813
814 /* Free helper structures */
815
816 fvm_writer_field_helper_destroy(&helper);
817 }
818
819 /*----------------------------------------------------------------------------
820 * Flush files associated with a given writer.
821 *
822 * parameters:
823 * writer <-- pointer to associated writer
824 *----------------------------------------------------------------------------*/
825
826 void
fvm_to_histogram_flush(void * writer)827 fvm_to_histogram_flush(void *writer)
828 {
829 fvm_to_histogram_writer_t *w = (fvm_to_histogram_writer_t *)writer;
830
831 if (w->f != NULL && w->buffer != NULL) {
832
833 if (fclose(w->f) != 0)
834 bft_error(__FILE__, __LINE__, errno,
835 _("Error closing file: \"%s\""), w->file_name);
836
837 w->f = NULL;
838
839 }
840
841 BFT_FREE(w->buffer);
842 }
843
844 /*----------------------------------------------------------------------------*/
845
846 END_C_DECLS
847