1 /*============================================================================
2 * \file Auxiliary structure used to read, write, and partition mesh data.
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 #include <limits.h>
38 #include <math.h>
39 #include <assert.h>
40
41 /*----------------------------------------------------------------------------
42 * Local headers
43 *----------------------------------------------------------------------------*/
44
45 #include "bft_mem.h"
46 #include "bft_printf.h"
47
48 #include "fvm_io_num.h"
49 #include "fvm_periodicity.h"
50 #include "fvm_selector.h"
51
52 #include "cs_base.h"
53 #include "cs_interface.h"
54
55 /*----------------------------------------------------------------------------
56 * Header for the current file
57 *----------------------------------------------------------------------------*/
58
59 #include "cs_mesh_builder.h"
60
61 /*----------------------------------------------------------------------------*/
62
63 BEGIN_C_DECLS
64
65 /*! \cond DOXYGEN_SHOULD_SKIP_THIS */
66
67 /*=============================================================================
68 * Local Macro Definitions
69 *============================================================================*/
70
71 /*============================================================================
72 * Local structure definitions
73 *============================================================================*/
74
75 /*=============================================================================
76 * Local Macro definitions
77 *============================================================================*/
78
79 /*============================================================================
80 * Static global variables
81 *============================================================================*/
82
83 /* Pointer on the temporary mesh used for building main mesh */
84
85 cs_mesh_builder_t *cs_glob_mesh_builder = NULL;
86
87 /*=============================================================================
88 * Private function definitions
89 *============================================================================*/
90
91 #if defined(HAVE_MPI)
92
93 /*----------------------------------------------------------------------------
94 * Compare periodic couples in global numbering form (qsort function).
95 *
96 * parameters:
97 * x <-> pointer to first couple
98 * y <-> pointer to second couple
99 *
100 * returns:
101 * lexicographical
102 *----------------------------------------------------------------------------*/
103
_compare_couples(const void * x,const void * y)104 static int _compare_couples(const void *x, const void *y)
105 {
106 int retval = 1;
107
108 const cs_gnum_t *c0 = x;
109 const cs_gnum_t *c1 = y;
110
111 if (c0[0] < c1[0])
112 retval = -1;
113
114 else if (c0[0] == c1[0]) {
115 if (c0[1] < c1[1])
116 retval = -1;
117 else if (c0[1] == c1[1])
118 retval = 0;
119 }
120
121 return retval;
122 }
123
124 #endif /* defined(HAVE_MPI) */
125
126 /*! (DOXYGEN_SHOULD_SKIP_THIS) \endcond */
127
128 /*=============================================================================
129 * Public function definitions
130 *============================================================================*/
131
132 /*----------------------------------------------------------------------------*/
133 /*!
134 * \brief Create an empty mesh builder structure.
135 *
136 * \return pointer to a mesh builder structure
137 */
138 /*----------------------------------------------------------------------------*/
139
140 cs_mesh_builder_t *
cs_mesh_builder_create(void)141 cs_mesh_builder_create(void)
142 {
143 cs_mesh_builder_t *mb = NULL;
144
145 BFT_MALLOC(mb, 1, cs_mesh_builder_t);
146
147 mb->n_g_faces = 0;
148 mb->n_g_face_connect_size = 0;
149
150 mb->n_perio = 0;
151
152 mb->have_cell_rank = false;
153
154 mb->have_face_r_gen = false;
155
156 /* Temporary mesh data */
157
158 mb->face_cells = NULL;
159 mb->face_vertices_idx = NULL;
160 mb->face_vertices = NULL;
161 mb->cell_gc_id = NULL;
162 mb->face_gc_id = NULL;
163 mb->vertex_coords = NULL;
164
165 /* Refinement features */
166
167 mb->face_r_gen = NULL;
168
169 /* Periodic features */
170
171 mb->periodicity_num = NULL;
172 mb->n_per_face_couples = NULL;
173 mb->n_per_face_couples = NULL;
174 mb->n_g_per_face_couples = NULL;
175 mb->per_face_couples = NULL;
176
177 /* Optional partitioning info */
178
179 mb->cell_rank = NULL;
180
181 /* Block ranges for parallel distribution */
182
183 mb->min_rank_step = 1;
184 memset(&(mb->cell_bi), 0, sizeof(cs_block_dist_info_t));
185 memset(&(mb->face_bi), 0, sizeof(cs_block_dist_info_t));
186 memset(&(mb->vertex_bi), 0, sizeof(cs_block_dist_info_t));
187 mb->per_face_bi = NULL;
188
189 return mb;
190 }
191
192 /*----------------------------------------------------------------------------*/
193 /*!
194 * \brief Destroy a cs_mesh_builder_t structure.
195 *
196 * \param[in, out] mb pointer to pointer of structure to destroy
197 */
198 /*----------------------------------------------------------------------------*/
199
200 void
cs_mesh_builder_destroy(cs_mesh_builder_t ** mb)201 cs_mesh_builder_destroy(cs_mesh_builder_t **mb)
202 {
203 if (mb != NULL) {
204
205 cs_mesh_builder_t *_mb = *mb;
206
207 if (_mb == NULL)
208 return;
209
210 /* Temporary mesh data */
211
212 BFT_FREE(_mb->face_cells);
213 BFT_FREE(_mb->face_vertices_idx);
214 BFT_FREE(_mb->face_vertices);
215 BFT_FREE(_mb->cell_gc_id);
216 BFT_FREE(_mb->face_gc_id);
217 BFT_FREE(_mb->vertex_coords);
218
219 /* Refinement features */
220
221 BFT_FREE(_mb->face_r_gen);
222
223 /* Periodic features */
224
225 BFT_FREE(_mb->periodicity_num);
226 BFT_FREE(_mb->n_per_face_couples);
227 BFT_FREE(_mb->n_g_per_face_couples);
228 if (_mb->per_face_couples != NULL) {
229 for (int i = 0; i < _mb->n_perio; i++)
230 BFT_FREE(_mb->per_face_couples[i]);
231 BFT_FREE(_mb->per_face_couples);
232 }
233
234 /* Optional partitioning info */
235
236 BFT_FREE(_mb->cell_rank);
237
238 /* Block ranges for parallel distribution */
239
240 BFT_FREE(_mb->per_face_bi);
241
242 BFT_FREE(*mb);
243 }
244 }
245
246 /*----------------------------------------------------------------------------*/
247 /*!
248 * \brief Define block distribution sizes for mesh builder.
249 *
250 * \param[in, out] mb pointer to mesh builder to update
251 * \param[in] rank_id id of local rank
252 * \param[in] n_ranks number of associated ranks
253 * \param[in] min_rank_step minimum rank step between blocks
254 * \param[in] min_block_size minimum number of entities per block
255 * \param[in] n_g_cells global number of cells
256 * \param[in] n_g_faces global number of faces
257 * \param[in] n_g_vertices global number of vertices
258 */
259 /*----------------------------------------------------------------------------*/
260
261 void
cs_mesh_builder_define_block_dist(cs_mesh_builder_t * mb,int rank_id,int n_ranks,int min_rank_step,int min_block_size,cs_gnum_t n_g_cells,cs_gnum_t n_g_faces,cs_gnum_t n_g_vertices)262 cs_mesh_builder_define_block_dist(cs_mesh_builder_t *mb,
263 int rank_id,
264 int n_ranks,
265 int min_rank_step,
266 int min_block_size,
267 cs_gnum_t n_g_cells,
268 cs_gnum_t n_g_faces,
269 cs_gnum_t n_g_vertices)
270 {
271 mb->min_rank_step = min_rank_step;
272
273 mb->cell_bi = cs_block_dist_compute_sizes(rank_id,
274 n_ranks,
275 min_rank_step,
276 min_block_size,
277 n_g_cells);
278
279 mb->face_bi = cs_block_dist_compute_sizes(rank_id,
280 n_ranks,
281 min_rank_step,
282 min_block_size,
283 n_g_faces);
284
285 mb->vertex_bi = cs_block_dist_compute_sizes(rank_id,
286 n_ranks,
287 min_rank_step,
288 min_block_size,
289 n_g_vertices);
290 }
291
292 #if defined(HAVE_MPI)
293
294 /*----------------------------------------------------------------------------*/
295 /*!
296 * \brief Extract periodic face connectivity information from faces interface
297 * for mesh builder when running in parallel mode.
298 *
299 * \param[in] n_init_perio number of initial periodicities
300 * \param[in, out] mb pointer to mesh builder structure
301 * \param[in] periodicity periodicity information
302 * \param[in] face_gnum global face numbers, or NULL
303 * \param[in] face_ifs parallel and periodic faces interfaces set
304 */
305 /*----------------------------------------------------------------------------*/
306
307 void
cs_mesh_builder_extract_periodic_faces_g(int n_init_perio,cs_mesh_builder_t * mb,fvm_periodicity_t * periodicity,const cs_gnum_t * face_gnum,const cs_interface_set_t * face_ifs)308 cs_mesh_builder_extract_periodic_faces_g(int n_init_perio,
309 cs_mesh_builder_t *mb,
310 fvm_periodicity_t *periodicity,
311 const cs_gnum_t *face_gnum,
312 const cs_interface_set_t *face_ifs)
313 {
314 int i, j;
315 cs_lnum_t k, l;
316
317 int perio_count = 0;
318 cs_lnum_t *send_index = NULL;
319 cs_gnum_t *recv_num = NULL;
320 int *tr_id = NULL;
321
322 cs_datatype_t gnum_type = CS_GNUM_TYPE;
323
324 const int n_perio = n_init_perio;
325 const int n_interfaces = cs_interface_set_size(face_ifs);
326
327 /* Free previous values if we are updating */
328
329 if (mb->n_perio > 0 && mb->n_per_face_couples != NULL) {
330 for (i = 0; i < n_perio; i++)
331 BFT_FREE(mb->per_face_couples[i]);
332 BFT_FREE(mb->n_per_face_couples);
333 BFT_FREE(mb->per_face_couples);
334 }
335
336 mb->n_perio = n_perio;
337
338 /* Allocate arrays in mesh builder (initializing per_face_idx) */
339
340 assert(periodicity != NULL);
341 assert(mb != NULL);
342 assert(mb->n_g_per_face_couples == 0);
343
344 BFT_MALLOC(mb->n_per_face_couples, n_perio, cs_lnum_t);
345 BFT_MALLOC(mb->per_face_couples, n_perio, cs_gnum_t *);
346
347 for (i = 0; i < n_perio; i++) {
348 mb->n_per_face_couples[i] = 0;
349 mb->per_face_couples[i] = NULL;
350 }
351
352 /* List direct and reverse transforms */
353
354 BFT_MALLOC(tr_id, n_perio*2, int);
355
356 for (i = 0; i < n_perio*2; i++) {
357 int rev_id = fvm_periodicity_get_reverse_id(periodicity, i);
358 if (i < rev_id) {
359 int parent_ids[2];
360 fvm_periodicity_get_parent_ids(periodicity, i, parent_ids);
361 if (parent_ids[0] < 0 && parent_ids[1] < 0) {
362 tr_id[perio_count*2] = i + 1;
363 tr_id[perio_count*2 + 1] = rev_id + 1;
364 perio_count++;
365 }
366 }
367 }
368 assert(perio_count == n_perio);
369
370 for (i = 0; i < n_interfaces; i++) {
371 const cs_interface_t *face_if = cs_interface_set_get(face_ifs, i);
372 const cs_lnum_t *tr_index = cs_interface_get_tr_index(face_if);
373 for (j = 0; j < n_perio; j++) {
374 const cs_lnum_t n_tr_faces = ( tr_index[tr_id[j*2] + 1]
375 - tr_index[tr_id[j*2]]);
376 mb->n_per_face_couples[j] += n_tr_faces;
377 }
378 }
379
380 BFT_MALLOC(recv_num, cs_interface_set_n_elts(face_ifs), cs_gnum_t);
381
382 cs_interface_set_copy_array(face_ifs,
383 gnum_type,
384 1,
385 true, /* src_on_parent */
386 face_gnum,
387 recv_num);
388
389 /* Prepare send buffer (send reverse transformation values) */
390
391 BFT_FREE(send_index);
392
393 for (i = 0; i < n_perio; i++)
394 BFT_MALLOC(mb->per_face_couples[i], mb->n_per_face_couples[i]*2, cs_gnum_t);
395
396 /* Reset couples count */
397
398 for (i = 0; i < n_perio; i++)
399 mb->n_per_face_couples[i] = 0;
400
401 /* Copy face couples to mesh builder */
402
403 for (i = 0, j = 0, l = 0; i < n_interfaces; i++) {
404
405 const cs_interface_t *face_if = cs_interface_set_get(face_ifs, i);
406 const cs_lnum_t *tr_index = cs_interface_get_tr_index(face_if);
407 const cs_lnum_t *elt_id = cs_interface_get_elt_ids(face_if);
408
409 l += tr_index[1];
410
411 for (j = 0; j < n_perio; j++) {
412
413 /* Count couples in direct periodicity */
414
415 cs_lnum_t nc = mb->n_per_face_couples[j]*2;
416 const cs_lnum_t start_id = tr_index[tr_id[j*2]];
417 const cs_lnum_t end_id = tr_index[tr_id[j*2] + 1];
418
419 for (k = start_id; k < end_id; k++) {
420 cs_lnum_t f_id = elt_id[k];
421 mb->per_face_couples[j][nc++] = face_gnum[f_id];
422 mb->per_face_couples[j][nc++] = recv_num[l++];
423 }
424 mb->n_per_face_couples[j] = nc/2;
425
426 /* Ignore couples in reverse periodicity */
427
428 l += tr_index[tr_id[j*2 + 1] + 1] - tr_index[tr_id[j*2 + 1]];
429
430 }
431
432 }
433
434 BFT_FREE(recv_num);
435 BFT_FREE(tr_id);
436
437 /* Now sort couples in place for future use (more for consistency
438 and ease of verification than absolutely necessary) */
439
440 for (i = 0; i < n_perio; i++) {
441 if (mb->n_per_face_couples[i] > 0)
442 qsort(mb->per_face_couples[i],
443 mb->n_per_face_couples[i],
444 sizeof(cs_gnum_t) * 2,
445 &_compare_couples);
446 }
447 }
448
449 #endif /* defined(HAVE_MPI) */
450
451 /*----------------------------------------------------------------------------*/
452
453 END_C_DECLS
454