1 /*============================================================================
2  * Set of structures and functions to handle probes and profiles
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 <stdio.h>
34 #include <stdlib.h>
35 #include <string.h>
36 #include <assert.h>
37 
38 /*----------------------------------------------------------------------------
39  * Local headers
40  *----------------------------------------------------------------------------*/
41 
42 #include "bft_mem.h"
43 #include "bft_error.h"
44 #include "bft_printf.h"
45 
46 #include "fvm_nodal.h"
47 #include "fvm_point_location.h"
48 
49 #include "cs_base.h"
50 #include "cs_map.h"
51 #include "cs_math.h"
52 #include "cs_mesh.h"
53 #include "cs_mesh_connect.h"
54 #include "cs_mesh_location.h"
55 #include "cs_mesh_quantities.h"
56 #include "cs_mesh_intersect.h"
57 #include "cs_order.h"
58 #include "cs_parall.h"
59 #include "cs_selector.h"
60 #include "cs_timer.h"
61 
62 /*----------------------------------------------------------------------------
63  *  Header for the current file
64  *----------------------------------------------------------------------------*/
65 
66 #include "cs_probe.h"
67 
68 /*----------------------------------------------------------------------------*/
69 
70 BEGIN_C_DECLS
71 
72 /*=============================================================================
73  * Additional doxygen documentation
74  *============================================================================*/
75 
76 /*!
77   \file cs_probe.c
78 
79   \brief Probes and profiles management.
80 */
81 
82 /*! \cond DOXYGEN_SHOULD_SKIP_THIS */
83 
84 /*=============================================================================
85  * Local Macro Definitions
86  *============================================================================*/
87 
88 /* Predefined masks to set a flag for a probe set structure */
89 
90 #define CS_PROBE_TRANSIENT   (1 << 0) //   1: locations of probes may change
91 #define CS_PROBE_BOUNDARY    (1 << 1) //   2: locations only on the border mesh
92 #define CS_PROBE_ON_CURVE    (1 << 2) //   4: locations are on a curve
93 #define CS_PROBE_OVERWRITE   (1 << 3) //   8: allow re-creation of probe set
94 #define CS_PROBE_AUTO_VAR    (1 << 4) //  16: automatic output of variables
95 #define CS_PROBE_AUTO_S      (1 << 5) //  32: output curvilinear coordinates
96 #define CS_PROBE_AUTO_COORD  (1 << 6) //  64: output cartesian coordinates
97 
98 /*=============================================================================
99  * Local Structure Definitions
100  *============================================================================*/
101 
102 /* Structure to handle a set of probes */
103 
104 struct _cs_probe_set_t {
105 
106   char            *name;          /* Associated name */
107 
108   int32_t          flags;         /* Metadata related to a set of probes */
109 
110   char            *sel_criter;    /* Selection criterion to filter entities
111                                      before the location step */
112   double           tolerance;     /* Criterion to define a threshold during
113                                      the location step. This is a relative
114                                      criterion. */
115   cs_probe_snap_t  snap_mode;     /* Indicate how the positions of the
116                                      probes are computed */
117 
118   cs_lnum_t     n_max_probes;   /* Number of probes initially requested */
119   cs_lnum_t     n_probes;       /* Number of probes really used */
120   cs_lnum_t     n_loc_probes;   /* Number of probes located on this rank */
121 
122   cs_real_3_t  *coords_ini;     /* Initial coordinates of the set of probes
123                                    (only when snap mode is on and
124                                    location is transient) */
125 
126   cs_real_3_t  *coords;         /* Coordinates of the set of probes
127                                    Initially allocated to n_max_probes. */
128   cs_real_t    *s_coords;       /* NULL or curvilinear coordinates
129                                    for a profile */
130 
131   char        **labels;         /* List of labels for each probe (optional)
132                                    By default, this is set to NULL */
133 
134   cs_probe_set_define_local_t  *p_define_func;   /* Advanced local definition
135                                                     function */
136   void                         *p_define_input;  /* Advanced local definition
137                                                     input */
138   void                        *_p_define_input;  /* Advanced local definition
139                                                     input, if owner */
140 
141   cs_lnum_t    *loc_id;         /* ids of probes located on local domain */
142 
143   cs_lnum_t    *elt_id;         /* Element ids where the probes have been
144                                    located (size: n_loc_probes); -1 for
145                                    unlocated probes assigned to local domain*/
146   cs_lnum_t    *vtx_id;         /* Vertex ids closest to probes; -1 for
147                                    unlocated probes assigned to local domain*/
148 
149   char         *located;        /* 1 for located probes, 0 for unlocated */
150 
151   int           interpolation;  /* 0: no interpolation;
152                                    1: local gradient-based interpolation */
153 
154   /* User-defined writers associated to this set of probes */
155 
156   int           n_writers;      /* Number of writers (-1 if unset) */
157   int          *writer_ids;     /* List of writer ids */
158 
159   /* User-defined fields associated to this set of probes */
160 
161   int           n_fields;       /* Number of associated fields in list */
162   int           n_max_fields;   /* Maximum number of associated fields in list */
163   int          *field_info;     /* List of associated field tuples:
164                                    (writer_id, field_id, component_id) */
165 
166 };
167 
168 /* List of available keys for setting a set of probes */
169 
170 typedef enum {
171 
172   PSETKEY_TRANSIENT_LOC,
173   PSETKEY_BOUNDARY,
174   PSETKEY_SELECT_CRIT,
175   PSETKEY_TOLERANCE,
176   PSETKEY_INTERPOLATION,
177   PSETKEY_ERROR
178 
179 } psetkey_t;
180 
181 /*============================================================================
182  *  Global static variables
183  *============================================================================*/
184 
185 /* Each structure stores a set of probes with some metadata */
186 static int  _n_probe_sets = 0;
187 static cs_probe_set_t  **_probe_set_array = NULL;
188 
189 static const char _err_empty_pset[]
190   = N_(" Stop execution since the given cs_probe_set_t structure is empty.\n"
191        " Please check your settings.\n");
192 static const char _err_truefalse_key[]
193   = N_(" Invalid value %s for setting key %s\n"
194        " Valid choices are true or false.\n"
195        " Please modify your setting.\n");
196 
197 /*============================================================================
198  * Private function definitions
199  *============================================================================*/
200 
201 /*----------------------------------------------------------------------------*/
202 /*!
203  * \brief  Print the name of the corresponding key
204  *
205  * \param[in] key        name of the key
206  *
207  * \return a string
208  */
209 /*----------------------------------------------------------------------------*/
210 
211 static const char *
_print_psetkey(psetkey_t key)212 _print_psetkey(psetkey_t  key)
213 {
214   switch (key) {
215 
216   case PSETKEY_TRANSIENT_LOC:
217     return "transient_location";
218   case PSETKEY_BOUNDARY:
219     return "boundary";
220   case PSETKEY_SELECT_CRIT:
221     return "selection_criteria";
222   case PSETKEY_TOLERANCE:
223     return "tolerance";
224   case PSETKEY_INTERPOLATION:
225     return "interpolation";
226   default:
227     assert(0);
228   }
229 
230   return NULL; // avoid a warning
231 }
232 
233 /*----------------------------------------------------------------------------*/
234 /*!
235  * \brief  Get the corresponding enum from the name of a key.
236  *         If not found, return a key error.
237  *
238  * \param[in] keyname    name of the key
239  *
240  * \return a psetkey_t
241  */
242 /*----------------------------------------------------------------------------*/
243 
244 static psetkey_t
_get_psetkey(const char * keyname)245 _get_psetkey(const char  *keyname)
246 {
247   psetkey_t  key = PSETKEY_ERROR;
248 
249   if (strcmp(keyname, "transient_location") == 0)
250     key = PSETKEY_TRANSIENT_LOC;
251   else if (strcmp(keyname, "boundary") == 0)
252     key = PSETKEY_BOUNDARY;
253   else if (strcmp(keyname, "selection_criteria") == 0)
254     key = PSETKEY_SELECT_CRIT;
255   else if (strcmp(keyname, "tolerance") == 0)
256     key = PSETKEY_TOLERANCE;
257   else if (strcmp(keyname, "interpolation") == 0)
258     key = PSETKEY_INTERPOLATION;
259 
260   return key;
261 }
262 
263 /*----------------------------------------------------------------------------*/
264 /*!
265  * \brief  Copy a label to a new string
266  *
267  * \param[in, out]  label    label to set
268  *
269  * \return   copy of label, or NULL
270  */
271 /*----------------------------------------------------------------------------*/
272 
273 inline static char *
_copy_label(const char * name)274 _copy_label(const char  *name)
275 {
276   char *label = NULL;
277 
278   if (name) {
279     size_t  len = strlen(name) + 1;
280     BFT_MALLOC(label, len, char);
281     strcpy(label, name);
282   }
283 
284   return label;
285 }
286 
287 /*----------------------------------------------------------------------------*/
288 /*!
289  * \brief  Free a cs_probe_set_t structure
290  *
291  * \param[in, out]  pset          pointer to a cs_probe_set_t structure
292  */
293 /*----------------------------------------------------------------------------*/
294 
295 static void
_probe_set_free(cs_probe_set_t * pset)296 _probe_set_free(cs_probe_set_t   *pset)
297 {
298   if (pset == NULL)
299     return;
300 
301   BFT_FREE(pset->name);
302   BFT_FREE(pset->coords_ini);
303   BFT_FREE(pset->coords);
304   BFT_FREE(pset->sel_criter);
305   BFT_FREE(pset->loc_id);
306   BFT_FREE(pset->elt_id);
307   BFT_FREE(pset->vtx_id);
308   BFT_FREE(pset->located);
309 
310   BFT_FREE(pset->_p_define_input);
311 
312   if (pset->labels != NULL) {
313     for (int i = 0; i < pset->n_probes; i++)
314       BFT_FREE(pset->labels[i]);
315     BFT_FREE(pset->labels);
316   }
317 
318   if (pset->s_coords != NULL)
319     BFT_FREE(pset->s_coords);
320 
321   if (pset->n_writers > 0)
322     BFT_FREE(pset->writer_ids);
323 
324   if (pset->n_fields > 0)
325     BFT_FREE(pset->field_info);
326 }
327 
328 /*----------------------------------------------------------------------------*/
329 /*!
330  * \brief  Allocate and initialize by default a new set of probes
331  *
332  * \param[in]  name          name of the probe set
333  * \param[in]  n_max_probes  maximal number of probes
334  *
335  * \return a pointer to a cs_probe_set_t structure
336  */
337 /*----------------------------------------------------------------------------*/
338 
339 static cs_probe_set_t *
_probe_set_create(const char * name,cs_lnum_t n_max_probes)340 _probe_set_create(const char    *name,
341                   cs_lnum_t      n_max_probes)
342 {
343   cs_probe_set_t *pset = cs_probe_set_get(name);
344 
345   /* Already defined */
346 
347   if (pset != NULL) {
348     if (pset->flags & CS_PROBE_OVERWRITE)
349       _probe_set_free(pset);
350     else
351       bft_error(__FILE__, __LINE__, 0,
352                 _(" Error adding a new set of probes.\n"
353                   " %s is already used as a name for a set of probes.\n"
354                   " Please check your settings."), name);
355   }
356 
357   /* Add a new set of probes */
358 
359   else {
360     int pset_id = _n_probe_sets;
361 
362     _n_probe_sets++;
363     BFT_REALLOC(_probe_set_array, _n_probe_sets, cs_probe_set_t *);
364     BFT_MALLOC(pset, 1, cs_probe_set_t);
365     _probe_set_array[pset_id] = pset;
366   }
367 
368   /* Copy name */
369   int  len = strlen(name) + 1;
370   BFT_MALLOC(pset->name, len, char);
371   strncpy(pset->name, name, len);
372 
373   pset->flags = CS_PROBE_AUTO_VAR;
374   pset->tolerance = 0.1;
375   pset->sel_criter = NULL;
376 
377   pset->snap_mode = CS_PROBE_SNAP_NONE;
378 
379   pset->n_max_probes = n_max_probes;
380   pset->n_probes = 0;
381   pset->n_loc_probes = 0;
382 
383   pset->coords_ini = NULL;
384 
385   BFT_MALLOC(pset->coords, n_max_probes, cs_real_3_t);
386   pset->s_coords = NULL;
387   pset->labels = NULL;
388 
389   pset->p_define_func = NULL;
390   pset->p_define_input = NULL;
391   pset->_p_define_input = NULL;
392 
393   pset->loc_id = NULL;
394   pset->elt_id = NULL;
395   pset->vtx_id = NULL;
396   pset->located = NULL;
397 
398   pset->interpolation = 0;
399 
400   pset->n_writers = -1;
401   pset->writer_ids = NULL;
402 
403   pset->n_fields = 0;
404   pset->n_max_fields = 0;
405   pset->field_info = NULL;
406 
407   return pset;
408 }
409 
410 /*----------------------------------------------------------------------------*/
411 /*!
412  * \brief  Check if the given name has the same name as the given set of probes
413  *
414  * \param[in]  pset        pointer to a cs_probe_set_t structure to test
415  * \param[in]  ref_name    name to check
416  *
417  * \return true if the name of the set is ref_name otherwise false
418  */
419 /*----------------------------------------------------------------------------*/
420 
421 static bool
_check_probe_set_name(const cs_probe_set_t * pset,const char * ref_name)422 _check_probe_set_name(const cs_probe_set_t   *pset,
423                       const char             *ref_name)
424 {
425   if (pset == NULL)
426     return false;
427 
428   int  reflen = strlen(ref_name);
429   int  len = strlen(pset->name);
430 
431   if (reflen == len) {
432     if (strcmp(ref_name, pset->name) == 0)
433       return true;
434     else
435       return false;
436   }
437   else
438     return false;
439 }
440 
441 /*----------------------------------------------------------------------------*/
442 /*!
443  * \brief Build probes based on a function-definition
444  *
445  * Note: if the p_define_input pointer is non-NULL, it must point to valid data
446  * when the selection function is called, so that value or structure should
447  * not be temporary (i.e. local);
448  *
449  * \param[in, out]  pset  pointer to a cs_probe_set_t structure to update
450  */
451 /*----------------------------------------------------------------------------*/
452 
453 static void
_build_local_probe_set(cs_probe_set_t * pset)454 _build_local_probe_set(cs_probe_set_t  *pset)
455 {
456   assert(pset->p_define_func != NULL);
457 
458   pset->n_max_probes = 0;
459   pset->n_probes = 0;
460   pset->n_loc_probes = 0;
461 
462   BFT_FREE(pset->coords);
463   BFT_FREE(pset->s_coords);
464 
465   cs_lnum_t     n_elts = 0;
466   cs_real_3_t  *coords = NULL;
467   cs_real_t    *s = NULL;
468 
469   pset->p_define_func(pset->p_define_input,
470                       &n_elts,
471                       &coords,
472                       &s);
473 
474   pset->n_probes = n_elts;
475 
476   /* Order by curvilinear coordinates to avoid issues with some
477      plot formats/tools in single-rank mode (not needed in parallel,
478      as this is already handled through IO numbering in that case) */
479 
480   if (cs_glob_n_ranks <= 1) {
481     BFT_MALLOC(pset->coords, n_elts, cs_real_3_t);
482     BFT_MALLOC(pset->s_coords, n_elts, cs_real_t);
483 
484     cs_lnum_t *order = NULL;
485     BFT_MALLOC(order, n_elts, cs_lnum_t);
486     cs_order_real_allocated(NULL, s, order, n_elts);
487 
488     for (cs_lnum_t i = 0; i < n_elts; i++) {
489       cs_lnum_t j = order[i];
490       for (cs_lnum_t k = 0; k < 3; k++)
491         pset->coords[i][k] = coords[j][k];
492       pset->s_coords[i] = s[j];
493     }
494 
495     BFT_FREE(order);
496 
497     BFT_FREE(coords);
498     BFT_FREE(s);
499   }
500 
501   else {
502     pset->coords = coords;
503     pset->s_coords = s;
504   }
505 }
506 
507 /*----------------------------------------------------------------------------*/
508 /*!
509  * \brief  Merge identical probes when snapped to center.
510  *
511  * \param[in, out]  pset    pointer to a cs_probe_set_t structure
512  * \param[in]       center  cell centers
513  */
514 /*----------------------------------------------------------------------------*/
515 
516 #if defined(__INTEL_COMPILER)
517 #if __INTEL_COMPILER < 1800
518 #pragma optimization_level 1 /* Bug with O2 or above with icc 17.0.0 20160721 */
519 #endif
520 #endif
521 
522 static void
_merge_snapped_to_center(cs_probe_set_t * pset,const cs_real_t centers[])523 _merge_snapped_to_center(cs_probe_set_t   *pset,
524                          const cs_real_t   centers[])
525 {
526   if (pset == NULL)
527     return;
528 
529   int *tag;
530   BFT_MALLOC(tag, pset->n_probes, int);
531 
532   for (int i = 0; i < pset->n_probes; i++)
533     tag[i] = 0;
534 
535   cs_lnum_t *order;
536   BFT_MALLOC(order, pset->n_loc_probes, cs_lnum_t);
537   cs_order_lnum_allocated(NULL, pset->elt_id, order, pset->n_loc_probes);
538 
539   cs_lnum_t e_id = 0;
540 
541   while (e_id < pset->n_loc_probes) {
542 
543     cs_lnum_t s_id = e_id;
544     cs_lnum_t j_min = order[s_id];
545     int l = pset->elt_id[j_min];
546     int k = pset->loc_id[j_min];
547     cs_real_t d_min = cs_math_3_distance(pset->coords[k], centers + l*3);
548     tag[pset->loc_id[j_min]] = 1;
549 
550     for (e_id = s_id+1; e_id < pset->n_loc_probes; e_id++) {
551       cs_lnum_t j = order[e_id];
552       if (pset->elt_id[j] != l)
553         break;
554       else {
555         cs_real_t d = cs_math_3_distance(pset->coords[k], centers + l*3);
556         if (d < d_min) {
557           tag[pset->loc_id[j_min]] = 0;
558           tag[pset->loc_id[j]] = 1;
559           j_min = j;
560           d_min = d;
561         }
562       }
563     }
564 
565   }
566 
567   BFT_FREE(order);
568 
569   /* Now all points to keep are tagged */
570 
571   cs_parall_max(pset->n_probes, CS_INT_TYPE, tag);
572 
573   /* Discard untagged points
574      ----------------------- */
575 
576   /* Update global arrays; transfer tag to renumbering array */
577 
578   int n_probes = 0;
579 
580   for (int i = 0; i < pset->n_probes; i++) {
581     if (tag[i] > 0) {
582       tag[i] = n_probes;
583       for (int j = 0; j < 3; j++)
584         pset->coords[n_probes][j] = pset->coords[i][j];
585       if (pset->s_coords != NULL)
586         pset->s_coords[n_probes] = pset->s_coords[i];
587       if (pset->labels != NULL)
588         pset->labels[n_probes] = pset->labels[i];
589       pset->located[n_probes] = pset->located[i];
590       n_probes++;
591     }
592     else {
593       if (pset->labels != NULL)
594         BFT_FREE(pset->labels[i]);
595       tag[i] = -1;
596     }
597   }
598   pset->n_probes = n_probes;
599 
600   BFT_REALLOC(pset->coords, n_probes, cs_real_3_t);
601   if (pset->s_coords != NULL)
602     BFT_REALLOC(pset->s_coords, n_probes, cs_real_t);
603   if (pset->labels != NULL)
604     BFT_REALLOC(pset->labels, n_probes, char *);
605   BFT_REALLOC(pset->located, n_probes, char);
606 
607   /* Update local arrays */
608 
609   int n_loc_probes = 0;
610   for (int i = 0; i < pset->n_loc_probes; i++) {
611     int j = pset->loc_id[i];
612     if (tag[j] > -1) {
613       pset->loc_id[n_loc_probes] = tag[j];
614       pset->elt_id[n_loc_probes] = pset->elt_id[i];
615       pset->vtx_id[n_loc_probes] = pset->vtx_id[i];
616       n_loc_probes += 1;
617     }
618   }
619   pset->n_loc_probes = n_loc_probes;
620 
621   BFT_REALLOC(pset->loc_id, n_loc_probes, cs_lnum_t);
622   BFT_REALLOC(pset->elt_id, n_loc_probes, cs_lnum_t);
623   BFT_REALLOC(pset->vtx_id, n_loc_probes, cs_lnum_t);
624 
625   BFT_FREE(tag);
626 }
627 
628 /*----------------------------------------------------------------------------*/
629 /*!
630  * \brief Define a profile based on centers of cells cut by a line segment
631  *
632  * In this case, the input points to a real array containing the segment's
633  * start and end coordinates.
634  *
635  * Note: the input pointer must point to valid data when this selection
636  * function is called, so either:
637  * - that value or structure should not be temporary (i.e. local);
638  * - post-processing output must be ensured using cs_post_write_meshes()
639  *   with a fixed-mesh writer before the data pointed to goes out of scope;
640  *
641  * \param[in]   input   pointer to segment start and end:
642  *                      [x0, y0, z0, x1, y1, z1]
643  * \param[out]  n_elts  number of selected coordinates
644  * \param[out]  coords  coordinates of selected elements.
645  * \param[out]  s       curvilinear coordinates of selected elements
646  *----------------------------------------------------------------------------*/
647 
648 static void
_cell_profile_probes_define(void * input,cs_lnum_t * n_elts,cs_real_3_t ** coords,cs_real_t ** s)649 _cell_profile_probes_define(void          *input,
650                             cs_lnum_t     *n_elts,
651                             cs_real_3_t  **coords,
652                             cs_real_t    **s)
653 {
654   /* Determine segment endpoints from input  */
655 
656   cs_real_t *sx = (cs_real_t *)input;
657 
658   const cs_real_t dx1[3] = {sx[3]-sx[0], sx[4]-sx[1], sx[5]-sx[2]};
659   const cs_real_t s_norm2 = cs_math_3_square_norm(dx1);
660 
661   const cs_real_3_t  *cell_cen
662     = (const cs_real_3_t *)(cs_glob_mesh_quantities->cell_cen);
663 
664   cs_lnum_t n_cells = 0;
665   cs_lnum_t *cell_ids = NULL;
666   cs_real_t *seg_c_len = NULL;
667 
668   /* Compared to cs_cell_segment_intersect_select, this function also
669      selects a cell if the segment is included in this cell */
670 
671   cs_mesh_intersect_polyline_cell_select(sx,
672                                          2, &n_cells, &cell_ids, &seg_c_len);
673 
674   cs_real_3_t *_coords;
675   cs_real_t *_s;
676   BFT_MALLOC(_coords, n_cells, cs_real_3_t);
677   BFT_MALLOC(_s, n_cells, cs_real_t);
678 
679   for (cs_lnum_t i = 0; i < n_cells; i++) {
680     cs_real_t dx[3], coo[3];
681     for (cs_lnum_t j = 0; j < 3; j++) {
682       coo[j] = cell_cen[cell_ids[i]][j];
683       dx[j] = coo[j] - sx[j];
684       _coords[i][j] = coo[j];
685     }
686     _s[i] = cs_math_3_dot_product(dx, dx1) / s_norm2;
687   }
688 
689   BFT_FREE(cell_ids);
690   BFT_FREE(seg_c_len);
691 
692   /* Set return values */
693 
694   *n_elts = n_cells;
695   *coords = _coords;
696   *s = _s;
697 }
698 
699 /*============================================================================
700  * Semi-private function definitions
701  *
702  * The following functions are intended to be used by the postprocessing layer
703  * (cs_post.c), not directly by the user.
704  *============================================================================*/
705 
706 /*----------------------------------------------------------------------------*/
707 /*!
708  * \brief  Free all structures related to a set of probes.
709  */
710 /*----------------------------------------------------------------------------*/
711 
712 void
cs_probe_finalize(void)713 cs_probe_finalize(void)
714 {
715   for (int i = 0; i < _n_probe_sets; i++) {
716     cs_probe_set_t *pset = _probe_set_array[i];
717     _probe_set_free(pset);
718     BFT_FREE(pset);
719   }
720 
721   _n_probe_sets = 0;
722   BFT_FREE(_probe_set_array);
723 }
724 
725 /*----------------------------------------------------------------------------*/
726 /*!
727  * \brief Transfer info on associated fields to the caller.
728  *
729  * This function transfers the property of the associated arrays to the caller
730  * and removes it from the probe set info.
731  *
732  * \param[in]   pset        pointer to a cs_probe_set_t structure
733  * \param[out]  n_fields    number of associated fields
734  * \param[out]  field_info  associated field info (array of tuples with
735  *                          writer id, field id, component id)
736  */
737 /*----------------------------------------------------------------------------*/
738 
739 void
cs_probe_set_transfer_associated_field_info(cs_probe_set_t * pset,int * n_fields,int ** field_info)740 cs_probe_set_transfer_associated_field_info(cs_probe_set_t   *pset,
741                                             int              *n_fields,
742                                             int             **field_info)
743 {
744   if (pset == NULL)
745     bft_error(__FILE__, __LINE__, 0, _(_err_empty_pset));
746 
747   *n_fields = pset->n_fields;
748   *field_info = pset->field_info;
749 
750   pset->n_fields = 0;
751   pset->n_max_fields = 0;
752   pset->field_info = NULL;
753 }
754 
755 /*! (DOXYGEN_SHOULD_SKIP_THIS) \endcond */
756 
757 /*============================================================================
758  * Public function definitions
759  *============================================================================*/
760 
761 /*----------------------------------------------------------------------------*/
762 /*!
763  * \brief  Retrieve the number of probe sets defined.
764  *
765  * \return the number of probe sets defined
766  */
767 /*----------------------------------------------------------------------------*/
768 
769 int
cs_probe_get_n_sets(void)770 cs_probe_get_n_sets(void)
771 {
772   return _n_probe_sets;
773 }
774 
775 /*----------------------------------------------------------------------------*/
776 /*!
777  * \brief  Retrieve a cs_probe_set_t structure.
778  *
779  * \param[in]   name        name of the set of probes to find
780  *
781  * \return a pointer to a cs_probes_t structure or NULL if not found
782  */
783 /*----------------------------------------------------------------------------*/
784 
785 cs_probe_set_t *
cs_probe_set_get(const char * name)786 cs_probe_set_get(const char  *name)
787 {
788   cs_probe_set_t  *pset = NULL;
789 
790   if (name == NULL)
791     bft_error(__FILE__, __LINE__, 0,
792               _(" The given name for this set of probes is empty."));
793 
794   /* Check if the given name is already used */
795   for (int pset_id = 0; pset_id < _n_probe_sets; pset_id++) {
796     if (_check_probe_set_name(_probe_set_array[pset_id], name)) {
797       pset = _probe_set_array[pset_id];
798       break;
799     }
800   }
801 
802   return pset;
803 }
804 
805 /*----------------------------------------------------------------------------*/
806 /*!
807  * \brief  Retrieve a cs_probe_set_t structure from its id.
808  *
809  * \param[in]   pset_id       id related to the set of probes to find
810  *
811  * \return a pointer to a cs_probes_t structure or NULL if not found
812  */
813 /*----------------------------------------------------------------------------*/
814 
815 cs_probe_set_t *
cs_probe_set_get_by_id(int pset_id)816 cs_probe_set_get_by_id(int   pset_id)
817 {
818   /* Check if the given name is already used */
819   if (pset_id < 0 || pset_id >= _n_probe_sets)
820     return  NULL;
821 
822   return _probe_set_array[pset_id];
823 }
824 
825 /*----------------------------------------------------------------------------*/
826 /*!
827  * \brief  Retrieve the name related to a cs_probe_set_t structure.
828  *
829  * \param[in]   pset       pointer to a cs_probe_set_t structure
830  *
831  * \return the name of the cs_probes_t structure or NULL if not found
832  */
833 /*----------------------------------------------------------------------------*/
834 
835 const char *
cs_probe_set_get_name(cs_probe_set_t * pset)836 cs_probe_set_get_name(cs_probe_set_t   *pset)
837 {
838   if (pset == NULL)
839     return  NULL;
840 
841   return pset->name;
842 }
843 
844 /*----------------------------------------------------------------------------*/
845 /*!
846  * \brief  Retrieve information useful for the postprocessing step.
847  *
848  * Output arguments may be set to NULL if we do not need to query them.
849  *
850  * \param[in]  pset            pointer to a cs_probe_set_t structure
851  * \param[out] time_varying    true if probe locations may change with time
852  * \param[out] on_boundary     true if probes are located on boundary
853  * \param[out] on_curve        true if the probe set has cuvilinear coordinates
854  * \param[out] auto_variables  true if set of variables to output is predefined
855  * \param[out] auto_curve_coo  true if curvilinear coordinates should be output
856  * \param[out] auto_cart_coo   true if cartesian coordinates should be output
857  * \param[out] n_writers       number of associated  user-defined writers,
858  *                             or -1 if default unchanged
859  * \param[out] writer_ids      pointer to a list of writer ids
860  */
861 /*----------------------------------------------------------------------------*/
862 
863 void
cs_probe_set_get_post_info(const cs_probe_set_t * pset,bool * time_varying,bool * on_boundary,bool * on_curve,bool * auto_variables,bool * auto_curve_coo,bool * auto_cart_coo,int * n_writers,int * writer_ids[])864 cs_probe_set_get_post_info(const cs_probe_set_t   *pset,
865                            bool                   *time_varying,
866                            bool                   *on_boundary,
867                            bool                   *on_curve,
868                            bool                   *auto_variables,
869                            bool                   *auto_curve_coo,
870                            bool                   *auto_cart_coo,
871                            int                    *n_writers,
872                            int                    *writer_ids[])
873 {
874   if (pset == NULL)
875     bft_error(__FILE__, __LINE__, 0, _err_empty_pset);
876 
877   if (time_varying != NULL)
878     *time_varying = (pset->flags & CS_PROBE_TRANSIENT) ? true : false;
879   if (auto_variables != NULL)
880     *auto_variables = (pset->flags & CS_PROBE_AUTO_VAR) ? true: false;
881   if (auto_curve_coo != NULL)
882     *auto_curve_coo = (pset->flags & CS_PROBE_AUTO_S) ? true: false;
883   if (auto_cart_coo != NULL)
884     *auto_cart_coo = (pset->flags & CS_PROBE_AUTO_COORD) ? true: false;
885   if (on_curve != NULL)
886     *on_curve = (pset->flags & CS_PROBE_ON_CURVE) ? true : false;
887   if (on_boundary != NULL)
888     *on_boundary = (pset->flags & CS_PROBE_BOUNDARY) ? true : false;
889 
890   if (n_writers != NULL)
891     *n_writers = pset->n_writers;
892   if (writer_ids != NULL)
893     *writer_ids = pset->writer_ids;
894 }
895 
896 /*----------------------------------------------------------------------------*/
897 /*!
898  * \brief  Return the location filter selection criteria string for a
899  *         given probe set
900  *
901  * \param[in]   pset       pointer to a cs_probe_set_t structure
902  *
903  * \return selection criteria string, or NULL if no filter defined
904  */
905 /*----------------------------------------------------------------------------*/
906 
907 const char *
cs_probe_set_get_location_criteria(cs_probe_set_t * pset)908 cs_probe_set_get_location_criteria(cs_probe_set_t   *pset)
909 {
910   if (pset == NULL)
911     return  NULL;
912 
913   return pset->sel_criter;
914 }
915 
916 /*----------------------------------------------------------------------------*/
917 /*!
918  * \brief  Return the interpolation option for a given probe set
919  *
920  * Interpolation will be applied only where possible.
921  *
922  * \param[in]   pset       pointer to a cs_probe_set_t structure
923  *
924  * \return interpolation option value
925  */
926 /*----------------------------------------------------------------------------*/
927 
928 int
cs_probe_set_get_interpolation(cs_probe_set_t * pset)929 cs_probe_set_get_interpolation(cs_probe_set_t   *pset) {
930   if (pset == NULL)
931     return  0;
932 
933   return pset->interpolation;
934 }
935 
936 /*----------------------------------------------------------------------------*/
937 /*!
938  * \brief  Create a new set of probes.
939  *
940  * \param[in]   name   name of the set of probes
941  *
942  * \return a pointer to a new allocated cs_probe_set_t structure
943  */
944 /*----------------------------------------------------------------------------*/
945 
946 cs_probe_set_t *
cs_probe_set_create(const char * name)947 cs_probe_set_create(const char  *name)
948 {
949   /* Default initialization of a set of probes (max number is set to 4 by
950      default but a realloc is available if the max. number of probes is
951      reached) */
952 
953   cs_probe_set_t *pset = _probe_set_create(name, 4);
954 
955   return  pset;
956 }
957 
958 /*----------------------------------------------------------------------------*/
959 /*!
960  * \brief  Add a new probe to an existing set of probes.
961  *
962  * \param[in, out]  pset    set of probes
963  * \param[in]       x       x coordinate  of the point to add
964  * \param[in]       y       y coordinate  of the point to add
965  * \param[in]       z       z coordinate  of the point to add
966  * \param[in]       label   NULL or the name of the point (optional)
967  */
968 /*----------------------------------------------------------------------------*/
969 
970 void
cs_probe_set_add_probe(cs_probe_set_t * pset,cs_real_t x,cs_real_t y,cs_real_t z,const char * label)971 cs_probe_set_add_probe(cs_probe_set_t     *pset,
972                        cs_real_t           x,
973                        cs_real_t           y,
974                        cs_real_t           z,
975                        const char         *label)
976 {
977   if (pset == NULL)
978     bft_error(__FILE__, __LINE__, 0, _err_empty_pset);
979 
980   cs_lnum_t  point_id = pset->n_probes;
981 
982   pset->n_probes++;
983 
984   if (point_id >= pset->n_max_probes) { /* Reallocate arrays */
985     pset->n_max_probes *= 2;
986     BFT_REALLOC(pset->coords, pset->n_max_probes, cs_real_3_t);
987     if (pset->labels != NULL)
988       BFT_REALLOC(pset->labels, pset->n_max_probes, char *);
989   }
990 
991   /* Set coordinates */
992   pset->coords[point_id][0] = x;
993   pset->coords[point_id][1] = y;
994   pset->coords[point_id][2] = z;
995 
996   if (label != NULL) { /* Manage the label */
997     if (pset->labels == NULL)
998       BFT_MALLOC(pset->labels, pset->n_max_probes, char *);
999 
1000     /* Copy label */
1001     pset->labels[point_id] = _copy_label(label);
1002   }
1003 }
1004 
1005 /*----------------------------------------------------------------------------*/
1006 /*!
1007  * \brief  Define a new set of probes from an array of coordinates.
1008  *
1009  * \param[in]   name      name of the set of probes
1010  * \param[in]   n_probes  number of probes in coords and labels
1011  * \param[in]   coords    list of coordinates related to each probe
1012  * \param[in]   labels    list of label related to each probe (optional)
1013  *
1014  * \return a pointer to a new allocated cs_probe_set_t structure
1015  */
1016 /*----------------------------------------------------------------------------*/
1017 
1018 cs_probe_set_t *
cs_probe_set_create_from_array(const char * name,int n_probes,const cs_real_3_t * coords,const char ** labels)1019 cs_probe_set_create_from_array(const char          *name,
1020                                int                  n_probes,
1021                                const cs_real_3_t   *coords,
1022                                const char         **labels)
1023 {
1024   cs_probe_set_t  *pset = _probe_set_create(name, n_probes);
1025 
1026   pset->n_probes = n_probes;
1027 
1028   /* Coordinates */
1029   for (int i = 0; i < n_probes; i++) {
1030     pset->coords[i][0] = coords[i][0];
1031     pset->coords[i][1] = coords[i][1];
1032     pset->coords[i][2] = coords[i][2];
1033   }
1034 
1035   /* Copy labels */
1036   if (labels != NULL) {
1037     BFT_MALLOC(pset->labels, n_probes, char *);
1038     for (int i = 0; i < n_probes; i++)
1039       pset->labels[i] = _copy_label(labels[i]);
1040   }
1041 
1042   return  pset;
1043 }
1044 
1045 /*----------------------------------------------------------------------------*/
1046 /*!
1047  * \brief Define a new set of probes from the segment spanned by two points.
1048  *
1049  * If n_probes > 0, the given number of probes will be used for sampling, and
1050  * will be evenly distributed along the segment.
1051  *
1052  * If n_probes <= 0, one probe will be defined for each cell intersected by
1053  * the line segment, and its position based on the projection of those cell
1054  * centers on the given segment.
1055  *
1056  * In both cases, using cs_probe_set_snap_mode can further modify this behavior.
1057  *
1058  * \param[in]  name          name of the set of probes
1059  * \param[in]  n_probes      number of probes, or <= 0 for mesh intersection
1060  * \param[in]  start_coords  coordinates of the starting point
1061  * \param[in]  end_coords    coordinates of the ending point
1062  *
1063  * \return a pointer to a new allocated cs_probe_set_t structure
1064  */
1065 /*----------------------------------------------------------------------------*/
1066 
1067 cs_probe_set_t *
cs_probe_set_create_from_segment(const char * name,int n_probes,const cs_real_t start_coords[3],const cs_real_t end_coords[3])1068 cs_probe_set_create_from_segment(const char        *name,
1069                                  int                n_probes,
1070                                  const cs_real_t    start_coords[3],
1071                                  const cs_real_t    end_coords[3])
1072 {
1073   cs_probe_set_t  *pset = NULL;
1074 
1075   if (n_probes > 0) {
1076 
1077     pset = _probe_set_create(name, n_probes);
1078 
1079     pset->n_probes = n_probes;
1080     pset->flags |= CS_PROBE_ON_CURVE;
1081     if (pset->flags & CS_PROBE_AUTO_VAR)
1082       pset->flags -= CS_PROBE_AUTO_VAR;
1083 
1084     BFT_MALLOC(pset->s_coords, n_probes, cs_real_t);
1085 
1086     /* 2 probes are already defined (the starting and ending points)
1087        Define the additional probes */
1088     cs_real_t  distance;
1089     cs_real_3_t  unitv, delta_vect;
1090 
1091     cs_math_3_length_unitv(start_coords, end_coords, &distance, unitv);
1092 
1093     const double  delta = distance / (n_probes - 1);
1094     for (int k = 0; k < 3; k++)
1095       delta_vect[k] = delta*unitv[k];
1096 
1097     /* Set the starting probe */
1098     pset->s_coords[0] = 0;
1099     for (int k = 0; k < 3; k++)
1100       pset->coords[0][k] = start_coords[k];
1101 
1102     /* Set additional probes */
1103     for (int i = 1; i < n_probes - 1; i++) {
1104 
1105       pset->s_coords[i] = pset->s_coords[i-1] + delta;
1106       for (int k = 0; k < 3; k++)
1107         pset->coords[i][k] = pset->coords[i-1][k] + delta_vect[k];
1108 
1109     }
1110 
1111     /* Set the ending probe */
1112     pset->s_coords[n_probes-1] = distance;
1113     for (int k = 0; k < 3; k++)
1114       pset->coords[n_probes-1][k] = end_coords[k];
1115   }
1116 
1117   else {
1118     cs_real_t *se_coords;
1119     BFT_MALLOC(se_coords, 6, cs_real_t);
1120 
1121     for (int i = 0; i < 3; i++) {
1122       se_coords[i] = start_coords[i];
1123       se_coords[i+3] = end_coords[i];
1124     }
1125 
1126     pset = cs_probe_set_create_from_local(name,
1127                                           _cell_profile_probes_define,
1128                                           (void *)se_coords);  /* input */
1129 
1130     pset->_p_define_input = se_coords;
1131   }
1132 
1133   return  pset;
1134 }
1135 
1136 /*----------------------------------------------------------------------------*/
1137 /*!
1138  * \brief Define a new set of probes from rank-local definition function.
1139  *
1140  * The local definition function given by the p_define_func pointer
1141  * is called just before locating probes on the parent mesh, so this allows
1142  * building probe sets based on subsets of the computational mesh.
1143  *
1144  * Note: if the p_define_input pointer is non-NULL, it must point to valid data
1145  * when the selection function is called, so that value or structure should
1146  * not be temporary (i.e. local);
1147  *
1148  * \param[in]  name           name of the set of probes
1149  * \param[in]  p_define_func  function used for local definition
1150  * \param[in]  p_define_input optional input for local definition function
1151  *
1152  * \return a pointer to a new allocated cs_probe_set_t structure
1153  */
1154 /*----------------------------------------------------------------------------*/
1155 
1156 cs_probe_set_t *
cs_probe_set_create_from_local(const char * name,cs_probe_set_define_local_t * p_define_func,void * p_define_input)1157 cs_probe_set_create_from_local(const char                   *name,
1158                                cs_probe_set_define_local_t  *p_define_func,
1159                                void                         *p_define_input)
1160 {
1161   cs_probe_set_t  *pset = _probe_set_create(name, 0);
1162 
1163   pset->flags |= CS_PROBE_ON_CURVE;
1164   if (pset->flags & CS_PROBE_AUTO_VAR)
1165     pset->flags -= CS_PROBE_AUTO_VAR;
1166 
1167   pset->p_define_func = p_define_func;
1168   pset->p_define_input = p_define_input;
1169 
1170   return  pset;
1171 }
1172 
1173 /*----------------------------------------------------------------------------*/
1174 /*!
1175  * \brief  allow overwriting the definition of a given probe set.
1176  *
1177  * If no a probe set of the given name exists, the operation is ignored.
1178  *
1179  * \param[in]  name  name of the probe set
1180  */
1181 /*----------------------------------------------------------------------------*/
1182 
1183 void
cs_probe_set_allow_overwrite(const char * name)1184 cs_probe_set_allow_overwrite(const char  *name)
1185 {
1186   cs_probe_set_t *pset = cs_probe_set_get(name);
1187 
1188   if (pset != NULL)
1189     pset->flags = pset->flags | CS_PROBE_OVERWRITE;
1190 }
1191 
1192 /*----------------------------------------------------------------------------*/
1193 /*!
1194  * \brief  Assign curvilinear abscissa for the given probe set
1195  *
1196  * This is effective only when using probe sets to which curvilinear
1197  * coordinates have not already been assigned.
1198  *
1199  * So this may be used following \ref cs_probe_set_create (and probe additions)
1200  * or \ref cs_probe_set_create_from_array, but will be ignored for probe sets
1201  * defined using \ref cs_probe_set_create_from_segment and
1202  * \ref cs_probe_set_create_from_local.
1203  *
1204  * If provided, the array of curvilinear absicssa must match the size of the
1205  * probe set. If set to NULL, it is assumed a uniform spacing is provided
1206  * between points.
1207  *
1208  * \param[in]  pset  pointer to a cs_probe_set_t structure
1209  * \param[in]  s     array of curvilinear abscissa, or NULL
1210  */
1211 /*----------------------------------------------------------------------------*/
1212 
1213 void
cs_probe_set_assign_curvilinear_abscissa(cs_probe_set_t * pset,const cs_real_t * s)1214 cs_probe_set_assign_curvilinear_abscissa(cs_probe_set_t   *pset,
1215                                          const cs_real_t  *s)
1216 {
1217   if (pset == NULL)
1218     return;
1219 
1220   if (pset->flags & CS_PROBE_ON_CURVE)
1221     return;
1222 
1223   pset->flags |= CS_PROBE_ON_CURVE;
1224 
1225   BFT_REALLOC(pset->s_coords, pset->n_probes, cs_real_t);
1226 
1227   if (s != NULL) {
1228     for (int i = 0; i < pset->n_probes; i++)
1229       pset->s_coords[i] = s[i];
1230   }
1231   else if (pset->n_probes > 0) {
1232     pset->s_coords[0] = 0.;
1233     if (pset->n_probes > 1) {
1234       int npm1 = pset->n_probes - 1;
1235       cs_real_t a = 1./(cs_real_t)npm1;
1236       for (int i = 1; i < npm1; i++)
1237         pset->s_coords[i] = i*a;;
1238       pset->s_coords[npm1] = 1.;
1239     }
1240   }
1241 }
1242 
1243 /*----------------------------------------------------------------------------*/
1244 /*!
1245  * \brief  Associate a list of writers to a probe set.
1246  *
1247  * \param[in, out] pset        pointer to a cs_probe_set_t structure to set
1248  * \param[in]      n_writers   number of writers assocuated to this probe set
1249  * \param[in]      writer_ids  list of writer ids
1250  */
1251 /*----------------------------------------------------------------------------*/
1252 
1253 void
cs_probe_set_associate_writers(cs_probe_set_t * pset,int n_writers,const int * writer_ids)1254 cs_probe_set_associate_writers(cs_probe_set_t   *pset,
1255                                int               n_writers,
1256                                const int        *writer_ids)
1257 {
1258   if (pset == NULL)
1259     bft_error(__FILE__, __LINE__, 0, _(_err_empty_pset));
1260 
1261   if (pset->n_writers < 0)
1262     pset->n_writers = 0;
1263 
1264   int  n_init_writers = pset->n_writers;
1265   pset->n_writers += n_writers;
1266   BFT_REALLOC(pset->writer_ids, pset->n_writers, int);
1267 
1268   for (int i = n_init_writers, j = 0; i < pset->n_writers; i++, j++)
1269     pset->writer_ids[i] = writer_ids[j];
1270 }
1271 
1272 /*----------------------------------------------------------------------------*/
1273 /*!
1274  * \brief Associate a field to a probe set.
1275  *
1276  * This function must be called during the postprocessing output definition
1277  * stage, before any output actually occurs.
1278  *
1279  * If the field should already be output automatically based on the mesh
1280  * category and field output keywords, it will be filtered at a later stage.
1281  *
1282  * \param[in]  pset       pointer to a cs_probe_set_t structure
1283  * \param[in]  writer_id  id of specified associated writer,
1284  *                        or 0 (\ref CS_POST_WRITER_ALL_ASSOCIATED) for all
1285  * \param[in]  field_id   id of field to attach
1286  * \param[in]  comp_id    id of field component (-1 for all)
1287  */
1288 /*----------------------------------------------------------------------------*/
1289 
1290 void
cs_probe_set_associate_field(cs_probe_set_t * pset,int writer_id,int field_id,int comp_id)1291 cs_probe_set_associate_field(cs_probe_set_t   *pset,
1292                              int               writer_id,
1293                              int               field_id,
1294                              int               comp_id)
1295 {
1296   if (pset == NULL)
1297     bft_error(__FILE__, __LINE__, 0, _(_err_empty_pset));
1298 
1299   if (pset->n_max_fields <= pset->n_fields) {
1300     if (pset->n_max_fields == 0)
1301       pset->n_max_fields = 8;
1302     else
1303       pset->n_max_fields *= 2;
1304     BFT_REALLOC(pset->field_info, 3*(pset->n_max_fields), int);
1305   }
1306 
1307   int *afi = pset->field_info + 3*pset->n_fields;
1308 
1309   afi[0] = writer_id;
1310   afi[1] = field_id;
1311   afi[2] = comp_id;
1312 
1313   pset->n_fields += 1;
1314 }
1315 
1316 /*----------------------------------------------------------------------------*/
1317 /*!
1318  * \brief  Set to true or false the automatic post-processing of variables
1319  *
1320  * \param[in, out] pset     pointer to a cs_probe_set_t structure
1321  * \param[in]      mode     true or false
1322  */
1323 /*----------------------------------------------------------------------------*/
1324 
1325 void
cs_probe_set_auto_var(cs_probe_set_t * pset,bool mode)1326 cs_probe_set_auto_var(cs_probe_set_t   *pset,
1327                       bool              mode)
1328 {
1329   if (pset == NULL)
1330     bft_error(__FILE__, __LINE__, 0, _(_err_empty_pset));
1331 
1332   if (mode == false) {
1333     if (pset->flags & CS_PROBE_AUTO_VAR)
1334       pset->flags -= CS_PROBE_AUTO_VAR;
1335   }
1336   else
1337     pset->flags |= CS_PROBE_AUTO_VAR;
1338 }
1339 
1340 /*----------------------------------------------------------------------------*/
1341 /*!
1342  * \brief  Set automatic output of curvilinear coordinates.
1343  *
1344  * \param[in, out] pset     pointer to a cs_probe_set_t structure
1345  * \param[in]      mode     true or false
1346  */
1347 /*----------------------------------------------------------------------------*/
1348 
1349 void
cs_probe_set_auto_curvilinear_coords(cs_probe_set_t * pset,bool mode)1350 cs_probe_set_auto_curvilinear_coords(cs_probe_set_t   *pset,
1351                                      bool              mode)
1352 {
1353   if (pset == NULL)
1354     bft_error(__FILE__, __LINE__, 0, _(_err_empty_pset));
1355 
1356   if (mode == false) {
1357     if (pset->flags & CS_PROBE_AUTO_S)
1358       pset->flags -= CS_PROBE_AUTO_S;
1359   }
1360   else
1361     pset->flags |= CS_PROBE_AUTO_S;
1362 }
1363 
1364 /*----------------------------------------------------------------------------*/
1365 /*!
1366  * \brief  Set automatic output of cartesian coordinates.
1367  *
1368  * \param[in, out] pset     pointer to a cs_probe_set_t structure
1369  * \param[in]      mode     true or false
1370  */
1371 /*----------------------------------------------------------------------------*/
1372 
1373 void
cs_probe_set_auto_cartesian_coords(cs_probe_set_t * pset,bool mode)1374 cs_probe_set_auto_cartesian_coords(cs_probe_set_t   *pset,
1375                                    bool              mode)
1376 {
1377   if (pset == NULL)
1378     bft_error(__FILE__, __LINE__, 0, _(_err_empty_pset));
1379 
1380   if (mode == false) {
1381     if (pset->flags & CS_PROBE_AUTO_COORD)
1382       pset->flags -= CS_PROBE_AUTO_COORD;
1383   }
1384   else
1385     pset->flags |= CS_PROBE_AUTO_COORD;
1386 }
1387 
1388 /*----------------------------------------------------------------------------*/
1389 /*!
1390  * \brief  Set snap mode related to the management of a set of probes.
1391  *
1392  * \param[in, out] pset        pointer to a cs_probe_set_t structure
1393  * \param[in]      snap_mode   snap mode to set
1394  */
1395 /*----------------------------------------------------------------------------*/
1396 
1397 void
cs_probe_set_snap_mode(cs_probe_set_t * pset,cs_probe_snap_t snap_mode)1398 cs_probe_set_snap_mode(cs_probe_set_t   *pset,
1399                        cs_probe_snap_t   snap_mode)
1400 {
1401   if (pset == NULL)
1402     bft_error(__FILE__, __LINE__, 0, _(_err_empty_pset));
1403 
1404   pset->snap_mode = snap_mode;
1405 }
1406 
1407 /*----------------------------------------------------------------------------*/
1408 /*!
1409  * \brief  Set optional parameters related to the management of a set of probes
1410  *
1411  * Available option key names accepting \c true or \c false:
1412  *
1413  * - \c \b transient_location if \c true, relocate probes relative to
1414  *         deforming or moving mesh (default: \c false)
1415  * - \c \b boundary  if \ c true, locate on boundary mesh; if \c false,
1416  *         locate on volume mesh (default)
1417  *
1418  * Other options:
1419  *
1420  * - \c \b selection_criteria where keyval is selection criteria string
1421  * - \c \b tolerance  where keyval is for instance "0.05" (default "0.10")
1422  * - \c \b interpolation if \ c 0, P0 interpolation (default); if \c 1,
1423  *         simple gradient-based interpolation for volume probes (ignored
1424  *         for boundaries). Other interpolation options might be added
1425  *         in the future.
1426  *
1427  * \param[in, out] pset     pointer to a cs_probe_set_t structure to set
1428  * \param[in]      keyname  name of the keyword related to the parameter to set
1429  * \param[in]      keyval   value of the keyword to set
1430  */
1431 /*----------------------------------------------------------------------------*/
1432 
1433 void
cs_probe_set_option(cs_probe_set_t * pset,const char * keyname,const char * keyval)1434 cs_probe_set_option(cs_probe_set_t   *pset,
1435                     const char       *keyname,
1436                     const char       *keyval)
1437 {
1438   if (pset == NULL)
1439     bft_error(__FILE__, __LINE__, 0, _(_err_empty_pset));
1440 
1441   psetkey_t  key = _get_psetkey(keyname);
1442 
1443   if (key == PSETKEY_ERROR) {
1444 
1445     bft_printf("\n\n Current key: %s\n", keyname);
1446     bft_printf(" Possible keys: ");
1447     for (int i = 0; i < PSETKEY_ERROR; i++) {
1448       bft_printf("%s ", _print_psetkey(i));
1449       if (i > 0 && i % 4 == 0) {
1450         bft_printf("\n");
1451         if (i + 1 < PSETKEY_ERROR)
1452           bft_printf("\t");
1453       }
1454     }
1455     bft_error(__FILE__, __LINE__, 0,
1456               _(" Invalid key for probe options %s.\n"
1457                 " Please read run_solver.log for more details and"
1458                 " modify your settings."), pset->name);
1459 
1460   } /* Error message */
1461 
1462   switch(key) {
1463 
1464   case PSETKEY_BOUNDARY:
1465     if (strcmp(keyval, "true") == 0)
1466       pset->flags |= CS_PROBE_BOUNDARY;
1467     else if (strcmp(keyval, "false") == 0) { // remove the flags if it is set
1468       if (pset->flags & CS_PROBE_BOUNDARY)
1469         pset->flags ^= CS_PROBE_BOUNDARY;
1470     }
1471     else
1472       bft_error(__FILE__, __LINE__, 0, _err_truefalse_key, keyval, keyname);
1473     break;
1474 
1475   case PSETKEY_SELECT_CRIT:
1476     {
1477       int  len = strlen(keyval) + 1;
1478       BFT_MALLOC(pset->sel_criter, len, char);
1479       strncpy(pset->sel_criter, keyval, len);
1480     }
1481     break;
1482 
1483   case PSETKEY_TRANSIENT_LOC:
1484     if (strcmp(keyval, "true") == 0)
1485       pset->flags |= CS_PROBE_TRANSIENT;
1486     else if (strcmp(keyval, "false") == 0) { // remove the flag if it is set
1487       if (pset->flags & CS_PROBE_TRANSIENT)
1488         pset->flags ^= CS_PROBE_TRANSIENT;
1489     }
1490     else
1491       bft_error(__FILE__, __LINE__, 0, _err_truefalse_key, keyval, keyname);
1492     break;
1493 
1494   case PSETKEY_TOLERANCE:
1495     pset->tolerance = atof(keyval);
1496     break;
1497 
1498   case PSETKEY_INTERPOLATION:
1499     pset->interpolation = atoi(keyval);
1500     break;
1501 
1502   default:
1503     bft_error(__FILE__, __LINE__, 0,
1504               _(" Key %s is not implemented yet."), keyname);
1505     break;
1506 
1507   } /* Switch on keys */
1508 
1509 }
1510 
1511 /*----------------------------------------------------------------------------*/
1512 /*!
1513  * \brief  Try to locate each probe and define the coordinate really used for
1514  *         the postprocessing step.
1515  *
1516  * For better performance when using multiple probe sets, a pointer to
1517  * an existing location mesh may be passed to this function. The caller is
1518  * responsible for ensuring this mesh matches selection criteria for the
1519  * probe set.
1520  *
1521  * \param[in, out]  pset           pointer to a cs_probe_set_t structure
1522  * \param[in]       location_mesh  optional pointer to mesh relative to which
1523  *                                 probe set should be located, or NULL
1524  */
1525 /*----------------------------------------------------------------------------*/
1526 
1527 void
cs_probe_set_locate(cs_probe_set_t * pset,const fvm_nodal_t * location_mesh)1528 cs_probe_set_locate(cs_probe_set_t     *pset,
1529                     const fvm_nodal_t  *location_mesh)
1530 {
1531   if (pset == NULL)
1532     return;
1533 
1534   bool first_location = false;
1535   const double  tolerance_base = 0.;
1536   const cs_mesh_t  *mesh = cs_glob_mesh;
1537   fvm_nodal_t *_location_mesh = NULL;
1538 
1539   const bool  on_boundary = (pset->flags & CS_PROBE_BOUNDARY) ? true : false;
1540 
1541   /* Build in local case, restore saved coordinates if snapped otherwise */
1542 
1543   if (pset->p_define_func != NULL)
1544     _build_local_probe_set(pset);
1545 
1546   else if (pset->coords_ini != NULL) {
1547     memcpy(pset->coords,
1548            pset->coords_ini,
1549            pset->n_probes*3*sizeof(cs_real_t));
1550   }
1551 
1552   /* Allocate on first pass */
1553 
1554   if (pset->located == NULL) {
1555     BFT_MALLOC(pset->located, pset->n_probes, char);
1556     first_location = true;
1557   }
1558 
1559   /* Reallocate on all passes, in case local sizes change */
1560 
1561   BFT_REALLOC(pset->loc_id, pset->n_probes, cs_lnum_t);
1562   BFT_REALLOC(pset->elt_id, pset->n_probes, cs_lnum_t);
1563   BFT_FREE(pset->vtx_id);
1564 
1565   if (location_mesh == NULL) {
1566 
1567     cs_lnum_t  n_select_elements = 0;
1568     cs_lnum_t  *selected_elements = NULL;
1569 
1570     if (on_boundary) { /* Deal with the surface mesh related to the boundary */
1571 
1572       n_select_elements = mesh->n_b_faces;
1573       if (pset->sel_criter != NULL) {
1574         if (strcmp(pset->sel_criter, "all[]")) {
1575           BFT_MALLOC(selected_elements, mesh->n_b_faces, cs_lnum_t);
1576           cs_selector_get_b_face_num_list(pset->sel_criter,
1577                                           &n_select_elements, selected_elements);
1578         }
1579       } /* Need to define a list of faces ? */
1580 
1581       _location_mesh = cs_mesh_connect_faces_to_nodal(mesh,
1582                                                       "probe_location_mesh",
1583                                                       false, // no family info
1584                                                       0,     // interior faces
1585                                                       n_select_elements,
1586                                                       NULL,  // interior faces
1587                                                       selected_elements);
1588 
1589     }
1590     else { /* Deal with the volumic mesh */
1591 
1592       n_select_elements = mesh->n_cells;
1593       if (pset->sel_criter != NULL) {
1594         if (strcmp(pset->sel_criter, "all[]")) {
1595           BFT_MALLOC(selected_elements, mesh->n_cells, cs_lnum_t);
1596           cs_selector_get_cell_num_list(pset->sel_criter,
1597                                         &n_select_elements, selected_elements);
1598         }
1599       } /* Need to define a list of cells ? */
1600 
1601       _location_mesh = cs_mesh_connect_cells_to_nodal(mesh,
1602                                                       "probe_location_mesh",
1603                                                       false, // no family info
1604                                                       n_select_elements,
1605                                                       selected_elements);
1606 
1607     } /* volumic or surfacic mesh */
1608 
1609     if (selected_elements != NULL)
1610       BFT_FREE(selected_elements);
1611 
1612     location_mesh = _location_mesh;
1613   }
1614 
1615   /* Locate probes on this location mesh */
1616 
1617   float *distance;
1618 
1619   BFT_MALLOC(distance, pset->n_probes, float);
1620 
1621   for (int i = 0; i < pset->n_probes; i++) {
1622     pset->elt_id[i] = -1;
1623     distance[i] = -1.0;
1624   }
1625 
1626   fvm_point_location_nodal(location_mesh,
1627                            tolerance_base,
1628                            pset->tolerance,
1629                            0, /* locate_on_parents */
1630                            pset->n_probes,
1631                            NULL, /* point_tag */
1632                            (const cs_coord_t *)(pset->coords),
1633                            pset->elt_id,
1634                            distance);
1635 
1636   for (int i = 0; i < pset->n_probes; i++) {
1637     if (pset->elt_id[i] < 0) /* Not found */
1638       distance[i] = HUGE_VAL;
1639   }
1640 
1641   /* Warning if points have not beeen located */
1642   cs_gnum_t  n_unlocated_probes = 0;
1643 
1644   int n_loc_probes = 0;
1645 
1646   if (cs_glob_n_ranks == 1 || pset->p_define_func != NULL) {
1647     for (int i = 0; i < pset->n_probes; i++) {
1648       if (distance[i] >= 0.5*HUGE_VAL) {
1649         pset->located[i] = 0;
1650         n_unlocated_probes++;
1651       }
1652       else {
1653         pset->loc_id[n_loc_probes] = i;
1654         pset->elt_id[n_loc_probes] = pset->elt_id[i];
1655         pset->located[i] = 1;
1656         n_loc_probes += 1;
1657       }
1658     }
1659   }
1660 
1661 #if defined(HAVE_MPI)
1662   if (cs_glob_n_ranks > 1 && pset->p_define_func == NULL) {
1663 
1664     cs_double_int_t  *gmin_loc = NULL, *loc = NULL;
1665 
1666     BFT_MALLOC(gmin_loc, pset->n_probes, cs_double_int_t);
1667     BFT_MALLOC(loc, pset->n_probes, cs_double_int_t);
1668 
1669     for (int i = 0; i < pset->n_probes; i++) {
1670       gmin_loc[i].id = loc[i].id = cs_glob_rank_id;
1671       gmin_loc[i].val = loc[i].val = distance[i];
1672     }
1673 
1674     MPI_Allreduce(loc, gmin_loc, pset->n_probes, MPI_DOUBLE_INT, MPI_MINLOC,
1675                   cs_glob_mpi_comm);
1676 
1677     for (int i = 0; i < pset->n_probes; i++) {
1678       if (gmin_loc[i].val >= 0.5*HUGE_VAL) {
1679         pset->located[i] = 0;
1680         n_unlocated_probes++;
1681       }
1682       else {
1683         pset->located[i] = 1;
1684         if (gmin_loc[i].id == cs_glob_rank_id) {
1685           pset->loc_id[n_loc_probes] = i;
1686           pset->elt_id[n_loc_probes] = pset->elt_id[i];
1687           n_loc_probes += 1;
1688         }
1689       }
1690     }
1691 
1692     BFT_FREE(gmin_loc);
1693     BFT_FREE(loc);
1694   }
1695 #endif
1696 
1697   BFT_FREE(distance);
1698 
1699   if (n_unlocated_probes > 0 && first_location) {
1700     bft_printf(_("\n Warning: probe set \"%s\"\n"
1701                  "   %lu (of %ld) probes are not located"
1702                  " on the associated mesh:\n"),
1703                pset->name, (unsigned long)n_unlocated_probes,
1704                (long)pset->n_probes);
1705     for (int i = 0; i < pset->n_probes; i++) {
1706       if (pset->located[i] == 0) {
1707         if (pset->labels == NULL)
1708           bft_printf(_("    %2d ([%8.3e, %8.3e, %8.3e])\n"),
1709                      i+1,
1710                      pset->coords[i][0],
1711                      pset->coords[i][1],
1712                      pset->coords[i][2]);
1713         else
1714           bft_printf(_("    %s ([%8.3e, %8.3e, %8.3e])\n"),
1715                      pset->labels[i],
1716                      pset->coords[i][0],
1717                      pset->coords[i][1],
1718                      pset->coords[i][2]);
1719       }
1720     }
1721   }
1722 
1723   pset->n_loc_probes = n_loc_probes;
1724 
1725   if (n_unlocated_probes) {
1726 
1727     /* For a transient mesh, non-located probes may vary, so we need to maintain
1728        them in the structure to avoid having a varying number of columns in a
1729        plot. We add those locations to the last rank */
1730 
1731     if (pset->flags & CS_PROBE_TRANSIENT) {
1732       if (cs_glob_rank_id == cs_glob_n_ranks - 1 || cs_glob_n_ranks == 1)
1733         pset->n_loc_probes += n_unlocated_probes;
1734     }
1735 
1736     /* For a fixed mesh, if we do not have labels, we add labels based on the
1737        probe global ids so as to maintain probe ids in the output metadata. */
1738 
1739     else if (   pset->labels == NULL
1740              && (pset->flags & CS_PROBE_ON_CURVE) == false) {
1741       BFT_MALLOC(pset->labels, pset->n_probes, char *);
1742       char label[16];
1743       for (int i = 0; i < pset->n_probes; i++) {
1744         snprintf(label, 15, "%d", i+1); label[15] = '\0';
1745         pset->labels[i] = _copy_label(label);
1746       }
1747     }
1748 
1749   }
1750 
1751   BFT_REALLOC(pset->loc_id, pset->n_loc_probes, cs_lnum_t);
1752   BFT_REALLOC(pset->elt_id, pset->n_loc_probes, cs_lnum_t);
1753   BFT_MALLOC(pset->vtx_id, pset->n_loc_probes, cs_lnum_t);
1754 
1755   /* Now also locate relative to closest vertices and update
1756      element num relative to parents */
1757 
1758   cs_coord_3_t  *probe_coords = NULL;
1759   BFT_MALLOC(probe_coords, pset->n_loc_probes, cs_coord_3_t);
1760   for (int i = 0; i < n_loc_probes; i++) {
1761     int j = pset->loc_id[i];
1762     for (int k = 0; k < 3; k++)
1763       probe_coords[i][k] = pset->coords[j][k];
1764   }
1765 
1766   fvm_point_location_closest_vertex(location_mesh,
1767                                     1, /* locate on parents */
1768                                     n_loc_probes,
1769                                     (const cs_coord_t *)(probe_coords),
1770                                     pset->elt_id,
1771                                     pset->vtx_id);
1772 
1773   BFT_FREE(probe_coords);
1774 
1775   /* Now switch to 0-based location */
1776 
1777   for (int i = 0; i < n_loc_probes; i++) {
1778     if (pset->elt_id[i] > -1) {
1779       pset->elt_id[i] -= 1;
1780       pset->vtx_id[i] -= 1;
1781     }
1782   }
1783 
1784   if (_location_mesh != NULL)
1785     _location_mesh = fvm_nodal_destroy(_location_mesh);
1786 
1787   /* Finally, effectively add non-located probe info when required */
1788 
1789   if (pset->n_loc_probes > n_loc_probes) {
1790     for (int i = 0; i < pset->n_probes; i++) {
1791       if (pset->located[i] == 0) {
1792         pset->loc_id[n_loc_probes] = i;
1793         pset->elt_id[n_loc_probes] = -1;
1794         pset->vtx_id[n_loc_probes] = -1;
1795         n_loc_probes++;
1796       }
1797     }
1798   }
1799 
1800   assert(n_loc_probes == pset->n_loc_probes);
1801 }
1802 
1803 /*----------------------------------------------------------------------------*/
1804 /*!
1805  * \brief  Define a fvm_nodal_t structure from the set of probes.
1806  *
1807  * \param[in, out]  pset        pointer to a cs_probe_set_t structure
1808  * \param[in]       mesh_name   name of the mesh to export
1809  *
1810  * \return a pointer to a fvm_nodal_t structure
1811  */
1812 /*----------------------------------------------------------------------------*/
1813 
1814 fvm_nodal_t *
cs_probe_set_export_mesh(cs_probe_set_t * pset,const char * mesh_name)1815 cs_probe_set_export_mesh(cs_probe_set_t   *pset,
1816                          const char       *mesh_name)
1817 {
1818   if (pset == NULL)
1819     return  NULL;
1820 
1821   cs_gnum_t  *global_num = NULL;
1822   cs_coord_3_t  *probe_coords = NULL;
1823   fvm_nodal_t  *exp_mesh = fvm_nodal_create(mesh_name, 3);
1824 
1825   const cs_mesh_t  *m = cs_glob_mesh;
1826   const cs_mesh_quantities_t  *mq = cs_glob_mesh_quantities;
1827 
1828   const cs_real_t *centers = mq->cell_cen;
1829 
1830   if (pset->flags & CS_PROBE_BOUNDARY) centers = mq->b_face_cog;
1831 
1832   /* In case of clip to cell center:
1833      - Remove redundant probes in case of fixed mesh
1834      - save original coordinates in case of variable mesh */
1835 
1836   if (   pset->snap_mode == CS_PROBE_SNAP_ELT_CENTER
1837       && pset->p_define_func == NULL) {
1838     if (! (pset->flags & CS_PROBE_TRANSIENT))
1839       _merge_snapped_to_center(pset, centers);
1840     else if (pset->coords_ini == NULL) {
1841       BFT_MALLOC(pset->coords_ini, pset->n_probes, cs_real_3_t);
1842       memcpy(pset->coords_ini,
1843              pset->coords,
1844              pset->n_probes*3*sizeof(cs_real_t));
1845     }
1846   }
1847 
1848   BFT_MALLOC(probe_coords, pset->n_loc_probes, cs_coord_3_t);
1849   BFT_MALLOC(global_num, pset->n_loc_probes, cs_gnum_t);
1850 
1851   /* Build the final list of probe coordinates */
1852 
1853   cs_real_t  max_distance = 0.;
1854 
1855   for (int i = 0; i < pset->n_loc_probes; i++) {
1856     int j = pset->loc_id[i];
1857 
1858     for (int k = 0; k < 3; k++)
1859       probe_coords[i][k] = pset->coords[j][k];
1860 
1861     global_num[i] = j+1;
1862 
1863     if (pset->elt_id[i] > -1) {
1864       const cs_real_t  *elt_coords = centers + pset->elt_id[i]*3;
1865       cs_real_t v[3];
1866       for (int k = 0; k < 3; k++)
1867         v[k] = elt_coords[k] - pset->coords[j][k];
1868       max_distance = fmax(max_distance, cs_math_3_square_norm(v));
1869     }
1870   }
1871 
1872   cs_real_t gmax_distance = max_distance;
1873 
1874   /* Handle snap mode if active */
1875 
1876   if (pset->snap_mode == CS_PROBE_SNAP_ELT_CENTER) {
1877     for (int i = 0; i < pset->n_loc_probes; i++) {
1878       if (pset->elt_id[i] > -1) {
1879         const int j = pset->loc_id[i];
1880         const cs_real_t  *elt_coords = centers + pset->elt_id[i]*3;
1881         for (int k = 0; k < 3; k++)
1882           pset->coords[j][k] = elt_coords[k];
1883       }
1884     }
1885   }
1886   else if (pset->snap_mode == CS_PROBE_SNAP_VERTEX) {
1887     for (int i = 0; i < pset->n_loc_probes; i++) {
1888       if (pset->vtx_id[i] > -1) {
1889         const int j = pset->loc_id[i];
1890         const cs_real_t  *vtx_coords = m->vtx_coord + pset->vtx_id[i]*3;
1891         for (int k = 0; k < 3; k++)
1892           pset->coords[j][k] = vtx_coords[k];
1893       }
1894     }
1895   }
1896 
1897   /* Update the probe set structure */
1898 
1899   fvm_nodal_define_vertex_list(exp_mesh, pset->n_loc_probes, NULL);
1900   fvm_nodal_transfer_vertices(exp_mesh, (cs_coord_t *)probe_coords);
1901 
1902   /* Set a global numbering if needed */
1903 
1904   if (pset->p_define_func != NULL) {
1905     cs_real_t *s = cs_probe_set_get_loc_curvilinear_abscissa(pset);
1906     fvm_io_num_t *vtx_io_num
1907       = fvm_io_num_create_from_real(s, pset->n_loc_probes);
1908     BFT_FREE(s);
1909     fvm_nodal_transfer_vertex_io_num(exp_mesh, &vtx_io_num);
1910   }
1911   else if (cs_glob_n_ranks > 1) {
1912     fvm_nodal_init_io_num(exp_mesh, global_num, 0); // 0 = vertices
1913 
1914 #if defined(HAVE_MPI)
1915     MPI_Reduce(&max_distance, &gmax_distance, 1, MPI_DOUBLE, MPI_MAX, 0,
1916                cs_glob_mpi_comm);
1917 #endif
1918   }
1919 
1920   if (! (   pset->flags & CS_PROBE_ON_CURVE
1921          || pset->flags & CS_PROBE_TRANSIENT))
1922     bft_printf(_("\n Probe set: \"%s\":\n"
1923                  "   maximum distance between cell centers and"
1924                  " requested coordinates:"
1925                  " %5.3e\n"), pset->name, gmax_distance);
1926 
1927   BFT_FREE(global_num);
1928 
1929   /* Add global labels */
1930 
1931   if (pset->labels != NULL) {
1932     char **g_labels;
1933     int ngl = fvm_nodal_get_n_g_vertices(exp_mesh);
1934     BFT_MALLOC(g_labels, ngl, char *);
1935 
1936     int j = 0;
1937     for (int i = 0; i < pset->n_probes; i++) {
1938       if (pset->located[i] != 0)
1939         g_labels[j++] = _copy_label(pset->labels[i]);
1940     }
1941     assert(j == ngl);
1942     fvm_nodal_transfer_global_vertex_labels(exp_mesh, g_labels);
1943   }
1944 
1945   /* probe_coords is managed by exp_mesh */
1946 
1947   return exp_mesh;
1948 }
1949 
1950 /*----------------------------------------------------------------------------*/
1951 /*!
1952  * \brief  Define a fvm_nodal_t structure from the set of unlocated probes.
1953  *
1954  * \param[in, out]  pset        pointer to a cs_probe_set_t structure
1955  * \param[in]       mesh_name   name of the mesh to export
1956  *
1957  * \return a pointer to a fvm_nodal_t structure
1958  */
1959 /*----------------------------------------------------------------------------*/
1960 
1961 fvm_nodal_t *
cs_probe_set_unlocated_export_mesh(cs_probe_set_t * pset,const char * mesh_name)1962 cs_probe_set_unlocated_export_mesh(cs_probe_set_t   *pset,
1963                                    const char       *mesh_name)
1964 {
1965   if (pset == NULL)
1966     return  NULL;
1967 
1968   int  n_exp_probes = 0;
1969   cs_gnum_t  *global_num = NULL;
1970   cs_coord_3_t  *probe_coords = NULL;
1971   fvm_nodal_t  *exp_mesh = fvm_nodal_create(mesh_name, 3);
1972 
1973   BFT_MALLOC(probe_coords, pset->n_probes, cs_coord_3_t);
1974   BFT_MALLOC(global_num, pset->n_loc_probes, cs_gnum_t);
1975 
1976   /* Build the final list of probe coordinates */
1977 
1978   for (int i = 0; i < pset->n_probes; i++) {
1979     if (pset->located[i] == 0) {
1980       for (int k = 0; k < 3; k++)
1981         probe_coords[n_exp_probes][k] = pset->coords[i][k];
1982       global_num[n_exp_probes] = i+1;
1983       n_exp_probes++;
1984     }
1985   }
1986 
1987   /* Update the probe set structure */
1988 
1989   fvm_nodal_define_vertex_list(exp_mesh, n_exp_probes, NULL);
1990   fvm_nodal_transfer_vertices(exp_mesh, (cs_coord_t *)probe_coords);
1991 
1992   /* Set a global numbering if needed */
1993 
1994   if (pset->p_define_func != NULL) {
1995     cs_real_t *s;
1996     BFT_MALLOC(s, pset->n_probes, cs_real_t);
1997     int j = 0;
1998     for (int i = 0; i < pset->n_probes; i++) {
1999       if (pset->located[i] == 0)
2000         s[j++] = pset->s_coords[i];
2001     }
2002     fvm_io_num_t *vtx_io_num
2003       = fvm_io_num_create_from_real(pset->s_coords, j);
2004     BFT_FREE(s);
2005     fvm_nodal_transfer_vertex_io_num(exp_mesh, &vtx_io_num);
2006   }
2007   else if (cs_glob_n_ranks > 1)
2008     fvm_nodal_init_io_num(exp_mesh, global_num, 0); // 0 = vertices
2009 
2010   BFT_FREE(global_num);
2011 
2012   /* Add global labels */
2013 
2014   if (pset->labels != NULL) {
2015     char **g_labels;
2016     int ngl = fvm_nodal_get_n_g_vertices(exp_mesh);
2017     BFT_MALLOC(g_labels, ngl, char *);
2018 
2019     int j = 0;
2020     for (int i = 0; i < pset->n_probes; i++) {
2021       if (pset->located[i] == 0)
2022         g_labels[j++] = _copy_label(pset->labels[i]);
2023     }
2024     assert(j == ngl);
2025     fvm_nodal_transfer_global_vertex_labels(exp_mesh, g_labels);
2026   }
2027 
2028   /* probe_coords is managed by exp_mesh */
2029 
2030   return exp_mesh;
2031 }
2032 
2033 /*----------------------------------------------------------------------------*/
2034 /*!
2035  * \brief  Dump a cs_probe_set_t structure.
2036  *
2037  * \param[in]  pset    pointer to a cs_probe_set_t structure
2038  */
2039 /*----------------------------------------------------------------------------*/
2040 
2041 void
cs_probe_set_dump(const cs_probe_set_t * pset)2042 cs_probe_set_dump(const cs_probe_set_t   *pset)
2043 {
2044   bft_printf("\n\n Dump cs_probe_set_t structure %p\n", (const void *)pset);
2045 
2046   if (pset == NULL)
2047     return;
2048 
2049   bft_printf(" name:                %s\n"
2050              " flags:               %d\n"
2051              " location criteria:   %s\n"
2052              " tolerance:           %5.3e\n",
2053              pset->name, pset->flags, pset->sel_criter,
2054              pset->tolerance);
2055 
2056   if (pset->sel_criter != NULL)
2057     bft_printf(" selection:  %s\n", pset->sel_criter);
2058 
2059   bft_printf(" n_probes:   %d; %d; %d (locally located; defined; max.)\n",
2060              (int)pset->n_loc_probes, (int)pset->n_probes,
2061              (int)pset->n_max_probes);
2062 
2063   for (int i = 0; i < pset->n_probes; i++) {
2064 
2065     bft_printf(" %4d | %-5.3e %-5.3e %-5.3e |", i,
2066                pset->coords[i][0], pset->coords[i][1], pset->coords[i][2]);
2067 
2068     if (pset->s_coords != NULL)
2069       bft_printf(" %5.3e |", pset->s_coords[i]);
2070     if (pset->elt_id != NULL && pset->located != NULL)
2071       bft_printf(" %6d | %c |", (int)pset->elt_id[i], pset->located[i]);
2072     if (pset->labels != NULL)
2073       if (pset->labels[i] != NULL)
2074         bft_printf(" %s", pset->labels[i]);
2075     bft_printf("\n");
2076 
2077   }
2078 }
2079 
2080 /*----------------------------------------------------------------------------*/
2081 /*!
2082  * \brief  Retrieve the main members of a cs_probe_set_t structure.
2083  *
2084  * \param[in]       pset       pointer to a cs_probe_set_t structure
2085  * \param[in, out]  snap_mode  mode of location
2086  * \param[in, out]  n_probes   number of probes
2087  * \param[in, out]  coords     probe coordinates
2088  */
2089 /*----------------------------------------------------------------------------*/
2090 
2091 void
cs_probe_set_get_members(const cs_probe_set_t * pset,cs_probe_snap_t * snap_mode,int * n_probes,cs_real_3_t * coords[])2092 cs_probe_set_get_members(const cs_probe_set_t   *pset,
2093                          cs_probe_snap_t        *snap_mode,
2094                          int                    *n_probes,
2095                          cs_real_3_t            *coords[])
2096 {
2097   if (pset == NULL)
2098     return;
2099 
2100   /* Return pointers */
2101 
2102   if (snap_mode != NULL)
2103     *snap_mode = pset->snap_mode;
2104 
2105   if (n_probes != NULL)
2106     *n_probes = pset->n_probes;
2107 
2108   if (coords != NULL)
2109     *coords = pset->coords;
2110 }
2111 
2112 /*----------------------------------------------------------------------------*/
2113 /*!
2114  * \brief  Return the number probes in the local domain.
2115  *
2116  * \param[in]       pset       pointer to a cs_probe_set_t structure
2117  *
2118  * \return  number of probes in local domain
2119  */
2120 /*----------------------------------------------------------------------------*/
2121 
2122 int
cs_probe_set_get_n_local(const cs_probe_set_t * pset)2123 cs_probe_set_get_n_local(const cs_probe_set_t   *pset)
2124 {
2125   int retval = 0;
2126 
2127   if (pset != NULL)
2128     retval = pset->n_loc_probes;
2129 
2130   return retval;
2131 }
2132 
2133 /*----------------------------------------------------------------------------*/
2134 /*!
2135  * \brief  Return the list of curvilinear abscissa for the given probe set
2136  *
2137  * \param[in]  pset              pointer to a cs_probe_set_t structure
2138  *
2139  * \return NULL or the pointer to the array of abscissa
2140  */
2141 /*----------------------------------------------------------------------------*/
2142 
2143 const cs_real_t *
cs_probe_set_get_curvilinear_abscissa(const cs_probe_set_t * pset)2144 cs_probe_set_get_curvilinear_abscissa(const cs_probe_set_t   *pset)
2145 {
2146   const cs_real_t *retval = NULL;
2147 
2148   if (pset == NULL)
2149     return retval;
2150   else
2151     return pset->s_coords;
2152 }
2153 
2154 /*----------------------------------------------------------------------------*/
2155 /*!
2156  * \brief  Return the list of curvilinear abscissa of probes located
2157  *         on the local ranks for the given probe set
2158  *
2159  * The caller is responsible for freeing the returned array.
2160  *
2161  * \param[in]  pset              pointer to a cs_probe_set_t structure
2162  *
2163  * \return NULL or the pointer to the array of abscissa
2164  */
2165 /*----------------------------------------------------------------------------*/
2166 
2167 cs_real_t *
cs_probe_set_get_loc_curvilinear_abscissa(const cs_probe_set_t * pset)2168 cs_probe_set_get_loc_curvilinear_abscissa(const cs_probe_set_t   *pset)
2169 {
2170   if (pset == NULL)
2171     return NULL;
2172 
2173   cs_real_t *s;
2174   BFT_MALLOC(s, pset->n_loc_probes, cs_real_t);
2175   for (int i = 0; i < pset->n_loc_probes; i++) {
2176     int j = pset->loc_id[i];
2177     s[i] = pset->s_coords[j];
2178   }
2179 
2180   return s;
2181 }
2182 
2183 /*----------------------------------------------------------------------------*/
2184 /*!
2185  * \brief  Return the ids of a probe set's local matching elements, relative
2186  *         to a given mesh location.
2187  *
2188  * The mesh_location id must match one of \ref CS_MESH_LOCATION_CELLS,
2189  * \ref CS_MESH_LOCATION_BOUNDARY_FACES, or \ref CS_MESH_LOCATION_VERTICES.
2190  *
2191  * \param[in]  pset              pointer to a cs_probe_set_t structure
2192  * \param[in]  mesh_location_id  id of parent mesh location
2193  */
2194 /*----------------------------------------------------------------------------*/
2195 
2196 const cs_lnum_t *
cs_probe_set_get_elt_ids(const cs_probe_set_t * pset,int mesh_location_id)2197 cs_probe_set_get_elt_ids(const cs_probe_set_t  *pset,
2198                          int                    mesh_location_id)
2199 {
2200   const cs_lnum_t *retval = NULL;
2201 
2202   if (pset == NULL)
2203     return retval;
2204 
2205   bool on_boundary = (pset->flags & CS_PROBE_BOUNDARY) ? true : false;
2206 
2207   if (mesh_location_id == CS_MESH_LOCATION_CELLS && on_boundary == false)
2208     retval = pset->elt_id;
2209   else if (mesh_location_id == CS_MESH_LOCATION_BOUNDARY_FACES && on_boundary)
2210     retval = pset->elt_id;
2211 
2212   else if (CS_MESH_LOCATION_VERTICES)
2213     retval = pset->vtx_id;
2214 
2215   return retval;
2216 }
2217 
2218 /*----------------------------------------------------------------------------*/
2219 
2220 END_C_DECLS
2221