1 /*============================================================================
2 * Main structure for a nodal representation associated with a mesh
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 <stdio.h>
35 #include <stdlib.h>
36 #include <string.h>
37
38 /*----------------------------------------------------------------------------
39 * Local headers
40 *----------------------------------------------------------------------------*/
41
42 #include "bft_mem.h"
43 #include "bft_printf.h"
44 #include "cs_mesh.h"
45
46 #include "fvm_defs.h"
47 #include "fvm_io_num.h"
48 #include "fvm_tesselation.h"
49
50 #include "cs_parall.h"
51
52 /*----------------------------------------------------------------------------
53 * Header for the current file
54 *----------------------------------------------------------------------------*/
55
56 #include "fvm_nodal.h"
57 #include "fvm_nodal_priv.h"
58
59 /*----------------------------------------------------------------------------*/
60
61 BEGIN_C_DECLS
62
63 /*=============================================================================
64 * Additional doxygen documentation
65 *============================================================================*/
66
67 /*!
68 \file fvm_nodal.c
69 Main structure for a nodal representation associated with a mesh.
70
71 The \c fvm_nodal_t structure is mostly used to handle post-processing
72 output in external formats, usually running in parallel using MPI.
73 This implies reconstructing a nodal connectivity (cells -> vertices) from
74 the main faces -> cells connectivity (using an intermediate
75 cells -> faces representation). It is also possible to directly
76 pass a nodal connectivity to such a structure.
77
78 So as to limit memory usage and avoid unecessary copies, the
79 \c fvm_nodal_t structure associated to a mesh defined by its nodal
80 connectivity is construed as a "view" on a mesh, and owns as little
81 data as possible. As such, most main structures associated with
82 this representation are defined by 2 arrays, one private, and one
83 shared. For example, an fvm_nodal_t structure has 2 coordinate
84 arrays:
85
86 <tt>const *cs_coord_t vertex_coords;</tt>
87
88 <tt>*cs_coord_t _vertex_coords;</tt>
89
90 If coordinates are shared with the calling code (and owned by that
91 code), we have <tt>_vertex_coords = NULL</tt>, and \c vertex_coords
92 points to the array passed from the calling code. If the
93 coordinates array belongs to FVM (either having been "given" by
94 the calling code, or generated by an operation wich invalidates
95 coorindates sharing with the parent mesh, such as mesh extrusion),
96 we have <tt>vertex_coords = _vertex_coords</tt>, which points to the
97 private array.
98
99 When an fvm_nodal_t object is destroyed, it destroys its private
100 data and frees the corresponding memory, but obviously leaves its
101 shared data untouched.
102
103 If a \c fvm_nodal_t structure B is built from a structure A with
104 which it shares its data, and a second \c fvm_nodal_t mesh C
105 is a view on B (for example a restriction to a part of the domain
106 that we whish to post-process more frequently), C must be destroyed
107 before B, which must be destroyed before A.
108 FVM does not use reference counters or a similar mechanism, so
109 good managment of object lifecycles is of the responsiblity of
110 the calling code. In practice, this logic is simple enough that
111 no issues have been encountered with this model so far in the
112 intended uses of the code.
113
114 Another associated concept is that of "parent_number": if a mesh
115 constitutes a "view" on another, sur un autre, a list of
116 parent entity numbers allows accessing variables associated
117 with the "parent" mesh, without needing to extract or duplicate
118 values at the level of the calling code.
119
120 Note that an FVM structure is global from an MPI point of
121 view, as it may participate in collective parallel operations.
122 Thus, if an FVM mesh represents a subset of a global domain,
123 it may very well be empty for some ranks, but it must still exist
124 on all ranks of the main communicator associated with FVM.
125
126 For parallelism, another important concept is that of "global numbering",
127 corresponding to entity numbers in a "serial" or "global" version
128 of an object: two entities (vertices, elements, ...) on different
129 processors having a same global number correspond to the same
130 "absolute" object. Except when explicitely mentioned, all other
131 data defining an object is based on local definitions and numberings.
132 Parallelism is thus made quite "transparent" to the calling code.
133 */
134
135 /*============================================================================
136 * Static global variables
137 *============================================================================*/
138
139 /* Number of vertices associated with each "nodal" element type */
140
141 const int fvm_nodal_n_vertices_element[] = {2, /* Edge */
142 3, /* Triangle */
143 4, /* Quadrangle */
144 0, /* Simple polygon */
145 4, /* Tetrahedron */
146 5, /* Pyramid */
147 6, /* Prism */
148 8, /* Hexahedron */
149 0}; /* Simple polyhedron */
150
151 /* Number of vertices associated with each "nodal" element type */
152
153 static int fvm_nodal_n_edges_element[] = {1, /* Edge */
154 3, /* Triangle */
155 4, /* Quadrangle */
156 0, /* Simple polygon */
157 6, /* Tetrahedron */
158 8, /* Pyramid */
159 9, /* Prism */
160 12, /* Hexahedron */
161 0}; /* Simple polyhedron */
162
163 /*! \cond DOXYGEN_SHOULD_SKIP_THIS */
164
165 /*============================================================================
166 * Private function definitions
167 *============================================================================*/
168
169 /*----------------------------------------------------------------------------
170 * Compare edges (qsort function).
171 *
172 * parameters:
173 * x <-> pointer to first edge definition
174 * y <-> pointer to second edge definition
175 *
176 * returns:
177 * result of strcmp() on group names
178 *----------------------------------------------------------------------------*/
179
_compare_edges(const void * x,const void * y)180 static int _compare_edges(const void *x, const void *y)
181 {
182 int retval = 1;
183
184 const cs_lnum_t *e0 = x;
185 const cs_lnum_t *e1 = y;
186
187 if (e0[0] < e1[0])
188 retval = -1;
189
190 else if (e0[0] == e1[0]) {
191 if (e0[1] < e1[1])
192 retval = -1;
193 else if (e0[1] == e1[1])
194 retval = 0;
195 }
196
197 return retval;
198 }
199
200 /*----------------------------------------------------------------------------
201 * Copy a nodal mesh section representation structure, sharing arrays
202 * with the original structure.
203 *
204 * parameters:
205 * this_section <-> pointer to structure that should be copied
206 *
207 * returns:
208 * pointer to created nodal mesh section representation structure
209 *----------------------------------------------------------------------------*/
210
211 static fvm_nodal_section_t *
_fvm_nodal_section_copy(const fvm_nodal_section_t * this_section)212 _fvm_nodal_section_copy(const fvm_nodal_section_t *this_section)
213 {
214 fvm_nodal_section_t *new_section;
215
216 BFT_MALLOC(new_section, 1, fvm_nodal_section_t);
217
218 /* Global information */
219
220 new_section->entity_dim = this_section->entity_dim;
221
222 new_section->n_elements = this_section->n_elements;
223 new_section->type = this_section->type;
224 new_section->boundary_flag = this_section->boundary_flag;
225
226 /* Connectivity */
227
228 new_section->connectivity_size = this_section->connectivity_size;
229 new_section->stride = this_section->stride;
230
231 new_section->n_faces = this_section->n_faces;
232
233 new_section->face_index = this_section->face_index;
234 new_section->face_num = this_section->face_num;
235 new_section->vertex_index = this_section->vertex_index;
236 new_section->vertex_num = this_section->vertex_num;
237
238 new_section->_face_index = NULL;
239 new_section->_face_num = NULL;
240 new_section->_vertex_index = NULL;
241 new_section->_vertex_num = NULL;
242
243 new_section->gc_id = NULL;
244 new_section->tag = NULL;
245
246 new_section->tesselation = NULL; /* TODO: copy tesselation */
247
248 /* Numbering */
249 /*-----------*/
250
251 new_section->parent_element_num = this_section->parent_element_num;
252 new_section->_parent_element_num = NULL;
253
254 if (this_section->global_element_num != NULL) {
255 cs_lnum_t n_ent
256 = fvm_io_num_get_local_count(this_section->global_element_num);
257 cs_gnum_t global_count
258 = fvm_io_num_get_global_count(this_section->global_element_num);
259 const cs_gnum_t *global_num
260 = fvm_io_num_get_global_num(this_section->global_element_num);
261
262 new_section->global_element_num
263 = fvm_io_num_create_shared(global_num, global_count, n_ent);
264 }
265 else
266 new_section->global_element_num = NULL;
267
268 return (new_section);
269 }
270
271 /*----------------------------------------------------------------------------
272 * Reduction of a nodal mesh representation section: only the associations
273 * (numberings) necessary to redistribution of fields for output are
274 * conserved, the full connectivity being no longer useful once it has been
275 * output.
276 *
277 * parameters:
278 * this_section <-> pointer to structure that should be reduced
279 *
280 * returns:
281 * true if connectivity has been reduced
282 *----------------------------------------------------------------------------*/
283
284 static bool
_fvm_nodal_section_reduce(fvm_nodal_section_t * this_section)285 _fvm_nodal_section_reduce(fvm_nodal_section_t * this_section)
286 {
287 bool retval = false;
288
289 /* If we have a tesselation of polyhedra (face index != NULL),
290 we may need to keep the connectivity information, to
291 interpolate nodal values to added vertices */
292
293 if ( this_section->tesselation == NULL
294 || this_section->_face_index == NULL) {
295
296 /* Connectivity */
297
298 this_section->connectivity_size = 0;
299
300 if (this_section->_face_index != NULL)
301 BFT_FREE(this_section->_face_index);
302 this_section->face_index = NULL;
303
304 if (this_section->_face_num != NULL)
305 BFT_FREE(this_section->_face_num);
306 this_section->face_num = NULL;
307
308 if (this_section->_vertex_index != NULL)
309 BFT_FREE(this_section->_vertex_index);
310 this_section->vertex_index = NULL;
311
312 if (this_section->_vertex_num != NULL)
313 BFT_FREE(this_section->_vertex_num);
314 this_section->vertex_num = NULL;
315
316 retval = true;
317 }
318
319 if (this_section->gc_id != NULL)
320 BFT_FREE(this_section->gc_id);
321
322 if (this_section->tag != NULL)
323 BFT_FREE(this_section->tag);
324
325 if (this_section->tesselation != NULL)
326 fvm_tesselation_reduce(this_section->tesselation);
327
328 return retval;
329 }
330
331 /*----------------------------------------------------------------------------
332 * Change entity parent numbering; this is useful when entities of the
333 * parent mesh have been renumbered after a nodal mesh representation
334 * structure's creation. As the parent_num[] array is defined only when
335 * non trivial (i.e. not 1, 2, ..., n), it may be allocated or freed
336 * by this function. The return argument corresponds to the new
337 * pointer which should replace the parent_num input argument.
338 *
339 * parameters:
340 * parent_num_size <-- size of local parent numbering array
341 * new_parent_num <-- pointer to local parent renumbering array
342 * ({1, ..., n} <-- {1, ..., n})
343 * parent_num <-> pointer to local parent numbering array
344 * _parent_num <-> pointer to local parent numbering array if
345 * owner, NULL otherwise
346 *
347 * returns:
348 * pointer to resulting parent_num[] array
349 *----------------------------------------------------------------------------*/
350
351 static cs_lnum_t *
_renumber_parent_num(cs_lnum_t parent_num_size,const cs_lnum_t new_parent_num[],const cs_lnum_t parent_num[],cs_lnum_t _parent_num[])352 _renumber_parent_num(cs_lnum_t parent_num_size,
353 const cs_lnum_t new_parent_num[],
354 const cs_lnum_t parent_num[],
355 cs_lnum_t _parent_num[])
356 {
357 int i;
358 cs_lnum_t old_num_id;
359 cs_lnum_t *parent_num_p = _parent_num;
360 bool trivial = true;
361
362 if (parent_num_size > 0 && new_parent_num != NULL) {
363
364 if (parent_num_p != NULL) {
365 for (i = 0; i < parent_num_size; i++) {
366 old_num_id = parent_num_p[i] - 1;
367 parent_num_p[i] = new_parent_num[old_num_id];
368 if (parent_num_p[i] != i+1)
369 trivial = false;
370 }
371 }
372 else {
373 BFT_MALLOC(parent_num_p, parent_num_size, cs_lnum_t);
374 if (parent_num != NULL) {
375 for (i = 0; i < parent_num_size; i++) {
376 old_num_id = parent_num[i] - 1;
377 parent_num_p[i] = new_parent_num[old_num_id];
378 if (parent_num_p[i] != i+1)
379 trivial = false;
380 }
381 }
382 else {
383 for (i = 0; i < parent_num_size; i++) {
384 parent_num_p[i] = new_parent_num[i];
385 if (parent_num_p[i] != i+1)
386 trivial = false;
387 }
388 }
389 }
390 }
391
392 if (trivial == true)
393 BFT_FREE(parent_num_p);
394
395 return parent_num_p;
396 }
397
398 /*----------------------------------------------------------------------------
399 * Renumber vertices based on those actually referenced, and update
400 * connectivity arrays and parent numbering in accordance.
401 *
402 * The number of vertices assigned to the nodal mesh (this_nodal->n_vertices)
403 * is computed and set by this function. If this number was previously
404 * non-zero (i.e. vertices have already been assigned to the structure),
405 * those vertices are considered as referenced. This is useful if we want
406 * to avoid discarding a given set of vertices, such as when building a
407 * nodal mesh representation containing only vertices.
408 *
409 * parameters:
410 * this_nodal <-> nodal mesh structure
411 *----------------------------------------------------------------------------*/
412
413 static void
_renumber_vertices(fvm_nodal_t * this_nodal)414 _renumber_vertices(fvm_nodal_t *this_nodal)
415 {
416 size_t i;
417 int section_id;
418 cs_lnum_t j;
419 cs_lnum_t vertex_id;
420 cs_lnum_t n_vertices;
421 fvm_nodal_section_t *section;
422
423 cs_lnum_t *loc_vertex_num = NULL;
424 cs_lnum_t max_vertex_num = 0;
425
426 /* Find maximum vertex reference */
427 /*-------------------------------*/
428
429 /* The mesh may already contain direct vertex references
430 (as in the case of a "mesh" only containing vertices) */
431
432 if (this_nodal->n_vertices > 0) {
433 if (this_nodal->parent_vertex_num != NULL) {
434 for (j = 0; j < this_nodal->n_vertices; j++) {
435 if (this_nodal->parent_vertex_num[j] > max_vertex_num)
436 max_vertex_num = this_nodal->parent_vertex_num[j];
437 }
438 }
439 else
440 max_vertex_num = this_nodal->n_vertices;
441 }
442
443 /* In most cases, the mesh will reference vertices through elements */
444
445 for (section_id = 0; section_id < this_nodal->n_sections; section_id++) {
446 section = this_nodal->sections[section_id];
447 if (this_nodal->parent_vertex_num != NULL) {
448 for (i = 0; i < section->connectivity_size; i++) {
449 cs_lnum_t vertex_num
450 = this_nodal->parent_vertex_num[section->vertex_num[i] - 1];
451 if (vertex_num > max_vertex_num)
452 max_vertex_num = vertex_num;
453 }
454 }
455 else {
456 for (i = 0; i < section->connectivity_size; i++) {
457 if (section->vertex_num[i] > max_vertex_num)
458 max_vertex_num = section->vertex_num[i];
459 }
460 }
461 }
462
463 /* Flag referenced vertices and compute size */
464 /*-------------------------------------------*/
465
466 BFT_MALLOC(loc_vertex_num, max_vertex_num, cs_lnum_t);
467
468 for (vertex_id = 0; vertex_id < max_vertex_num; vertex_id++)
469 loc_vertex_num[vertex_id] = 0;
470
471 if (this_nodal->n_vertices > 0) {
472 if (this_nodal->parent_vertex_num != NULL) {
473 for (j = 0; j < this_nodal->n_vertices; j++) {
474 vertex_id = this_nodal->parent_vertex_num[j] - 1;
475 if (loc_vertex_num[vertex_id] == 0)
476 loc_vertex_num[vertex_id] = 1;
477 }
478 }
479 else {
480 for (j = 0; j < this_nodal->n_vertices; j++) {
481 if (loc_vertex_num[j] == 0)
482 loc_vertex_num[j] = 1;
483 }
484 }
485 }
486
487 for (section_id = 0; section_id < this_nodal->n_sections; section_id++) {
488 section = this_nodal->sections[section_id];
489 if (this_nodal->parent_vertex_num != NULL) {
490 for (i = 0; i < section->connectivity_size; i++) {
491 vertex_id
492 = this_nodal->parent_vertex_num[section->vertex_num[i] - 1] - 1;
493 if (loc_vertex_num[vertex_id] == 0)
494 loc_vertex_num[vertex_id] = 1;
495 }
496 }
497 else {
498 for (i = 0; i < section->connectivity_size; i++) {
499 vertex_id = section->vertex_num[i] - 1;
500 if (loc_vertex_num[vertex_id] == 0)
501 loc_vertex_num[vertex_id] = 1;
502 }
503 }
504 }
505
506 /* Build vertices renumbering */
507 /*----------------------------*/
508
509 n_vertices = 0;
510
511 for (vertex_id = 0; vertex_id < max_vertex_num; vertex_id++) {
512 if (loc_vertex_num[vertex_id] == 1) {
513 n_vertices += 1;
514 loc_vertex_num[vertex_id] = n_vertices;
515 }
516 }
517 this_nodal->n_vertices = n_vertices;
518
519 /* Update connectivity and vertex parent numbering */
520 /*-------------------------------------------------*/
521
522 /* If all vertices are flagged, no need to renumber */
523
524 if (n_vertices == max_vertex_num)
525 BFT_FREE(loc_vertex_num);
526
527 else {
528
529 /* Update connectivity */
530
531 for (section_id = 0; section_id < this_nodal->n_sections; section_id++) {
532 section = this_nodal->sections[section_id];
533 if (section->_vertex_num == NULL)
534 fvm_nodal_section_copy_on_write(section, false, false, false, true);
535 if (this_nodal->parent_vertex_num != NULL) {
536 for (i = 0; i < section->connectivity_size; i++) {
537 vertex_id
538 = this_nodal->parent_vertex_num[section->vertex_num[i] - 1] - 1;
539 section->_vertex_num[i] = loc_vertex_num[vertex_id];
540 }
541 }
542 else {
543 for (i = 0; i < section->connectivity_size; i++) {
544 vertex_id = section->vertex_num[i] - 1;
545 section->_vertex_num[i] = loc_vertex_num[vertex_id];
546 }
547 }
548 }
549
550 /* Build or update vertex parent numbering */
551
552 this_nodal->parent_vertex_num = NULL;
553 if (this_nodal->_parent_vertex_num != NULL)
554 BFT_FREE(this_nodal->_parent_vertex_num);
555
556 if (loc_vertex_num != NULL) {
557 BFT_MALLOC(this_nodal->_parent_vertex_num, n_vertices, cs_lnum_t);
558 for (vertex_id = 0; vertex_id < max_vertex_num; vertex_id++) {
559 if (loc_vertex_num[vertex_id] > 0)
560 this_nodal->_parent_vertex_num[loc_vertex_num[vertex_id] - 1]
561 = vertex_id + 1;
562 }
563 this_nodal->parent_vertex_num = this_nodal->_parent_vertex_num;
564 }
565 }
566
567 /* Free renumbering array */
568
569 BFT_FREE(loc_vertex_num);
570 }
571
572 /*----------------------------------------------------------------------------
573 * Dump printout of a nodal representation structure section.
574 *
575 * parameters:
576 * this_section <-- pointer to structure that should be dumped
577 *----------------------------------------------------------------------------*/
578
579 static void
_fvm_nodal_section_dump(const fvm_nodal_section_t * this_section)580 _fvm_nodal_section_dump(const fvm_nodal_section_t *this_section)
581 {
582 cs_lnum_t n_elements, i, j;
583 const cs_lnum_t *idx, *num;
584
585 /* Global indicators */
586 /*--------------------*/
587
588 bft_printf("\n"
589 "Entity dimension: %d\n"
590 "Number of elements: %ld\n"
591 "Element type: %s\n"
592 "Boundary flag: %d\n",
593 this_section->entity_dim, (long)this_section->n_elements,
594 fvm_elements_type_name[this_section->type],
595 this_section->boundary_flag);
596
597 bft_printf("\n"
598 "Connectivity_size: %llu\n"
599 "Stride: %d\n"
600 "Number of faces: %ld\n",
601 (unsigned long long)(this_section->connectivity_size),
602 this_section->stride,
603 (long)(this_section->n_faces));
604
605 bft_printf("\n"
606 "Pointers to shareable arrays:\n"
607 " face_index: %p\n"
608 " face_num: %p\n"
609 " vertex_index: %p\n"
610 " vertex_num: %p\n"
611 " parent_element_num: %p\n",
612 (const void *)this_section->face_index,
613 (const void *)this_section->face_num,
614 (const void *)this_section->vertex_index,
615 (const void *)this_section->vertex_num,
616 (const void *)this_section->parent_element_num);
617
618 bft_printf("\n"
619 "Pointers to local arrays:\n"
620 " _face_index: %p\n"
621 " _face_num: %p\n"
622 " _vertex_index: %p\n"
623 " _vertex_num: %p\n"
624 " _parent_element_num: %p\n"
625 " gc_id: %p\n"
626 " tag: %p\n",
627 (const void *)this_section->_face_index,
628 (const void *)this_section->_face_num,
629 (const void *)this_section->_vertex_index,
630 (const void *)this_section->_vertex_num,
631 (const void *)this_section->_parent_element_num,
632 (const void *)this_section->gc_id,
633 (const void *)this_section->tag);
634
635 if (this_section->face_index != NULL) {
636 bft_printf("\nPolyhedra -> faces (polygons) connectivity:\n\n");
637 n_elements = this_section->n_elements;
638 idx = this_section->face_index;
639 num = this_section->face_num;
640 for (i = 0; i < n_elements; i++) {
641 bft_printf("%10ld (idx = %10ld) %10ld\n",
642 (long)i+1, (long)idx[i], (long)num[idx[i]]);
643 for (j = idx[i] + 1; j < idx[i + 1]; j++)
644 bft_printf(" %10ld\n", (long)num[j]);
645 }
646 bft_printf(" end (idx = %10ld)\n", (long)idx[n_elements]);
647 }
648
649 if (this_section->vertex_index != NULL) {
650 cs_lnum_t n_faces = (this_section->n_faces > 0) ?
651 this_section->n_faces : this_section->n_elements;
652 bft_printf("\nPolygons -> vertices connectivity:\n\n");
653 idx = this_section->vertex_index;
654 num = this_section->vertex_num;
655 for (i = 0; i < n_faces; i++) {
656 bft_printf("%10ld (idx = %10ld) %10ld\n",
657 (long)i + 1, (long)idx[i], (long)num[idx[i]]);
658 for (j = idx[i] + 1; j < idx[i + 1]; j++)
659 bft_printf(" %10ld\n", (long)num[j]);
660 }
661 bft_printf(" end (idx = %10ld)\n", (long)idx[n_faces]);
662 }
663
664 else {
665 bft_printf("\nElements -> vertices connectivity:\n\n");
666 n_elements = this_section->n_elements;
667 num = this_section->vertex_num;
668 switch(this_section->stride) {
669 case 2:
670 for (i = 0; i < n_elements; i++)
671 bft_printf("%10ld : %10ld %10ld\n",
672 (long)i+1, (long)num[i*2], (long)num[i*2+1]);
673 break;
674 case 3:
675 for (i = 0; i < n_elements; i++)
676 bft_printf("%10ld : %10ld %10ld %10ld\n",
677 (long)i+1, (long)num[i*3], (long)num[i*3+1],
678 (long)num[i*3+2]);
679 break;
680 case 4:
681 for (i = 0; i < n_elements; i++)
682 bft_printf("%10ld : %10ld %10ld %10ld %10ld\n",
683 (long)i+1, (long)num[i*4], (long)num[i*4+1],
684 (long)num[i*4+2], (long)num[i*4+3]);
685 break;
686 case 5:
687 for (i = 0; i < n_elements; i++)
688 bft_printf("%10ld : %10ld %10ld %10ld %10ld %10ld\n",
689 (long)i+1, (long)num[i*5], (long)num[i*5+1],
690 (long)num[i*5+2], (long)num[i*5+3], (long)num[i*5+4]);
691 break;
692 case 6:
693 for (i = 0; i < n_elements; i++)
694 bft_printf("%10ld : %10ld %10ld %10ld %10ld %10ld %10ld\n",
695 (long)i+1, (long)num[i*6], (long)num[i*6+1],
696 (long)num[i*6+2], (long)num[i*6+3], (long)num[i*6+4],
697 (long)num[i*6+5]);
698 break;
699 case 8:
700 for (i = 0; i < n_elements; i++)
701 bft_printf("%10ld : %10ld %10ld %10ld %10ld %10ld %10ld %10ld %10ld\n",
702 (long)i+1, (long)num[i*8], (long)num[i*8+1],
703 (long)num[i*8+2], (long)num[i*8+3], (long)num[i*8+4],
704 (long)num[i*8+5], (long)num[i*8+6], (long)num[i*8+7]);
705 break;
706 default:
707 for (i = 0; i < n_elements; i++) {
708 bft_printf("%10ld :", (long)i+1);
709 for (j = 0; j < this_section->stride; j++)
710 bft_printf(" %10ld", (long)num[i*this_section->stride + j]);
711 bft_printf("\n");
712 }
713 }
714 }
715
716 if (this_section->gc_id != NULL) {
717 bft_printf("\nGroup class ids:\n\n");
718 for (i = 0; i < this_section->n_elements; i++)
719 bft_printf("%10ld : %10d\n", (long)i+1, this_section->gc_id[i]);
720 bft_printf("\n");
721 }
722
723 if (this_section->tag != NULL) {
724 bft_printf("\nTags:\n\n");
725 for (i = 0; i < this_section->n_elements; i++)
726 bft_printf("%10ld : %10d\n", (long)i+1, this_section->tag[i]);
727 bft_printf("\n");
728 }
729
730 /* Faces tesselation */
731
732 if (this_section->tesselation != NULL)
733 fvm_tesselation_dump(this_section->tesselation);
734
735 /* Numbers of associated elements in the parent mesh */
736
737 bft_printf("\nLocal element numbers in parent mesh:\n");
738 if (this_section->parent_element_num == NULL)
739 bft_printf("\n Nil\n\n");
740 else {
741 for (i = 0; i < this_section->n_elements; i++)
742 bft_printf(" %10ld %10ld\n", (long)i+1,
743 (long)this_section->parent_element_num[i]);
744 }
745
746 /* Global element numbers (only for parallel execution) */
747
748 if (this_section->global_element_num != NULL) {
749 bft_printf("\nGlobal element numbers:\n");
750 fvm_io_num_dump(this_section->global_element_num);
751 }
752 }
753
754 /*----------------------------------------------------------------------------
755 * Remove global vertex labels from a nodal mesh.
756 *
757 * parameters:
758 * this_nodal <-> nodal mesh structure
759 *----------------------------------------------------------------------------*/
760
761 static void
_remove_global_vertex_labels(fvm_nodal_t * this_nodal)762 _remove_global_vertex_labels(fvm_nodal_t *this_nodal)
763 {
764 assert(this_nodal != NULL);
765
766 if (this_nodal->global_vertex_labels != NULL) {
767
768 cs_gnum_t n_g_vertices = fvm_nodal_n_g_vertices(this_nodal);
769
770 for (cs_gnum_t i = 0; i < n_g_vertices; i++)
771 BFT_FREE(this_nodal->global_vertex_labels[i]);
772
773 BFT_FREE(this_nodal->global_vertex_labels);
774
775 }
776 }
777
778 /*============================================================================
779 * Semi-private function definitions (prototypes in fvm_nodal_priv.h)
780 *============================================================================*/
781
782 /*----------------------------------------------------------------------------
783 * Creation of a nodal mesh section representation structure.
784 *
785 * parameters:
786 * type <-- type of element defined by this section
787 *
788 * returns:
789 * pointer to created nodal mesh section representation structure
790 *----------------------------------------------------------------------------*/
791
792 fvm_nodal_section_t *
fvm_nodal_section_create(const fvm_element_t type)793 fvm_nodal_section_create(const fvm_element_t type)
794 {
795 fvm_nodal_section_t *this_section;
796
797 BFT_MALLOC(this_section, 1, fvm_nodal_section_t);
798
799 /* Global information */
800
801 if (type == FVM_EDGE)
802 this_section->entity_dim = 1;
803 else if (type >= FVM_FACE_TRIA && type <= FVM_FACE_POLY)
804 this_section->entity_dim = 2;
805 else
806 this_section->entity_dim = 3;
807
808 this_section->n_elements = 0;
809 this_section->type = type;
810
811 this_section->boundary_flag = -1;
812
813 /* Connectivity */
814
815 this_section->connectivity_size = 0;
816
817 if (type != FVM_FACE_POLY && type != FVM_CELL_POLY)
818 this_section->stride = fvm_nodal_n_vertices_element[type];
819 else
820 this_section->stride = 0;
821
822 this_section->n_faces = 0;
823 this_section->face_index = NULL;
824 this_section->face_num = NULL;
825 this_section->vertex_index = NULL;
826 this_section->vertex_num = NULL;
827
828 this_section->_face_index = NULL;
829 this_section->_face_num = NULL;
830 this_section->_vertex_index = NULL;
831 this_section->_vertex_num = NULL;
832
833 this_section->gc_id = NULL;
834 this_section->tag = NULL;
835
836 this_section->tesselation = NULL;
837
838 /* Numbering */
839 /*-----------*/
840
841 this_section->parent_element_num = NULL;
842 this_section->_parent_element_num = NULL;
843
844 this_section->global_element_num = NULL;
845
846 return (this_section);
847 }
848
849 /*----------------------------------------------------------------------------
850 * Destruction of a nodal mesh section representation structure.
851 *
852 * parameters:
853 * this_section <-> pointer to structure that should be destroyed
854 *
855 * returns:
856 * NULL pointer
857 *----------------------------------------------------------------------------*/
858
859 fvm_nodal_section_t *
fvm_nodal_section_destroy(fvm_nodal_section_t * this_section)860 fvm_nodal_section_destroy(fvm_nodal_section_t * this_section)
861 {
862 /* Connectivity */
863
864 if (this_section->_face_index != NULL)
865 BFT_FREE(this_section->_face_index);
866 if (this_section->_face_num != NULL)
867 BFT_FREE(this_section->_face_num);
868
869 if (this_section->_vertex_index != NULL)
870 BFT_FREE(this_section->_vertex_index);
871 if (this_section->_vertex_num != NULL)
872 BFT_FREE(this_section->_vertex_num);
873
874 if (this_section->gc_id != NULL)
875 BFT_FREE(this_section->gc_id);
876
877 if (this_section->tag != NULL)
878 BFT_FREE(this_section->tag);
879
880 if (this_section->tesselation != NULL)
881 fvm_tesselation_destroy(this_section->tesselation);
882
883 /* Numbering */
884 /*-----------*/
885
886 if (this_section->parent_element_num != NULL) {
887 this_section->parent_element_num = NULL;
888 BFT_FREE(this_section->_parent_element_num);
889 }
890
891 if (this_section->global_element_num != NULL)
892 fvm_io_num_destroy(this_section->global_element_num);
893
894 /* Main structure destroyed and NULL returned */
895
896 BFT_FREE(this_section);
897
898 return (this_section);
899 }
900
901 /*----------------------------------------------------------------------------
902 * Copy selected shared connectivity information to private connectivity
903 * for a nodal mesh section.
904 *
905 * parameters:
906 * this_section <-> pointer to section structure
907 * copy_face_index <-- copy face index (polyhedra only) ?
908 * copy_face_num <-- copy face numbers (polyhedra only) ?
909 * copy_vertex_index <-- copy vertex index (polyhedra/polygons only) ?
910 * copy_vertex_num <-- copy vertex numbers ?
911 *----------------------------------------------------------------------------*/
912
913 void
fvm_nodal_section_copy_on_write(fvm_nodal_section_t * this_section,bool copy_face_index,bool copy_face_num,bool copy_vertex_index,bool copy_vertex_num)914 fvm_nodal_section_copy_on_write(fvm_nodal_section_t *this_section,
915 bool copy_face_index,
916 bool copy_face_num,
917 bool copy_vertex_index,
918 bool copy_vertex_num)
919 {
920 cs_lnum_t n_faces;
921 size_t i;
922
923 if (copy_face_index == true
924 && this_section->face_index != NULL && this_section->_face_index == NULL) {
925 BFT_MALLOC(this_section->_face_index, this_section->n_elements + 1, cs_lnum_t);
926 for (i = 0; i < (size_t)(this_section->n_elements + 1); i++) {
927 this_section->_face_index[i] = this_section->face_index[i];
928 }
929 this_section->face_index = this_section->_face_index;
930 }
931
932 if (copy_face_num == true
933 && this_section->face_num != NULL && this_section->_face_num == NULL) {
934 n_faces = this_section->face_index[this_section->n_elements];
935 BFT_MALLOC(this_section->_face_num, n_faces, cs_lnum_t);
936 for (i = 0; i < (size_t)n_faces; i++) {
937 this_section->_face_num[i] = this_section->face_num[i];
938 }
939 this_section->face_num = this_section->_face_num;
940 }
941
942 if ( copy_vertex_index == true
943 && this_section->vertex_index != NULL
944 && this_section->_vertex_index == NULL) {
945 if (this_section->n_faces != 0)
946 n_faces = this_section->n_faces;
947 else
948 n_faces = this_section->n_elements;
949 BFT_MALLOC(this_section->_vertex_index, n_faces + 1, cs_lnum_t);
950 for (i = 0; i < (size_t)n_faces + 1; i++) {
951 this_section->_vertex_index[i] = this_section->vertex_index[i];
952 }
953 this_section->vertex_index = this_section->_vertex_index;
954 }
955
956 if (copy_vertex_num == true && this_section->_vertex_num == NULL) {
957 BFT_MALLOC(this_section->_vertex_num,
958 this_section->connectivity_size, cs_lnum_t);
959 for (i = 0; i < this_section->connectivity_size; i++) {
960 this_section->_vertex_num[i] = this_section->vertex_num[i];
961 }
962 this_section->vertex_num = this_section->_vertex_num;
963 }
964
965 }
966
967 /*----------------------------------------------------------------------------
968 * Return global number of elements associated with section.
969 *
970 * parameters:
971 * this_section <-- pointer to section structure
972 *
973 * returns:
974 * global number of elements associated with section
975 *----------------------------------------------------------------------------*/
976
977 cs_gnum_t
fvm_nodal_section_n_g_elements(const fvm_nodal_section_t * this_section)978 fvm_nodal_section_n_g_elements(const fvm_nodal_section_t *this_section)
979 {
980 if (this_section->global_element_num != NULL)
981 return fvm_io_num_get_global_count(this_section->global_element_num);
982 else
983 return this_section->n_elements;
984 }
985
986 /*----------------------------------------------------------------------------
987 * Return global number of vertices associated with nodal mesh.
988 *
989 * parameters:
990 * this_nodal <-- pointer to nodal mesh structure
991 *
992 * returns:
993 * global number of vertices associated with nodal mesh
994 *----------------------------------------------------------------------------*/
995
996 cs_gnum_t
fvm_nodal_n_g_vertices(const fvm_nodal_t * this_nodal)997 fvm_nodal_n_g_vertices(const fvm_nodal_t *this_nodal)
998 {
999 cs_gnum_t n_g_vertices;
1000
1001 if (this_nodal->global_vertex_num != NULL)
1002 n_g_vertices = fvm_io_num_get_global_count(this_nodal->global_vertex_num);
1003 else
1004 n_g_vertices = this_nodal->n_vertices;
1005
1006 return n_g_vertices;
1007 }
1008
1009 /*----------------------------------------------------------------------------
1010 * Define cell->face connectivity for strided cell types.
1011 *
1012 * parameters:
1013 * element_type <-- type of strided element
1014 * n_faces --> number of element faces
1015 * n_face_vertices --> number of vertices of each face
1016 * face_vertices --> face -> vertex base connectivity (0 to n-1)
1017 *----------------------------------------------------------------------------*/
1018
1019 void
fvm_nodal_cell_face_connect(fvm_element_t element_type,int * n_faces,int n_face_vertices[6],int face_vertices[6][4])1020 fvm_nodal_cell_face_connect(fvm_element_t element_type,
1021 int *n_faces,
1022 int n_face_vertices[6],
1023 int face_vertices[6][4])
1024 {
1025 int i, j;
1026
1027 /* Initialization */
1028
1029 *n_faces = 0;
1030
1031 for (i = 0; i < 6; i++) {
1032 n_face_vertices[i] = 0;
1033 for (j = 0; j < 4; j++)
1034 face_vertices[i][j] = 0;
1035 }
1036
1037 /* Define connectivity based on element type */
1038
1039 switch(element_type) {
1040
1041 case FVM_CELL_TETRA:
1042 {
1043 cs_lnum_t _face_vertices[4][3] = {{1, 3, 2}, /* x 4 */
1044 {1, 2, 4}, /* /|\ */
1045 {1, 4, 3}, /* / | \ */
1046 {2, 3, 4}}; /* / | \ */
1047 /* 1 x- -|- -x 3 */
1048 for (i = 0; i < 4; i++) { /* \ | / */
1049 n_face_vertices[i] = 3; /* \ | / */
1050 for (j = 0; j < 3; j++) /* \|/ */
1051 face_vertices[i][j] = _face_vertices[i][j]; /* x 2 */
1052 }
1053 *n_faces = 4;
1054 }
1055 break;
1056
1057 case FVM_CELL_PYRAM:
1058 {
1059 cs_lnum_t _n_face_vertices[5] = {3, 3, 3, 3, 4};
1060 cs_lnum_t _face_vertices[5][4] = {{1, 2, 5, 0}, /* 5 x */
1061 {1, 5, 4, 0}, /* /|\ */
1062 {2, 3, 5, 0}, /* //| \ */
1063 {3, 4, 5, 0}, /* // | \ */
1064 {1, 4, 3, 2}}; /* 4 x/--|---x 3 */
1065 /* // | / */
1066 for (i = 0; i < 5; i++) { /* // | / */
1067 n_face_vertices[i] = _n_face_vertices[i]; /* 1 x-------x 2 */
1068 for (j = 0; j < 4; j++)
1069 face_vertices[i][j] = _face_vertices[i][j];
1070 }
1071 *n_faces = 5;
1072 }
1073 break;
1074
1075 case FVM_CELL_PRISM:
1076 {
1077 cs_lnum_t _n_face_vertices[5] = {3, 3, 4, 4, 4};
1078 cs_lnum_t _face_vertices[5][4] = {{1, 3, 2, 0}, /* 4 x-------x 6 */
1079 {4, 5, 6, 0}, /* |\ /| */
1080 {1, 2, 5, 4}, /* | \ / | */
1081 {1, 4, 6, 3}, /* 1 x- \-/ -x 3 */
1082 {2, 3, 6, 5}}; /* \ 5x / */
1083 /* \ | / */
1084 for (i = 0; i < 5; i++) { /* \|/ */
1085 n_face_vertices[i] = _n_face_vertices[i]; /* x 2 */
1086 for (j = 0; j < 4; j++)
1087 face_vertices[i][j] = _face_vertices[i][j];
1088 }
1089 *n_faces = 5;
1090 }
1091 break;
1092
1093 case FVM_CELL_HEXA:
1094 {
1095 cs_lnum_t _n_face_vertices[6] = {4, 4, 4, 4, 4, 4};
1096 cs_lnum_t _face_vertices[6][4] = {{1, 4, 3, 2}, /* 8 x-------x 7 */
1097 {1, 2, 6, 5}, /* /| /| */
1098 {1, 5, 8, 4}, /* / | / | */
1099 {2, 3, 7, 6}, /* 5 x-------x6 | */
1100 {3, 4, 8, 7}, /* | 4x----|--x 3 */
1101 {5, 6, 7, 8}}; /* | / | / */
1102 for (i = 0; i < 6; i++) { /* |/ |/ */
1103 n_face_vertices[i] = _n_face_vertices[i]; /* 1 x-------x 2 */
1104 for (j = 0; j < 4; j++)
1105 face_vertices[i][j] = _face_vertices[i][j];
1106 }
1107 *n_faces = 6;
1108 }
1109 break;
1110
1111 default:
1112 assert(0);
1113 }
1114
1115 /* Switch from (1, n) to (0, n-1) numbering */
1116
1117 for (i = 0; i < 6; i++) {
1118 for (j = 0; j < 4; j++)
1119 face_vertices[i][j] -= 1;
1120 }
1121 }
1122
1123 /*! (DOXYGEN_SHOULD_SKIP_THIS) \endcond */
1124
1125 /*============================================================================
1126 * Public function definitions
1127 *============================================================================*/
1128
1129 /*----------------------------------------------------------------------------
1130 * Creation of a nodal mesh representation structure.
1131 *
1132 * parameters:
1133 * name <-- name that should be assigned to the nodal mesh
1134 * dim <-- spatial dimension
1135 *
1136 * returns:
1137 * pointer to created nodal mesh representation structure
1138 *----------------------------------------------------------------------------*/
1139
1140 fvm_nodal_t *
fvm_nodal_create(const char * name,int dim)1141 fvm_nodal_create(const char *name,
1142 int dim)
1143 {
1144 fvm_nodal_t *this_nodal;
1145
1146 BFT_MALLOC(this_nodal, 1, fvm_nodal_t);
1147
1148 /* Global indicators */
1149
1150 if (name != NULL) {
1151 BFT_MALLOC(this_nodal->name, strlen(name) + 1, char);
1152 strcpy(this_nodal->name, name);
1153 }
1154 else
1155 this_nodal->name = NULL;
1156
1157 this_nodal->dim = dim;
1158 this_nodal->num_dom = (cs_glob_rank_id >= 0) ? cs_glob_rank_id + 1 : 1;
1159 this_nodal->n_doms = cs_glob_n_ranks;
1160 this_nodal->n_sections = 0;
1161
1162 /* Local dimensions */
1163
1164 this_nodal->n_cells = 0;
1165 this_nodal->n_faces = 0;
1166 this_nodal->n_edges = 0;
1167 this_nodal->n_vertices = 0;
1168
1169 /* Local structures */
1170
1171 this_nodal->vertex_coords = NULL;
1172 this_nodal->_vertex_coords = NULL;
1173
1174 this_nodal->parent_vertex_num = NULL;
1175 this_nodal->_parent_vertex_num = NULL;
1176
1177 this_nodal->global_vertex_num = NULL;
1178
1179 this_nodal->sections = NULL;
1180
1181 this_nodal->gc_set = NULL;
1182
1183 this_nodal->global_vertex_labels = NULL;
1184
1185 this_nodal->parent = NULL;
1186
1187 return (this_nodal);
1188 }
1189
1190 /*----------------------------------------------------------------------------
1191 * Destruction of a nodal mesh representation structure.
1192 *
1193 * parameters:
1194 * this_nodal <-> pointer to structure that should be destroyed
1195 *
1196 * returns:
1197 * NULL pointer
1198 *----------------------------------------------------------------------------*/
1199
1200 fvm_nodal_t *
fvm_nodal_destroy(fvm_nodal_t * this_nodal)1201 fvm_nodal_destroy(fvm_nodal_t *this_nodal)
1202 {
1203 int i;
1204
1205 if (this_nodal == NULL)
1206 return NULL;
1207
1208 /* Local structures */
1209
1210 _remove_global_vertex_labels(this_nodal);
1211
1212 if (this_nodal->name != NULL)
1213 BFT_FREE(this_nodal->name);
1214
1215 if (this_nodal->_vertex_coords != NULL)
1216 BFT_FREE(this_nodal->_vertex_coords);
1217
1218 if (this_nodal->parent_vertex_num != NULL) {
1219 this_nodal->parent_vertex_num = NULL;
1220 BFT_FREE(this_nodal->_parent_vertex_num);
1221 }
1222
1223 if (this_nodal->global_vertex_num != NULL)
1224 fvm_io_num_destroy(this_nodal->global_vertex_num);
1225
1226 for (i = 0; i < this_nodal->n_sections; i++)
1227 fvm_nodal_section_destroy(this_nodal->sections[i]);
1228
1229 if (this_nodal->sections != NULL)
1230 BFT_FREE(this_nodal->sections);
1231
1232 if (this_nodal->gc_set != NULL)
1233 this_nodal->gc_set = fvm_group_class_set_destroy(this_nodal->gc_set);
1234
1235 /* Main structure destroyed and NULL returned */
1236
1237 BFT_FREE(this_nodal);
1238
1239 return (this_nodal);
1240 }
1241
1242 /*----------------------------------------------------------------------------
1243 * Copy a nodal mesh representation structure, sharing arrays with the
1244 * original structure.
1245 *
1246 * Element group classes and mesh group class descriptions are not currently
1247 * copied.
1248 *
1249 * parameters:
1250 * this_nodal <-> pointer to structure that should be copied
1251 *
1252 * returns:
1253 * pointer to created nodal mesh representation structure
1254 *----------------------------------------------------------------------------*/
1255
1256 fvm_nodal_t *
fvm_nodal_copy(const fvm_nodal_t * this_nodal)1257 fvm_nodal_copy(const fvm_nodal_t *this_nodal)
1258 {
1259 int i;
1260 fvm_nodal_t *new_nodal;
1261
1262 BFT_MALLOC(new_nodal, 1, fvm_nodal_t);
1263
1264 /* Global indicators */
1265
1266 if (this_nodal->name != NULL) {
1267 BFT_MALLOC(new_nodal->name, strlen(this_nodal->name) + 1, char);
1268 strcpy(new_nodal->name, this_nodal->name);
1269 }
1270 else
1271 new_nodal->name = NULL;
1272
1273 new_nodal->dim = this_nodal->dim;
1274 new_nodal->num_dom = this_nodal->num_dom;
1275 new_nodal->n_doms = this_nodal->n_doms;
1276 new_nodal->n_sections = this_nodal->n_sections;
1277
1278 /* Local dimensions */
1279
1280 new_nodal->n_cells = this_nodal->n_cells;
1281 new_nodal->n_faces = this_nodal->n_faces;
1282 new_nodal->n_edges = this_nodal->n_edges;
1283 new_nodal->n_vertices = this_nodal->n_vertices;
1284
1285 /* Local structures */
1286
1287 new_nodal->vertex_coords = this_nodal->vertex_coords;
1288 new_nodal->_vertex_coords = NULL;
1289
1290 new_nodal->parent_vertex_num = this_nodal->parent_vertex_num;
1291 new_nodal->_parent_vertex_num = NULL;
1292
1293 if (this_nodal->global_vertex_num != NULL) {
1294 cs_lnum_t n_ent
1295 = fvm_io_num_get_local_count(this_nodal->global_vertex_num);
1296 cs_gnum_t global_count
1297 = fvm_io_num_get_global_count(this_nodal->global_vertex_num);
1298 const cs_gnum_t *global_num
1299 = fvm_io_num_get_global_num(this_nodal->global_vertex_num);
1300
1301 new_nodal->global_vertex_num
1302 = fvm_io_num_create_shared(global_num, global_count, n_ent);
1303 }
1304 else
1305 new_nodal->global_vertex_num = NULL;
1306
1307 BFT_MALLOC(new_nodal->sections,
1308 new_nodal->n_sections,
1309 fvm_nodal_section_t *);
1310 for (i = 0; i < new_nodal->n_sections; i++)
1311 new_nodal->sections[i] = _fvm_nodal_section_copy(this_nodal->sections[i]);
1312
1313 new_nodal->gc_set = NULL;
1314 new_nodal->global_vertex_labels = NULL;
1315
1316 return (new_nodal);
1317 }
1318
1319 /*----------------------------------------------------------------------------
1320 * Reduction of a nodal mesh representation structure: only the associations
1321 * (numberings) necessary to redistribution of fields for output are
1322 * conserved, the full connectivity being in many cases no longer useful
1323 * once it has been output. If the del_vertex_num value is set
1324 * to true, vertex-based values may not be output in parallel mode
1325 * after this function is called.
1326 *
1327 * parameters:
1328 * this_nodal <-> pointer to structure that should be reduced
1329 * del_vertex_num <-- indicates if vertex parent indirection and
1330 * I/O numbering are destroyed (1) or not (0)
1331 *----------------------------------------------------------------------------*/
1332
1333 void
fvm_nodal_reduce(fvm_nodal_t * this_nodal,int del_vertex_num)1334 fvm_nodal_reduce(fvm_nodal_t *this_nodal,
1335 int del_vertex_num)
1336 {
1337 int i;
1338 bool reduce_vertices = true;
1339
1340 /* Connectivity */
1341
1342 for (i = 0; i < this_nodal->n_sections; i++) {
1343 if (_fvm_nodal_section_reduce(this_nodal->sections[i]) == false)
1344 reduce_vertices = false;
1345 }
1346
1347 /* Vertices */
1348
1349 if (reduce_vertices == true) {
1350
1351 if (this_nodal->_vertex_coords != NULL)
1352 BFT_FREE(this_nodal->_vertex_coords);
1353 this_nodal->vertex_coords = NULL;
1354
1355 }
1356
1357 /* Depending on this option, output on vertices may not remain possible */
1358
1359 if (del_vertex_num > 0) {
1360
1361 if (this_nodal->parent_vertex_num != NULL) {
1362 this_nodal->parent_vertex_num = NULL;
1363 BFT_FREE(this_nodal->_parent_vertex_num);
1364 }
1365
1366 if (this_nodal->global_vertex_num != NULL)
1367 this_nodal->global_vertex_num
1368 = fvm_io_num_destroy(this_nodal->global_vertex_num);
1369
1370 }
1371
1372 if (this_nodal->gc_set != NULL)
1373 this_nodal->gc_set = fvm_group_class_set_destroy(this_nodal->gc_set);
1374 }
1375
1376 /*----------------------------------------------------------------------------
1377 * Change entity parent numbering; this is useful when entities of the
1378 * parent mesh have been renumbered after a nodal mesh representation
1379 * structure's creation.
1380 *
1381 * parameters:
1382 * this_nodal <-- nodal mesh structure
1383 * new_parent_num <-- pointer to local parent renumbering array
1384 * ({1, ..., n} <-- {1, ..., n})
1385 * entity_dim <-- 3 for cells, 2 for faces, 1 for edges,
1386 * and 0 for vertices
1387 *----------------------------------------------------------------------------*/
1388
1389 void
fvm_nodal_change_parent_num(fvm_nodal_t * this_nodal,const cs_lnum_t new_parent_num[],int entity_dim)1390 fvm_nodal_change_parent_num(fvm_nodal_t *this_nodal,
1391 const cs_lnum_t new_parent_num[],
1392 int entity_dim)
1393 {
1394 /* Vertices */
1395
1396 if (entity_dim == 0) {
1397
1398 this_nodal->_parent_vertex_num
1399 = _renumber_parent_num(this_nodal->n_vertices,
1400 new_parent_num,
1401 this_nodal->parent_vertex_num,
1402 this_nodal->_parent_vertex_num);
1403 this_nodal->parent_vertex_num = this_nodal->_parent_vertex_num;
1404
1405 }
1406
1407 /* Other elements */
1408
1409 else {
1410
1411 int i = 0;
1412 fvm_nodal_section_t *section = NULL;
1413
1414 for (i = 0; i < this_nodal->n_sections; i++) {
1415 section = this_nodal->sections[i];
1416 if (section->entity_dim == entity_dim) {
1417 section->_parent_element_num
1418 = _renumber_parent_num(section->n_elements,
1419 new_parent_num,
1420 section->parent_element_num,
1421 section->_parent_element_num);
1422 section->parent_element_num = section->_parent_element_num;
1423 }
1424 }
1425
1426 }
1427
1428 }
1429
1430 /*----------------------------------------------------------------------------
1431 * Remove entity parent numbering; this is useful for example when we
1432 * want to assign coordinates or fields to an extracted mesh using
1433 * arrays relative to the mesh, and not to its parent.
1434 *
1435 * This is equivalent to calling fvm_nodal_change_parent_num(), with
1436 * 'trivial' (1 o n) new_parent_num[] values.
1437 *
1438 * parameters:
1439 * this_nodal <-- nodal mesh structure
1440 * entity_dim <-- 3 for cells, 2 for faces, 1 for edges,
1441 * and 0 for vertices
1442 *----------------------------------------------------------------------------*/
1443
1444 void
fvm_nodal_remove_parent_num(fvm_nodal_t * this_nodal,int entity_dim)1445 fvm_nodal_remove_parent_num(fvm_nodal_t *this_nodal,
1446 int entity_dim)
1447 {
1448 /* Vertices */
1449
1450 if (entity_dim == 0) {
1451 this_nodal->parent_vertex_num = NULL;
1452 if (this_nodal->_parent_vertex_num != NULL)
1453 BFT_FREE(this_nodal->_parent_vertex_num);
1454 }
1455
1456 /* Other elements */
1457
1458 else {
1459
1460 int i = 0;
1461 fvm_nodal_section_t *section = NULL;
1462
1463 for (i = 0; i < this_nodal->n_sections; i++) {
1464 section = this_nodal->sections[i];
1465 if (section->entity_dim == entity_dim) {
1466 section->parent_element_num = NULL;
1467 if (section->_parent_element_num != NULL)
1468 BFT_FREE(section->_parent_element_num);
1469 }
1470 }
1471
1472 }
1473
1474 }
1475
1476 /*----------------------------------------------------------------------------
1477 * Build external numbering for entities based on global numbers.
1478 *
1479 * parameters:
1480 * this_nodal <-- nodal mesh structure
1481 * parent_global_number <-- pointer to list of global (i.e. domain splitting
1482 * independent) parent entity numbers
1483 * entity_dim <-- 3 for cells, 2 for faces, 1 for edges,
1484 * and 0 for vertices
1485 *----------------------------------------------------------------------------*/
1486
1487 void
fvm_nodal_init_io_num(fvm_nodal_t * this_nodal,const cs_gnum_t parent_global_numbers[],int entity_dim)1488 fvm_nodal_init_io_num(fvm_nodal_t *this_nodal,
1489 const cs_gnum_t parent_global_numbers[],
1490 int entity_dim)
1491 {
1492 int i;
1493 fvm_nodal_section_t *section;
1494
1495 if (entity_dim == 0) {
1496 this_nodal->global_vertex_num
1497 = fvm_io_num_create(this_nodal->parent_vertex_num,
1498 parent_global_numbers,
1499 this_nodal->n_vertices,
1500 0);
1501 _remove_global_vertex_labels(this_nodal);
1502 }
1503
1504 else {
1505 for (i = 0; i < this_nodal->n_sections; i++) {
1506 section = this_nodal->sections[i];
1507 if (section->entity_dim == entity_dim) {
1508 section->global_element_num
1509 = fvm_io_num_create(section->parent_element_num,
1510 parent_global_numbers,
1511 section->n_elements,
1512 0);
1513 }
1514 }
1515 }
1516
1517 }
1518
1519 /*----------------------------------------------------------------------------
1520 * Transfer an existing global numbering structure to a nodal mesh.
1521 *
1522 * This assumes the local number of vertices is identical to that of
1523 * the numberign structure.
1524 *
1525 * The argument pointer to the global numbering structure is set to NULL
1526 * so the caller does not spuriously modify the structure after this call.
1527 *
1528 * A call to this function may be used instead of fvm_nodal_init_io_num
1529 * (for vertices only).
1530 *
1531 * parameters:
1532 * this_nodal <-- nodal mesh structure
1533 * io_num <-> pointer to external numbering structure
1534 * whose property is transferred to the mesh
1535 * independent) parent entity numbers
1536 *----------------------------------------------------------------------------*/
1537
1538 void
fvm_nodal_transfer_vertex_io_num(fvm_nodal_t * this_nodal,fvm_io_num_t ** io_num)1539 fvm_nodal_transfer_vertex_io_num(fvm_nodal_t *this_nodal,
1540 fvm_io_num_t **io_num)
1541 {
1542 this_nodal->global_vertex_num = *io_num;
1543 *io_num = NULL;
1544 _remove_global_vertex_labels(this_nodal);
1545 }
1546
1547 /*----------------------------------------------------------------------------
1548 * Set entity tags (for non-vertex entities).
1549 *
1550 * The number of entities of the given dimension may be obtained
1551 * through fvm_nodal_get_n_entities(), the tag[] array is populated
1552 * in local section order, section by section).
1553 *
1554 * parameters:
1555 * this_nodal <-- nodal mesh structure
1556 * tag <-- tag values to assign
1557 * entity_dim <-- 3 for cells, 2 for faces, 1 for edges
1558 *----------------------------------------------------------------------------*/
1559
1560 void
fvm_nodal_set_tag(fvm_nodal_t * this_nodal,const int tag[],int entity_dim)1561 fvm_nodal_set_tag(fvm_nodal_t *this_nodal,
1562 const int tag[],
1563 int entity_dim)
1564 {
1565 cs_lnum_t entitiy_count = 0;
1566 for (int i = 0; i < this_nodal->n_sections; i++) {
1567 fvm_nodal_section_t *section = this_nodal->sections[i];
1568 if (section->entity_dim == entity_dim) {
1569 BFT_REALLOC(section->tag, section->n_elements, int);
1570 for (cs_lnum_t j = 0; j < section->n_elements; j++)
1571 section->tag[j] = tag[entitiy_count + j];
1572 entitiy_count += section->n_elements;
1573 }
1574 }
1575 }
1576
1577 /*----------------------------------------------------------------------------
1578 * Remove entity tags.
1579 *
1580 * parameters:
1581 * this_nodal <-- nodal mesh structure
1582 * entity_dim <-- 3 for cells, 2 for faces, 1 for edges
1583 *----------------------------------------------------------------------------*/
1584
1585 void
fvm_nodal_remove_tag(fvm_nodal_t * this_nodal,int entity_dim)1586 fvm_nodal_remove_tag(fvm_nodal_t *this_nodal,
1587 int entity_dim)
1588 {
1589 for (int i = 0; i < this_nodal->n_sections; i++) {
1590 fvm_nodal_section_t *section = this_nodal->sections[i];
1591 if (section->entity_dim == entity_dim)
1592 BFT_FREE(section->tag);
1593 }
1594 }
1595
1596 /*----------------------------------------------------------------------------
1597 * Preset number and list of vertices to assign to a nodal mesh.
1598 *
1599 * If the parent_vertex_num argument is NULL, the list is assumed to
1600 * be {1, 2, ..., n}. If parent_vertex_num is given, it specifies a
1601 * list of n vertices from a larger set (1 to n numbering).
1602 *
1603 * Ownership of the given parent vertex numbering array is
1604 * transferred to the nodal mesh representation structure.
1605 *
1606 * This function should be called before fvm_nodal_set_shared_vertices()
1607 * or fvm_nodal_transfer_vertices() if we want to force certain
1608 * vertices to appear in the mesh (especially if we want to define
1609 * a mesh containing only vertices).
1610 *
1611 * parameters:
1612 * this_nodal <-> nodal mesh structure
1613 * n_vertices <-- number of vertices to assign
1614 * parent_vertex_num <-- parent numbers of vertices to assign
1615 *----------------------------------------------------------------------------*/
1616
1617 void
fvm_nodal_define_vertex_list(fvm_nodal_t * this_nodal,cs_lnum_t n_vertices,cs_lnum_t parent_vertex_num[])1618 fvm_nodal_define_vertex_list(fvm_nodal_t *this_nodal,
1619 cs_lnum_t n_vertices,
1620 cs_lnum_t parent_vertex_num[])
1621 {
1622 assert(this_nodal != NULL);
1623
1624 this_nodal->n_vertices = n_vertices;
1625
1626 this_nodal->parent_vertex_num = NULL;
1627 if (this_nodal->_parent_vertex_num != NULL)
1628 BFT_FREE(this_nodal->_parent_vertex_num);
1629
1630 if (parent_vertex_num != NULL) {
1631 this_nodal->_parent_vertex_num = parent_vertex_num;
1632 this_nodal->parent_vertex_num = parent_vertex_num;
1633 }
1634
1635 _remove_global_vertex_labels(this_nodal);
1636 }
1637
1638 /*----------------------------------------------------------------------------
1639 * Assign shared vertex coordinates to an extracted nodal mesh,
1640 * renumbering vertex numbers based on those really referenced,
1641 * and updating connectivity arrays in accordance.
1642 *
1643 * This function should only be called once all element sections
1644 * have been added to a nodal mesh representation.
1645 *
1646 * parameters:
1647 * this_nodal <-> nodal mesh structure
1648 * vertex_coords <-- coordinates of parent vertices (interlaced)
1649 *----------------------------------------------------------------------------*/
1650
1651 void
fvm_nodal_set_shared_vertices(fvm_nodal_t * this_nodal,const cs_coord_t vertex_coords[])1652 fvm_nodal_set_shared_vertices(fvm_nodal_t *this_nodal,
1653 const cs_coord_t vertex_coords[])
1654 {
1655 assert(this_nodal != NULL);
1656
1657 /* Map vertex coordinates to array passed as argument
1658 (this_nodal->_vertex_coords remains NULL, so only
1659 the const pointer may be used for a shared array) */
1660
1661 this_nodal->vertex_coords = vertex_coords;
1662
1663 /* If the mesh contains only vertices, its n_vertices and
1664 parent_vertex_num must already have been set, and do not
1665 require updating */
1666
1667 if (this_nodal->n_sections == 0)
1668 return;
1669
1670 /* Renumber vertices based on those really referenced */
1671
1672 _renumber_vertices(this_nodal);
1673
1674 _remove_global_vertex_labels(this_nodal);
1675 }
1676
1677 /*----------------------------------------------------------------------------
1678 * Assign private vertex coordinates to a nodal mesh,
1679 * renumbering vertex numbers based on those really referenced,
1680 * and updating connectivity arrays in accordance.
1681 *
1682 * Ownership of the given coordinates array is transferred to
1683 * the nodal mesh representation structure.
1684 *
1685 * This function should only be called once all element sections
1686 * have been added to a nodal mesh representation.
1687 *
1688 * parameters:
1689 * this_nodal <-> nodal mesh structure
1690 * vertex_coords <-- coordinates of parent vertices (interlaced)
1691 *
1692 * returns:
1693 * updated pointer to vertex_coords (may be different from initial
1694 * argument if vertices were renumbered).
1695 *----------------------------------------------------------------------------*/
1696
1697 cs_coord_t *
fvm_nodal_transfer_vertices(fvm_nodal_t * this_nodal,cs_coord_t vertex_coords[])1698 fvm_nodal_transfer_vertices(fvm_nodal_t *this_nodal,
1699 cs_coord_t vertex_coords[])
1700 {
1701 cs_lnum_t i;
1702 int j;
1703
1704 cs_coord_t *_vertex_coords = vertex_coords;
1705
1706 assert(this_nodal != NULL);
1707
1708 /* Renumber vertices based on those really referenced, and
1709 update connectivity arrays in accordance. */
1710
1711 _renumber_vertices(this_nodal);
1712
1713 /* If renumbering is necessary, update connectivity */
1714
1715 if (this_nodal->parent_vertex_num != NULL) {
1716
1717 int dim = this_nodal->dim;
1718 const cs_lnum_t *parent_vertex_num = this_nodal->parent_vertex_num;
1719
1720 BFT_MALLOC(_vertex_coords, this_nodal->n_vertices * dim, cs_coord_t);
1721
1722 for (i = 0; i < this_nodal->n_vertices; i++) {
1723 for (j = 0; j < dim; j++)
1724 _vertex_coords[i*dim + j]
1725 = vertex_coords[(parent_vertex_num[i]-1)*dim + j];
1726 }
1727
1728 BFT_FREE(vertex_coords);
1729
1730 this_nodal->parent_vertex_num = NULL;
1731 if (this_nodal->_parent_vertex_num != NULL)
1732 BFT_FREE(this_nodal->_parent_vertex_num);
1733 }
1734
1735 this_nodal->_vertex_coords = _vertex_coords;
1736 this_nodal->vertex_coords = _vertex_coords;
1737
1738 _remove_global_vertex_labels(this_nodal);
1739
1740 return _vertex_coords;
1741 }
1742
1743 /*----------------------------------------------------------------------------
1744 * Make vertex coordinates of a nodal mesh private.
1745 *
1746 * If vertex coordinates were previously shared, those coordinates that
1747 * are actually references are copied, and the relation to parent vertices
1748 * is discarded.
1749 *
1750 * If vertices were already private, the mesh is not modified.
1751 *
1752 * parameters:
1753 * this_nodal <-> nodal mesh structure
1754 *----------------------------------------------------------------------------*/
1755
1756 void
fvm_nodal_make_vertices_private(fvm_nodal_t * this_nodal)1757 fvm_nodal_make_vertices_private(fvm_nodal_t *this_nodal)
1758 {
1759 assert(this_nodal != NULL);
1760
1761 if (this_nodal->_vertex_coords == NULL) {
1762
1763 cs_coord_t *_vertex_coords = NULL;
1764 const cs_coord_t *vertex_coords = this_nodal->vertex_coords;
1765 const cs_lnum_t n_vertices = this_nodal->n_vertices;
1766 const int dim = this_nodal->dim;
1767
1768 BFT_MALLOC(_vertex_coords, n_vertices * dim, cs_coord_t);
1769
1770 /* If renumbering is necessary, update connectivity */
1771
1772 if (this_nodal->parent_vertex_num != NULL) {
1773
1774 cs_lnum_t i;
1775 int j;
1776 const cs_lnum_t *parent_vertex_num = this_nodal->parent_vertex_num;
1777
1778 for (i = 0; i < n_vertices; i++) {
1779 for (j = 0; j < dim; j++)
1780 _vertex_coords[i*dim + j]
1781 = vertex_coords[(parent_vertex_num[i]-1)*dim + j];
1782 }
1783
1784 this_nodal->parent_vertex_num = NULL;
1785 if (this_nodal->_parent_vertex_num != NULL)
1786 BFT_FREE(this_nodal->_parent_vertex_num);
1787 }
1788 else
1789 memcpy(_vertex_coords, vertex_coords, n_vertices*dim*sizeof(cs_coord_t));
1790
1791 /* Assign new array to structure */
1792
1793 this_nodal->_vertex_coords = _vertex_coords;
1794 this_nodal->vertex_coords = _vertex_coords;
1795 }
1796 }
1797
1798 /*----------------------------------------------------------------------------
1799 * Assign group class set descriptions to a nodal mesh.
1800 *
1801 * The structure builds its own copy of the group class sets,
1802 * renumbering them so as to discard those not referenced.
1803 * Empty group classes are also renumbered to zero.
1804 *
1805 * This function should only be called once all element sections
1806 * have been added to a nodal mesh representation.
1807 *
1808 * parameters:
1809 * this_nodal <-> nodal mesh structure
1810 * gc_set <-- group class set descriptions
1811 *----------------------------------------------------------------------------*/
1812
1813 void
fvm_nodal_set_group_class_set(fvm_nodal_t * this_nodal,const fvm_group_class_set_t * gc_set)1814 fvm_nodal_set_group_class_set(fvm_nodal_t *this_nodal,
1815 const fvm_group_class_set_t *gc_set)
1816 {
1817 int gc_id, section_id;
1818 int n_gc = fvm_group_class_set_size(gc_set);
1819 int n_gc_new = 0;
1820 int *gc_renum = NULL;
1821
1822 assert(this_nodal != NULL);
1823
1824 if (this_nodal->gc_set != NULL)
1825 this_nodal->gc_set = fvm_group_class_set_destroy(this_nodal->gc_set);
1826
1827 if (gc_set == NULL)
1828 return;
1829
1830 /* Mark referenced group classes */
1831
1832 BFT_MALLOC(gc_renum, n_gc, int);
1833
1834 for (gc_id = 0; gc_id < n_gc; gc_id++)
1835 gc_renum[gc_id] = 0;
1836
1837 for (section_id = 0; section_id < this_nodal->n_sections; section_id++) {
1838 cs_lnum_t i;
1839 const fvm_nodal_section_t *section = this_nodal->sections[section_id];
1840 if (section->gc_id == NULL)
1841 continue;
1842 for (i = 0; i < section->n_elements; i++) {
1843 if (section->gc_id[i] != 0)
1844 gc_renum[section->gc_id[i] - 1] = 1;
1845 }
1846 }
1847
1848 cs_parall_max(n_gc, CS_INT_TYPE, gc_renum);
1849
1850 /* Renumber group classes if necessary */
1851
1852 for (gc_id = 0; gc_id < n_gc; gc_id++) {
1853 if (gc_renum[gc_id] != 0) {
1854 gc_renum[gc_id] = n_gc_new + 1;
1855 n_gc_new++;
1856 }
1857 }
1858
1859 if (n_gc_new < n_gc) {
1860 for (section_id = 0; section_id < this_nodal->n_sections; section_id++) {
1861 cs_lnum_t i;
1862 const fvm_nodal_section_t *section = this_nodal->sections[section_id];
1863 if (section->gc_id == NULL)
1864 continue;
1865 for (i = 0; i < section->n_elements; i++) {
1866 if (section->gc_id[i] != 0)
1867 section->gc_id[i] = gc_renum[section->gc_id[i] - 1];
1868 }
1869 }
1870 }
1871
1872 /* Transform renumbering array to list */
1873
1874 n_gc_new = 0;
1875 for (gc_id = 0; gc_id < n_gc; gc_id++) {
1876 if (gc_renum[gc_id] != 0)
1877 gc_renum[n_gc_new++] = gc_id;
1878 }
1879
1880 if (n_gc_new > 0)
1881 this_nodal->gc_set = fvm_group_class_set_copy(gc_set,
1882 n_gc_new,
1883 gc_renum);
1884
1885 BFT_FREE(gc_renum);
1886 }
1887
1888 /*----------------------------------------------------------------------------
1889 * Assign global vertex labels to a nodal mesh.
1890 *
1891 * As these are expected to be used only for small sets (i.e. probes)
1892 * where the point set is built from a global definition and data movement
1893 * would add complexity and overhead, the labels refer to a global view
1894 * on rank 0.
1895 *
1896 * The size of the labels pointers array should be the same as that
1897 * returned by fvm_nodal_n_g_vertices();
1898 *
1899 * This function should only be called once the nodal mesh representation
1900 * has been completed, as most functions modifying its vertex definitions
1901 * will remove these labels.
1902 *
1903 * parameters:
1904 * this_nodal <-> nodal mesh structure
1905 * g_labels <-- global vertex labels, or NULL
1906 *----------------------------------------------------------------------------*/
1907
1908 void
fvm_nodal_transfer_global_vertex_labels(fvm_nodal_t * this_nodal,char * g_labels[])1909 fvm_nodal_transfer_global_vertex_labels(fvm_nodal_t *this_nodal,
1910 char *g_labels[])
1911 {
1912 assert(this_nodal != NULL);
1913
1914 _remove_global_vertex_labels(this_nodal);
1915
1916 this_nodal->global_vertex_labels = g_labels;
1917 }
1918
1919 /*----------------------------------------------------------------------------
1920 * Obtain the name of a nodal mesh.
1921 *
1922 * parameters:
1923 * this_nodal <-- pointer to nodal mesh structure
1924 *
1925 * returns:
1926 * pointer to constant string containing the mesh name
1927 *----------------------------------------------------------------------------*/
1928
1929 const char *
fvm_nodal_get_name(const fvm_nodal_t * this_nodal)1930 fvm_nodal_get_name(const fvm_nodal_t *this_nodal)
1931 {
1932 assert(this_nodal != NULL);
1933
1934 return this_nodal->name;
1935 }
1936
1937 /*----------------------------------------------------------------------------
1938 * Return spatial dimension of the nodal mesh.
1939 *
1940 * parameters:
1941 * this_nodal <-- pointer to nodal mesh structure
1942 *
1943 * returns:
1944 * spatial dimension.
1945 *----------------------------------------------------------------------------*/
1946
1947 int
fvm_nodal_get_dim(const fvm_nodal_t * this_nodal)1948 fvm_nodal_get_dim(const fvm_nodal_t *this_nodal)
1949 {
1950 return this_nodal->dim;
1951 }
1952
1953 /*----------------------------------------------------------------------------
1954 * Return maximum dimension of entities in a nodal mesh.
1955 *
1956 * parameters:
1957 * this_nodal <-- pointer to nodal mesh structure
1958 *
1959 * returns:
1960 * maximum dimension of entities in mesh (0 to 3)
1961 *----------------------------------------------------------------------------*/
1962
1963 int
fvm_nodal_get_max_entity_dim(const fvm_nodal_t * this_nodal)1964 fvm_nodal_get_max_entity_dim(const fvm_nodal_t *this_nodal)
1965 {
1966 int section_id;
1967 int max_entity_dim = 0;
1968
1969 assert(this_nodal != NULL);
1970
1971 for (section_id = 0; section_id < this_nodal->n_sections; section_id++) {
1972 const fvm_nodal_section_t *section = this_nodal->sections[section_id];
1973 if (section->entity_dim > max_entity_dim)
1974 max_entity_dim = section->entity_dim;
1975 }
1976
1977 return max_entity_dim;
1978 }
1979
1980 /*----------------------------------------------------------------------------
1981 * Return number of entities of a given dimension in a nodal mesh.
1982 *
1983 * parameters:
1984 * this_nodal <-- pointer to nodal mesh structure
1985 * entity_dim <-- dimension of entities we want to count (0 to 3)
1986 *
1987 * returns:
1988 * number of entities of given dimension in mesh
1989 *----------------------------------------------------------------------------*/
1990
1991 cs_lnum_t
fvm_nodal_get_n_entities(const fvm_nodal_t * this_nodal,int entity_dim)1992 fvm_nodal_get_n_entities(const fvm_nodal_t *this_nodal,
1993 int entity_dim)
1994 {
1995 cs_lnum_t n_entities;
1996
1997 assert(this_nodal != NULL);
1998
1999 switch(entity_dim) {
2000 case 0:
2001 n_entities = this_nodal->n_vertices;
2002 break;
2003 case 1:
2004 n_entities = this_nodal->n_edges;
2005 break;
2006 case 2:
2007 n_entities = this_nodal->n_faces;
2008 break;
2009 case 3:
2010 n_entities = this_nodal->n_cells;
2011 break;
2012 default:
2013 n_entities = 0;
2014 }
2015
2016 return n_entities;
2017 }
2018
2019 /*----------------------------------------------------------------------------
2020 * Return global number of vertices associated with nodal mesh.
2021 *
2022 * parameters:
2023 * this_nodal <-- pointer to nodal mesh structure
2024 *
2025 * returns:
2026 * global number of vertices associated with nodal mesh
2027 *----------------------------------------------------------------------------*/
2028
2029 cs_gnum_t
fvm_nodal_get_n_g_vertices(const fvm_nodal_t * this_nodal)2030 fvm_nodal_get_n_g_vertices(const fvm_nodal_t *this_nodal)
2031 {
2032 return fvm_nodal_n_g_vertices(this_nodal);
2033 }
2034
2035 /*----------------------------------------------------------------------------
2036 * Return global number of elements of a given type associated with nodal mesh.
2037 *
2038 * parameters:
2039 * this_nodal <-- pointer to nodal mesh structure
2040 * element_type <-- type of elements for query
2041 *
2042 * returns:
2043 * global number of elements of the given type associated with nodal mesh
2044 *----------------------------------------------------------------------------*/
2045
2046 cs_gnum_t
fvm_nodal_get_n_g_elements(const fvm_nodal_t * this_nodal,fvm_element_t element_type)2047 fvm_nodal_get_n_g_elements(const fvm_nodal_t *this_nodal,
2048 fvm_element_t element_type)
2049 {
2050 int i;
2051 cs_gnum_t n_g_elements = 0;
2052
2053 assert(this_nodal != NULL);
2054
2055 for (i = 0; i < this_nodal->n_sections; i++) {
2056 fvm_nodal_section_t *section = this_nodal->sections[i];
2057 if (section->type == element_type)
2058 n_g_elements += fvm_nodal_section_n_g_elements(section);
2059 }
2060
2061 return n_g_elements;
2062 }
2063
2064 /*----------------------------------------------------------------------------
2065 * Return local number of elements of a given type associated with nodal mesh.
2066 *
2067 * parameters:
2068 * this_nodal <-- pointer to nodal mesh structure
2069 * element_type <-- type of elements for query
2070 *
2071 * returns:
2072 * local number of elements of the given type associated with nodal mesh
2073 *----------------------------------------------------------------------------*/
2074
2075 cs_lnum_t
fvm_nodal_get_n_elements(const fvm_nodal_t * this_nodal,fvm_element_t element_type)2076 fvm_nodal_get_n_elements(const fvm_nodal_t *this_nodal,
2077 fvm_element_t element_type)
2078 {
2079 int i;
2080 cs_lnum_t n_elements = 0;
2081
2082 assert(this_nodal != NULL);
2083
2084 for (i = 0; i < this_nodal->n_sections; i++) {
2085 fvm_nodal_section_t *section = this_nodal->sections[i];
2086 if (section->type == element_type)
2087 n_elements += section->n_elements;
2088 }
2089
2090 return n_elements;
2091 }
2092
2093 /*----------------------------------------------------------------------------
2094 * Return local parent numbering array for all entities of a given
2095 * dimension in a nodal mesh.
2096 *
2097 * The number of entities of the given dimension may be obtained
2098 * through fvm_nodal_get_n_entities(), the parent_num[] array is populated
2099 * with the parent entity numbers of those entities, in order (i.e. in
2100 * local section order, section by section).
2101 *
2102 * This function is similar to fvm_nodal_get_parent_num(), but returns
2103 * numbers (1 to n) instead of ids (0 to n-1).
2104 *
2105 * parameters:
2106 * this_nodal <-- pointer to nodal mesh structure
2107 * entity_dim <-- dimension of entities we are interested in (0 to 3)
2108 * parent_num --> entity parent numbering (array must be pre-allocated)
2109 *----------------------------------------------------------------------------*/
2110
2111 void
fvm_nodal_get_parent_num(const fvm_nodal_t * this_nodal,int entity_dim,cs_lnum_t parent_num[])2112 fvm_nodal_get_parent_num(const fvm_nodal_t *this_nodal,
2113 int entity_dim,
2114 cs_lnum_t parent_num[])
2115 {
2116 int section_id;
2117 cs_lnum_t i;
2118
2119 cs_lnum_t entity_count = 0;
2120
2121 assert(this_nodal != NULL);
2122
2123 /* Entity dimension 0: vertices */
2124
2125 if (entity_dim == 0) {
2126 if (this_nodal->parent_vertex_num != NULL) {
2127 for (i = 0; i < this_nodal->n_vertices; i++)
2128 parent_num[entity_count++] = this_nodal->parent_vertex_num[i];
2129 }
2130 else {
2131 for (i = 0; i < this_nodal->n_vertices; i++)
2132 parent_num[entity_count++] = i + 1;
2133 }
2134 }
2135
2136 /* Entity dimension > 0: edges, faces, or cells */
2137
2138 else {
2139
2140 for (section_id = 0; section_id < this_nodal->n_sections; section_id++) {
2141
2142 const fvm_nodal_section_t *section = this_nodal->sections[section_id];
2143
2144 if (section->entity_dim == entity_dim) {
2145 if (section->parent_element_num != NULL) {
2146 for (i = 0; i < section->n_elements; i++)
2147 parent_num[entity_count++] = section->parent_element_num[i];
2148 }
2149 else {
2150 for (i = 0; i < section->n_elements; i++)
2151 parent_num[entity_count++] = i + 1;
2152 }
2153 }
2154
2155 } /* end loop on sections */
2156
2157 }
2158 }
2159
2160 /*----------------------------------------------------------------------------
2161 * Return local parent id array for all entities of a given
2162 * dimension in a nodal mesh.
2163 *
2164 * The number of entities of the given dimension may be obtained
2165 * through fvm_nodal_get_n_entities(), the parent_num[] array is populated
2166 * with the parent entity numbers of those entities, in order (i.e. in
2167 * local section order, section by section).
2168 *
2169 * parameters:
2170 * this_nodal <-- pointer to nodal mesh structure
2171 * entity_dim <-- dimension of entities we are interested in (0 to 3)
2172 * parent_id --> entity parent id (array must be pre-allocated)
2173 *----------------------------------------------------------------------------*/
2174
2175 void
fvm_nodal_get_parent_id(const fvm_nodal_t * this_nodal,int entity_dim,cs_lnum_t parent_id[])2176 fvm_nodal_get_parent_id(const fvm_nodal_t *this_nodal,
2177 int entity_dim,
2178 cs_lnum_t parent_id[])
2179 {
2180 int section_id;
2181 cs_lnum_t i;
2182
2183 cs_lnum_t entity_count = 0;
2184
2185 assert(this_nodal != NULL);
2186
2187 /* Entity dimension 0: vertices */
2188
2189 if (entity_dim == 0) {
2190 if (this_nodal->parent_vertex_num != NULL) {
2191 for (i = 0; i < this_nodal->n_vertices; i++)
2192 parent_id[entity_count++] = this_nodal->parent_vertex_num[i] - 1;
2193 }
2194 else {
2195 for (i = 0; i < this_nodal->n_vertices; i++)
2196 parent_id[entity_count++] = i;
2197 }
2198 }
2199
2200 /* Entity dimension > 0: edges, faces, or cells */
2201
2202 else {
2203
2204 for (section_id = 0; section_id < this_nodal->n_sections; section_id++) {
2205
2206 const fvm_nodal_section_t *section = this_nodal->sections[section_id];
2207
2208 if (section->entity_dim == entity_dim) {
2209 if (section->parent_element_num != NULL) {
2210 for (i = 0; i < section->n_elements; i++)
2211 parent_id[entity_count++] = section->parent_element_num[i] - 1;
2212 }
2213 else {
2214 for (i = 0; i < section->n_elements; i++)
2215 parent_id[entity_count++] = i;
2216 }
2217 }
2218
2219 } /* end loop on sections */
2220
2221 }
2222 }
2223
2224 /*----------------------------------------------------------------------------
2225 * Return pointer to global vertex labels of a nodal mesh.
2226 *
2227 * As these are expected to be used only for small sets (i.e. probes)
2228 * where the point set is built from a global definition and data movement
2229 * would adds complexity and overhead, the labels refer to a global view
2230 * on rank 0; for the same reason, only shared labels are needed.
2231 *
2232 * The size of the labels pointers array returned should be the same
2233 * as that returned by fvm_nodal_n_g_vertices();
2234 *
2235 * parameters:
2236 * this_nodal <-> nodal mesh structure
2237 *
2238 * returns:
2239 * pointer to global vertex labels, or NULL
2240 *----------------------------------------------------------------------------*/
2241
2242 const char **
fvm_nodal_get_global_vertex_labels(const fvm_nodal_t * this_nodal)2243 fvm_nodal_get_global_vertex_labels(const fvm_nodal_t *this_nodal)
2244 {
2245 assert(this_nodal != NULL);
2246
2247 return (const char **)(this_nodal->global_vertex_labels);
2248 }
2249
2250 /*----------------------------------------------------------------------------
2251 * Return a const pointer to a nodal mesh representation's parent mesh.
2252 *
2253 * parameters:
2254 * this_nodal <-> pointer to structure that should be reduced
2255 *
2256 * return:
2257 * const pointer to parent mesh
2258 *----------------------------------------------------------------------------*/
2259
2260 const cs_mesh_t *
fvm_nodal_get_parent(const fvm_nodal_t * this_nodal)2261 fvm_nodal_get_parent(const fvm_nodal_t *this_nodal)
2262 {
2263 return this_nodal->parent;
2264 }
2265
2266 /*----------------------------------------------------------------------------
2267 * Associate a parent mesh to a nodal mesh representation structure.
2268 *
2269 * parameters:
2270 * this_nodal <-> pointer to structure that should be reduced
2271 * parent <-- const pointer to parent mesh
2272 *----------------------------------------------------------------------------*/
2273
2274 void
fvm_nodal_set_parent(fvm_nodal_t * this_nodal,const cs_mesh_t * parent)2275 fvm_nodal_set_parent(fvm_nodal_t *this_nodal,
2276 const cs_mesh_t *parent)
2277 {
2278 this_nodal->parent = parent;
2279 }
2280
2281 /*----------------------------------------------------------------------------
2282 * Compute tesselation a a nodal mesh's sections of a given type, and add the
2283 * corresponding structure to the mesh representation.
2284 *
2285 * If global element numbers are used (i.e. in parallel mode), this function
2286 * should be only be used after calling fvm_nodal_init_io_num().
2287 *
2288 * If some mesh sections have already been tesselated, their tesselation
2289 * is unchanged.
2290 *
2291 * parameters:
2292 * this_nodal <-> pointer to nodal mesh structure
2293 * type <-> element type that should be tesselated
2294 * error_count --> number of elements with a tesselation error
2295 * counter (optional)
2296 *----------------------------------------------------------------------------*/
2297
2298 void
fvm_nodal_tesselate(fvm_nodal_t * this_nodal,fvm_element_t type,cs_lnum_t * error_count)2299 fvm_nodal_tesselate(fvm_nodal_t *this_nodal,
2300 fvm_element_t type,
2301 cs_lnum_t *error_count)
2302 {
2303 int section_id;
2304 cs_lnum_t section_error_count;
2305
2306 assert(this_nodal != NULL);
2307
2308 if (error_count != NULL)
2309 *error_count = 0;
2310
2311 for (section_id = 0; section_id < this_nodal->n_sections; section_id++) {
2312
2313 fvm_nodal_section_t *section = this_nodal->sections[section_id];
2314
2315 if (section->type == type && section->tesselation == NULL) {
2316
2317 section->tesselation = fvm_tesselation_create(type,
2318 section->n_elements,
2319 section->face_index,
2320 section->face_num,
2321 section->vertex_index,
2322 section->vertex_num,
2323 section->global_element_num);
2324
2325 fvm_tesselation_init(section->tesselation,
2326 this_nodal->dim,
2327 this_nodal->vertex_coords,
2328 this_nodal->parent_vertex_num,
2329 §ion_error_count);
2330
2331 if (error_count != NULL)
2332 *error_count += section_error_count;
2333 }
2334
2335 }
2336 }
2337
2338 /*----------------------------------------------------------------------------
2339 * Build a nodal representation structure based on extraction of a
2340 * mesh's edges.
2341 *
2342 * parameters:
2343 * name <-- name to assign to extracted mesh
2344 * this_nodal <-> pointer to nodal mesh structure
2345 *----------------------------------------------------------------------------*/
2346
2347 fvm_nodal_t *
fvm_nodal_copy_edges(const char * name,const fvm_nodal_t * this_nodal)2348 fvm_nodal_copy_edges(const char *name,
2349 const fvm_nodal_t *this_nodal)
2350 {
2351 int i;
2352 cs_lnum_t j, k;
2353 cs_lnum_t n_edges = 0, n_max_edges = 0;
2354 fvm_nodal_t *new_nodal = NULL;
2355 fvm_nodal_section_t *new_section = NULL;
2356
2357 BFT_MALLOC(new_nodal, 1, fvm_nodal_t);
2358
2359 /* Global indicators */
2360
2361 if (name != NULL) {
2362 BFT_MALLOC(new_nodal->name, strlen(name) + 1, char);
2363 strcpy(new_nodal->name, name);
2364 }
2365 else
2366 new_nodal->name = NULL;
2367
2368 new_nodal->dim = this_nodal->dim;
2369 new_nodal->num_dom = this_nodal->num_dom;
2370 new_nodal->n_doms = this_nodal->n_doms;
2371 new_nodal->n_sections = 1;
2372
2373 /* Local dimensions */
2374
2375 new_nodal->n_cells = 0;
2376 new_nodal->n_faces = 0;
2377 new_nodal->n_edges = 0;
2378 new_nodal->n_vertices = this_nodal->n_vertices;
2379
2380 /* Local structures */
2381
2382 new_nodal->vertex_coords = this_nodal->vertex_coords;
2383 new_nodal->_vertex_coords = NULL;
2384
2385 new_nodal->parent_vertex_num = this_nodal->parent_vertex_num;
2386 new_nodal->_parent_vertex_num = NULL;
2387
2388 if (this_nodal->global_vertex_num != NULL) {
2389 cs_lnum_t n_ent
2390 = fvm_io_num_get_local_count(this_nodal->global_vertex_num);
2391 cs_gnum_t global_count
2392 = fvm_io_num_get_global_count(this_nodal->global_vertex_num);
2393 const cs_gnum_t *global_num
2394 = fvm_io_num_get_global_num(this_nodal->global_vertex_num);
2395
2396 new_nodal->global_vertex_num
2397 = fvm_io_num_create_shared(global_num, global_count, n_ent);
2398 }
2399 else
2400 new_nodal->global_vertex_num = NULL;
2401
2402 /* Counting step */
2403
2404 for (i = 0; i < this_nodal->n_sections; i++) {
2405 const fvm_nodal_section_t *this_section = this_nodal->sections[i];
2406 if (this_section->vertex_index == NULL)
2407 n_max_edges += ( fvm_nodal_n_edges_element[this_section->type]
2408 * this_section->n_elements);
2409 else if (this_section->type == FVM_FACE_POLY)
2410 n_max_edges += this_section->vertex_index[this_section->n_elements];
2411 else if (this_section->type == FVM_CELL_POLY)
2412 n_max_edges += this_section->vertex_index[this_section->n_faces];
2413 }
2414
2415 BFT_MALLOC(new_nodal->sections, 1, fvm_nodal_section_t *);
2416
2417 new_section = fvm_nodal_section_create(FVM_EDGE);
2418 new_nodal->sections[0] = new_section;
2419
2420 BFT_MALLOC(new_section->_vertex_num, n_max_edges*2, cs_lnum_t);
2421
2422 /* Add edges */
2423
2424 for (i = 0; i < this_nodal->n_sections; i++) {
2425
2426 const fvm_nodal_section_t *this_section = this_nodal->sections[i];
2427
2428 if ( this_section->type == FVM_FACE_POLY
2429 || this_section->type == FVM_CELL_POLY) {
2430
2431 cs_lnum_t n_faces = this_section->type == FVM_FACE_POLY ?
2432 this_section->n_elements : this_section->n_faces;
2433
2434 for (j = 0; j < n_faces; j++) {
2435 const cs_lnum_t face_start_id = this_section->vertex_index[j];
2436 const cs_lnum_t n_face_edges
2437 = this_section->vertex_index[j+1] - this_section->vertex_index[j];
2438 for (k = 0; k < n_face_edges; k++) {
2439 new_section->_vertex_num[n_edges*2]
2440 = this_section->vertex_num[face_start_id + k];
2441 new_section->_vertex_num[n_edges*2 + 1]
2442 = this_section->vertex_num[face_start_id + (k + 1)%n_face_edges];
2443 n_edges += 1;
2444 }
2445 }
2446
2447 }
2448 else {
2449
2450 cs_lnum_t edges[2][12];
2451
2452 cs_lnum_t n_elt_edges = fvm_nodal_n_edges_element[this_section->type];
2453 cs_lnum_t n_elts = this_section->n_elements;
2454 cs_lnum_t stride = this_section->stride;
2455
2456 switch (this_section->type) {
2457
2458 case FVM_EDGE:
2459 case FVM_FACE_TRIA:
2460 case FVM_FACE_QUAD:
2461 for (j = 0; j < n_elt_edges; j++) {
2462 edges[0][j] = j;
2463 edges[1][j] = (j+1)%n_elt_edges;
2464 }
2465 break;
2466
2467 case FVM_CELL_TETRA:
2468 edges[0][0] = 0; edges[1][0] = 1;
2469 edges[0][1] = 1; edges[1][1] = 2;
2470 edges[0][2] = 2; edges[1][2] = 0;
2471 edges[0][3] = 0; edges[1][3] = 3;
2472 edges[0][4] = 1; edges[1][4] = 3;
2473 edges[0][5] = 2; edges[1][5] = 3;
2474 break;
2475
2476 case FVM_CELL_PYRAM:
2477 edges[0][0] = 0; edges[1][0] = 1;
2478 edges[0][1] = 1; edges[1][1] = 2;
2479 edges[0][2] = 2; edges[1][2] = 3;
2480 edges[0][3] = 3; edges[1][3] = 0;
2481 edges[0][4] = 0; edges[1][4] = 4;
2482 edges[0][5] = 1; edges[1][5] = 4;
2483 edges[0][6] = 2; edges[1][6] = 4;
2484 edges[0][7] = 3; edges[1][7] = 4;
2485 break;
2486
2487 case FVM_CELL_PRISM:
2488 edges[0][0] = 0; edges[1][0] = 1;
2489 edges[0][1] = 1; edges[1][1] = 2;
2490 edges[0][2] = 2; edges[1][2] = 0;
2491 edges[0][3] = 0; edges[1][3] = 3;
2492 edges[0][4] = 1; edges[1][4] = 4;
2493 edges[0][5] = 2; edges[1][5] = 5;
2494 edges[0][6] = 3; edges[1][6] = 4;
2495 edges[0][7] = 4; edges[1][7] = 5;
2496 edges[0][8] = 5; edges[1][8] = 3;
2497 break;
2498
2499 case FVM_CELL_HEXA:
2500 edges[0][0] = 0; edges[1][0] = 1;
2501 edges[0][1] = 1; edges[1][1] = 2;
2502 edges[0][2] = 2; edges[1][2] = 3;
2503 edges[0][3] = 3; edges[1][3] = 0;
2504 edges[0][4] = 0; edges[1][4] = 4;
2505 edges[0][5] = 1; edges[1][5] = 5;
2506 edges[0][6] = 2; edges[1][6] = 6;
2507 edges[0][7] = 3; edges[1][7] = 7;
2508 edges[0][8] = 4; edges[1][8] = 5;
2509 edges[0][9] = 5; edges[1][9] = 6;
2510 edges[0][10] = 6; edges[1][10] = 7;
2511 edges[0][11] = 7; edges[1][11] = 4;
2512 break;
2513
2514 default:
2515 assert(0);
2516 edges[0][0] = -1; /* For nonempty default clause */
2517 }
2518
2519 for (j = 0; j < n_elts; j++) {
2520 const cs_lnum_t *_vertex_num = this_section->vertex_num + (j*stride);
2521 for (k = 0; k < n_elt_edges; k++) {
2522 new_section->_vertex_num[n_edges*2] = _vertex_num[edges[0][k]];
2523 new_section->_vertex_num[n_edges*2 + 1] = _vertex_num[edges[1][k]];
2524 n_edges += 1;
2525 }
2526 }
2527 }
2528 } /* End of loop on sections */
2529
2530 assert(n_edges == n_max_edges);
2531
2532 /* Ensure edges are oriented in the same direction */
2533
2534 if (this_nodal->global_vertex_num != NULL) {
2535
2536 const cs_gnum_t *v_num_g
2537 = fvm_io_num_get_global_num(this_nodal->global_vertex_num);
2538
2539 for (j = 0; j < n_max_edges; j++) {
2540 cs_lnum_t vnum_1 = new_section->_vertex_num[j*2];
2541 cs_lnum_t vnum_2 = new_section->_vertex_num[j*2 + 1];
2542 if (v_num_g[vnum_1 - 1] > v_num_g[vnum_2 - 1]) {
2543 new_section->_vertex_num[j*2] = vnum_2;
2544 new_section->_vertex_num[j*2 + 1] = vnum_1;
2545 }
2546
2547 }
2548
2549 }
2550 else {
2551
2552 for (j = 0; j < n_max_edges; j++) {
2553 cs_lnum_t vnum_1 = new_section->_vertex_num[j*2];
2554 cs_lnum_t vnum_2 = new_section->_vertex_num[j*2 + 1];
2555 if (vnum_1 > vnum_2) {
2556 new_section->_vertex_num[j*2] = vnum_2;
2557 new_section->_vertex_num[j*2 + 1] = vnum_1;
2558 }
2559 }
2560
2561 }
2562
2563 /* Sort and remove duplicates
2564 (use qsort rather than cs_order_..._s() so as to sort in place) */
2565
2566 qsort(new_section->_vertex_num,
2567 n_max_edges,
2568 sizeof(cs_lnum_t) * 2,
2569 &_compare_edges);
2570
2571 {
2572 cs_lnum_t vn_1_p = -1;
2573 cs_lnum_t vn_2_p = -1;
2574
2575 n_edges = 0;
2576
2577 for (j = 0; j < n_max_edges; j++) {
2578
2579 cs_lnum_t vn_1 = new_section->_vertex_num[j*2];
2580 cs_lnum_t vn_2 = new_section->_vertex_num[j*2 + 1];
2581
2582 if (vn_1 != vn_1_p || vn_2 != vn_2_p) {
2583 new_section->_vertex_num[n_edges*2] = vn_1;
2584 new_section->_vertex_num[n_edges*2 + 1] = vn_2;
2585 vn_1_p = vn_1;
2586 vn_2_p = vn_2;
2587 n_edges += 1;
2588 }
2589 }
2590 }
2591
2592 /* Resize edge connectivity to adjust to final size */
2593
2594 BFT_REALLOC(new_section->_vertex_num, n_edges*2, cs_lnum_t);
2595 new_section->vertex_num = new_section->_vertex_num;
2596
2597 new_section->n_elements = n_edges;
2598 new_nodal->n_edges = n_edges;
2599
2600 /* Build global edge numbering if necessary */
2601
2602 if (new_nodal->n_doms > 1) {
2603
2604 cs_gnum_t *edge_vertices_g; /* edges -> global vertices */
2605
2606 BFT_MALLOC(edge_vertices_g, n_edges*2, cs_gnum_t);
2607
2608 if (this_nodal->global_vertex_num != NULL) {
2609 const cs_gnum_t *v_num_g
2610 = fvm_io_num_get_global_num(this_nodal->global_vertex_num);
2611 for (j = 0; j < n_edges; j++) {
2612 edge_vertices_g[j*2] = v_num_g[new_section->_vertex_num[j*2] - 1];
2613 edge_vertices_g[j*2+1] = v_num_g[new_section->_vertex_num[j*2+1] - 1];
2614 }
2615 }
2616 else {
2617 for (j = 0; j < n_edges; j++) {
2618 edge_vertices_g[j*2] = new_section->_vertex_num[j*2];
2619 edge_vertices_g[j*2 + 1] = new_section->_vertex_num[j*2 + 1];
2620 }
2621 }
2622
2623 new_section->global_element_num
2624 = fvm_io_num_create_from_adj_s(NULL, edge_vertices_g, n_edges, 2);
2625
2626
2627 BFT_FREE(edge_vertices_g);
2628 };
2629
2630 new_nodal->gc_set = NULL;
2631
2632 return (new_nodal);
2633 }
2634
2635 /*----------------------------------------------------------------------------
2636 * Dump printout of a nodal representation structure.
2637 *
2638 * parameters:
2639 * this_nodal <-- pointer to structure that should be dumped
2640 *----------------------------------------------------------------------------*/
2641
2642 void
fvm_nodal_dump(const fvm_nodal_t * this_nodal)2643 fvm_nodal_dump(const fvm_nodal_t *this_nodal)
2644 {
2645 cs_lnum_t i;
2646 cs_lnum_t num_vertex = 1;
2647 const cs_coord_t *coord = this_nodal->vertex_coords;
2648
2649 /* Global indicators */
2650 /*--------------------*/
2651
2652 bft_printf("\n"
2653 "Mesh name:\"%s\"\n",
2654 this_nodal->name);
2655
2656 bft_printf("\n"
2657 "Mesh dimension: %d\n"
2658 "Domain number: %d\n"
2659 "Number of domains: %d\n"
2660 "Number of sections: %d\n",
2661 this_nodal->dim, this_nodal->num_dom, this_nodal->n_doms,
2662 this_nodal->n_sections);
2663
2664 bft_printf("\n"
2665 "Number of cells: %ld\n"
2666 "Number of faces: %ld\n"
2667 "Number of edges: %ld\n"
2668 "Number of vertices: %ld\n",
2669 (long)this_nodal->n_cells,
2670 (long)this_nodal->n_faces,
2671 (long)this_nodal->n_edges,
2672 (long)this_nodal->n_vertices);
2673
2674 if (this_nodal->n_vertices > 0) {
2675
2676 bft_printf("\n"
2677 "Pointers to shareable arrays:\n"
2678 " vertex_coords: %p\n"
2679 " parent_vertex_num: %p\n",
2680 (const void *)this_nodal->vertex_coords,
2681 (const void *)this_nodal->parent_vertex_num);
2682
2683 bft_printf("\n"
2684 "Pointers to local arrays:\n"
2685 " _vertex_coords: %p\n"
2686 " _parent_vertex_num: %p\n",
2687 (const void *)this_nodal->_vertex_coords,
2688 (const void *)this_nodal->_parent_vertex_num);
2689
2690 /* Output coordinates depending on parent numbering */
2691
2692 if (this_nodal->parent_vertex_num == NULL) {
2693
2694 bft_printf("\nVertex coordinates:\n\n");
2695 switch(this_nodal->dim) {
2696 case 1:
2697 for (i = 0; i < this_nodal->n_vertices; i++)
2698 bft_printf("%10ld : %12.5f\n",
2699 (long)num_vertex++, (double)(coord[i]));
2700 break;
2701 case 2:
2702 for (i = 0; i < this_nodal->n_vertices; i++)
2703 bft_printf("%10ld : %12.5f %12.5f\n",
2704 (long)num_vertex++, (double)(coord[i*2]),
2705 (double)(coord[i*2+1]));
2706 break;
2707 case 3:
2708 for (i = 0; i < this_nodal->n_vertices; i++)
2709 bft_printf("%10ld : %12.5f %12.5f %12.5f\n",
2710 (long)num_vertex++, (double)(coord[i*3]),
2711 (double)(coord[i*3+1]), (double)(coord[i*3+2]));
2712 break;
2713 default:
2714 bft_printf("coordinates not output\n"
2715 "dimension = %d unsupported\n", this_nodal->dim);
2716 }
2717
2718 }
2719 else { /* if (this_nodal->parent_vertex_num != NULL) */
2720
2721 bft_printf("\nVertex parent and coordinates:\n\n");
2722
2723 switch(this_nodal->dim) {
2724 case 1:
2725 for (i = 0; i < this_nodal->n_vertices; i++) {
2726 coord = this_nodal->vertex_coords
2727 + (this_nodal->parent_vertex_num[i]-1);
2728 bft_printf("%10ld : %12.5f\n",
2729 (long)num_vertex++, (double)(coord[0]));
2730 }
2731 break;
2732 case 2:
2733 for (i = 0; i < this_nodal->n_vertices; i++) {
2734 coord = this_nodal->vertex_coords
2735 + ((this_nodal->parent_vertex_num[i]-1)*2);
2736 bft_printf("%10ld : %12.5f %12.5f\n",
2737 (long)num_vertex++, (double)(coord[0]), (double)(coord[1]));
2738 }
2739 break;
2740 case 3:
2741 for (i = 0; i < this_nodal->n_vertices; i++) {
2742 coord = this_nodal->vertex_coords
2743 + ((this_nodal->parent_vertex_num[i]-1)*3);
2744 bft_printf("%10ld : %12.5f %12.5f %12.5f\n",
2745 (long)num_vertex++, (double)(coord[0]), (double)(coord[1]),
2746 (double)(coord[2]));
2747 }
2748 break;
2749 default:
2750 bft_printf("coordinates not output\n"
2751 "dimension = %d unsupported\n", this_nodal->dim);
2752 }
2753
2754 }
2755
2756 }
2757
2758 /* Global vertex numbers (only for parallel execution) */
2759 if (this_nodal->global_vertex_num != NULL) {
2760 bft_printf("\nGlobal vertex numbers:\n\n");
2761 fvm_io_num_dump(this_nodal->global_vertex_num);
2762 }
2763
2764 /* Dump element sections */
2765 /*-----------------------*/
2766
2767 for (i = 0; i < this_nodal->n_sections; i++)
2768 _fvm_nodal_section_dump(this_nodal->sections[i]);
2769
2770 /* Dump metadata */
2771 /*---------------*/
2772
2773 /* Dump group class set (NULL allowed) */
2774
2775 fvm_group_class_set_dump(this_nodal->gc_set);
2776
2777 /* Dump vertex labels */
2778
2779 if (this_nodal->global_vertex_labels != NULL) {
2780 bft_printf("\nGlobal vertex labels:\n");
2781 cs_lnum_t n_labels = (cs_lnum_t)fvm_nodal_n_g_vertices(this_nodal);
2782 for (i = 0; i < n_labels; i++) {
2783 if (this_nodal->global_vertex_labels[i] != NULL)
2784 bft_printf("%10ld : %s\n",
2785 (long)i+1, this_nodal->global_vertex_labels[i]);
2786 else
2787 bft_printf("%10ld : (null)\n", (long)i+1);
2788 }
2789 }
2790 }
2791
2792 /*----------------------------------------------------------------------------*/
2793
2794 END_C_DECLS
2795