7.2
general documentation
Examples of mesh modifications

Introduction

C user functions for optional modification of the mesh. These subroutines are called in all cases.

Several functions are present in the file, each specific to different modification types.

General mesh modifications

Mesh modifications not available through specialized functions should be defined in cs_user_mesh_modify.

Coordinates modification

For example, to modify coordinates, the following code can be added:

{
cs_lnum_t vtx_id;
const double coo_mult = 1. / 1000.;
for (vtx_id = 0; vtx_id < mesh->n_vertices; vtx_id++) {
mesh->vtx_coord[vtx_id*3] *= coo_mult;
mesh->vtx_coord[vtx_id*3 + 1] *= coo_mult;
mesh->vtx_coord[vtx_id*3 + 2] *= coo_mult;
}
/* Set mesh modification flag if it should be saved for future re-use. */
mesh->modified |= CS_MESH_MODIFIED;
}

Boundary mesh patch extrusion

It is possible to add cells by extruding selected boundary faces. A simplified usage is available through the cs_mesh_extrude_constant, while extruding a mesh with vertex-local extrusion parameters is available through cs_mesh_extrude.

The example which follows illustrates the use of the simplified function.

{
int n_layers = 2;
double thickness = 1.0;
double expansion_factor = 1.5;
const char criteria[] = "outlet";
/* Select boundary faces */
cs_lnum_t n_selected_faces = 0;
cs_lnum_t *selected_faces = NULL;
BFT_MALLOC(selected_faces, mesh->n_b_faces, cs_lnum_t);
cs_selector_get_b_face_list(criteria,
&n_selected_faces,
selected_faces);
/* Extrude selected boundary */
cs_mesh_extrude_constant(mesh,
false, /* Maintain groups on interior faces? */
n_layers,
thickness,
expansion_factor,
n_selected_faces,
selected_faces);
/* Free temporary memory */
BFT_FREE(selected_faces);
}

The example which follows illustrates the use of the advanced function to impose the vector of extrusion.

{
/* Define extrusion parameters for each face */
int n_zones = 2;
const char *sel_criteria[] = {"wall_1", "wall_2"};
const int zone_layers[] = {2, 4};
const double zone_thickness[] = {0.1, 0.3};
const float zone_expansion[] = {0.8, 0.7};
cs_mesh_extrude_face_info_t *efi = cs_mesh_extrude_face_info_create(mesh);
cs_lnum_t n_faces;
cs_lnum_t *face_list;
BFT_MALLOC(face_list, mesh->n_b_faces, cs_lnum_t);
for (int z_id = 0; z_id < n_zones; z_id++) {
cs_selector_get_b_face_list(sel_criteria[z_id], &n_faces, face_list);
cs_mesh_extrude_set_info_by_zone(efi,
zone_layers[z_id],
zone_thickness[z_id],
zone_expansion[z_id],
n_faces,
face_list);
}
BFT_FREE(face_list);
/* Determine vertex values for extrusion */
cs_mesh_extrude_vectors_t *e = cs_mesh_extrude_vectors_create(efi);
/* Overwrite the total coord_shift */
cs_real_3_t *vtx_coords = (cs_real_3_t *)mesh->vtx_coord;
for (cs_lnum_t id = 0; id < e->n_vertices; id++) {
cs_lnum_t v_id = e->vertex_ids[id];
e->coord_shift[id][0] = 0.;
e->coord_shift[id][1] = 0.;
/* Percentage of the original z coordinate */
e->coord_shift[id][2] = 0.01 * vtx_coords[v_id][2];
}
/* Extrude mesh with this */
cs_mesh_extrude_face_info_destroy(&efi);
cs_mesh_extrude(mesh,
e,
true); /* Maintain group classes of interior
faces previously on boundary */
cs_mesh_extrude_vectors_destroy(&e);
}

The example which follows illustrates how to tag the extruded cells with a user defined criteria.

{
int n_layers = 2;
double thickness = 1.0;
double expansion_factor = 1.5;
const char criteria[] = "walls";
/* Save the initial number of cells */
cs_lnum_t n_prev_cells = mesh->n_cells ;
/* Select boundary faces */
cs_lnum_t n_selected_faces = 0;
cs_lnum_t *selected_faces = NULL;
BFT_MALLOC(selected_faces, mesh->n_b_faces, cs_lnum_t);
cs_selector_get_b_face_list(criteria,
&n_selected_faces,
selected_faces);
/* Extrude selected boundary */
cs_mesh_extrude_constant(mesh,
false, /* Maintain groups on interior faces? */
n_layers,
thickness,
expansion_factor,
n_selected_faces,
selected_faces);
/* Free temporary memory */
BFT_FREE(selected_faces);
/* Compute the number of extruded cells */
cs_lnum_t n_selected_elts = mesh->n_cells - n_prev_cells ;
/* Among all the cells, only select the cells above
* the initial number of cells (before extrusion). */
cs_lnum_t *selected_elts = NULL;
BFT_MALLOC(selected_elts, mesh->n_cells, cs_lnum_t);
for(int i=0; i<n_selected_elts; i++)
selected_elts[i] = n_prev_cells + i;
/* Add selected cells to a new group called "solid" */
cs_mesh_group_cells_add(mesh,
"solid",
n_selected_elts,
selected_elts);
BFT_FREE(selected_elts);
}

Boundary layer insertion

Boundary faces extrusion can also be used to insert boundary layer cells, by first shrinking the mesh around the selected zones so as to accommodate for added cells. The following example shows how this can be done for 2 different zones, using different parameters for each zone. Note that adjacent zones with a different number of inserted layers are possible. Also a positive zone thickness implies a fixed thicknesss, while a negative thickness is interpreted as a ratio relative to the mean adjacent cell size, allowing for automatic and adapted local thickness definition.

{
/* Define extrusion parameters for each face */
int n_zones = 2;
const char *sel_criteria[] = {"wall_1", "wall_2"};
const int zone_layers[] = {2, 4};
const double zone_thickness[] = {0.1, 0.3};
const float zone_expansion[] = {0.8, 0.7};
cs_mesh_extrude_face_info_t *efi = cs_mesh_extrude_face_info_create(mesh);
cs_lnum_t n_faces;
cs_lnum_t *face_list;
BFT_MALLOC(face_list, mesh->n_b_faces, cs_lnum_t);
for (int z_id = 0; z_id < n_zones; z_id++) {
cs_selector_get_b_face_list(sel_criteria[z_id], &n_faces, face_list);
cs_mesh_extrude_set_info_by_zone(efi,
zone_layers[z_id],
zone_thickness[z_id],
zone_expansion[z_id],
n_faces,
face_list);
}
BFT_FREE(face_list);
/* Determine vertex values for extrusion */
cs_mesh_extrude_vectors_t *e = cs_mesh_extrude_vectors_create(efi);
/* Insert boundary layer */
cs_mesh_extrude_face_info_destroy(&efi);
cs_mesh_boundary_layer_insert(mesh, e, 0.2, false, 0, NULL);
cs_mesh_extrude_vectors_destroy(&e);
}

Groups of cells, interior, and boundary faces may be created or modified, using the cs_mesh_group_cells_set, cs_mesh_group_i_faces_set, cs_mesh_group_b_faces_set functions to assign a group to selected elements, removing all previous group attributions for those elements, and cs_mesh_group_cells_add, cs_mesh_group_i_faces_add, cs_mesh_group_b_faces_add may be used to add those elements to a group while also keeping previous group information for those elements.

The mesh is not marked as modified by default for this "light" modification, so the user may force this using a modification flag, as in the example here:

{
cs_lnum_t n_selected_elts = 0;
cs_lnum_t *selected_elts = NULL;
const char criteria[] = "box[0.5, 0.5, 0, 1, 1, 0.05]";
BFT_MALLOC(selected_elts, mesh->n_cells, cs_lnum_t);
cs_selector_get_cell_list(criteria,
&n_selected_elts,
selected_elts);
cs_mesh_group_cells_add(mesh,
"source_region",
n_selected_elts,
selected_elts);
BFT_FREE(selected_elts);
/* Mark mesh as modified to save it */
mesh->modified |= CS_MESH_MODIFIED;
}

The user can also add groups for boundary faces:

{
cs_lnum_t n_selected_elts = 0;
cs_lnum_t *selected_elts = NULL;
const char criteria[] = "box[0.5, 0.5, 0, 1, 1, 0.05]";
BFT_MALLOC(selected_elts, mesh->n_b_faces, cs_lnum_t);
cs_selector_get_b_face_list(criteria,
&n_selected_elts,
selected_elts);
cs_mesh_group_b_faces_add(mesh,
"source_region",
n_selected_elts,
selected_elts);
BFT_FREE(selected_elts);
/* Mark mesh as modified to save it */
mesh->modified |= CS_MESH_MODIFIED;
}

Mesh refinement

Cells may be refined automatically, using templates adapted to each cell type (see cs_mesh_refine.c").

This simply required selecting the cells to refine and whether the refinement should be partially propagated to neighboring cells (conforming) or not.

The following code shows an example of mesh refinement for a given region.

{
const char criteria[] = "box[0, 0, 0, 0.5, 0.5, 0.5]";
cs_lnum_t n_selected_cells = 0;
cs_lnum_t *selected_cells = NULL;
BFT_MALLOC(selected_cells, mesh->n_cells, cs_lnum_t);
cs_selector_get_cell_list(criteria,
&n_selected_cells,
selected_cells);
cs_mesh_refine_simple_selected(mesh,
true, /* conforming or not */
n_selected_cells,
selected_cells);
BFT_FREE(selected_cells);
}

Mesh reading and modification

The user function cs_user_mesh_input allows a detailed selection of imported meshes, allows reading files multiple times, applying geometric transformations, and renaming groups.

The following code shows an example of mesh reading with no transformation.

{
/* Determine list of files to add */
/*--------------------------------*/
/* Read input mesh with no modification */
cs_preprocessor_data_add_file("mesh_input/mesh_01", 0, NULL, NULL);
}

Multiple inputs read through this function are automatically concatenated, and may be later joined using the mesh joining options.

A mesh can also be read while its groups are renamed, and its geometry transformed.

{
/* Add same mesh with transformations */
const char *renames[] = {"Inlet", "Injection_2",
"Group_to_remove", NULL};
const double transf_matrix[3][4] = {{1., 0., 0., 5.},
{0., 1., 0., 0.},
{0., 0., 1., 0.}};
cs_preprocessor_data_add_file("mesh_input/mesh_02",
2, renames,
transf_matrix);
}

Geometric transformations are defined using a homogeneous coordinates transformation matrix. Such a matrix has 3 lines and 4 columns, with the first 3 columns describing a rotation/scaling factor, and the last column describing a translation. A 4th line is implicit, containing zeroes off-diagonal, and 1 on the diagonal. The advantages of this representation is that any rotation/translation/scaling combination may be expressed by matrix multiplication, while simple rotations or translations may still be defined easily.

Cells from a selection can also be removed.

{
cs_lnum_t n_selected_elts = 0;
cs_lnum_t *selected_elts = NULL;
const char criteria[] = "box[-250, -250, 0, 250, 250, 100]";
BFT_MALLOC(selected_elts, mesh->n_cells, cs_lnum_t);
cs_selector_get_cell_list(criteria,
&n_selected_elts,
selected_elts);
char *flag;
BFT_MALLOC(flag, mesh->n_cells, char);
for (cs_lnum_t i = 0; i < mesh->n_cells; i++) {
flag[i] = 0;
}
for (cs_lnum_t i = 0; i < n_selected_elts; i++) {
flag[selected_elts[i]] = 1;
}
cs_mesh_remove_cells(mesh, flag, "[Building]");
BFT_FREE(selected_elts);
BFT_FREE(flag);
/* Mark for re-partitioning */
mesh->modified |= CS_MESH_MODIFIED_BALANCE;
}

Mesh saving

The user function cs_user_mesh_save can enable or disable mesh saving. By default, mesh is saved when modified. The following code shows an example of disabled saving.

{
/* Disable saving of modified mesh by setting flag to 0;
Force saving of unmodified mesh by setting it to 2. */
mesh->save_if_modified = 0;
}

Mesh quality modifications

Mesh warping

The cs_user_mesh_warping function allows the user to cut the warped faces of his mesh using the cs_mesh_warping_set_defaults function to define the maximum warped angle.

double max_warp_angle = 3; /* bounded between 0 and 90 degrees */
int postprocess = 0;
cs_mesh_warping_set_defaults(max_warp_angle,
postprocess);

Mesh smoothing

The smoothing utilities may be useful when the calculation mesh has local defects. The principle of smoothers is to mitigate the local defects by averaging the mesh quality. This procedure can help for calculation robustness or/and results quality. The user function cs_user_mesh_smoothe allows to use different smoothing functions detailed below.

The following code shows an example of use of the cs_mesh_smoother functions, cs_mesh_smoother_fix_by_feature which fixes all boundary vertices that have one of their feature angles less than the maximum feature angle defined by the user and cs_mesh_smoother_unwarp which reduces face warping in the calculation mesh.

double feature_angle = 25; /* bounded between 0 and 90 degrees */
int *vtx_is_fixed = NULL;
BFT_MALLOC(vtx_is_fixed, mesh->n_vertices, int);
/* Get fixed boundary vertices flag */
cs_mesh_smoother_fix_by_feature(mesh,
feature_angle,
vtx_is_fixed);
/* Call unwarping smoother */
cs_mesh_smoother_unwarp(mesh, vtx_is_fixed);
/* Free memory */
BFT_FREE(vtx_is_fixed);

Bad cells tagging

Bad cells of a mesh can be tagged based on user-defined geometric criteria. The following example shows how to tag cells that have a volume below a certain value and then post-process the tagged cells. This is done using the cs_user_mesh_bad_cells_tag function.

const cs_lnum_t n_cells = mesh->n_cells;
const cs_lnum_t n_cells_wghosts = mesh->n_cells_with_ghosts;
unsigned *bad_cell_flag = mesh_quantities->bad_cell_flag;
double *volume = mesh_quantities->cell_vol;
cs_lnum_t cell_id;
cs_gnum_t n_cells_tot, iwarning, ibad;
/*---------------------------------------*/
/* User condition: check for cell volume */
/*---------------------------------------*/
cs_lnum_t *bad_vol_cells = NULL;
BFT_MALLOC(bad_vol_cells, n_cells_wghosts, cs_lnum_t);
for (cell_id = 0; cell_id < n_cells_wghosts; cell_id++)
bad_vol_cells[cell_id] = 0;
/* Get the total number of cells within the domain */
n_cells_tot = n_cells;
cs_parall_counter(&n_cells_tot, 1);
for (cell_id = 0; cell_id < n_cells; cell_id++) {
/* Compare the cell volume to the user condition */
if (volume[cell_id] < 0.01) {
/* Local array used to post-process results
--> user tagged cells are set to 1 */
bad_vol_cells[cell_id] = 1;
/* Array used to store bad cells flag
--> user tagged cells are flagged using the mask */
bad_cell_flag[cell_id] = bad_cell_flag[cell_id] | CS_BAD_CELL_USER;
}
}
ibad = 0;
iwarning = 0;
for (cell_id = 0; cell_id < n_cells; cell_id++) {
if (bad_cell_flag[cell_id] & CS_BAD_CELL_USER) {
ibad++;
iwarning++;
}
}
/* Parallel processing */
if (cs_glob_rank_id >= 0) {
cs_parall_counter(&ibad, 1);
cs_parall_counter(&iwarning, 1);
}
/* Display log output */
bft_printf("\n Criteria 6: User Specific Tag:\n");
bft_printf(" Number of bad cells detected: %llu --> %3.0f %%\n",
(unsigned long long)ibad,
(double)ibad / (double)n_cells_tot * 100.0);
if (iwarning > 0)
bft_printf
(" Warning:\n"
" --------\n"
" Mesh quality issue based on user criteria has been detected\n\n"
" The mesh should be re-considered using the listed criteria.\n\n");
/* Post processing is activated automatically in mesh quality mode
for the CS_BAD_CELL_USER tag, as for all tags defined in
cs_user_bad_cells.h.
For a different tag, postprocessing could be done with a user
postprocessing function, or calling directly cs_post_write_var(). */
BFT_FREE(bad_vol_cells);

Mesh joining

Simple mesh joining

Conforming joining of possibly non-conforming meshes may be done by the cs_user_join user function. For a simple mesh joining, the cs_join_add subroutine is sufficient.

int join_num;
/* Add a joining operation */
/* ----------------------- */
int verbosity = 1; /* per-task dump if > 1, debug level if >= 3 */
int visualization = 1; /* debug level if >= 3 */
float fraction = 0.10, plane = 25.;
join_num = cs_join_add("98 or 99",
fraction,
plane,
verbosity,
visualization);

Advanced mesh joining

For lower quality meshes or curved and anisotropically refined meshes leading to joining difficulties, or to reduce memory usage, using the the cs_join_set_advanced_param function refering to a defined joning allows finer control (see the function reference for parameter details).

cs_join_set_advanced_param(join_num,
1.00, /* merge tolerance factor */
0.10, /* pre-merge factor */
1, /* tolerance computation mode */
1, /* intersection computation mode */
500, /* max. nb. of equivalence breaks */
100, /* maximum number of sub-faces */
30, /* tree max level */
25, /* tree max boxes per node */
5.0, /* tree max ratio */
2.0); /* distribution tree max ratio */

Mesh periodicity

Handling of periodicity can be performed with the cs_user_periodicity function.

Periodicity of translation

The following example illustrates the periodicity of translation case using the cs_join_perio_add_translation subroutine.

{
int join_num;
int verbosity = 1; /* per-task dump if > 1, debug level if >= 3 */
int visualization = 1; /* debug level if >= 3 */
float fraction = 0.10, plane = 25.;
const double translation[3] = {1.0, 0.0, 0.0}; /* Translation vector */
join_num = cs_join_perio_add_translation("98 or 99",
fraction,
plane,
verbosity,
visualization,
translation);
}

Periodicity of rotation

The following example illustrates the periodicity of rotation case using the cs_join_perio_add_rotation subroutine.

{
int join_num;
int verbosity = 1; /* per-task dump if > 1, debug level if >= 3 */
int visualization = 1; /* debug level if >= 3 */
float fraction = 0.10, plane = 25.;
double theta = 20; /* angle in degrees */
double axis[3] = {1.0, 0, 0}; /* axis of rotation */
double invariant[3] = {0, 0, 0}; /* invariant point */
/* change default values */
fraction = 0.2;
verbosity = 2;
join_num = cs_join_perio_add_rotation("3",
fraction,
plane,
verbosity,
visualization,
theta,
axis,
invariant);
}

General periodicity

The following example illustrates a more general case of periodicity which combines different kinds of transformation. The function cs_join_perio_add_mixed is used.

{
/* We define a general transformation using a homogeneous matrix:
*
* We define the first 3 rows of a 4x4 matrix:
* _ _
* | r11 r12 r13 tx | t(x,y,z) : translation vector
* | r21 r22 r23 ty | r(i,j) : rotation matrix
* | r31 r32 r33 tz |
* |_ 0 0 0 1 _|
*
* Transformations may be combined using matrix multiplication,
* so this be used for helecoidal transformations for instance. */
int join_num;
int verbosity = 1; /* per-task dump if > 1, debug level if >= 3 */
int visualization = 1; /* debug level if >= 3 */
float fraction = 0.10, plane = 25.;
double matrix[3][4] = {{1., 0., 0., 0.5},
{0., 1., 0., 0.},
{0., 0., 1., 0.}};
join_num = cs_join_perio_add_mixed("all[]",
fraction,
plane,
verbosity,
visualization,
matrix);
}

Periodicity advanced parameters

As with the Advanced mesh joining subsection, a more complex periodicity can be defined using the cs_join_set_advanced_param subroutine.

{
/* Merge tolerance factor:
* used to locally modify the tolerance associated to each
* vertex BEFORE the merge step.
* = 0 => no vertex merge;
* < 1 => vertex merge is more strict. It may increase the number
* of tolerance reduction and so define smaller subset of
* vertices to merge together but it can drive to loose
* intersections;
* = 1 => no change;
* > 1 => vertex merge is less strict. The subset of vertices able
* to be merged together is greater. */
double mtf = 1.00;
/* Pre-merge factor:
* this parameter is used to define a bound under which two vertices
* are merged before the merge step.
* Tolerance limit for the pre-merge = pmf * fraction. */
double pmf = 0.10;
/* Tolerance computation mode:
*
* 1: (default) tol = min. edge length related to a vertex * fraction
* 2: tolerance is computed like in mode 1 with in addition, the
* multiplication by a coefficient equal to the max sin(e1, e2)
* where e1 and e2 are two edges sharing the same vertex V for which
* we want to compute the tolerance.
* 11: as 1 but taking into account only the selected faces
* 12: as 2 but taking into account only the selected faces */
int tcm = 1;
/* Intersection computation mode:
* 1: (default) Original algorithm.
* Try to snap intersection to extremity.
* 2: New intersection algorithm.
* Avoid snapping intersection to extremity. */
int icm = 1;
/* Maximum number of equivalence breaks which is
* enabled during the merge step. */
int max_break = 500;
/* Maximum number of sub-faces when splitting a selected face. */
int max_sub_face = 100;
/* tml, tmb and tmr are parameters of the searching algorithm for
* face intersections between selected faces (octree structure).
* Useful to adjust speed vs. memory consumption. */
/* Tree Max Level:
* deepest level reachable when building the tree */
int tml = 30;
/* Tree Max Boxes:
* max. number of bounding boxes which can be linked to a leaf of the tree
* (not necessary true for the deepest level) */
int tmb = 25;
/* Tree Max. Ratio:
* stop refining the tree structure when
* number of bounding boxes > tmr * number of faces to locate.
* In parallel, a separate (usually lower) value may be set for
* the initial coarse tree used to determine distribution.
* Reducing this will help reduce memory consumption. */
double tmr = 5.0;
double tmr_distrib = 2.0;
/* Set advanced parameters */
cs_join_set_advanced_param(join_num,
mtf, pmf, tcm, icm,
max_break, max_sub_face,
tml, tmb, tmr, tmr_distrib);
}

Mesh boundary insertion

The user function cs_user_mesh_boundary allows insertion of boundaries in the calculation mesh. This function transforms the selected interior faces into boundary faces, on which boundary conditions can (and must) be applied. Vertices are also split so as to be distinct on each side of the boundary.

Boundaries can be directly inserted based on a selection of interior faces, such as shown here:

{
cs_lnum_t n_selected_faces = 0;
cs_lnum_t *selected_faces = NULL;
/* example of multi-line character string */
const char criteria[] = "plane[0, -1, 0, 0.5, epsilon = 0.0001]"
" or plane[-1, 0, 0, 0.5, epsilon = 0.0001]";
BFT_MALLOC(selected_faces, mesh->n_i_faces, cs_lnum_t);
cs_selector_get_i_face_list(criteria,
&n_selected_faces,
selected_faces);
cs_mesh_boundary_insert(mesh,
n_selected_faces,
selected_faces);
BFT_FREE(selected_faces);
}

Boundaries can also be inserted between a set of selected cells and the rest of the mesh. In this case, a mesh group name can be assigned to the added boundary faces.

{
cs_lnum_t n_selected_cells = 0;
cs_lnum_t *selected_cells = NULL;
const char criteria[] = "box[0.5, 0.5, 0, 1, 1, 0.05]";
BFT_MALLOC(selected_cells, mesh->n_cells, cs_lnum_t);
cs_selector_get_cell_list(criteria,
&n_selected_cells,
selected_cells);
cs_mesh_boundary_insert_separating_cells(mesh,
"zone_interface",
n_selected_cells,
selected_cells);
BFT_FREE(selected_cells);
}