@@ -146,6 +146,10 @@ _transfer_bl_faces_selection(void *input,
146146/*!
147147 * \brief Flag vertices for limiter.
148148 *
149+ * We mark cells using cell_vol_cmp = -3 for negative volumes, -2 for volumes
150+ * reduced below the required threshold, and -1 for cells marked through
151+ * adjacency with one of the above.
152+ *
149153 * \param[in] m mesh
150154 * \param[in] cell_vol_cmp comparative cell volume (< 0 for limit)
151155 * \param[out] vtx_flag vertex flag (0 for unlimited, 1 for limited)
@@ -282,6 +286,77 @@ _extrude_vector_limit(const char *vtx_flag,
282286 return n_limited ;
283287}
284288
289+ /*----------------------------------------------------------------------------*/
290+ /*!
291+ * \brief Expand limiter to neighboring cells.
292+ *
293+ * We mark cells using cell_vol_cmp = -3 for negative volumes, -2 for volumes
294+ * reduced below the required threshold, and -1 for cells marked through
295+ * adjacency with one of the above.
296+ *
297+ * \param[in] m mesh
298+ * \param[in] vtx_flag vertex flag (0 for unlimited, 1 for limited)
299+ * \param[in, out] cell_vol_cmp comparative cell volume (< 0 for limit)
300+ *
301+ * \return:
302+ * number of cells marked by adjacency
303+ */
304+ /*----------------------------------------------------------------------------*/
305+
306+ static cs_lnum_t
307+ _expand_limit (const cs_mesh_t * m ,
308+ cs_real_t * cell_vol_cmp ,
309+ char * vtx_flag )
310+ {
311+ const cs_lnum_t n_cells = m -> n_cells ;
312+ const cs_lnum_t n_i_faces = m -> n_i_faces ;
313+ const cs_lnum_t n_b_faces = m -> n_b_faces ;
314+
315+ /* Use vertices flag to mark adjacent cells with bad volumes */
316+
317+ for (cs_lnum_t f_id = 0 ; f_id < n_i_faces ; f_id ++ ) {
318+ bool flag_cells = false;
319+ cs_lnum_t s_id = m -> i_face_vtx_idx [f_id ];
320+ cs_lnum_t e_id = m -> i_face_vtx_idx [f_id + 1 ];
321+ for (cs_lnum_t i = s_id ; i < e_id ; i ++ ) {
322+ if (vtx_flag [m -> i_face_vtx_lst [i ]] != 0 )
323+ flag_cells = true;
324+ }
325+ if (flag_cells ) {
326+ cs_lnum_t c_id0 = m -> i_face_cells [f_id ][0 ];
327+ cs_lnum_t c_id1 = m -> i_face_cells [f_id ][0 ];
328+ if (c_id0 > -1 && c_id0 < n_cells )
329+ cell_vol_cmp [c_id0 ] = CS_MIN (cell_vol_cmp [c_id0 ], -1 );
330+ if (c_id1 > -1 && c_id1 < n_cells )
331+ cell_vol_cmp [c_id1 ] = CS_MIN (cell_vol_cmp [c_id0 ], -1 );
332+ }
333+ }
334+
335+ for (cs_lnum_t f_id = 0 ; f_id < n_b_faces ; f_id ++ ) {
336+ bool flag_cells = false;
337+ cs_lnum_t s_id = m -> b_face_vtx_idx [f_id ];
338+ cs_lnum_t e_id = m -> b_face_vtx_idx [f_id + 1 ];
339+ for (cs_lnum_t i = s_id ; i < e_id ; i ++ ) {
340+ if (vtx_flag [m -> b_face_vtx_lst [i ]] != 0 )
341+ flag_cells = true;
342+ }
343+ if (flag_cells ) {
344+ cs_lnum_t c_id0 = m -> b_face_cells [f_id ];
345+ if (c_id0 > -1 && c_id0 < n_cells )
346+ cell_vol_cmp [c_id0 ] = CS_MIN (cell_vol_cmp [c_id0 ], -1 );
347+ }
348+
349+ }
350+
351+ cs_lnum_t count = 0 ;
352+ for (cs_lnum_t c_id = 0 ; c_id < n_cells ; c_id ++ ) {
353+ if (fabs (cell_vol_cmp [c_id ] + 1 ) < 0.1 )
354+ count ++ ;
355+ }
356+
357+ return count ;
358+ }
359+
285360/*! (DOXYGEN_SHOULD_SKIP_THIS) \endcond */
286361
287362/*============================================================================
@@ -410,23 +485,26 @@ cs_mesh_boundary_layer_insert(cs_mesh_t *m,
410485 m -> vtx_coord [i * 3 + 2 ] += vd [i ][2 ];
411486 }
412487
413- /* Check deformation is acceptable */
488+ /* Check deformation is acceptable;
489+ * We mark cells using cell_vol_cmp = -3 for negative volumes,
490+ * -2 for volumes reduced below the required threshold, and -1 for cells
491+ * marked through adjacency with one of the above. */
414492
415493 compute_displacement = false;
416494
417495 if (min_volume_factor > 0 && min_volume_factor < 1 ) {
418496
419- cs_gnum_t counts [3 ] = {0 , 0 , 0 };
497+ cs_gnum_t counts [4 ] = {0 , 0 , 0 , 0 };
420498
421499 cs_real_t * cell_vol_cmp = cs_mesh_quantities_cell_volume (m );
422500
423501 for (cs_lnum_t i = 0 ; i < n_cells_ini ; i ++ ) {
424502 if (cell_vol_cmp [i ] <= 0 ) {
425- cell_vol_cmp [i ] = -1 ;
503+ cell_vol_cmp [i ] = -3 ;
426504 counts [0 ] += 1 ;
427505 }
428506 else if (cell_vol_cmp [i ] < cell_vol_ref [i ]* min_volume_factor ) {
429- cell_vol_cmp [i ] = -1 ;
507+ cell_vol_cmp [i ] = -2 ;
430508 counts [1 ] += 1 ;
431509 }
432510 }
@@ -442,36 +520,48 @@ cs_mesh_boundary_layer_insert(cs_mesh_t *m,
442520 cell_vol_cmp ,
443521 vtx_flag );
444522
445- BFT_FREE (cell_vol_cmp );
446-
447523 /* Now adjust extrusion vectors structure,
448524 removing a layer at flagged vertices */
449525
450526 counts [2 ] = _extrude_vector_limit (vtx_flag , e );
451527
452- BFT_FREE (vtx_flag );
453-
454528 cs_parall_sum (3 , CS_GNUM_TYPE , counts );
455529
530+ /* If bad volumes are present but not adjacent to a boundary
531+ layer insertion zone, propagate to neighboring volumes */
532+
533+ const int max_propagation_iter = 30 ;
534+ for (int p_iter = 0 ;
535+ counts [0 ] > 0 && counts [2 ] == 0 && p_iter < max_propagation_iter ;
536+ p_iter ++ ) {
537+ counts [3 ] = _expand_limit (m , cell_vol_cmp , vtx_flag );
538+ _flag_vertices_for_limiter (m , cell_vol_cmp , vtx_flag );
539+ counts [2 ] = _extrude_vector_limit (vtx_flag , e );
540+ cs_parall_sum (2 , CS_GNUM_TYPE , counts + 2 ); /* do not change initial cell counts */
541+ }
542+
543+ BFT_FREE (vtx_flag );
544+ BFT_FREE (cell_vol_cmp );
545+
456546 if (counts [2 ] > 0 ) {
457547
458548 bft_printf
459549 (_ ("\nBoundary layer insertion:\n"
460550 " %llu cells would have a negative volume\n"
461551 " %llu cells would have a volume reduced by more than %g\n"
462552 " (which is the user-defined threshold)\n"
463- " reducing insertion at adjacent boundary vertices.\n" ),
553+ " reducing insertion at nearby boundary vertices.\n" ),
464554 (unsigned long long )counts [0 ], (unsigned long long )counts [1 ],
465555 min_volume_factor );
466556
467557 compute_displacement = true;
468558
469559 }
560+
470561 else if (counts [0 ] > 0 ) {
471- cs_base_warn (__FILE__ , __LINE__ );
472562 bft_printf
473563 (_ ("%llu cells would have a negative volume after boundary insertion\n"
474- "but none of these are adjacent to an inserted boundary.\n"
564+ "but none of these are near to an inserted boundary.\n"
475565 "Unable to detemine appropriate insertion limitation." ),
476566 (unsigned long long )counts [0 ]);
477567 }
0 commit comments