1
/*============================================================================
2
* Ordering of nodal mesh entity lists and connectivity
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_nodal.h"
54
#include "fvm_nodal_priv.h"
55
#include "fvm_order.h"
56
#include "fvm_parall.h"
58
/*----------------------------------------------------------------------------
59
* Header for the current file
60
*----------------------------------------------------------------------------*/
62
#include "fvm_nodal_order.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
/*============================================================================
78
* Private function definitions
79
*============================================================================*/
81
/*----------------------------------------------------------------------------
82
* Create ordered parent entity list or order existing list.
85
* _list <-> pointer to optional list (1 to n numbering) of selected
86
* entities. An existing list is ordered, otherwise one is
88
* list <-> pointer tooptional list (1 to n numbering) of selected
89
* entities. A shared list is copied to _list and ordered,
90
* a private list (pointed to by _list) is simply ordered.
91
* order <-- ordering of entities (0 to n-1).
92
* nb_ent <-- number of entities considered.
93
*----------------------------------------------------------------------------*/
96
_fvm_nodal_order_parent_list(fvm_lnum_t * _list[],
97
const fvm_lnum_t * list[],
98
const fvm_lnum_t order[],
103
fvm_lnum_t *ordered_list = NULL;
105
BFT_MALLOC(ordered_list, nb_ent, fvm_lnum_t);
109
for (i = 0 ; i < nb_ent ; i++)
110
ordered_list[i] = (*list)[order[i]];
112
if (*_list != NULL) {
113
for (i = 0 ; i < nb_ent ; i++)
114
(*_list)[i] = ordered_list[i];
115
BFT_FREE(ordered_list);
118
*_list = ordered_list;
123
assert(*list == NULL);
125
for (i = 0 ; i < nb_ent ; i++)
126
ordered_list[i] = order[i] + 1;
127
*_list = ordered_list;
134
/*----------------------------------------------------------------------------
135
* Reorder strided connectivity array.
138
* connect <-> connectivity array (nb_ent * stride) to be ordered.
139
* order <-- ordering of entities (0 to n-1).
140
* nb_ent <-- number of entities considered.
141
*----------------------------------------------------------------------------*/
144
_fvm_nodal_order_strided_connect(fvm_lnum_t connect[],
145
const fvm_lnum_t order[],
152
fvm_lnum_t *tmp_connect = NULL;
154
BFT_MALLOC(tmp_connect, nb_ent * stride, fvm_lnum_t);
156
/* Temporary ordered copy */
158
for (i = 0 ; i < nb_ent ; i++) {
159
p1 = tmp_connect + i*stride;
160
p2 = connect + (order[i] * stride);
161
for (j = 0 ; j < stride ; j++)
165
/* Now put back in initial location */
167
memcpy(connect, tmp_connect, stride * nb_ent * sizeof(fvm_lnum_t));
169
BFT_FREE(tmp_connect);
172
/*----------------------------------------------------------------------------
173
* Reorder indexed connectivity array.
176
* connect_idx <-> connectivity index array (0 to n -1) to be ordered.
177
* connect_num <-> connectivity numbers array to be ordered.
178
* order <-- ordering of entities (0 to n-1).
179
* nb_ent <-- number of entities considered.
180
*----------------------------------------------------------------------------*/
183
_fvm_nodal_order_indexed_connect(fvm_lnum_t connect_idx[],
184
fvm_lnum_t connect_num[],
185
const fvm_lnum_t order[],
188
size_t i, j, nb_ent_max, nb_loc;
191
fvm_lnum_t *tmp_connect = NULL;
193
nb_ent_max = connect_idx[nb_ent] ; /* size of connect_num */
194
if (nb_ent > nb_ent_max) /* only if some entities have no connectivity */
197
BFT_MALLOC(tmp_connect, nb_ent_max, fvm_lnum_t);
199
/* Temporary ordered copy of values */
202
for (i = 0 ; i < nb_ent ; i++) {
203
nb_loc = connect_idx[order[i]+1] - connect_idx[order[i]];
204
p2 = connect_num + connect_idx[order[i]];
205
for (j = 0 ; j < nb_loc ; j++)
209
/* Now put back in initial location */
211
memcpy(connect_num, tmp_connect,
212
(size_t)(connect_idx[nb_ent]) * sizeof(fvm_lnum_t));
214
/* Index to size : size associated with entity i in position i+1 */
216
for (i = nb_ent ; i > 0 ; i--)
217
connect_idx[i] = connect_idx[i] - connect_idx[i-1];
219
/* Temporary ordered copy of transformed index */
223
for (i = 0 ; i < nb_ent ; i++)
224
*p1++ = connect_idx[order[i] + 1];
226
/* Put back in initial location and re-convert to index*/
228
memcpy(connect_idx, tmp_connect, (size_t)(nb_ent + 1) * sizeof(fvm_lnum_t));
230
for (i = 0 ; i < nb_ent ; i++)
231
connect_idx[i+1] = connect_idx[i+1] + connect_idx[i];
233
BFT_FREE(tmp_connect);
237
/*============================================================================
238
* Public function definitions
239
*============================================================================*/
241
/*----------------------------------------------------------------------------
242
* Locally order cells and associated connectivity for a nodal mesh
245
* this_nodal <-- pointer to nodal mesh structure.
246
* parent_global_number <-- global numbers of parent cells (if NULL, a
247
* default 1 to n numbering is considered).
248
*----------------------------------------------------------------------------*/
251
fvm_nodal_order_cells(fvm_nodal_t *this_nodal,
252
const fvm_gnum_t parent_global_number[])
255
fvm_lnum_t *order = NULL;
256
fvm_nodal_section_t *section = NULL;
258
if (this_nodal == NULL)
261
/* Order locally if necessary */
263
for (i = 0 ; i < this_nodal->n_sections ; i++) {
265
section = this_nodal->sections[i];
267
if (section->entity_dim == 3) {
269
assert(section->global_element_num == NULL);
271
if (fvm_order_local_test(section->parent_element_num,
272
parent_global_number,
273
section->n_elements) == false) {
275
order = fvm_order_local(section->parent_element_num,
276
parent_global_number,
277
section->n_elements);
279
_fvm_nodal_order_parent_list(&(section->_parent_element_num),
280
&(section->parent_element_num),
282
section->n_elements);
284
if (section->type != FVM_CELL_POLY) {
285
fvm_nodal_section_copy_on_write(section, false, false, false, true);
286
_fvm_nodal_order_strided_connect(section->_vertex_num,
288
(size_t)(section->stride),
289
section->n_elements);
292
fvm_nodal_section_copy_on_write(section, true, true, false, false);
293
_fvm_nodal_order_indexed_connect(section->_face_index,
296
section->n_elements);
309
/*----------------------------------------------------------------------------
310
* Locally order faces and associated connectivity for a nodal mesh
313
* this_nodal <-- pointer to nodal mesh structure.
314
* parent_global_number <-- global numbers of parent faces (if NULL, a
315
* default 1 to n numbering is considered).
316
*----------------------------------------------------------------------------*/
319
fvm_nodal_order_faces(fvm_nodal_t *this_nodal,
320
const fvm_gnum_t parent_global_number[])
323
fvm_lnum_t *order = NULL;
324
fvm_nodal_section_t *section = NULL;
326
if (this_nodal == NULL)
329
/* Order locally if necessary */
331
for (i = 0 ; i < this_nodal->n_sections ; i++) {
333
section = this_nodal->sections[i];
335
if (section->entity_dim == 2) {
337
assert(section->global_element_num == NULL);
339
if (fvm_order_local_test(section->parent_element_num,
340
parent_global_number,
341
section->n_elements) == false) {
343
order = fvm_order_local(section->parent_element_num,
344
parent_global_number,
345
section->n_elements);
347
_fvm_nodal_order_parent_list(&(section->_parent_element_num),
348
&(section->parent_element_num),
350
section->n_elements);
352
if (section->type != FVM_FACE_POLY) {
353
fvm_nodal_section_copy_on_write(section, false, false, false, true);
354
_fvm_nodal_order_strided_connect(section->_vertex_num,
356
(size_t)(section->stride),
357
section->n_elements);
360
fvm_nodal_section_copy_on_write(section, false, false, true, true);
361
_fvm_nodal_order_indexed_connect(section->_vertex_index,
362
section->_vertex_num,
364
section->n_elements);
377
/*----------------------------------------------------------------------------
378
* Locally order vertices and update connectivity for a nodal mesh
381
* this_nodal <-- pointer to nodal mesh structure.
382
* parent_global_number <-- global numbers of parent vertices (if NULL, a
383
* default 1 to n numbering is considered).
384
*----------------------------------------------------------------------------*/
387
fvm_nodal_order_vertices(fvm_nodal_t *this_nodal,
388
const fvm_gnum_t parent_global_number[])
393
fvm_lnum_t *order = NULL;
394
fvm_lnum_t *renumber = NULL;
395
fvm_nodal_section_t *section = NULL;
397
/* Do nothing for trivial cases */
399
if (this_nodal == NULL)
402
else if (this_nodal->n_vertices < 2)
405
/* Return if already ordered */
407
if (fvm_order_local_test(this_nodal->parent_vertex_num,
408
parent_global_number,
409
this_nodal->n_vertices) == true)
412
/* Else, we must re-order vertices and update connectivity */
414
order = fvm_order_local(this_nodal->parent_vertex_num,
415
parent_global_number,
416
this_nodal->n_vertices);
418
/* Re-order parent list */
420
_fvm_nodal_order_parent_list(&(this_nodal->_parent_vertex_num),
421
&(this_nodal->parent_vertex_num),
423
this_nodal->n_vertices);
425
/* Calculate renumbering table for associated connectivity
426
and free ordering array, no longer needed after that */
428
renumber = fvm_order_local_renumbering(order,
429
this_nodal->n_vertices);
433
/* Update element connectivities */
435
for (i = 0 ; i < this_nodal->n_sections ; i++) {
437
section = this_nodal->sections[i];
439
fvm_nodal_section_copy_on_write(section, false, false, false, true);
441
for (j = 0 ; j < section->connectivity_size ; j++)
442
section->_vertex_num[j] = renumber[section->_vertex_num[j] - 1] + 1;
446
/* Free renumbering table */
452
/*----------------------------------------------------------------------------*/
456
#endif /* __cplusplus */