1 /*============================================================================
2 * Manage the (generic) evaluation of extended definitions
3 *============================================================================*/
4
5 /*
6 This file is part of Code_Saturne, a general-purpose CFD tool.
7
8 Copyright (C) 1998-2021 EDF S.A.
9
10 This program is free software; you can redistribute it and/or modify it under
11 the terms of the GNU General Public License as published by the Free Software
12 Foundation; either version 2 of the License, or (at your option) any later
13 version.
14
15 This program is distributed in the hope that it will be useful, but WITHOUT
16 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
17 FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
18 details.
19
20 You should have received a copy of the GNU General Public License along with
21 this program; if not, write to the Free Software Foundation, Inc., 51 Franklin
22 Street, Fifth Floor, Boston, MA 02110-1301, USA.
23 */
24
25 /*----------------------------------------------------------------------------*/
26
27 #include "cs_defs.h"
28
29 /*----------------------------------------------------------------------------
30 * Standard C library headers
31 *----------------------------------------------------------------------------*/
32
33 #include <assert.h>
34 #include <ctype.h>
35 #include <stdio.h>
36 #include <stdlib.h>
37 #include <string.h>
38 #include <math.h>
39
40 /*----------------------------------------------------------------------------
41 * Local headers
42 *----------------------------------------------------------------------------*/
43
44 #include <bft_mem.h>
45
46 #include "cs_defs.h"
47 #include "cs_field.h"
48 #include "cs_mesh_location.h"
49 #include "cs_reco.h"
50
51 /*----------------------------------------------------------------------------
52 * Header for the current file
53 *----------------------------------------------------------------------------*/
54
55 #include "cs_xdef_eval.h"
56
57 /*----------------------------------------------------------------------------*/
58
59 BEGIN_C_DECLS
60
61 /*=============================================================================
62 * Local Macro definitions and structure definitions
63 *============================================================================*/
64
65 /* Redefined the name of functions from cs_math to get shorter names */
66 #define _dp3 cs_math_3_dot_product
67
68 /*============================================================================
69 * Private function prototypes
70 *============================================================================*/
71
72 /*============================================================================
73 * Public function prototypes
74 *============================================================================*/
75
76 /*----------------------------------------------------------------------------*/
77 /*!
78 * \brief Evaluate a scalar-valued quantity for a list of elements
79 * This function complies with the generic function type defined as
80 * cs_xdef_eval_t
81 *
82 * \param[in] n_elts number of elements to consider
83 * \param[in] elt_ids list of element ids
84 * \param[in] dense_output perform an indirection for output (true/false)
85 * \param[in] mesh pointer to a cs_mesh_t structure
86 * \param[in] connect pointer to a cs_cdo_connect_t structure
87 * \param[in] quant pointer to a cs_cdo_quantities_t structure
88 * \param[in] time_eval physical time at which one evaluates the term
89 * \param[in] context NULL or pointer to a context structure
90 * \param[in, out] eval array storing the result (must be allocated)
91 */
92 /*----------------------------------------------------------------------------*/
93
94 void
cs_xdef_eval_scalar_by_val(cs_lnum_t n_elts,const cs_lnum_t * elt_ids,bool dense_output,const cs_mesh_t * mesh,const cs_cdo_connect_t * connect,const cs_cdo_quantities_t * quant,cs_real_t time_eval,void * context,cs_real_t * eval)95 cs_xdef_eval_scalar_by_val(cs_lnum_t n_elts,
96 const cs_lnum_t *elt_ids,
97 bool dense_output,
98 const cs_mesh_t *mesh,
99 const cs_cdo_connect_t *connect,
100 const cs_cdo_quantities_t *quant,
101 cs_real_t time_eval,
102 void *context,
103 cs_real_t *eval)
104 {
105 CS_UNUSED(mesh);
106 CS_UNUSED(quant);
107 CS_UNUSED(connect);
108 CS_UNUSED(time_eval);
109
110 if (n_elts == 0)
111 return;
112
113 const cs_real_t *constant_val = (cs_real_t *)context;
114
115 assert(eval != NULL && constant_val != NULL);
116
117 if (elt_ids != NULL && !dense_output) {
118
119 # pragma omp parallel for if (n_elts > CS_THR_MIN)
120 for (cs_lnum_t i = 0; i < n_elts; i++)
121 eval[elt_ids[i]] = constant_val[0];
122
123 }
124 else {
125
126 # pragma omp parallel for if (n_elts > CS_THR_MIN)
127 for (cs_lnum_t i = 0; i < n_elts; i++)
128 eval[i] = constant_val[0];
129
130 }
131 }
132
133 /*----------------------------------------------------------------------------*/
134 /*!
135 * \brief Evaluate a vector-valued quantity for a list of elements
136 * This function complies with the generic function type defined as
137 * cs_xdef_eval_t
138 *
139 * \param[in] n_elts number of elements to consider
140 * \param[in] elt_ids list of element ids
141 * \param[in] dense_output perform an indirection for output (true/false)
142 * \param[in] mesh pointer to a cs_mesh_t structure
143 * \param[in] connect pointer to a cs_cdo_connect_t structure
144 * \param[in] quant pointer to a cs_cdo_quantities_t structure
145 * \param[in] time_eval physical time at which one evaluates the term
146 * \param[in] context NULL or pointer to a context structure
147 * \param[in, out] eval array storing the result (must be allocated)
148 */
149 /*----------------------------------------------------------------------------*/
150
151 void
cs_xdef_eval_vector_by_val(cs_lnum_t n_elts,const cs_lnum_t * elt_ids,bool dense_output,const cs_mesh_t * mesh,const cs_cdo_connect_t * connect,const cs_cdo_quantities_t * quant,cs_real_t time_eval,void * context,cs_real_t * eval)152 cs_xdef_eval_vector_by_val(cs_lnum_t n_elts,
153 const cs_lnum_t *elt_ids,
154 bool dense_output,
155 const cs_mesh_t *mesh,
156 const cs_cdo_connect_t *connect,
157 const cs_cdo_quantities_t *quant,
158 cs_real_t time_eval,
159 void *context,
160 cs_real_t *eval)
161 {
162 CS_UNUSED(mesh);
163 CS_UNUSED(quant);
164 CS_UNUSED(connect);
165 CS_UNUSED(time_eval);
166
167 if (n_elts == 0)
168 return;
169
170 const cs_real_t *constant_val = (cs_real_t *)context;
171
172 /* Sanity checks */
173 assert(eval != NULL && constant_val != NULL);
174
175 if (elt_ids != NULL && !dense_output) {
176
177 # pragma omp parallel for if (n_elts > CS_THR_MIN)
178 for (cs_lnum_t i = 0; i < n_elts; i++) {
179
180 cs_real_t *_res = eval + 3*elt_ids[i];
181
182 _res[0] = constant_val[0];
183 _res[1] = constant_val[1];
184 _res[2] = constant_val[2];
185
186 }
187
188 }
189 else {
190
191 # pragma omp parallel for if (n_elts > CS_THR_MIN)
192 for (cs_lnum_t i = 0; i < n_elts; i++) {
193
194 cs_real_t *_res = eval + 3*i;
195
196 _res[0] = constant_val[0];
197 _res[1] = constant_val[1];
198 _res[2] = constant_val[2];
199
200 }
201
202 }
203 }
204
205 /*----------------------------------------------------------------------------*/
206 /*!
207 * \brief Evaluate a tensor-valued quantity for a list of elements with
208 * symmetric storage.
209 * This function complies with the generic function type defined as
210 * cs_xdef_eval_t
211 *
212 * \param[in] n_elts number of elements to consider
213 * \param[in] elt_ids list of element ids
214 * \param[in] dense_output perform an indirection for output (true/false)
215 * \param[in] mesh pointer to a cs_mesh_t structure
216 * \param[in] connect pointer to a cs_cdo_connect_t structure
217 * \param[in] quant pointer to a cs_cdo_quantities_t structure
218 * \param[in] time_eval physical time at which one evaluates the term
219 * \param[in] context NULL or pointer to a context structure
220 * \param[in, out] eval array storing the result (must be allocated)
221 */
222 /*----------------------------------------------------------------------------*/
223
224 void
cs_xdef_eval_symtens_by_val(cs_lnum_t n_elts,const cs_lnum_t * elt_ids,bool dense_output,const cs_mesh_t * mesh,const cs_cdo_connect_t * connect,const cs_cdo_quantities_t * quant,cs_real_t time_eval,void * context,cs_real_t * eval)225 cs_xdef_eval_symtens_by_val(cs_lnum_t n_elts,
226 const cs_lnum_t *elt_ids,
227 bool dense_output,
228 const cs_mesh_t *mesh,
229 const cs_cdo_connect_t *connect,
230 const cs_cdo_quantities_t *quant,
231 cs_real_t time_eval,
232 void *context,
233 cs_real_t *eval)
234 {
235 CS_UNUSED(quant);
236 CS_UNUSED(mesh);
237 CS_UNUSED(connect);
238 CS_UNUSED(time_eval);
239
240 if (n_elts == 0)
241 return;
242
243 const cs_real_t *constant_val = (const cs_real_t *)context;
244
245 /* Sanity checks */
246 assert(eval != NULL && constant_val != NULL);
247
248 if (elt_ids != NULL && !dense_output) {
249
250 # pragma omp parallel for if (n_elts > CS_THR_MIN)
251 for (cs_lnum_t i = 0; i < n_elts; i++) {
252
253 cs_real_t *shift_eval = eval + 6*elt_ids[i];
254 for (int k = 0; k < 6; k++)
255 shift_eval[k] = constant_val[k];
256
257 }
258
259 }
260 else {
261
262 # pragma omp parallel for if (n_elts > CS_THR_MIN)
263 for (cs_lnum_t i = 0; i < n_elts; i++) {
264
265 cs_real_t *shift_eval = eval + 6*i;
266 for (int k = 0; k < 6; k++)
267 shift_eval[k] = constant_val[k];
268
269 }
270
271 }
272 }
273
274 /*----------------------------------------------------------------------------*/
275 /*!
276 * \brief Evaluate a tensor-valued quantity for a list of elements
277 * This function complies with the generic function type defined as
278 * cs_xdef_eval_t
279 *
280 * \param[in] n_elts number of elements to consider
281 * \param[in] elt_ids list of element ids
282 * \param[in] dense_output perform an indirection for output (true/false)
283 * \param[in] mesh pointer to a cs_mesh_t structure
284 * \param[in] connect pointer to a cs_cdo_connect_t structure
285 * \param[in] quant pointer to a cs_cdo_quantities_t structure
286 * \param[in] time_eval physical time at which one evaluates the term
287 * \param[in] context NULL or pointer to a context structure
288 * \param[in, out] eval array storing the result (must be allocated)
289 */
290 /*----------------------------------------------------------------------------*/
291
292 void
cs_xdef_eval_tensor_by_val(cs_lnum_t n_elts,const cs_lnum_t * elt_ids,bool dense_output,const cs_mesh_t * mesh,const cs_cdo_connect_t * connect,const cs_cdo_quantities_t * quant,cs_real_t time_eval,void * context,cs_real_t * eval)293 cs_xdef_eval_tensor_by_val(cs_lnum_t n_elts,
294 const cs_lnum_t *elt_ids,
295 bool dense_output,
296 const cs_mesh_t *mesh,
297 const cs_cdo_connect_t *connect,
298 const cs_cdo_quantities_t *quant,
299 cs_real_t time_eval,
300 void *context,
301 cs_real_t *eval)
302 {
303 CS_UNUSED(quant);
304 CS_UNUSED(mesh);
305 CS_UNUSED(connect);
306 CS_UNUSED(time_eval);
307
308 if (n_elts == 0)
309 return;
310
311 const cs_real_3_t *constant_val = (const cs_real_3_t *)context;
312
313 /* Sanity checks */
314 assert(eval != NULL && constant_val != NULL);
315
316 if (elt_ids != NULL && !dense_output) {
317
318 # pragma omp parallel for if (n_elts > CS_THR_MIN)
319 for (cs_lnum_t i = 0; i < n_elts; i++) {
320
321 cs_real_t *shift_eval = eval + 9*elt_ids[i];
322 for (int ki = 0; ki < 3; ki++)
323 for (int kj = 0; kj < 3; kj++)
324 shift_eval[3*ki+kj] = constant_val[ki][kj];
325
326 }
327
328 }
329 else {
330
331 # pragma omp parallel for if (n_elts > CS_THR_MIN)
332 for (cs_lnum_t i = 0; i < n_elts; i++) {
333
334 cs_real_t *shift_eval = eval + 9*i;
335 for (int ki = 0; ki < 3; ki++)
336 for (int kj = 0; kj < 3; kj++)
337 shift_eval[3*ki+kj] = constant_val[ki][kj];
338
339 }
340
341 }
342 }
343
344 /*----------------------------------------------------------------------------*/
345 /*!
346 * \brief Evaluate a scalar-valued quantity with only a time-dependent
347 * variation for a list of elements
348 * This function complies with the generic function type defined as
349 * cs_xdef_eval_t
350 *
351 * \param[in] n_elts number of elements to consider
352 * \param[in] elt_ids list of element ids
353 * \param[in] dense_output perform an indirection for output (true/false)
354 * \param[in] mesh pointer to a cs_mesh_t structure
355 * \param[in] connect pointer to a cs_cdo_connect_t structure
356 * \param[in] quant pointer to a cs_cdo_quantities_t structure
357 * \param[in] time_eval physical time at which one evaluates the term
358 * \param[in] context NULL or pointer to a context structure
359 * \param[in, out] eval array storing the result (must be allocated)
360 */
361 /*----------------------------------------------------------------------------*/
362
363 void
cs_xdef_eval_scalar_at_cells_by_time_func(cs_lnum_t n_elts,const cs_lnum_t * elt_ids,bool dense_output,const cs_mesh_t * mesh,const cs_cdo_connect_t * connect,const cs_cdo_quantities_t * quant,cs_real_t time_eval,void * context,cs_real_t * eval)364 cs_xdef_eval_scalar_at_cells_by_time_func(cs_lnum_t n_elts,
365 const cs_lnum_t *elt_ids,
366 bool dense_output,
367 const cs_mesh_t *mesh,
368 const cs_cdo_connect_t *connect,
369 const cs_cdo_quantities_t *quant,
370 cs_real_t time_eval,
371 void *context,
372 cs_real_t *eval)
373 {
374 CS_UNUSED(mesh);
375 CS_UNUSED(quant);
376 CS_UNUSED(connect);
377
378 cs_xdef_time_func_context_t *tfc = (cs_xdef_time_func_context_t *)context;
379
380 /* Sanity checks */
381 assert(tfc != NULL);
382
383 /* Evaluate the quantity only once */
384 cs_real_t _eval;
385 tfc->func(time_eval, tfc->input, &_eval);
386
387 if (elt_ids != NULL && !dense_output) {
388
389 # pragma omp parallel for if (n_elts > CS_THR_MIN)
390 for (cs_lnum_t i = 0; i < n_elts; i++)
391 eval[elt_ids[i]] = _eval;
392
393 }
394 else {
395
396 for (cs_lnum_t i = 0; i < n_elts; i++)
397 eval[i] = _eval;
398
399 }
400 }
401
402 /*----------------------------------------------------------------------------*/
403 /*!
404 * \brief Evaluate a vector-valued quantity with only a time-dependent
405 * variation for a list of elements
406 * This function complies with the generic function type defined as
407 * cs_xdef_eval_t
408 *
409 * \param[in] n_elts number of elements to consider
410 * \param[in] elt_ids list of element ids
411 * \param[in] dense_output perform an indirection for output (true/false)
412 * \param[in] mesh pointer to a cs_mesh_t structure
413 * \param[in] connect pointer to a cs_cdo_connect_t structure
414 * \param[in] quant pointer to a cs_cdo_quantities_t structure
415 * \param[in] time_eval physical time at which one evaluates the term
416 * \param[in] context NULL or pointer to a context structure
417 * \param[in, out] eval array storing the result (must be allocated)
418 */
419 /*----------------------------------------------------------------------------*/
420
421 void
cs_xdef_eval_vector_at_cells_by_time_func(cs_lnum_t n_elts,const cs_lnum_t * elt_ids,bool dense_output,const cs_mesh_t * mesh,const cs_cdo_connect_t * connect,const cs_cdo_quantities_t * quant,cs_real_t time_eval,void * context,cs_real_t * eval)422 cs_xdef_eval_vector_at_cells_by_time_func(cs_lnum_t n_elts,
423 const cs_lnum_t *elt_ids,
424 bool dense_output,
425 const cs_mesh_t *mesh,
426 const cs_cdo_connect_t *connect,
427 const cs_cdo_quantities_t *quant,
428 cs_real_t time_eval,
429 void *context,
430 cs_real_t *eval)
431 {
432 CS_UNUSED(mesh);
433 CS_UNUSED(quant);
434 CS_UNUSED(connect);
435
436 cs_xdef_time_func_context_t *tfc = (cs_xdef_time_func_context_t *)context;
437
438 /* Sanity checks */
439 assert(tfc != NULL);
440
441 /* Evaluate the quantity */
442 cs_real_t _eval[3];
443 tfc->func(time_eval, tfc->input, _eval);
444
445 if (elt_ids != NULL && !dense_output) {
446
447 # pragma omp parallel for if (n_elts > CS_THR_MIN)
448 for (cs_lnum_t i = 0; i < n_elts; i++)
449 for (int k = 0; k < 3; k++)
450 eval[3*elt_ids[i] + k] = _eval[k];
451
452 }
453 else {
454
455 # pragma omp parallel for if (n_elts > CS_THR_MIN)
456 for (cs_lnum_t i = 0; i < n_elts; i++)
457 for (int k = 0; k < 3; k++)
458 eval[3*i+k] = _eval[k];
459
460 }
461 }
462
463 /*----------------------------------------------------------------------------*/
464 /*!
465 * \brief Evaluate a tensor-valued quantity with a symmetric storage and with
466 * only a time-dependent variation for a list of elements
467 * This function complies with the generic function type defined as
468 * cs_xdef_eval_t
469 *
470 * \param[in] n_elts number of elements to consider
471 * \param[in] elt_ids list of element ids
472 * \param[in] dense_output perform an indirection for output (true/false)
473 * \param[in] mesh pointer to a cs_mesh_t structure
474 * \param[in] connect pointer to a cs_cdo_connect_t structure
475 * \param[in] quant pointer to a cs_cdo_quantities_t structure
476 * \param[in] time_eval physical time at which one evaluates the term
477 * \param[in] context NULL or pointer to a context structure
478 * \param[in, out] eval array storing the result (must be allocated)
479 */
480 /*----------------------------------------------------------------------------*/
481
482 void
cs_xdef_eval_symtens_at_cells_by_time_func(cs_lnum_t n_elts,const cs_lnum_t * elt_ids,bool dense_output,const cs_mesh_t * mesh,const cs_cdo_connect_t * connect,const cs_cdo_quantities_t * quant,cs_real_t time_eval,void * context,cs_real_t * eval)483 cs_xdef_eval_symtens_at_cells_by_time_func(cs_lnum_t n_elts,
484 const cs_lnum_t *elt_ids,
485 bool dense_output,
486 const cs_mesh_t *mesh,
487 const cs_cdo_connect_t *connect,
488 const cs_cdo_quantities_t *quant,
489 cs_real_t time_eval,
490 void *context,
491 cs_real_t *eval)
492 {
493 CS_UNUSED(mesh);
494 CS_UNUSED(quant);
495 CS_UNUSED(connect);
496
497 cs_xdef_time_func_context_t *tfc = (cs_xdef_time_func_context_t *)context;
498
499 /* Sanity checks */
500 assert(tfc != NULL);
501
502 /* Evaluate the quantity */
503 cs_real_t _eval[6];
504 tfc->func(time_eval, tfc->input, _eval);
505
506 if (elt_ids != NULL && !dense_output) {
507
508 # pragma omp parallel for if (n_elts > CS_THR_MIN)
509 for (cs_lnum_t i = 0; i < n_elts; i++)
510 for (int k = 0; k < 6; k++)
511 eval[6*elt_ids[i] + k] = _eval[k];
512
513 }
514 else {
515
516 # pragma omp parallel for if (n_elts > CS_THR_MIN)
517 for (cs_lnum_t i = 0; i < n_elts; i++)
518 for (int k = 0; k < 6; k++)
519 eval[6*i+k] = _eval[k];
520
521 }
522 }
523
524 /*----------------------------------------------------------------------------*/
525 /*!
526 * \brief Evaluate a tensor-valued quantity with only a time-dependent
527 * variation for a list of elements
528 * This function complies with the generic function type defined as
529 * cs_xdef_eval_t
530 *
531 * \param[in] n_elts number of elements to consider
532 * \param[in] elt_ids list of element ids
533 * \param[in] dense_output perform an indirection for output (true/false)
534 * \param[in] mesh pointer to a cs_mesh_t structure
535 * \param[in] connect pointer to a cs_cdo_connect_t structure
536 * \param[in] quant pointer to a cs_cdo_quantities_t structure
537 * \param[in] time_eval physical time at which one evaluates the term
538 * \param[in] context NULL or pointer to a context structure
539 * \param[in, out] eval array storing the result (must be allocated)
540 */
541 /*----------------------------------------------------------------------------*/
542
543 void
cs_xdef_eval_tensor_at_cells_by_time_func(cs_lnum_t n_elts,const cs_lnum_t * elt_ids,bool dense_output,const cs_mesh_t * mesh,const cs_cdo_connect_t * connect,const cs_cdo_quantities_t * quant,cs_real_t time_eval,void * context,cs_real_t * eval)544 cs_xdef_eval_tensor_at_cells_by_time_func(cs_lnum_t n_elts,
545 const cs_lnum_t *elt_ids,
546 bool dense_output,
547 const cs_mesh_t *mesh,
548 const cs_cdo_connect_t *connect,
549 const cs_cdo_quantities_t *quant,
550 cs_real_t time_eval,
551 void *context,
552 cs_real_t *eval)
553 {
554 CS_UNUSED(mesh);
555 CS_UNUSED(quant);
556 CS_UNUSED(connect);
557
558 cs_xdef_time_func_context_t *tfc = (cs_xdef_time_func_context_t *)context;
559
560 /* Sanity checks */
561 assert(tfc != NULL);
562
563 /* Evaluate the quantity */
564 cs_real_t _eval[9];
565 tfc->func(time_eval, tfc->input, _eval);
566
567 if (elt_ids != NULL && !dense_output) {
568
569 # pragma omp parallel for if (n_elts > CS_THR_MIN)
570 for (cs_lnum_t i = 0; i < n_elts; i++)
571 for (int k = 0; k < 9; k++)
572 eval[9*elt_ids[i] + k] = _eval[k];
573
574 }
575 else {
576
577 # pragma omp parallel for if (n_elts > CS_THR_MIN)
578 for (cs_lnum_t i = 0; i < n_elts; i++)
579 for (int k = 0; k < 9; k++)
580 eval[9*i+k] = _eval[k];
581
582 }
583 }
584
585 /*----------------------------------------------------------------------------*/
586 /*!
587 * \brief Evaluate a quantity defined at cells using an analytic function
588 * This function complies with the generic function type defined as
589 * cs_xdef_eval_t
590 *
591 * \param[in] n_elts number of elements to consider
592 * \param[in] elt_ids list of element ids
593 * \param[in] dense_output perform an indirection for output (true/false)
594 * \param[in] mesh pointer to a cs_mesh_t structure
595 * \param[in] connect pointer to a cs_cdo_connect_t structure
596 * \param[in] quant pointer to a cs_cdo_quantities_t structure
597 * \param[in] time_eval physical time at which one evaluates the term
598 * \param[in] context NULL or pointer to a context structure
599 * \param[in, out] eval array storing the result (must be allocated)
600 */
601 /*----------------------------------------------------------------------------*/
602
603 void
cs_xdef_eval_at_cells_by_analytic(cs_lnum_t n_elts,const cs_lnum_t * elt_ids,bool dense_output,const cs_mesh_t * mesh,const cs_cdo_connect_t * connect,const cs_cdo_quantities_t * quant,cs_real_t time_eval,void * context,cs_real_t * eval)604 cs_xdef_eval_at_cells_by_analytic(cs_lnum_t n_elts,
605 const cs_lnum_t *elt_ids,
606 bool dense_output,
607 const cs_mesh_t *mesh,
608 const cs_cdo_connect_t *connect,
609 const cs_cdo_quantities_t *quant,
610 cs_real_t time_eval,
611 void *context,
612 cs_real_t *eval)
613 {
614 CS_UNUSED(mesh);
615 CS_UNUSED(connect);
616
617 const cs_real_t *cell_centers = (quant != NULL) ? quant->cell_centers : NULL;
618
619 cs_xdef_analytic_context_t *cx = (cs_xdef_analytic_context_t *)context;
620
621 /* Sanity checks */
622 assert(cx != NULL);
623
624 /* Evaluate the function for this time at the cell center */
625 cx->func(time_eval, n_elts, elt_ids, cell_centers, dense_output, cx->input,
626 eval);
627 }
628
629 /*----------------------------------------------------------------------------*/
630 /*!
631 * \brief Evaluate a quantity defined at border faces using an analytic
632 * function
633 * This function complies with the generic function type defined as
634 * cs_xdef_eval_t
635 *
636 * \param[in] n_elts number of elements to consider
637 * \param[in] elt_ids list of element ids
638 * \param[in] dense_output perform an indirection for output (true/false)
639 * \param[in] mesh pointer to a cs_mesh_t structure
640 * \param[in] connect pointer to a cs_cdo_connect_t structure
641 * \param[in] quant pointer to a cs_cdo_quantities_t structure
642 * \param[in] time_eval physical time at which one evaluates the term
643 * \param[in] context NULL or pointer to a context structure
644 * \param[in, out] eval array storing the result (must be allocated)
645 */
646 /*----------------------------------------------------------------------------*/
647
648 void
cs_xdef_eval_at_b_faces_by_analytic(cs_lnum_t n_elts,const cs_lnum_t * elt_ids,bool dense_output,const cs_mesh_t * mesh,const cs_cdo_connect_t * connect,const cs_cdo_quantities_t * quant,cs_real_t time_eval,void * context,cs_real_t * eval)649 cs_xdef_eval_at_b_faces_by_analytic(cs_lnum_t n_elts,
650 const cs_lnum_t *elt_ids,
651 bool dense_output,
652 const cs_mesh_t *mesh,
653 const cs_cdo_connect_t *connect,
654 const cs_cdo_quantities_t *quant,
655 cs_real_t time_eval,
656 void *context,
657 cs_real_t *eval)
658 {
659 CS_UNUSED(mesh);
660 CS_UNUSED(connect);
661
662 const cs_real_t *bf_centers = (quant != NULL) ? quant->b_face_center : NULL;
663
664 cs_xdef_analytic_context_t *cx = (cs_xdef_analytic_context_t *)context;
665
666 /* Sanity checks */
667 assert(cx != NULL);
668
669 /* Evaluate the function for this time at the border face center */
670 cx->func(time_eval, n_elts, elt_ids, bf_centers, dense_output, cx->input,
671 eval);
672 }
673
674 /*----------------------------------------------------------------------------*/
675 /*!
676 * \brief Evaluate a quantity defined at vertices using an analytic function
677 * This function complies with the generic function type defined as
678 * cs_xdef_eval_t
679 *
680 * \param[in] n_elts number of elements to consider
681 * \param[in] elt_ids list of element ids
682 * \param[in] dense_output perform an indirection for output (true/false)
683 * \param[in] mesh pointer to a cs_mesh_t structure
684 * \param[in] connect pointer to a cs_cdo_connect_t structure
685 * \param[in] quant pointer to a cs_cdo_quantities_t structure
686 * \param[in] time_eval physical time at which one evaluates the term
687 * \param[in] context NULL or pointer to a context structure
688 * \param[in, out] eval array storing the result (must be allocated)
689 */
690 /*----------------------------------------------------------------------------*/
691
692 void
cs_xdef_eval_at_vertices_by_analytic(cs_lnum_t n_elts,const cs_lnum_t * elt_ids,bool dense_output,const cs_mesh_t * mesh,const cs_cdo_connect_t * connect,const cs_cdo_quantities_t * quant,cs_real_t time_eval,void * context,cs_real_t * eval)693 cs_xdef_eval_at_vertices_by_analytic(cs_lnum_t n_elts,
694 const cs_lnum_t *elt_ids,
695 bool dense_output,
696 const cs_mesh_t *mesh,
697 const cs_cdo_connect_t *connect,
698 const cs_cdo_quantities_t *quant,
699 cs_real_t time_eval,
700 void *context,
701 cs_real_t *eval)
702 {
703 CS_UNUSED(connect);
704
705 if (n_elts == 0)
706 return;
707
708 cs_xdef_analytic_context_t *cx = (cs_xdef_analytic_context_t *)context;
709
710 /* Sanity checks */
711 assert(eval != NULL || cx != NULL);
712
713 const cs_real_t *v_coords;
714 if (quant != NULL)
715 v_coords = quant->vtx_coord;
716 else if (mesh != NULL)
717 v_coords = mesh->vtx_coord;
718 else {
719 v_coords = NULL;/* avoid a compilation warning */
720 bft_error(__FILE__, __LINE__, 0, "%s: No vertex coordinates available.",
721 __func__);
722 }
723
724 /* Evaluate the function for this time at the cell center */
725 cx->func(time_eval, n_elts, elt_ids, v_coords, dense_output, cx->input,
726 eval);
727 }
728
729 /*----------------------------------------------------------------------------*/
730 /*!
731 * \brief Evaluate a quantity defined at primal cells using a function
732 * associated to dof (dof = degrees of freedom).
733 * This function complies with the generic function type defined as
734 * cs_xdef_eval_t
735 *
736 * \param[in] n_elts number of elements to consider
737 * \param[in] elt_ids list of element ids
738 * \param[in] dense_output perform an indirection for output (true/false)
739 * \param[in] mesh pointer to a cs_mesh_t structure
740 * \param[in] connect pointer to a cs_cdo_connect_t structure
741 * \param[in] quant pointer to a cs_cdo_quantities_t structure
742 * \param[in] time_eval physical time at which one evaluates the term
743 * \param[in] context NULL or pointer to a context structure
744 * \param[in, out] eval array storing the result (must be allocated)
745 */
746 /*----------------------------------------------------------------------------*/
747
748 void
cs_xdef_eval_at_cells_by_dof_func(cs_lnum_t n_elts,const cs_lnum_t * elt_ids,bool dense_output,const cs_mesh_t * mesh,const cs_cdo_connect_t * connect,const cs_cdo_quantities_t * quant,cs_real_t time_eval,void * context,cs_real_t * eval)749 cs_xdef_eval_at_cells_by_dof_func(cs_lnum_t n_elts,
750 const cs_lnum_t *elt_ids,
751 bool dense_output,
752 const cs_mesh_t *mesh,
753 const cs_cdo_connect_t *connect,
754 const cs_cdo_quantities_t *quant,
755 cs_real_t time_eval,
756 void *context,
757 cs_real_t *eval)
758 {
759 CS_UNUSED(mesh);
760 CS_UNUSED(connect);
761 CS_UNUSED(quant);
762 CS_UNUSED(time_eval);
763
764 cs_xdef_dof_context_t *cx = (cs_xdef_dof_context_t *)context;
765
766 /* Sanity check */
767 assert(cx != NULL);
768
769 /* Values of the function are defined at the cells */
770 if (cs_flag_test(cx->loc, cs_flag_primal_cell))
771 cx->func(n_elts, elt_ids, dense_output, cx->input,
772 eval);
773 else
774 bft_error(__FILE__, __LINE__, 0, "%s: Invalid location.\n", __func__);
775 }
776
777 /*----------------------------------------------------------------------------*/
778 /*!
779 * \brief Evaluate a quantity defined at vertices using a function
780 * associated to dof (dof = degrees of freedom).
781 * This function complies with the generic function type defined as
782 * cs_xdef_eval_t
783 *
784 * \param[in] n_elts number of elements to consider
785 * \param[in] elt_ids list of element ids
786 * \param[in] dense_output perform an indirection for output (true/false)
787 * \param[in] mesh pointer to a cs_mesh_t structure
788 * \param[in] connect pointer to a cs_cdo_connect_t structure
789 * \param[in] quant pointer to a cs_cdo_quantities_t structure
790 * \param[in] time_eval physical time at which one evaluates the term
791 * \param[in] context NULL or pointer to a context structure
792 * \param[in, out] eval array storing the result (must be allocated)
793 */
794 /*----------------------------------------------------------------------------*/
795
796 void
cs_xdef_eval_at_vertices_by_dof_func(cs_lnum_t n_elts,const cs_lnum_t * elt_ids,bool dense_output,const cs_mesh_t * mesh,const cs_cdo_connect_t * connect,const cs_cdo_quantities_t * quant,cs_real_t time_eval,void * context,cs_real_t * eval)797 cs_xdef_eval_at_vertices_by_dof_func(cs_lnum_t n_elts,
798 const cs_lnum_t *elt_ids,
799 bool dense_output,
800 const cs_mesh_t *mesh,
801 const cs_cdo_connect_t *connect,
802 const cs_cdo_quantities_t *quant,
803 cs_real_t time_eval,
804 void *context,
805 cs_real_t *eval)
806 {
807 CS_UNUSED(mesh);
808 CS_UNUSED(connect);
809 CS_UNUSED(quant);
810 CS_UNUSED(time_eval);
811
812 cs_xdef_dof_context_t *cx = (cs_xdef_dof_context_t *)context;
813
814 /* Sanity check */
815 assert(cx != NULL);
816
817 /* Values of the function are defined at vertices */
818 if (cs_flag_test(cx->loc, cs_flag_primal_vtx))
819 cx->func(n_elts, elt_ids, dense_output, cx->input,
820 eval);
821 else
822 bft_error(__FILE__, __LINE__, 0, "%s: Invalid location.\n", __func__);
823 }
824
825 /*----------------------------------------------------------------------------*/
826 /*!
827 * \brief Evaluate a quantity defined at boundary faces using a function
828 * associated to dof (dof = degrees of freedom).
829 * This function complies with the generic function type defined as
830 * cs_xdef_eval_t
831 *
832 * \param[in] n_elts number of elements to consider
833 * \param[in] elt_ids list of element ids
834 * \param[in] dense_output perform an indirection for output (true/false)
835 * \param[in] mesh pointer to a cs_mesh_t structure
836 * \param[in] connect pointer to a cs_cdo_connect_t structure
837 * \param[in] quant pointer to a cs_cdo_quantities_t structure
838 * \param[in] time_eval physical time at which one evaluates the term
839 * \param[in] context NULL or pointer to a context structure
840 * \param[in, out] eval array storing the result (must be allocated)
841 */
842 /*----------------------------------------------------------------------------*/
843
844 void
cs_xdef_eval_at_b_faces_by_dof_func(cs_lnum_t n_elts,const cs_lnum_t * elt_ids,bool dense_output,const cs_mesh_t * mesh,const cs_cdo_connect_t * connect,const cs_cdo_quantities_t * quant,cs_real_t time_eval,void * context,cs_real_t * eval)845 cs_xdef_eval_at_b_faces_by_dof_func(cs_lnum_t n_elts,
846 const cs_lnum_t *elt_ids,
847 bool dense_output,
848 const cs_mesh_t *mesh,
849 const cs_cdo_connect_t *connect,
850 const cs_cdo_quantities_t *quant,
851 cs_real_t time_eval,
852 void *context,
853 cs_real_t *eval)
854 {
855 CS_UNUSED(mesh);
856 CS_UNUSED(connect);
857 CS_UNUSED(quant);
858 CS_UNUSED(time_eval);
859
860 cs_xdef_dof_context_t *cx = (cs_xdef_dof_context_t *)context;
861
862 /* Sanity check */
863 assert(cx != NULL);
864
865 /* Values of the function are defined at the boundary faces */
866 if (cs_flag_test(cx->loc, cs_flag_boundary_face))
867 cx->func(n_elts, elt_ids, dense_output, cx->input, eval);
868 else
869 bft_error(__FILE__, __LINE__, 0, "%s: Invalid location.\n", __func__);
870 }
871
872 /*----------------------------------------------------------------------------*/
873 /*!
874 * \brief Evaluate a scalar-valued quantity at cells defined by an array.
875 * Array is assumed to be interlaced.
876 * This function complies with the generic function type defined as
877 * cs_xdef_eval_t
878 *
879 * \param[in] n_elts number of elements to consider
880 * \param[in] elt_ids list of element ids
881 * \param[in] dense_output perform an indirection for output (true/false)
882 * \param[in] mesh pointer to a cs_mesh_t structure
883 * \param[in] connect pointer to a cs_cdo_connect_t structure
884 * \param[in] quant pointer to a cs_cdo_quantities_t structure
885 * \param[in] time_eval physical time at which one evaluates the term
886 * \param[in] context NULL or pointer to a context structure
887 * \param[in, out] eval array storing the result (must be allocated)
888 */
889 /*----------------------------------------------------------------------------*/
890
891 void
cs_xdef_eval_scalar_at_cells_by_array(cs_lnum_t n_elts,const cs_lnum_t * elt_ids,bool dense_output,const cs_mesh_t * mesh,const cs_cdo_connect_t * connect,const cs_cdo_quantities_t * quant,cs_real_t time_eval,void * context,cs_real_t * eval)892 cs_xdef_eval_scalar_at_cells_by_array(cs_lnum_t n_elts,
893 const cs_lnum_t *elt_ids,
894 bool dense_output,
895 const cs_mesh_t *mesh,
896 const cs_cdo_connect_t *connect,
897 const cs_cdo_quantities_t *quant,
898 cs_real_t time_eval,
899 void *context,
900 cs_real_t *eval)
901 {
902 CS_UNUSED(mesh);
903 CS_UNUSED(time_eval);
904
905 if (n_elts == 0)
906 return;
907
908 cs_xdef_array_context_t *cx = (cs_xdef_array_context_t *)context;
909
910 /* Sanity checks */
911 assert(eval != NULL || cx != NULL);
912 assert(cx->stride == 1);
913
914 if (cs_flag_test(cx->loc, cs_flag_primal_cell)) {
915
916 if (elt_ids != NULL && !dense_output) {
917
918 for (cs_lnum_t i = 0; i < n_elts; i++) {
919 const cs_lnum_t c_id = elt_ids[i];
920 eval[c_id] = cx->values[c_id];
921 }
922
923 }
924 else if (elt_ids != NULL && dense_output) {
925
926 for (cs_lnum_t i = 0; i < n_elts; i++)
927 eval[i] = cx->values[elt_ids[i]];
928
929 }
930 else {
931
932 assert(elt_ids == NULL);
933 memcpy(eval, cx->values, n_elts * sizeof(cs_real_t));
934
935 }
936
937 }
938 else if (cs_flag_test(cx->loc, cs_flag_primal_vtx)) {
939
940 assert(connect != NULL && quant != NULL);
941 if (elt_ids != NULL && !dense_output) {
942
943 for (cs_lnum_t i = 0; i < n_elts; i++) {
944 const cs_lnum_t c_id = elt_ids[i];
945 cs_reco_pv_at_cell_center(c_id, connect->c2v, quant, cx->values,
946 eval + c_id);
947 }
948
949 }
950 else if (elt_ids != NULL && dense_output) {
951
952 for (cs_lnum_t i = 0; i < n_elts; i++)
953 cs_reco_pv_at_cell_center(elt_ids[i], connect->c2v, quant, cx->values,
954 eval + i);
955
956 }
957 else {
958
959 assert(elt_ids == NULL);
960 for (cs_lnum_t i = 0; i < n_elts; i++)
961 cs_reco_pv_at_cell_center(i, connect->c2v, quant, cx->values,
962 eval + i);
963
964 }
965
966 }
967 else
968 bft_error(__FILE__, __LINE__, 0,
969 " %s: Invalid support for the input array", __func__);
970
971 }
972
973 /*----------------------------------------------------------------------------*/
974 /*!
975 * \brief Evaluate a nd-valued quantity at cells defined by an array.
976 * Array is assumed to be interlaced.
977 * This function complies with the generic function type defined as
978 * cs_xdef_eval_t
979 *
980 * \param[in] n_elts number of elements to consider
981 * \param[in] elt_ids list of element ids
982 * \param[in] dense_output perform an indirection for output (true/false)
983 * \param[in] mesh pointer to a cs_mesh_t structure
984 * \param[in] connect pointer to a cs_cdo_connect_t structure
985 * \param[in] quant pointer to a cs_cdo_quantities_t structure
986 * \param[in] time_eval physical time at which one evaluates the term
987 * \param[in] context NULL or pointer to a context structure
988 * \param[in, out] eval array storing the result (must be allocated)
989 */
990 /*----------------------------------------------------------------------------*/
991
992 void
cs_xdef_eval_nd_at_cells_by_array(cs_lnum_t n_elts,const cs_lnum_t * elt_ids,bool dense_output,const cs_mesh_t * mesh,const cs_cdo_connect_t * connect,const cs_cdo_quantities_t * quant,cs_real_t time_eval,void * context,cs_real_t * eval)993 cs_xdef_eval_nd_at_cells_by_array(cs_lnum_t n_elts,
994 const cs_lnum_t *elt_ids,
995 bool dense_output,
996 const cs_mesh_t *mesh,
997 const cs_cdo_connect_t *connect,
998 const cs_cdo_quantities_t *quant,
999 cs_real_t time_eval,
1000 void *context,
1001 cs_real_t *eval)
1002 {
1003 CS_UNUSED(mesh);
1004 CS_UNUSED(time_eval);
1005
1006 if (n_elts == 0)
1007 return;
1008
1009 cs_xdef_array_context_t *cx = (cs_xdef_array_context_t *)context;
1010
1011 /* Sanity checks */
1012 assert(eval != NULL || cx != NULL);
1013
1014 const int stride = cx->stride;
1015
1016 if (cs_flag_test(cx->loc, cs_flag_primal_cell)) {
1017
1018 assert(stride > 1);
1019 if (elt_ids != NULL && !dense_output) {
1020
1021 for (cs_lnum_t i = 0; i < n_elts; i++) {
1022 const cs_lnum_t c_id = elt_ids[i];
1023 for (int k = 0; k < stride; k++)
1024 eval[stride*c_id + k] = cx->values[stride*c_id + k];
1025 }
1026
1027 }
1028 else if (elt_ids != NULL && dense_output) {
1029
1030 for (cs_lnum_t i = 0; i < n_elts; i++) {
1031 const cs_lnum_t c_id = elt_ids[i];
1032 for (int k = 0; k < stride; k++)
1033 eval[stride*i + k] = cx->values[stride*c_id + k];
1034 }
1035
1036 }
1037 else { /* All elements are considered */
1038
1039 assert(elt_ids == NULL);
1040 memcpy(eval, cx->values, stride*n_elts * sizeof(cs_real_t));
1041
1042 }
1043
1044 }
1045 else if (cs_flag_test(cx->loc, cs_flag_dual_face_byc)) {
1046
1047 /* Sanity checks */
1048 assert(stride == 3);
1049 assert(connect!= NULL && quant != NULL);
1050 assert(cx->index == connect->c2e->idx);
1051
1052 if (elt_ids != NULL && !dense_output) {
1053
1054 for (cs_lnum_t i = 0; i < n_elts; i++) {
1055 const cs_lnum_t c_id = elt_ids[i];
1056 cs_reco_dfbyc_at_cell_center(c_id, connect->c2e, quant, cx->values,
1057 eval + c_id*stride);
1058 }
1059
1060 }
1061 else if (elt_ids != NULL && dense_output) {
1062
1063 for (cs_lnum_t i = 0; i < n_elts; i++)
1064 cs_reco_dfbyc_at_cell_center(elt_ids[i], connect->c2e, quant,
1065 cx->values,
1066 eval + i*stride);
1067
1068 }
1069 else {
1070
1071 for (cs_lnum_t i = 0; i < n_elts; i++)
1072 cs_reco_dfbyc_at_cell_center(i, connect->c2e, quant, cx->values,
1073 eval + i*stride);
1074
1075 }
1076
1077 }
1078 else
1079 bft_error(__FILE__, __LINE__, 0,
1080 " %s: Invalid case for the input array", __func__);
1081
1082 }
1083
1084 /*----------------------------------------------------------------------------*/
1085 /*!
1086 * \brief Evaluate a quantity defined at vertices using an array
1087 * This function complies with the generic function type defined as
1088 * cs_xdef_eval_t
1089 *
1090 * \param[in] n_elts number of elements to consider
1091 * \param[in] elt_ids list of element ids
1092 * \param[in] dense_output perform an indirection for output (true/false)
1093 * \param[in] mesh pointer to a cs_mesh_t structure
1094 * \param[in] connect pointer to a cs_cdo_connect_t structure
1095 * \param[in] quant pointer to a cs_cdo_quantities_t structure
1096 * \param[in] time_eval physical time at which one evaluates the term
1097 * \param[in] context NULL or pointer to a context structure
1098 * \param[in, out] eval array storing the result (must be allocated)
1099 */
1100 /*----------------------------------------------------------------------------*/
1101
1102 void
cs_xdef_eval_at_vertices_by_array(cs_lnum_t n_elts,const cs_lnum_t * elt_ids,bool dense_output,const cs_mesh_t * mesh,const cs_cdo_connect_t * connect,const cs_cdo_quantities_t * quant,cs_real_t time_eval,void * context,cs_real_t * eval)1103 cs_xdef_eval_at_vertices_by_array(cs_lnum_t n_elts,
1104 const cs_lnum_t *elt_ids,
1105 bool dense_output,
1106 const cs_mesh_t *mesh,
1107 const cs_cdo_connect_t *connect,
1108 const cs_cdo_quantities_t *quant,
1109 cs_real_t time_eval,
1110 void *context,
1111 cs_real_t *eval)
1112 {
1113 CS_UNUSED(mesh);
1114 CS_UNUSED(connect);
1115 CS_UNUSED(quant);
1116 CS_UNUSED(time_eval);
1117
1118 if (n_elts == 0)
1119 return;
1120
1121 cs_xdef_array_context_t *cx = (cs_xdef_array_context_t *)context;
1122
1123 /* Sanity checks */
1124 assert(eval != NULL || cx != NULL);
1125
1126 const int stride = cx->stride;
1127
1128 if (cs_flag_test(cx->loc, cs_flag_primal_vtx)) {
1129
1130 if (elt_ids != NULL && !dense_output) {
1131
1132 switch (stride) {
1133
1134 case 1: /* Scalar-valued */
1135 for (cs_lnum_t i = 0; i < n_elts; i++) {
1136 const cs_lnum_t v_id = elt_ids[i];
1137 eval[v_id] = cx->values[v_id];
1138 }
1139 break;
1140
1141 default:
1142 for (cs_lnum_t i = 0; i < n_elts; i++) {
1143 const cs_lnum_t v_id = elt_ids[i];
1144 for (int j = 0; j < stride; j++)
1145 eval[stride*v_id + j] = cx->values[stride*v_id+j];
1146 }
1147 break;
1148
1149 } /* End of switch */
1150
1151 }
1152 else if (elt_ids != NULL && dense_output) {
1153
1154 switch (stride) {
1155
1156 case 1: /* Scalar-valued */
1157 for (cs_lnum_t i = 0; i < n_elts; i++)
1158 eval[i] = cx->values[elt_ids[i]];
1159 break;
1160
1161 default:
1162 for (cs_lnum_t i = 0; i < n_elts; i++) {
1163 for (int j = 0; j < stride; j++)
1164 eval[stride*i + j] = cx->values[stride*elt_ids[i] + j];
1165 }
1166 break;
1167
1168 } /* End of switch */
1169
1170 }
1171 else {
1172
1173 assert(elt_ids == NULL);
1174 memcpy(eval, cx->values, n_elts*stride * sizeof(cs_real_t));
1175
1176 }
1177
1178 }
1179 else
1180 bft_error(__FILE__, __LINE__, 0,
1181 " %s: Invalid support for the input array", __func__);
1182 }
1183
1184 /*----------------------------------------------------------------------------*/
1185 /*!
1186 * \brief Evaluate a quantity inside a cell defined using a field
1187 * This function complies with the generic function type defined as
1188 * cs_xdef_eval_t
1189 *
1190 * \param[in] n_elts number of elements to consider
1191 * \param[in] elt_ids list of element ids
1192 * \param[in] dense_output perform an indirection for output (true/false)
1193 * \param[in] mesh pointer to a cs_mesh_t structure
1194 * \param[in] connect pointer to a cs_cdo_connect_t structure
1195 * \param[in] quant pointer to a cs_cdo_quantities_t structure
1196 * \param[in] time_eval physical time at which one evaluates the term
1197 * \param[in] context NULL or pointer to a context structure
1198 * \param[in, out] eval array storing the result (must be allocated)
1199 */
1200 /*----------------------------------------------------------------------------*/
1201
1202 void
cs_xdef_eval_cell_by_field(cs_lnum_t n_elts,const cs_lnum_t * elt_ids,bool dense_output,const cs_mesh_t * mesh,const cs_cdo_connect_t * connect,const cs_cdo_quantities_t * quant,cs_real_t time_eval,void * context,cs_real_t * eval)1203 cs_xdef_eval_cell_by_field(cs_lnum_t n_elts,
1204 const cs_lnum_t *elt_ids,
1205 bool dense_output,
1206 const cs_mesh_t *mesh,
1207 const cs_cdo_connect_t *connect,
1208 const cs_cdo_quantities_t *quant,
1209 cs_real_t time_eval,
1210 void *context,
1211 cs_real_t *eval)
1212 {
1213 CS_UNUSED(mesh);
1214 CS_UNUSED(time_eval);
1215
1216 if (n_elts == 0)
1217 return;
1218
1219 cs_field_t *field = (cs_field_t *)context;
1220
1221 /* Sanity checks */
1222 assert(eval != NULL || field != NULL);
1223
1224 cs_real_t *values = field->val;
1225
1226 const int c_ml_id = cs_mesh_location_get_id_by_name(N_("cells"));
1227 const int v_ml_id = cs_mesh_location_get_id_by_name(N_("vertices"));
1228
1229 if (field->location_id == c_ml_id) {
1230
1231 if (elt_ids != NULL && !dense_output) {
1232 for (cs_lnum_t i = 0; i < n_elts; i++) {
1233 const cs_lnum_t c_id = elt_ids[i];
1234 for (int k = 0; k < field->dim; k++)
1235 eval[field->dim*c_id + k] = values[field->dim*c_id + k];
1236 }
1237 }
1238 else if (elt_ids != NULL && dense_output) {
1239
1240 for (cs_lnum_t i = 0; i < n_elts; i++) {
1241 const cs_lnum_t c_id = elt_ids[i];
1242 for (int k = 0; k < field->dim; k++)
1243 eval[field->dim*i + k] = values[field->dim*c_id + k];
1244 }
1245
1246 }
1247 else {
1248
1249 assert(elt_ids == NULL);
1250 memcpy(eval, values, field->dim*n_elts * sizeof(cs_real_t));
1251
1252 }
1253
1254 }
1255 else if (field->location_id == v_ml_id) {
1256
1257 assert(connect != NULL);
1258 if (field->dim > 1)
1259 bft_error(__FILE__, __LINE__, 0,
1260 " %s: Invalid field dimension.", __func__);
1261
1262 if (elt_ids != NULL && !dense_output) {
1263 for (cs_lnum_t i = 0; i < n_elts; i++) {
1264
1265 const cs_lnum_t c_id = elt_ids[i];
1266 cs_reco_pv_at_cell_center(c_id,
1267 connect->c2v,
1268 quant,
1269 values,
1270 eval + c_id);
1271
1272 }
1273 }
1274 else if (elt_ids != NULL && dense_output) {
1275
1276 for (cs_lnum_t i = 0; i < n_elts; i++) {
1277
1278 const cs_lnum_t c_id = elt_ids[i];
1279 cs_reco_pv_at_cell_center(c_id,
1280 connect->c2v,
1281 quant,
1282 values,
1283 eval + i);
1284
1285 }
1286
1287 }
1288 else {
1289
1290 assert(elt_ids == NULL);
1291 for (cs_lnum_t c_id = 0; c_id < n_elts; c_id++)
1292 cs_reco_pv_at_cell_center(c_id,
1293 connect->c2v,
1294 quant,
1295 values,
1296 eval + c_id);
1297
1298 }
1299
1300 }
1301 else
1302 bft_error(__FILE__, __LINE__, 0,
1303 " %s: Invalid case for the input field", __func__);
1304 }
1305
1306 /*----------------------------------------------------------------------------*/
1307 /*!
1308 * \brief Evaluate a quantity defined at border faces using an analytic
1309 * function
1310 *
1311 * \param[in] n_elts number of elements to consider
1312 * \param[in] elt_ids list of element ids
1313 * \param[in] dense_output perform an indirection for output (true/false)
1314 * \param[in] mesh pointer to a cs_mesh_t structure
1315 * \param[in] connect pointer to a cs_cdo_connect_t structure
1316 * \param[in] quant pointer to a cs_cdo_quantities_t structure
1317 * \param[in] time_eval physical time at which one evaluates the term
1318 * \param[in] context NULL or pointer to a context structure
1319 * \param[in] qtype quadrature type
1320 * \param[in] dim dimension of the analytic function return
1321 * \param[in, out] eval array storing the result (must be allocated)
1322 */
1323 /*----------------------------------------------------------------------------*/
1324
1325 void
cs_xdef_eval_avg_at_b_faces_by_analytic(cs_lnum_t n_elts,const cs_lnum_t * elt_ids,bool dense_output,const cs_mesh_t * mesh,const cs_cdo_connect_t * connect,const cs_cdo_quantities_t * quant,cs_real_t time_eval,void * context,cs_quadrature_type_t qtype,int dim,cs_real_t * eval)1326 cs_xdef_eval_avg_at_b_faces_by_analytic(cs_lnum_t n_elts,
1327 const cs_lnum_t *elt_ids,
1328 bool dense_output,
1329 const cs_mesh_t *mesh,
1330 const cs_cdo_connect_t *connect,
1331 const cs_cdo_quantities_t *quant,
1332 cs_real_t time_eval,
1333 void *context,
1334 cs_quadrature_type_t qtype,
1335 int dim,
1336 cs_real_t *eval)
1337 {
1338 CS_UNUSED(mesh);
1339
1340 cs_xdef_analytic_context_t *cx = (cs_xdef_analytic_context_t *)context;
1341
1342 /* Sanity checks */
1343 assert(cx != NULL);
1344 assert(connect != NULL && quant != NULL);
1345
1346 cs_quadrature_tria_integral_t
1347 *qfunc = cs_quadrature_get_tria_integral(dim, qtype);
1348
1349 const cs_adjacency_t *f2e = connect->f2e;
1350 const cs_adjacency_t *e2v = connect->e2v;
1351 const cs_real_t *xv = quant->vtx_coord;
1352
1353 if (elt_ids == NULL) { /* All boundary faces are selected */
1354
1355 # pragma omp parallel for if (quant->n_b_faces > CS_THR_MIN)
1356 for (cs_lnum_t bf_id = 0; bf_id < quant->n_b_faces; bf_id++) {
1357
1358 const cs_lnum_t f_id = quant->n_i_faces + bf_id;
1359 const cs_quant_t pfq = cs_quant_set_face(f_id, quant);
1360 const cs_lnum_t start = f2e->idx[f_id], end = f2e->idx[f_id+1];
1361 double *val_i = eval + dim*bf_id;
1362
1363 /* Resetting */
1364 memset(val_i, 0, dim*sizeof(double));
1365
1366 switch (end - start) {
1367
1368 case CS_TRIANGLE_CASE:
1369 {
1370 cs_lnum_t v1, v2, v3;
1371 cs_connect_get_next_3_vertices(f2e->ids, e2v->ids, start,
1372 &v1, &v2, &v3);
1373 qfunc(time_eval, xv + 3*v1, xv + 3*v2, xv + 3*v3, pfq.meas,
1374 cx->func, cx->input, val_i);
1375 }
1376 break;
1377
1378 default:
1379 for (cs_lnum_t j = start; j < end; j++) {
1380
1381 const cs_lnum_t _2e = 2*f2e->ids[j];
1382 const cs_lnum_t v1 = e2v->ids[_2e];
1383 const cs_lnum_t v2 = e2v->ids[_2e+1];
1384
1385 qfunc(time_eval, xv + 3*v1, xv + 3*v2, pfq.center,
1386 cs_math_surftri(xv + 3*v1, xv + 3*v2, pfq.center),
1387 cx->func, cx->input, val_i);
1388
1389 } /* Loop on edges */
1390
1391 } /* Switch on the type of face. Special case for triangles */
1392
1393 /* Compute the average */
1394 const double _os = 1./pfq.meas;
1395 for (int k = 0; k < dim; k++)
1396 val_i[k] *= _os;
1397
1398 } /* Loop on faces */
1399
1400 }
1401 else { /* There is an indirection list */
1402
1403 # pragma omp parallel for if (n_elts > CS_THR_MIN)
1404 for (cs_lnum_t i = 0; i < n_elts; i++) { /* Loop on selected faces */
1405
1406 const cs_lnum_t bf_id = elt_ids[i];
1407 const cs_lnum_t f_id = quant->n_i_faces + bf_id;
1408 const cs_quant_t pfq = cs_quant_set_face(f_id, quant);
1409 const cs_lnum_t start = f2e->idx[f_id], end = f2e->idx[f_id+1];
1410
1411 double *val_i = dense_output ? eval + dim*i : eval + dim*bf_id;
1412
1413 /* Resetting */
1414 memset(val_i, 0, dim*sizeof(double));
1415
1416 switch (end - start) {
1417
1418 case CS_TRIANGLE_CASE:
1419 {
1420 /* If triangle, one-shot computation */
1421 cs_lnum_t v1, v2, v3;
1422 cs_connect_get_next_3_vertices(f2e->ids, e2v->ids, start,
1423 &v1, &v2, &v3);
1424 qfunc(time_eval, xv + 3*v1, xv + 3*v2, xv + 3*v3,
1425 pfq.meas, cx->func, cx->input, val_i);
1426 }
1427 break;
1428
1429 default:
1430 for (cs_lnum_t j = start; j < end; j++) {
1431
1432 const cs_lnum_t _2e = 2*f2e->ids[j];
1433 const cs_lnum_t v1 = e2v->ids[_2e];
1434 const cs_lnum_t v2 = e2v->ids[_2e+1];
1435
1436 qfunc(time_eval, xv + 3*v1, xv + 3*v2, pfq.center,
1437 cs_math_surftri(xv + 3*v1, xv + 3*v2, pfq.center),
1438 cx->func, cx->input, val_i);
1439
1440 } /* Loop on edges */
1441
1442 } /* Switch on the type of face. Special case for triangles */
1443
1444 /* Compute the average */
1445 const double _os = 1./pfq.meas;
1446 for (int k = 0; k < dim; k++)
1447 val_i[k] *= _os;
1448
1449 } /* Loop on selected faces */
1450
1451 }
1452 }
1453
1454 /*----------------------------------------------------------------------------*/
1455
1456 #undef _dp3
1457
1458 END_C_DECLS
1459