1
/*============================================================================
2
* Main structure for a nodal representation associated with a mesh
3
*============================================================================*/
6
This file is part of Code_Saturne, a general-purpose CFD tool.
8
Copyright (C) 1998-2011 EDF S.A.
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
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
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.
25
/*----------------------------------------------------------------------------*/
27
#if defined(HAVE_CONFIG_H)
28
#include "cs_config.h"
31
/*----------------------------------------------------------------------------
32
* Standard C library headers
33
*----------------------------------------------------------------------------*/
40
/*----------------------------------------------------------------------------
42
*----------------------------------------------------------------------------*/
45
#include <bft_printf.h>
47
/*----------------------------------------------------------------------------
49
*----------------------------------------------------------------------------*/
51
#include "fvm_config_defs.h"
53
#include "fvm_io_num.h"
54
#include "fvm_parall.h"
55
#include "fvm_tesselation.h"
57
/*----------------------------------------------------------------------------
58
* Header for the current file
59
*----------------------------------------------------------------------------*/
61
#include "fvm_nodal.h"
62
#include "fvm_nodal_priv.h"
64
/*----------------------------------------------------------------------------*/
69
} /* Fake brace to force back Emacs auto-indentation back to column 0 */
71
#endif /* __cplusplus */
73
/*============================================================================
74
* Static global variables
75
*============================================================================*/
77
/* Number of vertices associated with each "nodal" element type */
79
const int fvm_nodal_n_vertices_element[] = {2, /* Edge */
82
0, /* Simple polygon */
87
0}; /* Simple polyhedron */
89
/* Number of vertices associated with each "nodal" element type */
91
static int fvm_nodal_n_edges_element[] = {1, /* Edge */
94
0, /* Simple polygon */
99
0}; /* Simple polyhedron */
101
/*============================================================================
102
* Private function definitions
103
*============================================================================*/
105
/*----------------------------------------------------------------------------
106
* Compare edges (qsort function).
109
* x <-> pointer to first edge definition
110
* y <-> pointer to second edge definition
113
* result of strcmp() on group names
114
*----------------------------------------------------------------------------*/
116
static int _compare_edges(const void *x, const void *y)
120
const fvm_lnum_t *e0 = x;
121
const fvm_lnum_t *e1 = y;
126
else if (e0[0] == e1[0]) {
129
else if (e0[1] == e1[1])
136
/*----------------------------------------------------------------------------
137
* Copy a nodal mesh section representation structure, sharing arrays
138
* with the original structure.
141
* this_section <-> pointer to structure that should be copied
144
* pointer to created nodal mesh section representation structure
145
*----------------------------------------------------------------------------*/
147
static fvm_nodal_section_t *
148
_fvm_nodal_section_copy(const fvm_nodal_section_t *this_section)
150
fvm_nodal_section_t *new_section;
152
BFT_MALLOC(new_section, 1, fvm_nodal_section_t);
154
/* Global information */
156
new_section->entity_dim = this_section->entity_dim;
158
new_section->n_elements = this_section->n_elements;
159
new_section->type = this_section->type;
163
new_section->connectivity_size = this_section->connectivity_size;
164
new_section->stride = this_section->stride;
166
new_section->n_faces = this_section->n_faces;
168
new_section->face_index = this_section->face_index;
169
new_section->face_num = this_section->face_num;
170
new_section->vertex_index = this_section->vertex_index;
171
new_section->vertex_num = this_section->vertex_num;
173
new_section->_face_index = NULL;
174
new_section->_face_num = NULL;
175
new_section->_vertex_index = NULL;
176
new_section->_vertex_num = NULL;
178
new_section->gc_id = NULL;
180
new_section->tesselation = NULL; /* TODO: copy tesselation */
185
new_section->parent_element_num = this_section->parent_element_num;
186
new_section->_parent_element_num = NULL;
188
if (this_section->global_element_num != NULL) {
190
= fvm_io_num_get_local_count(this_section->global_element_num);
191
fvm_gnum_t global_count
192
= fvm_io_num_get_global_count(this_section->global_element_num);
193
const fvm_gnum_t *global_num
194
= fvm_io_num_get_global_num(this_section->global_element_num);
196
new_section->global_element_num
197
= fvm_io_num_create_shared(global_num, global_count, n_ent);
200
new_section->global_element_num = NULL;
202
return (new_section);
205
/*----------------------------------------------------------------------------
206
* Reduction of a nodal mesh representation section: only the associations
207
* (numberings) necessary to redistribution of fields for output are
208
* conserved, the full connectivity being no longer useful once it has been
212
* this_section <-> pointer to structure that should be reduced
215
* true if connectivity has been reduced
216
*----------------------------------------------------------------------------*/
219
_fvm_nodal_section_reduce(fvm_nodal_section_t * this_section)
221
_Bool retval = false;
223
/* If we have a tesselation of polyhedra (face index != NULL),
224
we may need to keep the connectivity information, to
225
interpolate nodal values to added vertices */
227
if ( this_section->tesselation == NULL
228
|| this_section->_face_index == NULL) {
232
this_section->connectivity_size = 0;
234
if (this_section->_face_index != NULL)
235
BFT_FREE(this_section->_face_index);
236
this_section->face_index = NULL;
238
if (this_section->_face_num != NULL)
239
BFT_FREE(this_section->_face_num);
240
this_section->face_num = NULL;
242
if (this_section->_vertex_index != NULL)
243
BFT_FREE(this_section->_vertex_index);
244
this_section->vertex_index = NULL;
246
if (this_section->_vertex_num != NULL)
247
BFT_FREE(this_section->_vertex_num);
248
this_section->vertex_num = NULL;
253
if (this_section->gc_id != NULL)
254
BFT_FREE(this_section->gc_id);
256
if (this_section->tesselation != NULL)
257
fvm_tesselation_reduce(this_section->tesselation);
262
/*----------------------------------------------------------------------------
263
* Change entity parent numbering; this is useful when entities of the
264
* parent mesh have been renumbered after a nodal mesh representation
265
* structure's creation. As the parent_num[] array is defined only when
266
* non trivial (i.e. not 1, 2, ..., n), it may be allocated or freed
267
* by this function. The return argument corresponds to the new
268
* pointer which should replace the parent_num input argument.
271
* parent_num_size <-- size of local parent numbering array
272
* new_parent_num <-- pointer to local parent renumbering array
273
* ({1, ..., n} <-- {1, ..., n})
274
* parent_num <-> pointer to local parent numbering array
275
* _parent_num <-> pointer to local parent numbering array if
276
* owner, NULL otherwise
279
* pointer to resulting parent_num[] array
280
*----------------------------------------------------------------------------*/
283
_renumber_parent_num(fvm_lnum_t parent_num_size,
284
const fvm_lnum_t new_parent_num[],
285
const fvm_lnum_t parent_num[],
286
fvm_lnum_t _parent_num[])
289
fvm_lnum_t old_num_id;
290
fvm_lnum_t *parent_num_p = _parent_num;
291
_Bool trivial = true;
293
if (parent_num_size > 0 && new_parent_num != NULL) {
295
if (parent_num_p != NULL) {
296
for (i = 0; i < parent_num_size; i++) {
297
old_num_id = parent_num_p[i] - 1;
298
parent_num_p[i] = new_parent_num[old_num_id];
299
if (parent_num_p[i] != i+1)
304
BFT_MALLOC(parent_num_p, parent_num_size, fvm_lnum_t);
305
if (parent_num != NULL) {
306
for (i = 0; i < parent_num_size; i++) {
307
old_num_id = parent_num[i] - 1;
308
parent_num_p[i] = new_parent_num[old_num_id];
309
if (parent_num_p[i] != i+1)
314
for (i = 0; i < parent_num_size; i++) {
315
parent_num_p[i] = new_parent_num[i];
316
if (parent_num_p[i] != i+1)
324
BFT_FREE(parent_num_p);
329
/*----------------------------------------------------------------------------
330
* Renumber vertices based on those actually referenced, and update
331
* connectivity arrays and parent numbering in accordance.
333
* The number of vertices assigned to the nodal mesh (this_nodal->n_vertices)
334
* is computed and set by this function. If this number was previously
335
* non-zero (i.e. vertices have already been assigned to the structure),
336
* those vertices are considered as referenced. This is useful if we want
337
* to avoid discarding a given set of vertices, such as when building a
338
* nodal mesh representation containing only vertices.
341
* this_nodal <-> nodal mesh structure
342
*----------------------------------------------------------------------------*/
345
_renumber_vertices(fvm_nodal_t *this_nodal)
350
fvm_lnum_t vertex_id;
351
fvm_lnum_t n_vertices;
352
fvm_nodal_section_t *section;
354
fvm_lnum_t *loc_vertex_num = NULL;
355
fvm_lnum_t max_vertex_num = 0;
357
/* Find maximum vertex reference */
358
/*-------------------------------*/
360
/* The mesh may already contain direct vertex references
361
(as in the case of a "mesh" only containing vertices) */
363
if (this_nodal->n_vertices > 0) {
364
if (this_nodal->parent_vertex_num != NULL) {
365
for (j = 0; j < this_nodal->n_vertices; j++) {
366
if (this_nodal->parent_vertex_num[j] > max_vertex_num)
367
max_vertex_num = this_nodal->parent_vertex_num[j];
371
max_vertex_num = this_nodal->n_vertices;
374
/* In most cases, the mesh will reference vertices through elements */
376
for (section_id = 0; section_id < this_nodal->n_sections; section_id++) {
377
section = this_nodal->sections[section_id];
378
if (this_nodal->parent_vertex_num != NULL) {
379
for (i = 0; i < section->connectivity_size; i++) {
380
fvm_lnum_t vertex_num
381
= this_nodal->parent_vertex_num[section->vertex_num[i] - 1];
382
if (vertex_num > max_vertex_num)
383
max_vertex_num = vertex_num;
387
for (i = 0; i < section->connectivity_size; i++) {
388
if (section->vertex_num[i] > max_vertex_num)
389
max_vertex_num = section->vertex_num[i];
394
/* Flag referenced vertices and compute size */
395
/*-------------------------------------------*/
397
BFT_MALLOC(loc_vertex_num, max_vertex_num, fvm_lnum_t);
399
for (vertex_id = 0; vertex_id < max_vertex_num; vertex_id++)
400
loc_vertex_num[vertex_id] = 0;
402
if (this_nodal->n_vertices > 0) {
403
if (this_nodal->parent_vertex_num != NULL) {
404
for (j = 0; j < this_nodal->n_vertices; j++) {
405
vertex_id = this_nodal->parent_vertex_num[j] - 1;
406
if (loc_vertex_num[vertex_id] == 0)
407
loc_vertex_num[vertex_id] = 1;
411
for (j = 0; j < this_nodal->n_vertices; j++) {
412
if (loc_vertex_num[j] == 0)
413
loc_vertex_num[j] = 1;
418
for (section_id = 0; section_id < this_nodal->n_sections; section_id++) {
419
section = this_nodal->sections[section_id];
420
if (this_nodal->parent_vertex_num != NULL) {
421
for (i = 0; i < section->connectivity_size; i++) {
423
= this_nodal->parent_vertex_num[section->vertex_num[i] - 1] - 1;
424
if (loc_vertex_num[vertex_id] == 0)
425
loc_vertex_num[vertex_id] = 1;
429
for (i = 0; i < section->connectivity_size; i++) {
430
vertex_id = section->vertex_num[i] - 1;
431
if (loc_vertex_num[vertex_id] == 0)
432
loc_vertex_num[vertex_id] = 1;
437
/* Build vertices renumbering */
438
/*----------------------------*/
442
for (vertex_id = 0; vertex_id < max_vertex_num; vertex_id++) {
443
if (loc_vertex_num[vertex_id] == 1) {
445
loc_vertex_num[vertex_id] = n_vertices;
448
this_nodal->n_vertices = n_vertices;
450
/* Update connectivity and vertex parent numbering */
451
/*-------------------------------------------------*/
453
/* If all vertices are flagged, no need to renumber */
455
if (n_vertices == max_vertex_num)
456
BFT_FREE(loc_vertex_num);
460
/* Update connectivity */
462
for (section_id = 0; section_id < this_nodal->n_sections; section_id++) {
463
section = this_nodal->sections[section_id];
464
if (section->_vertex_num == NULL)
465
fvm_nodal_section_copy_on_write(section, false, false, false, true);
466
if (this_nodal->parent_vertex_num != NULL) {
467
for (i = 0; i < section->connectivity_size; i++) {
469
= this_nodal->parent_vertex_num[section->vertex_num[i] - 1] - 1;
470
section->_vertex_num[i] = loc_vertex_num[vertex_id];
474
for (i = 0; i < section->connectivity_size; i++) {
475
vertex_id = section->vertex_num[i] - 1;
476
section->_vertex_num[i] = loc_vertex_num[vertex_id];
481
/* Build or update vertex parent numbering */
483
this_nodal->parent_vertex_num = NULL;
484
if (this_nodal->_parent_vertex_num != NULL)
485
BFT_FREE(this_nodal->_parent_vertex_num);
487
if (loc_vertex_num != NULL) {
488
BFT_MALLOC(this_nodal->_parent_vertex_num, n_vertices, fvm_lnum_t);
489
for (vertex_id = 0; vertex_id < max_vertex_num; vertex_id++) {
490
if (loc_vertex_num[vertex_id] > 0)
491
this_nodal->_parent_vertex_num[loc_vertex_num[vertex_id] - 1]
494
this_nodal->parent_vertex_num = this_nodal->_parent_vertex_num;
498
/* Free renumbering array */
500
BFT_FREE(loc_vertex_num);
503
/*----------------------------------------------------------------------------
504
* Dump printout of a nodal representation structure section.
507
* this_section <-- pointer to structure that should be dumped
508
*----------------------------------------------------------------------------*/
511
_fvm_nodal_section_dump(const fvm_nodal_section_t *this_section)
513
fvm_lnum_t n_elements, i, j;
514
const fvm_lnum_t *idx, *num;
516
/* Global indicators */
517
/*--------------------*/
520
"Entity dimension: %d\n"
521
"Number of elements: %ld\n"
522
"Element type: %s\n",
523
this_section->entity_dim, (long)this_section->n_elements,
524
fvm_elements_type_name[this_section->type]);
527
"Connectivity_size: %llu\n"
529
"Number of faces: %d\n",
530
(unsigned long long)(this_section->connectivity_size),
531
this_section->stride,
532
(long)(this_section->n_faces));
535
"Pointers to shareable arrays:\n"
538
" vertex_index: %p\n"
540
" parent_element_num: %p\n",
541
this_section->face_index, this_section->face_num,
542
this_section->vertex_index, this_section->vertex_num,
543
this_section->parent_element_num);
546
"Pointers to local arrays:\n"
549
" _vertex_index: %p\n"
551
" _parent_element_num: %p\n",
553
this_section->_face_index, this_section->_face_num,
554
this_section->_vertex_index, this_section->_vertex_num,
555
this_section->_parent_element_num, this_section->gc_id);
557
if (this_section->face_index != NULL) {
558
bft_printf("\nPolyhedra -> faces (polygons) connectivity:\n\n");
559
n_elements = this_section->n_elements;
560
idx = this_section->face_index;
561
num = this_section->face_num;
562
for (i = 0; i < n_elements; i++) {
563
bft_printf("%10d (idx = %10d) %10d\n",
564
i+1, idx[i], num[idx[i]]);
565
for (j = idx[i] + 1; j < idx[i + 1]; j++)
566
bft_printf(" %10d\n", num[j]);
568
bft_printf(" end (idx = %10d)\n", idx[n_elements]);
571
if (this_section->vertex_index != NULL) {
572
fvm_lnum_t n_faces = (this_section->n_faces > 0) ?
573
this_section->n_faces : this_section->n_elements;
574
bft_printf("\nPolygons -> vertices connectivity:\n\n");
575
idx = this_section->vertex_index;
576
num = this_section->vertex_num;
577
for (i = 0; i < n_faces; i++) {
578
bft_printf("%10d (idx = %10d) %10d\n",
579
i + 1, idx[i], num[idx[i]]);
580
for (j = idx[i] + 1; j < idx[i + 1]; j++)
581
bft_printf(" %10d\n", num[j]);
583
bft_printf(" end (idx = %10d)\n", idx[n_faces]);
587
bft_printf("\nElements -> vertices connectivity:\n\n");
588
n_elements = this_section->n_elements;
589
num = this_section->vertex_num;
590
switch(this_section->stride) {
592
for (i = 0; i < n_elements; i++)
593
bft_printf("%10d : %10d %10d\n",
594
i+1, num[i*2], num[i*2+1]);
597
for (i = 0; i < n_elements; i++)
598
bft_printf("%10d : %10d %10d %10d\n",
599
i+1, num[i*3], num[i*3+1], num[i*3+2]);
602
for (i = 0; i < n_elements; i++)
603
bft_printf("%10d : %10d %10d %10d %10d\n",
604
i+1, num[i*4], num[i*4+1], num[i*4+2],
608
for (i = 0; i < n_elements; i++)
609
bft_printf("%10d : %10d %10d %10d %10d %10d\n",
610
i+1, num[i*5], num[i*5+1], num[i*5+2],
611
num[i*5+3], num[i*5+4]);
614
for (i = 0; i < n_elements; i++)
615
bft_printf("%10d : %10d %10d %10d %10d %10d %10d\n",
616
i+1, num[i*6], num[i*6+1], num[i*6+2],
617
num[i*6+3], num[i*6+4], num[i*6+5]);
620
for (i = 0; i < n_elements; i++)
621
bft_printf("%10d : %10d %10d %10d %10d %10d %10d %10d %10d\n",
622
i+1, num[i*8], num[i*8+1], num[i*8+2], num[i*8+3],
623
num[i*8+4], num[i*8+5], num[i*8+6], num[i*8+7]);
626
for (i = 0; i < n_elements; i++) {
627
bft_printf("%10d :", i+1);
628
for (j = 0; j < this_section->stride; j++)
629
bft_printf(" %10d", num[i*this_section->stride + j]);
635
if (this_section->gc_id != NULL) {
636
bft_printf("\nGroup class ids:\n\n");
637
for (i = 0; i < this_section->n_elements; i++)
638
bft_printf("%10d : %10d\n", i + 1, this_section->gc_id[i]);
642
/* Faces tesselation */
644
if (this_section->tesselation != NULL)
645
fvm_tesselation_dump(this_section->tesselation);
647
/* Numbers of associated elements in the parent mesh */
649
bft_printf("\nLocal element numbers in parent mesh:\n");
650
if (this_section->parent_element_num == NULL)
651
bft_printf("\n Nil\n\n");
653
for (i = 0; i < this_section->n_elements; i++)
654
bft_printf(" %10d %10d\n", i+1, this_section->parent_element_num[i]);
657
/* Global element numbers (only for parallel execution) */
659
if (this_section->global_element_num != NULL) {
660
bft_printf("\nGlobal element numbers:\n");
661
fvm_io_num_dump(this_section->global_element_num);
665
/*============================================================================
666
* Semi-private function definitions (prototypes in fvm_nodal_priv.h)
667
*============================================================================*/
669
/*----------------------------------------------------------------------------
670
* Creation of a nodal mesh section representation structure.
673
* type <-- type of element defined by this section
676
* pointer to created nodal mesh section representation structure
677
*----------------------------------------------------------------------------*/
679
fvm_nodal_section_t *
680
fvm_nodal_section_create(const fvm_element_t type)
682
fvm_nodal_section_t *this_section;
684
BFT_MALLOC(this_section, 1, fvm_nodal_section_t);
686
/* Global information */
688
if (type == FVM_EDGE)
689
this_section->entity_dim = 1;
690
else if (type >= FVM_FACE_TRIA && type <= FVM_FACE_POLY)
691
this_section->entity_dim = 2;
693
this_section->entity_dim = 3;
695
this_section->n_elements = 0;
696
this_section->type = type;
700
this_section->connectivity_size = 0;
702
if (type != FVM_FACE_POLY && type != FVM_CELL_POLY)
703
this_section->stride = fvm_nodal_n_vertices_element[type];
705
this_section->stride = 0;
707
this_section->n_faces = 0;
708
this_section->face_index = NULL;
709
this_section->face_num = NULL;
710
this_section->vertex_index = NULL;
711
this_section->vertex_num = NULL;
713
this_section->_face_index = NULL;
714
this_section->_face_num = NULL;
715
this_section->_vertex_index = NULL;
716
this_section->_vertex_num = NULL;
718
this_section->gc_id = NULL;
720
this_section->tesselation = NULL;
725
this_section->parent_element_num = NULL;
726
this_section->_parent_element_num = NULL;
728
this_section->global_element_num = NULL;
730
return (this_section);
733
/*----------------------------------------------------------------------------
734
* Destruction of a nodal mesh section representation structure.
737
* this_section <-> pointer to structure that should be destroyed
741
*----------------------------------------------------------------------------*/
743
fvm_nodal_section_t *
744
fvm_nodal_section_destroy(fvm_nodal_section_t * this_section)
748
if (this_section->_face_index != NULL)
749
BFT_FREE(this_section->_face_index);
750
if (this_section->_face_num != NULL)
751
BFT_FREE(this_section->_face_num);
753
if (this_section->_vertex_index != NULL)
754
BFT_FREE(this_section->_vertex_index);
755
if (this_section->_vertex_num != NULL)
756
BFT_FREE(this_section->_vertex_num);
758
if (this_section->gc_id != NULL)
759
BFT_FREE(this_section->gc_id);
761
if (this_section->tesselation != NULL)
762
fvm_tesselation_destroy(this_section->tesselation);
767
if (this_section->parent_element_num != NULL) {
768
this_section->parent_element_num = NULL;
769
BFT_FREE(this_section->_parent_element_num);
772
if (this_section->global_element_num != NULL)
773
fvm_io_num_destroy(this_section->global_element_num);
775
/* Main structure destroyed and NULL returned */
777
BFT_FREE(this_section);
779
return (this_section);
782
/*----------------------------------------------------------------------------
783
* Copy selected shared connectivity information to private connectivity
784
* for a nodal mesh section.
787
* this_section <-> pointer to section structure
788
* copy_face_index <-- copy face index (polyhedra only) ?
789
* copy_face_num <-- copy face numbers (polyhedra only) ?
790
* copy_vertex_index <-- copy vertex index (polyhedra/polygons only) ?
791
* copy_vertex_num <-- copy vertex numbers ?
792
*----------------------------------------------------------------------------*/
795
fvm_nodal_section_copy_on_write(fvm_nodal_section_t *this_section,
796
_Bool copy_face_index,
798
_Bool copy_vertex_index,
799
_Bool copy_vertex_num)
804
if (copy_face_index == true
805
&& this_section->face_index != NULL && this_section->_face_index == NULL) {
806
BFT_MALLOC(this_section->_face_index, this_section->n_elements + 1, fvm_lnum_t);
807
for (i = 0; i < (size_t)(this_section->n_elements + 1); i++) {
808
this_section->_face_index[i] = this_section->face_index[i];
810
this_section->face_index = this_section->_face_index;
813
if (copy_face_num == true
814
&& this_section->face_num != NULL && this_section->_face_num == NULL) {
815
n_faces = this_section->face_index[this_section->n_elements];
816
BFT_MALLOC(this_section->_face_num, n_faces, fvm_lnum_t);
817
for (i = 0; i < (size_t)n_faces; i++) {
818
this_section->_face_num[i] = this_section->face_num[i];
820
this_section->face_num = this_section->_face_num;
823
if ( copy_vertex_index == true
824
&& this_section->vertex_index != NULL
825
&& this_section->_vertex_index == NULL) {
826
if (this_section->n_faces != 0)
827
n_faces = this_section->n_faces;
829
n_faces = this_section->n_elements;
830
BFT_MALLOC(this_section->_vertex_index, n_faces + 1, fvm_lnum_t);
831
for (i = 0; i < (size_t)n_faces + 1; i++) {
832
this_section->_vertex_index[i] = this_section->vertex_index[i];
834
this_section->vertex_index = this_section->_vertex_index;
837
if (copy_vertex_num == true && this_section->_vertex_num == NULL) {
838
BFT_MALLOC(this_section->_vertex_num,
839
this_section->connectivity_size, fvm_lnum_t);
840
for (i = 0; i < this_section->connectivity_size; i++) {
841
this_section->_vertex_num[i] = this_section->vertex_num[i];
843
this_section->vertex_num = this_section->_vertex_num;
848
/*----------------------------------------------------------------------------
849
* Return global number of elements associated with section.
852
* this_section <-- pointer to section structure
855
* global number of elements associated with section
856
*----------------------------------------------------------------------------*/
859
fvm_nodal_section_n_g_elements(const fvm_nodal_section_t *this_section)
861
if (this_section->global_element_num != NULL)
862
return fvm_io_num_get_global_count(this_section->global_element_num);
864
return this_section->n_elements;
867
/*----------------------------------------------------------------------------
868
* Return global number of vertices associated with nodal mesh.
871
* this_nodal <-- pointer to nodal mesh structure
874
* global number of vertices associated with nodal mesh
875
*----------------------------------------------------------------------------*/
878
fvm_nodal_n_g_vertices(const fvm_nodal_t *this_nodal)
880
fvm_gnum_t n_g_vertices;
882
if (this_nodal->global_vertex_num != NULL)
883
n_g_vertices = fvm_io_num_get_global_count(this_nodal->global_vertex_num);
885
n_g_vertices = this_nodal->n_vertices;
890
/*----------------------------------------------------------------------------
891
* Define cell->face connectivity for strided cell types.
894
* element_type <-- type of strided element
895
* n_faces --> number of element faces
896
* n_face_vertices --> number of vertices of each face
897
* face_vertices --> face -> vertex base connectivity (0 to n-1)
898
*----------------------------------------------------------------------------*/
901
fvm_nodal_cell_face_connect(fvm_element_t element_type,
903
int n_face_vertices[6],
904
int face_vertices[6][4])
912
for (i = 0; i < 6; i++) {
913
n_face_vertices[i] = 0;
914
for (j = 0; j < 4; j++)
915
face_vertices[i][j] = 0;
918
/* Define connectivity based on element type */
920
switch(element_type) {
924
fvm_lnum_t _face_vertices[4][3] = {{1, 3, 2}, /* x 4 */
926
{1, 4, 3}, /* / | \ */
927
{2, 3, 4}}; /* / | \ */
929
for (i = 0; i < 4; i++) { /* \ | / */
930
n_face_vertices[i] = 3; /* \ | / */
931
for (j = 0; j < 3; j++) /* \|/ */
932
face_vertices[i][j] = _face_vertices[i][j]; /* x 2 */
940
fvm_lnum_t _n_face_vertices[5] = {3, 3, 3, 3, 4};
941
fvm_lnum_t _face_vertices[5][4] = {{1, 2, 5, 0}, /* 5 x */
942
{1, 5, 4, 0}, /* /|\ */
943
{2, 3, 5, 0}, /* //| \ */
944
{3, 4, 5, 0}, /* // | \ */
945
{1, 4, 3, 2}}; /* 4 x/--|---x 3 */
947
for (i = 0; i < 5; i++) { /* // | / */
948
n_face_vertices[i] = _n_face_vertices[i]; /* 1 x-------x 2 */
949
for (j = 0; j < 4; j++)
950
face_vertices[i][j] = _face_vertices[i][j];
958
fvm_lnum_t _n_face_vertices[5] = {3, 3, 4, 4, 4};
959
fvm_lnum_t _face_vertices[5][4] = {{1, 3, 2, 0}, /* 4 x-------x 6 */
960
{4, 5, 6, 0}, /* |\ /| */
961
{1, 2, 5, 4}, /* | \ / | */
962
{1, 4, 6, 3}, /* 1 x- \-/ -x 3 */
963
{2, 3, 6, 5}}; /* \ 5x / */
965
for (i = 0; i < 5; i++) { /* \|/ */
966
n_face_vertices[i] = _n_face_vertices[i]; /* x 2 */
967
for (j = 0; j < 4; j++)
968
face_vertices[i][j] = _face_vertices[i][j];
976
fvm_lnum_t _n_face_vertices[6] = {4, 4, 4, 4, 4, 4};
977
fvm_lnum_t _face_vertices[6][4] = {{1, 4, 3, 2}, /* 8 x-------x 7 */
978
{1, 2, 6, 5}, /* /| /| */
979
{1, 5, 8, 4}, /* / | / | */
980
{2, 3, 7, 6}, /* 5 x-------x6 | */
981
{3, 4, 8, 7}, /* | 4x----|--x 3 */
982
{5, 6, 7, 8}}; /* | / | / */
983
for (i = 0; i < 6; i++) { /* |/ |/ */
984
n_face_vertices[i] = _n_face_vertices[i]; /* 1 x-------x 2 */
985
for (j = 0; j < 4; j++)
986
face_vertices[i][j] = _face_vertices[i][j];
996
/* Switch from (1, n) to (0, n-1) numbering */
998
for (i = 0; i < 6; i++) {
999
for (j = 0; j < 4; j++)
1000
face_vertices[i][j] -= 1;
1004
/*============================================================================
1005
* Public function definitions
1006
*============================================================================*/
1008
/*----------------------------------------------------------------------------
1009
* Creation of a nodal mesh representation structure.
1012
* name <-- name that should be assigned to the nodal mesh
1013
* dim <-- spatial dimension
1016
* pointer to created nodal mesh representation structure
1017
*----------------------------------------------------------------------------*/
1020
fvm_nodal_create(const char *name,
1023
fvm_nodal_t *this_nodal;
1025
BFT_MALLOC(this_nodal, 1, fvm_nodal_t);
1027
/* Global indicators */
1030
BFT_MALLOC(this_nodal->name, strlen(name) + 1, char);
1031
strcpy(this_nodal->name, name);
1034
this_nodal->name = NULL;
1036
this_nodal->dim = dim;
1037
this_nodal->num_dom = fvm_parall_get_rank() + 1;
1038
this_nodal->n_doms = fvm_parall_get_size();
1039
this_nodal->n_sections = 0;
1041
/* Local dimensions */
1043
this_nodal->n_cells = 0;
1044
this_nodal->n_faces = 0;
1045
this_nodal->n_edges = 0;
1046
this_nodal->n_vertices = 0;
1048
/* Local structures */
1050
this_nodal->vertex_coords = NULL;
1051
this_nodal->_vertex_coords = NULL;
1053
this_nodal->parent_vertex_num = NULL;
1054
this_nodal->_parent_vertex_num = NULL;
1056
this_nodal->global_vertex_num = NULL;
1058
this_nodal->sections = NULL;
1060
this_nodal->gc_set = NULL;
1062
return (this_nodal);
1065
/*----------------------------------------------------------------------------
1066
* Destruction of a nodal mesh representation structure.
1069
* this_nodal <-> pointer to structure that should be destroyed
1073
*----------------------------------------------------------------------------*/
1076
fvm_nodal_destroy(fvm_nodal_t * this_nodal)
1080
/* Local structures */
1082
if (this_nodal->name != NULL)
1083
BFT_FREE(this_nodal->name);
1085
if (this_nodal->_vertex_coords != NULL)
1086
BFT_FREE(this_nodal->_vertex_coords);
1088
if (this_nodal->parent_vertex_num != NULL) {
1089
this_nodal->parent_vertex_num = NULL;
1090
BFT_FREE(this_nodal->_parent_vertex_num);
1093
if (this_nodal->global_vertex_num != NULL)
1094
fvm_io_num_destroy(this_nodal->global_vertex_num);
1096
for (i = 0; i < this_nodal->n_sections; i++)
1097
fvm_nodal_section_destroy(this_nodal->sections[i]);
1099
if (this_nodal->sections != NULL)
1100
BFT_FREE(this_nodal->sections);
1102
if (this_nodal->gc_set != NULL)
1103
this_nodal->gc_set = fvm_group_class_set_destroy(this_nodal->gc_set);
1105
/* Main structure destroyed and NULL returned */
1107
BFT_FREE(this_nodal);
1109
return (this_nodal);
1112
/*----------------------------------------------------------------------------
1113
* Copy a nodal mesh representation structure, sharing arrays with the
1114
* original structure.
1116
* Element group classes and mesh group class descriptions are not currently
1120
* this_nodal <-> pointer to structure that should be copied
1123
* pointer to created nodal mesh representation structure
1124
*----------------------------------------------------------------------------*/
1127
fvm_nodal_copy(const fvm_nodal_t *this_nodal)
1130
fvm_nodal_t *new_nodal;
1132
BFT_MALLOC(new_nodal, 1, fvm_nodal_t);
1134
/* Global indicators */
1136
if (this_nodal->name != NULL) {
1137
BFT_MALLOC(new_nodal->name, strlen(this_nodal->name) + 1, char);
1138
strcpy(new_nodal->name, this_nodal->name);
1141
new_nodal->name = NULL;
1143
new_nodal->dim = this_nodal->dim;
1144
new_nodal->num_dom = this_nodal->num_dom;
1145
new_nodal->n_doms = this_nodal->n_doms;
1146
new_nodal->n_sections = this_nodal->n_sections;
1148
/* Local dimensions */
1150
new_nodal->n_cells = this_nodal->n_cells;
1151
new_nodal->n_faces = this_nodal->n_faces;
1152
new_nodal->n_edges = this_nodal->n_edges;
1153
new_nodal->n_vertices = this_nodal->n_vertices;
1155
/* Local structures */
1157
new_nodal->vertex_coords = this_nodal->vertex_coords;
1158
new_nodal->_vertex_coords = NULL;
1160
new_nodal->parent_vertex_num = this_nodal->parent_vertex_num;
1161
new_nodal->_parent_vertex_num = NULL;
1163
if (this_nodal->global_vertex_num != NULL) {
1165
= fvm_io_num_get_local_count(this_nodal->global_vertex_num);
1166
fvm_gnum_t global_count
1167
= fvm_io_num_get_global_count(this_nodal->global_vertex_num);
1168
const fvm_gnum_t *global_num
1169
= fvm_io_num_get_global_num(this_nodal->global_vertex_num);
1171
new_nodal->global_vertex_num
1172
= fvm_io_num_create_shared(global_num, global_count, n_ent);
1175
new_nodal->global_vertex_num = NULL;
1177
BFT_MALLOC(new_nodal->sections,
1178
new_nodal->n_sections,
1179
fvm_nodal_section_t *);
1180
for (i = 0; i < new_nodal->n_sections; i++)
1181
new_nodal->sections[i] = _fvm_nodal_section_copy(this_nodal->sections[i]);
1183
new_nodal->gc_set = NULL;
1188
/*----------------------------------------------------------------------------
1189
* Reduction of a nodal mesh representation structure: only the associations
1190
* (numberings) necessary to redistribution of fields for output are
1191
* conserved, the full connectivity being in many cases no longer useful
1192
* once it has been output. If the del_vertex_num value is set
1193
* to true, vertex-based values may not be output in parallel mode
1194
* after this function is called.
1197
* this_nodal <-> pointer to structure that should be reduced
1198
* del_vertex_num <-- indicates if vertex parent indirection and
1199
* I/O numbering are destroyed (1) or not (0)
1200
*----------------------------------------------------------------------------*/
1203
fvm_nodal_reduce(fvm_nodal_t *this_nodal,
1207
_Bool reduce_vertices = true;
1211
for (i = 0; i < this_nodal->n_sections; i++) {
1212
if (_fvm_nodal_section_reduce(this_nodal->sections[i]) == false)
1213
reduce_vertices = false;
1218
if (reduce_vertices == true) {
1220
if (this_nodal->_vertex_coords != NULL)
1221
BFT_FREE(this_nodal->_vertex_coords);
1222
this_nodal->vertex_coords = NULL;
1226
/* Depending on this option, output on vertices may not remain possible */
1228
if (del_vertex_num > 0) {
1230
if (this_nodal->parent_vertex_num != NULL) {
1231
this_nodal->parent_vertex_num = NULL;
1232
BFT_FREE(this_nodal->_parent_vertex_num);
1235
if (this_nodal->global_vertex_num != NULL)
1236
this_nodal->global_vertex_num
1237
= fvm_io_num_destroy(this_nodal->global_vertex_num);
1241
if (this_nodal->gc_set != NULL)
1242
this_nodal->gc_set = fvm_group_class_set_destroy(this_nodal->gc_set);
1245
/*----------------------------------------------------------------------------
1246
* Change entity parent numbering; this is useful when entities of the
1247
* parent mesh have been renumbered after a nodal mesh representation
1248
* structure's creation.
1251
* this_nodal <-- nodal mesh structure
1252
* new_parent_num <-- pointer to local parent renumbering array
1253
* ({1, ..., n} <-- {1, ..., n})
1254
* entity_dim <-- 3 for cells, 2 for faces, 1 for edges,
1255
* and 0 for vertices
1256
*----------------------------------------------------------------------------*/
1259
fvm_nodal_change_parent_num(fvm_nodal_t *this_nodal,
1260
const fvm_lnum_t new_parent_num[],
1265
if (entity_dim == 0) {
1267
this_nodal->_parent_vertex_num
1268
= _renumber_parent_num(this_nodal->n_vertices,
1270
this_nodal->parent_vertex_num,
1271
this_nodal->_parent_vertex_num);
1272
this_nodal->parent_vertex_num = this_nodal->_parent_vertex_num;
1276
/* Other elements */
1281
fvm_nodal_section_t *section = NULL;
1283
for (i = 0; i < this_nodal->n_sections; i++) {
1284
section = this_nodal->sections[i];
1285
if (section->entity_dim == entity_dim) {
1286
section->_parent_element_num
1287
= _renumber_parent_num(section->n_elements,
1289
section->parent_element_num,
1290
section->_parent_element_num);
1291
section->parent_element_num = section->_parent_element_num;
1299
/*----------------------------------------------------------------------------
1300
* Remove entity parent numbering; this is useful for example when we
1301
* want to assign coordinates or fields to an extracted mesh using
1302
* arrays relative to the mesh, and not to its parent.
1304
* This is equivalent to calling fvm_nodal_change_parent_num(), with
1305
* 'trivial' (1 o n) new_parent_num[] values.
1308
* this_nodal <-- nodal mesh structure
1309
* entity_dim <-- 3 for cells, 2 for faces, 1 for edges,
1310
* and 0 for vertices
1311
*----------------------------------------------------------------------------*/
1314
fvm_nodal_remove_parent_num(fvm_nodal_t *this_nodal,
1319
if (entity_dim == 0) {
1320
this_nodal->parent_vertex_num = NULL;
1321
if (this_nodal->_parent_vertex_num != NULL)
1322
BFT_FREE(this_nodal->_parent_vertex_num);
1325
/* Other elements */
1330
fvm_nodal_section_t *section = NULL;
1332
for (i = 0; i < this_nodal->n_sections; i++) {
1333
section = this_nodal->sections[i];
1334
if (section->entity_dim == entity_dim) {
1335
section->parent_element_num = NULL;
1336
if (section->_parent_element_num != NULL)
1337
BFT_FREE(section->_parent_element_num);
1345
/*----------------------------------------------------------------------------
1346
* Build external numbering for entities based on global numbers.
1349
* this_nodal <-- nodal mesh structure
1350
* parent_global_number <-- pointer to list of global (i.e. domain splitting
1351
* independent) parent entity numbers
1352
* entity_dim <-- 3 for cells, 2 for faces, 1 for edges,
1353
* and 0 for vertices
1354
*----------------------------------------------------------------------------*/
1357
fvm_nodal_init_io_num(fvm_nodal_t *this_nodal,
1358
const fvm_gnum_t parent_global_numbers[],
1362
fvm_nodal_section_t *section;
1364
if (entity_dim == 0)
1365
this_nodal->global_vertex_num
1366
= fvm_io_num_create(this_nodal->parent_vertex_num,
1367
parent_global_numbers,
1368
this_nodal->n_vertices,
1372
for (i = 0; i < this_nodal->n_sections; i++) {
1373
section = this_nodal->sections[i];
1374
if (section->entity_dim == entity_dim) {
1375
section->global_element_num
1376
= fvm_io_num_create(section->parent_element_num,
1377
parent_global_numbers,
1378
section->n_elements,
1386
/*----------------------------------------------------------------------------
1387
* Preset number and list of vertices to assign to a nodal mesh.
1389
* If the parent_vertex_num argument is NULL, the list is assumed to
1390
* be {1, 2, ..., n}. If parent_vertex_num is given, it specifies a
1391
* list of n vertices from a larger set (1 to n numbering).
1393
* Ownership of the given parent vertex numbering array is
1394
* transferred to the nodal mesh representation structure.
1396
* This function should be called before fvm_nodal_set_shared_vertices()
1397
* or fvm_nodal_transfer_vertices() if we want to force certain
1398
* vertices to appear in the mesh (especially if we want to define
1399
* a mesh containing only vertices).
1402
* this_nodal <-> nodal mesh structure
1403
* n_vertices <-- number of vertices to assign
1404
* parent_vertex_num <-- parent numbers of vertices to assign
1405
*----------------------------------------------------------------------------*/
1408
fvm_nodal_define_vertex_list(fvm_nodal_t *this_nodal,
1409
fvm_lnum_t n_vertices,
1410
fvm_lnum_t parent_vertex_num[])
1412
assert(this_nodal != NULL);
1414
this_nodal->n_vertices = n_vertices;
1416
this_nodal->parent_vertex_num = NULL;
1417
if (this_nodal->_parent_vertex_num != NULL)
1418
BFT_FREE(this_nodal->_parent_vertex_num);
1420
if (parent_vertex_num != NULL) {
1421
this_nodal->_parent_vertex_num = parent_vertex_num;
1422
this_nodal->parent_vertex_num = parent_vertex_num;
1426
/*----------------------------------------------------------------------------
1427
* Assign shared vertex coordinates to an extracted nodal mesh,
1428
* renumbering vertex numbers based on those really referenced,
1429
* and updating connectivity arrays in accordance.
1431
* This function should only be called once all element sections
1432
* have been added to a nodal mesh representation.
1435
* this_nodal <-> nodal mesh structure
1436
* vertex_coords <-- coordinates of parent vertices (interlaced)
1437
*----------------------------------------------------------------------------*/
1440
fvm_nodal_set_shared_vertices(fvm_nodal_t *this_nodal,
1441
const fvm_coord_t vertex_coords[])
1443
assert(this_nodal != NULL);
1445
/* Map vertex coordinates to array passed as argument
1446
(this_nodal->_vertex_coords remains NULL, so only
1447
the const pointer may be used for a shared array) */
1449
this_nodal->vertex_coords = vertex_coords;
1451
/* If the mesh contains only vertices, its n_vertices and
1452
parent_vertex_num must already have been set, and do not
1455
if (this_nodal->n_sections == 0)
1458
/* Renumber vertices based on those really referenced */
1460
_renumber_vertices(this_nodal);
1464
/*----------------------------------------------------------------------------
1465
* Assign private vertex coordinates to a nodal mesh,
1466
* renumbering vertex numbers based on those really referenced,
1467
* and updating connectivity arrays in accordance.
1469
* Ownership of the given coordinates array is transferred to
1470
* the nodal mesh representation structure.
1472
* This function should only be called once all element sections
1473
* have been added to a nodal mesh representation.
1476
* this_nodal <-> nodal mesh structure
1477
* vertex_coords <-- coordinates of parent vertices (interlaced)
1480
* updated pointer to vertex_coords (may be different from initial
1481
* argument if vertices were renumbered).
1482
*----------------------------------------------------------------------------*/
1485
fvm_nodal_transfer_vertices(fvm_nodal_t *this_nodal,
1486
fvm_coord_t vertex_coords[])
1491
fvm_coord_t *_vertex_coords = vertex_coords;
1493
assert(this_nodal != NULL);
1495
/* Renumber vertices based on those really referenced, and
1496
update connectivity arrays in accordance. */
1498
_renumber_vertices(this_nodal);
1500
/* If renumbering is necessary, update connectivity */
1502
if (this_nodal->parent_vertex_num != NULL) {
1504
int dim = this_nodal->dim;
1505
const fvm_lnum_t *parent_vertex_num = this_nodal->parent_vertex_num;
1507
BFT_MALLOC(_vertex_coords, this_nodal->n_vertices * dim, fvm_coord_t);
1509
for (i = 0; i < this_nodal->n_vertices; i++) {
1510
for (j = 0; j < dim; j++)
1511
_vertex_coords[i*dim + j]
1512
= vertex_coords[(parent_vertex_num[i]-1)*dim + j];
1515
BFT_FREE(vertex_coords);
1517
this_nodal->parent_vertex_num = NULL;
1518
if (this_nodal->_parent_vertex_num != NULL)
1519
BFT_FREE(this_nodal->_parent_vertex_num);
1522
this_nodal->_vertex_coords = _vertex_coords;
1523
this_nodal->vertex_coords = _vertex_coords;
1525
return _vertex_coords;
1528
/*----------------------------------------------------------------------------
1529
* Make vertex coordinates of a nodal mesh private.
1531
* If vertex coordinates were previously shared, those coordinates that
1532
* are actually refernces are copied, and the relation to parent vertices
1535
* If vertices were already private, the mesh is not modified.
1538
* this_nodal <-> nodal mesh structure
1539
*----------------------------------------------------------------------------*/
1542
fvm_nodal_make_vertices_private(fvm_nodal_t *this_nodal)
1544
assert(this_nodal != NULL);
1546
if (this_nodal->_vertex_coords == NULL) {
1548
fvm_coord_t *_vertex_coords = NULL;
1549
const fvm_coord_t *vertex_coords = this_nodal->vertex_coords;
1550
const fvm_lnum_t n_vertices = this_nodal->n_vertices;
1551
const int dim = this_nodal->dim;
1553
BFT_MALLOC(vertex_coords, n_vertices * dim, fvm_coord_t);
1555
/* If renumbering is necessary, update connectivity */
1557
if (this_nodal->parent_vertex_num != NULL) {
1561
const fvm_lnum_t *parent_vertex_num = this_nodal->parent_vertex_num;
1563
for (i = 0; i < n_vertices; i++) {
1564
for (j = 0; j < dim; j++)
1565
_vertex_coords[i*dim + j]
1566
= vertex_coords[(parent_vertex_num[i]-1)*dim + j];
1569
this_nodal->parent_vertex_num = NULL;
1570
if (this_nodal->_parent_vertex_num != NULL)
1571
BFT_FREE(this_nodal->_parent_vertex_num);
1574
memcpy(_vertex_coords, vertex_coords, n_vertices*dim*sizeof(fvm_coord_t));
1576
/* Assign new array to structure */
1578
this_nodal->_vertex_coords = _vertex_coords;
1579
this_nodal->vertex_coords = _vertex_coords;
1583
/*----------------------------------------------------------------------------
1584
* Assign group class set descriptions to a nodal mesh.
1586
* The structure builds its own copy of the group class sets,
1587
* renumbering them so as to discard those not referenced.
1588
* Empty group classes are also renumbered to zero.
1590
* This function should only be called once all element sections
1591
* have been added to a nodal mesh representation.
1594
* this_nodal <-> nodal mesh structure
1595
* gc_set <-- group class set descriptions
1596
*----------------------------------------------------------------------------*/
1599
fvm_nodal_set_group_class_set(fvm_nodal_t *this_nodal,
1600
const fvm_group_class_set_t *gc_set)
1602
int gc_id, section_id;
1603
int n_gc = fvm_group_class_set_size(gc_set);
1605
fvm_lnum_t *gc_renum = NULL;
1607
assert(this_nodal != NULL);
1609
if (this_nodal->gc_set != NULL)
1610
this_nodal->gc_set = fvm_group_class_set_destroy(this_nodal->gc_set);
1615
/* Mark referenced group classes */
1617
BFT_MALLOC(gc_renum, n_gc, fvm_lnum_t);
1619
for (gc_id = 0; gc_id < n_gc; gc_id++)
1620
gc_renum[gc_id] = 0;
1622
for (section_id = 0; section_id < this_nodal->n_sections; section_id++) {
1624
const fvm_nodal_section_t *section = this_nodal->sections[section_id];
1625
if (section->gc_id == NULL)
1627
for (i = 0; i < section->n_elements; i++) {
1628
if (section->gc_id[i] != 0)
1629
gc_renum[section->gc_id[i] - 1] = 1;
1633
fvm_parall_counter_max(gc_renum, n_gc);
1635
/* Renumber group classes if necessary */
1637
for (gc_id = 0; gc_id < n_gc; gc_id++) {
1638
if (gc_renum[gc_id] != 0) {
1639
gc_renum[gc_id] = n_gc_new + 1;
1644
if (n_gc_new < n_gc) {
1645
for (section_id = 0; section_id < this_nodal->n_sections; section_id++) {
1647
const fvm_nodal_section_t *section = this_nodal->sections[section_id];
1648
if (section->gc_id == NULL)
1650
for (i = 0; i < section->n_elements; i++) {
1651
if (section->gc_id[i] != 0)
1652
section->gc_id[i] = gc_renum[section->gc_id[i] - 1];
1657
/* Transform renumbering array to list */
1660
for (gc_id = 0; gc_id < n_gc; gc_id++) {
1661
if (gc_renum[gc_id] != 0)
1662
gc_renum[n_gc_new++] = gc_id;
1666
this_nodal->gc_set = fvm_group_class_set_copy(gc_set,
1673
/*----------------------------------------------------------------------------
1674
* Obtain the name of a nodal mesh.
1677
* this_nodal <-- pointer to nodal mesh structure
1680
* pointer to constant string containing the mesh name
1681
*----------------------------------------------------------------------------*/
1684
fvm_nodal_get_name(const fvm_nodal_t *this_nodal)
1686
assert(this_nodal != NULL);
1688
return this_nodal->name;
1691
/*----------------------------------------------------------------------------
1692
* Return spatial dimension of the nodal mesh.
1695
* this_nodal <-- pointer to nodal mesh structure
1698
* spatial dimension.
1699
*----------------------------------------------------------------------------*/
1702
fvm_nodal_get_dim(const fvm_nodal_t *this_nodal)
1704
return this_nodal->dim;
1707
/*----------------------------------------------------------------------------
1708
* Return maximum dimension of entities in a nodal mesh.
1711
* this_nodal <-- pointer to nodal mesh structure
1714
* maximum dimension of entities in mesh (0 to 3)
1715
*----------------------------------------------------------------------------*/
1718
fvm_nodal_get_max_entity_dim(const fvm_nodal_t *this_nodal)
1721
int max_entity_dim = 0;
1723
assert(this_nodal != NULL);
1725
for (section_id = 0; section_id < this_nodal->n_sections; section_id++) {
1726
const fvm_nodal_section_t *section = this_nodal->sections[section_id];
1727
if (section->entity_dim > max_entity_dim)
1728
max_entity_dim = section->entity_dim;
1731
return max_entity_dim;
1734
/*----------------------------------------------------------------------------
1735
* Return number of entities of a given dimension in a nodal mesh.
1738
* this_nodal <-- pointer to nodal mesh structure
1739
* entity_dim <-- dimension of entities we want to count (0 to 3)
1742
* number of entities of given dimension in mesh
1743
*----------------------------------------------------------------------------*/
1746
fvm_nodal_get_n_entities(const fvm_nodal_t *this_nodal,
1749
fvm_lnum_t n_entities;
1751
assert(this_nodal != NULL);
1753
switch(entity_dim) {
1755
n_entities = this_nodal->n_vertices;
1758
n_entities = this_nodal->n_edges;
1761
n_entities = this_nodal->n_faces;
1764
n_entities = this_nodal->n_cells;
1773
/*----------------------------------------------------------------------------
1774
* Return global number of vertices associated with nodal mesh.
1777
* this_nodal <-- pointer to nodal mesh structure
1780
* global number of vertices associated with nodal mesh
1781
*----------------------------------------------------------------------------*/
1784
fvm_nodal_get_n_g_vertices(const fvm_nodal_t *this_nodal)
1786
return fvm_nodal_n_g_vertices(this_nodal);
1789
/*----------------------------------------------------------------------------
1790
* Return global number of elements of a given type associated with nodal mesh.
1793
* this_nodal <-- pointer to nodal mesh structure
1794
* element_type <-- type of elements for query
1797
* global number of elements of the given type associated with nodal mesh
1798
*----------------------------------------------------------------------------*/
1801
fvm_nodal_get_n_g_elements(const fvm_nodal_t *this_nodal,
1802
fvm_element_t element_type)
1805
fvm_gnum_t n_g_elements = 0;
1807
assert(this_nodal != NULL);
1809
for (i = 0; i < this_nodal->n_sections; i++) {
1810
fvm_nodal_section_t *section = this_nodal->sections[i];
1811
if (section->type == element_type)
1812
n_g_elements += fvm_nodal_section_n_g_elements(section);
1815
return n_g_elements;
1818
/*----------------------------------------------------------------------------
1819
* Return local number of elements of a given type associated with nodal mesh.
1822
* this_nodal <-- pointer to nodal mesh structure
1823
* element_type <-- type of elements for query
1826
* local number of elements of the given type associated with nodal mesh
1827
*----------------------------------------------------------------------------*/
1830
fvm_nodal_get_n_elements(const fvm_nodal_t *this_nodal,
1831
fvm_element_t element_type)
1834
fvm_lnum_t n_elements = 0;
1836
assert(this_nodal != NULL);
1838
for (i = 0; i < this_nodal->n_sections; i++) {
1839
fvm_nodal_section_t *section = this_nodal->sections[i];
1840
if (section->type == element_type)
1841
n_elements += section->n_elements;
1847
/*----------------------------------------------------------------------------
1848
* Return local parent numbering array for all entities of a given
1849
* dimension in a nodal mesh.
1851
* The number of entities of the given dimension may be obtained
1852
* through fvm_nodal_get_n_entities(), the parent_num[] array is populated
1853
* with the parent entity numbers of those entities, in order (i.e. in
1854
* local section order, section by section).
1857
* this_nodal <-- pointer to nodal mesh structure
1858
* entity_dim <-- dimension of entities we are interested in (0 to 3)
1859
* parent_num --> entity parent numbering (array must be pre-allocated)
1860
*----------------------------------------------------------------------------*/
1863
fvm_nodal_get_parent_num(const fvm_nodal_t *this_nodal,
1865
fvm_lnum_t parent_num[])
1870
fvm_lnum_t entity_count = 0;
1872
assert(this_nodal != NULL);
1874
/* Entity dimension 0: vertices */
1876
if (entity_dim == 0) {
1877
if (this_nodal->parent_vertex_num != NULL) {
1878
for (i = 0; i < this_nodal->n_vertices; i++)
1879
parent_num[entity_count++] = this_nodal->parent_vertex_num[i];
1882
for (i = 0; i < this_nodal->n_vertices; i++)
1883
parent_num[entity_count++] = i + 1;
1887
/* Entity dimension > 0: edges, faces, or cells */
1891
for (section_id = 0; section_id < this_nodal->n_sections; section_id++) {
1893
const fvm_nodal_section_t *section = this_nodal->sections[section_id];
1895
if (section->entity_dim == entity_dim) {
1896
if (section->parent_element_num != NULL) {
1897
for (i = 0; i < section->n_elements; i++)
1898
parent_num[entity_count++] = section->parent_element_num[i];
1901
for (i = 0; i < section->n_elements; i++)
1902
parent_num[entity_count++] = i + 1;
1906
} /* end loop on sections */
1911
/*----------------------------------------------------------------------------
1912
* Compute tesselation a a nodal mesh's sections of a given type, and add the
1913
* corresponding structure to the mesh representation.
1915
* If global element numbers are used (i.e. in parallel mode), this function
1916
* should be only be used after calling fvm_nodal_init_io_num().
1918
* If some mesh sections have already been tesselated, their tesselation
1922
* this_nodal <-> pointer to nodal mesh structure
1923
* type <-> element type that should be tesselated
1924
* error_count --> number of elements with a tesselation error
1925
* counter (optional)
1926
*----------------------------------------------------------------------------*/
1929
fvm_nodal_tesselate(fvm_nodal_t *this_nodal,
1931
fvm_lnum_t *error_count)
1934
fvm_lnum_t section_error_count;
1936
assert(this_nodal != NULL);
1938
if (error_count != NULL)
1941
for (section_id = 0; section_id < this_nodal->n_sections; section_id++) {
1943
fvm_nodal_section_t *section = this_nodal->sections[section_id];
1945
if (section->type == type && section->tesselation == NULL) {
1947
section->tesselation = fvm_tesselation_create(type,
1948
section->n_elements,
1949
section->face_index,
1951
section->vertex_index,
1952
section->vertex_num,
1953
section->global_element_num);
1955
fvm_tesselation_init(section->tesselation,
1957
this_nodal->vertex_coords,
1958
this_nodal->parent_vertex_num,
1959
§ion_error_count);
1961
if (error_count != NULL)
1962
*error_count += section_error_count;
1968
/*----------------------------------------------------------------------------
1969
* Build a nodal representation structure based on extraction of a
1973
* name <-- name to assign to extracted mesh
1974
* this_nodal <-> pointer to nodal mesh structure
1975
*----------------------------------------------------------------------------*/
1978
fvm_nodal_copy_edges(const char *name,
1979
const fvm_nodal_t *this_nodal)
1983
fvm_lnum_t n_edges = 0, n_max_edges = 0;
1984
fvm_nodal_t *new_nodal = NULL;
1985
fvm_nodal_section_t *new_section = NULL;
1987
BFT_MALLOC(new_nodal, 1, fvm_nodal_t);
1989
/* Global indicators */
1992
BFT_MALLOC(new_nodal->name, strlen(name) + 1, char);
1993
strcpy(new_nodal->name, name);
1996
new_nodal->name = NULL;
1998
new_nodal->dim = this_nodal->dim;
1999
new_nodal->num_dom = this_nodal->num_dom;
2000
new_nodal->n_doms = this_nodal->n_doms;
2001
new_nodal->n_sections = 1;
2003
/* Local dimensions */
2005
new_nodal->n_cells = 0;
2006
new_nodal->n_faces = 0;
2007
new_nodal->n_edges = 0;
2008
new_nodal->n_vertices = this_nodal->n_vertices;
2010
/* Local structures */
2012
new_nodal->vertex_coords = this_nodal->vertex_coords;
2013
new_nodal->_vertex_coords = NULL;
2015
new_nodal->parent_vertex_num = this_nodal->parent_vertex_num;
2016
new_nodal->_parent_vertex_num = NULL;
2018
if (this_nodal->global_vertex_num != NULL) {
2020
= fvm_io_num_get_local_count(this_nodal->global_vertex_num);
2021
fvm_gnum_t global_count
2022
= fvm_io_num_get_global_count(this_nodal->global_vertex_num);
2023
const fvm_gnum_t *global_num
2024
= fvm_io_num_get_global_num(this_nodal->global_vertex_num);
2026
new_nodal->global_vertex_num
2027
= fvm_io_num_create_shared(global_num, global_count, n_ent);
2030
new_nodal->global_vertex_num = NULL;
2034
for (i = 0; i < this_nodal->n_sections; i++) {
2035
const fvm_nodal_section_t *this_section = this_nodal->sections[i];
2036
if (this_section->vertex_index == NULL)
2037
n_max_edges += ( fvm_nodal_n_edges_element[this_section->type]
2038
* this_section->n_elements);
2039
else if (this_section->type == FVM_FACE_POLY)
2040
n_max_edges += this_section->vertex_index[this_section->n_elements];
2041
else if (this_section->type == FVM_CELL_POLY)
2042
n_max_edges += this_section->vertex_index[this_section->n_faces];
2045
BFT_MALLOC(new_nodal->sections, 1, fvm_nodal_section_t *);
2047
new_section = fvm_nodal_section_create(FVM_EDGE);
2048
new_nodal->sections[0] = new_section;
2050
BFT_MALLOC(new_section->_vertex_num, n_max_edges*2, fvm_lnum_t);
2054
for (i = 0; i < this_nodal->n_sections; i++) {
2056
const fvm_nodal_section_t *this_section = this_nodal->sections[i];
2058
if ( this_section->type == FVM_FACE_POLY
2059
|| this_section->type == FVM_CELL_POLY) {
2061
fvm_lnum_t n_faces = this_section->type == FVM_FACE_POLY ?
2062
this_section->n_elements : this_section->n_faces;
2064
for (j = 0; j < n_faces; j++) {
2065
const fvm_lnum_t face_start_id = this_section->vertex_index[j];
2066
const fvm_lnum_t n_face_edges
2067
= this_section->vertex_index[j+1] - this_section->vertex_index[j];
2068
for (k = 0; k < n_face_edges; k++) {
2069
new_section->_vertex_num[n_edges*2]
2070
= this_section->vertex_num[face_start_id + k];
2071
new_section->_vertex_num[n_edges*2 + 1]
2072
= this_section->vertex_num[face_start_id + (k + 1)%n_face_edges];
2080
fvm_lnum_t edges[2][12];
2082
fvm_lnum_t n_elt_edges = fvm_nodal_n_edges_element[this_section->type];
2083
fvm_lnum_t n_elts = this_section->n_elements;
2084
fvm_lnum_t stride = this_section->stride;
2086
switch (this_section->type) {
2091
for (j = 0; j < n_elt_edges; j++) {
2093
edges[1][j] = (j+1)%n_elt_edges;
2097
case FVM_CELL_TETRA:
2098
edges[0][0] = 0; edges[1][0] = 1;
2099
edges[0][1] = 1; edges[1][1] = 2;
2100
edges[0][2] = 2; edges[1][2] = 0;
2101
edges[0][3] = 0; edges[1][3] = 3;
2102
edges[0][4] = 1; edges[1][4] = 3;
2103
edges[0][5] = 2; edges[1][5] = 3;
2106
case FVM_CELL_PYRAM:
2107
edges[0][0] = 0; edges[1][0] = 1;
2108
edges[0][1] = 1; edges[1][1] = 2;
2109
edges[0][2] = 2; edges[1][2] = 3;
2110
edges[0][3] = 3; edges[1][3] = 0;
2111
edges[0][4] = 0; edges[1][4] = 4;
2112
edges[0][5] = 1; edges[1][5] = 4;
2113
edges[0][6] = 2; edges[1][6] = 4;
2114
edges[0][7] = 3; edges[1][7] = 4;
2117
case FVM_CELL_PRISM:
2118
edges[0][0] = 0; edges[1][0] = 1;
2119
edges[0][1] = 1; edges[1][1] = 2;
2120
edges[0][2] = 2; edges[1][2] = 0;
2121
edges[0][3] = 0; edges[1][3] = 3;
2122
edges[0][4] = 1; edges[1][4] = 4;
2123
edges[0][5] = 2; edges[1][5] = 5;
2124
edges[0][6] = 3; edges[1][6] = 4;
2125
edges[0][7] = 4; edges[1][7] = 5;
2126
edges[0][8] = 5; edges[1][8] = 3;
2130
edges[0][0] = 0; edges[1][0] = 1;
2131
edges[0][1] = 1; edges[1][1] = 2;
2132
edges[0][2] = 2; edges[1][2] = 3;
2133
edges[0][3] = 3; edges[1][3] = 0;
2134
edges[0][4] = 0; edges[1][4] = 4;
2135
edges[0][5] = 1; edges[1][5] = 5;
2136
edges[0][6] = 2; edges[1][6] = 6;
2137
edges[0][7] = 3; edges[1][7] = 7;
2138
edges[0][8] = 4; edges[1][8] = 5;
2139
edges[0][9] = 5; edges[1][9] = 6;
2140
edges[0][10] = 6; edges[1][10] = 7;
2141
edges[0][11] = 7; edges[1][11] = 4;
2146
edges[0][0] = -1; /* For nonempty default clause */
2149
for (j = 0; j < n_elts; j++) {
2150
const fvm_lnum_t *_vertex_num = this_section->vertex_num + (j*stride);
2151
for (k = 0; k < n_elt_edges; k++) {
2152
new_section->_vertex_num[n_edges*2] = _vertex_num[edges[0][k]];
2153
new_section->_vertex_num[n_edges*2 + 1] = _vertex_num[edges[1][k]];
2158
} /* End of loop on sections */
2160
assert(n_edges == n_max_edges);
2162
/* Ensure edges are oriented in the same direction */
2164
if (this_nodal->global_vertex_num != NULL) {
2166
const fvm_gnum_t *v_num_g
2167
= fvm_io_num_get_global_num(this_nodal->global_vertex_num);
2169
for (j = 0; j < n_max_edges; j++) {
2170
fvm_lnum_t vnum_1 = new_section->_vertex_num[j*2];
2171
fvm_lnum_t vnum_2 = new_section->_vertex_num[j*2 + 1];
2172
if (v_num_g[vnum_1 - 1] > v_num_g[vnum_2 - 1]) {
2173
new_section->_vertex_num[j*2] = vnum_2;
2174
new_section->_vertex_num[j*2 + 1] = vnum_1;
2182
for (j = 0; j < n_max_edges; j++) {
2183
fvm_lnum_t vnum_1 = new_section->_vertex_num[j*2];
2184
fvm_lnum_t vnum_2 = new_section->_vertex_num[j*2 + 1];
2185
if (vnum_1 > vnum_2) {
2186
new_section->_vertex_num[j*2] = vnum_2;
2187
new_section->_vertex_num[j*2 + 1] = vnum_1;
2193
/* Sort and remove duplicates
2194
(use qsort rather than fvm_order_local_s() so as to sort in place) */
2196
qsort(new_section->_vertex_num,
2198
sizeof(fvm_lnum_t) * 2,
2202
fvm_lnum_t vn_1_p = -1;
2203
fvm_lnum_t vn_2_p = -1;
2207
for (j = 0; j < n_max_edges; j++) {
2209
fvm_lnum_t vn_1 = new_section->_vertex_num[j*2];
2210
fvm_lnum_t vn_2 = new_section->_vertex_num[j*2 + 1];
2212
if (vn_1 != vn_1_p || vn_2 != vn_2_p) {
2213
new_section->_vertex_num[n_edges*2] = vn_1;
2214
new_section->_vertex_num[n_edges*2 + 1] = vn_2;
2222
/* Resize edge connectivity to adjust to final size */
2224
BFT_REALLOC(new_section->_vertex_num, n_edges*2, fvm_lnum_t);
2225
new_section->vertex_num = new_section->_vertex_num;
2227
new_section->n_elements = n_edges;
2228
new_nodal->n_edges = n_edges;
2230
/* Build global edge numbering if necessary */
2232
if (new_nodal->n_doms > 1) {
2234
fvm_gnum_t *edge_vertices_g; /* edges -> global vertices */
2236
BFT_MALLOC(edge_vertices_g, n_edges*2, fvm_gnum_t);
2238
if (this_nodal->global_vertex_num != NULL) {
2239
const fvm_gnum_t *v_num_g
2240
= fvm_io_num_get_global_num(this_nodal->global_vertex_num);
2241
for (j = 0; j < n_edges; j++) {
2242
edge_vertices_g[j*2] = v_num_g[new_section->_vertex_num[j*2] - 1];
2243
edge_vertices_g[j*2+1] = v_num_g[new_section->_vertex_num[j*2+1] - 1];
2247
for (j = 0; j < n_edges; j++) {
2248
edge_vertices_g[j*2] = new_section->_vertex_num[j*2];
2249
edge_vertices_g[j*2 + 1] = new_section->_vertex_num[j*2 + 1];
2253
new_section->global_element_num
2254
= fvm_io_num_create_from_adj_s(NULL, edge_vertices_g, n_edges, 2);
2257
BFT_FREE(edge_vertices_g);
2260
new_nodal->gc_set = NULL;
2265
/*----------------------------------------------------------------------------
2266
* Dump printout of a nodal representation structure.
2269
* this_nodal <-- pointer to structure that should be dumped
2270
*----------------------------------------------------------------------------*/
2273
fvm_nodal_dump(const fvm_nodal_t *this_nodal)
2276
fvm_lnum_t num_vertex = 1;
2277
const fvm_coord_t *coord = this_nodal->vertex_coords;
2279
/* Global indicators */
2280
/*--------------------*/
2283
"Mesh name:\"%s\"\n",
2287
"Mesh dimension: %d\n"
2288
"Domain number: %d\n"
2289
"Number of domains: %d\n"
2290
"Number of sections: %d\n",
2291
this_nodal->dim, this_nodal->num_dom, this_nodal->n_doms,
2292
this_nodal->n_sections);
2295
"Number of cells: %d\n"
2296
"Number of faces: %d\n"
2297
"Number of edges: %d\n"
2298
"Number of vertices: %d\n",
2299
this_nodal->n_cells,
2300
this_nodal->n_faces,
2301
this_nodal->n_edges,
2302
this_nodal->n_vertices);
2304
if (this_nodal->n_vertices > 0) {
2307
"Pointers to shareable arrays:\n"
2308
" vertex_coords: %p\n"
2309
" parent_vertex_num: %p\n",
2310
this_nodal->vertex_coords,
2311
this_nodal->parent_vertex_num);
2314
"Pointers to local arrays:\n"
2315
" _vertex_coords: %p\n"
2316
" _parent_vertex_num: %p\n",
2317
this_nodal->_vertex_coords,
2318
this_nodal->_parent_vertex_num);
2320
/* Output coordinates depending on parent numbering */
2322
if (this_nodal->parent_vertex_num == NULL) {
2324
bft_printf("\nVertex coordinates:\n\n");
2325
switch(this_nodal->dim) {
2327
for (i = 0; i < this_nodal->n_vertices; i++)
2328
bft_printf("%10d : %12.5f\n",
2329
num_vertex++, (double)(coord[i]));
2332
for (i = 0; i < this_nodal->n_vertices; i++)
2333
bft_printf("%10d : %12.5f %12.5f\n",
2334
num_vertex++, (double)(coord[i*2]),
2335
(double)(coord[i*2+1]));
2338
for (i = 0; i < this_nodal->n_vertices; i++)
2339
bft_printf("%10d : %12.5f %12.5f %12.5f\n",
2340
num_vertex++, (double)(coord[i*3]),
2341
(double)(coord[i*3+1]), (double)(coord[i*3+2]));
2344
bft_printf("coordinates not output\n"
2345
"dimension = %d unsupported\n", this_nodal->dim);
2349
else { /* if (this_nodal->parent_vertex_num != NULL) */
2351
bft_printf("\nVertex parent and coordinates:\n\n");
2353
switch(this_nodal->dim) {
2355
for (i = 0; i < this_nodal->n_vertices; i++) {
2356
coord = this_nodal->vertex_coords
2357
+ (this_nodal->parent_vertex_num[i]-1);
2358
bft_printf("%10d : %12.5f\n",
2359
num_vertex++, (double)(coord[0]));
2363
for (i = 0; i < this_nodal->n_vertices; i++) {
2364
coord = this_nodal->vertex_coords
2365
+ ((this_nodal->parent_vertex_num[i]-1)*2);
2366
bft_printf("%10d : %12.5f %12.5f\n",
2367
num_vertex++, (double)(coord[0]), (double)(coord[1]));
2371
for (i = 0; i < this_nodal->n_vertices; i++) {
2372
coord = this_nodal->vertex_coords
2373
+ ((this_nodal->parent_vertex_num[i]-1)*3);
2374
bft_printf("%10d : %12.5f %12.5f %12.5f\n",
2375
num_vertex++, (double)(coord[0]), (double)(coord[1]),
2376
(double)(coord[2]));
2380
bft_printf("coordinates not output\n"
2381
"dimension = %d unsupported\n", this_nodal->dim);
2388
/* Global vertex numbers (only for parallel execution) */
2389
if (this_nodal->global_vertex_num != NULL) {
2390
bft_printf("\nGlobal vertex numbers:\n\n");
2391
fvm_io_num_dump(this_nodal->global_vertex_num);
2394
/* Dump element sections */
2395
/*-----------------------*/
2397
for (i = 0; i < this_nodal->n_sections; i++)
2398
_fvm_nodal_section_dump(this_nodal->sections[i]);
2400
/* Dump group class set (NULL allowed) */
2402
fvm_group_class_set_dump(this_nodal->gc_set);
2405
/*----------------------------------------------------------------------------*/
2409
#endif /* __cplusplus */
added
removed
src/fvm/fvm_nodal.c
