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:
{
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;
}
}
int cs_lnum_t
local mesh entity id
Definition: cs_defs.h:313
#define CS_MESH_MODIFIED
Definition: cs_mesh.h:62
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";
&n_selected_faces,
selected_faces);
false,
n_layers,
thickness,
expansion_factor,
n_selected_faces,
selected_faces);
}
#define BFT_MALLOC(_ptr, _ni, _type)
Allocate memory for _ni elements of type _type.
Definition: bft_mem.h:62
#define BFT_FREE(_ptr)
Free allocated memory.
Definition: bft_mem.h:101
void cs_mesh_extrude_constant(cs_mesh_t *m, bool interior_gc, cs_lnum_t n_layers, double thickness, double expansion_factor, cs_lnum_t n_faces, const cs_lnum_t faces[])
Extrude mesh boundary faces in the normal direction by a constant thickness.
Definition: cs_mesh_extrude.c:2206
void cs_selector_get_b_face_list(const char *criteria, cs_lnum_t *n_b_faces, cs_lnum_t b_face_list[])
Fill a list of boundary faces verifying a given selection criteria.
Definition: cs_selector.c:213
The example which follows illustrates the use of the advanced function to impose the vector of extrusion.
{
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};
for (int z_id = 0; z_id < n_zones; z_id++) {
zone_layers[z_id],
zone_thickness[z_id],
zone_expansion[z_id],
n_faces,
face_list);
}
}
e,
true);
}
cs_real_t cs_real_3_t[3]
vector of 3 floating-point values
Definition: cs_defs.h:332
void cs_mesh_extrude_face_info_destroy(cs_mesh_extrude_face_info_t **efi)
Destroy a mesh extrusion face information structure.
Definition: cs_mesh_extrude.c:2283
void cs_mesh_extrude_vectors_destroy(cs_mesh_extrude_vectors_t **e)
Destroy a mesh extrusion vectors definition.
Definition: cs_mesh_extrude.c:2392
cs_mesh_extrude_face_info_t * cs_mesh_extrude_face_info_create(const cs_mesh_t *m)
Create a mesh extrusion face information structure.
Definition: cs_mesh_extrude.c:2249
void cs_mesh_extrude(cs_mesh_t *m, const cs_mesh_extrude_vectors_t *e, bool interior_gc)
Extrude mesh boundary faces in the normal direction.
Definition: cs_mesh_extrude.c:1975
cs_mesh_extrude_vectors_t * cs_mesh_extrude_vectors_create(const cs_mesh_extrude_face_info_t *efi)
Create and build a mesh extrusion vectors definition.
Definition: cs_mesh_extrude.c:2361
void cs_mesh_extrude_set_info_by_zone(cs_mesh_extrude_face_info_t *efi, int n_layers, double distance, float expansion_factor, const cs_lnum_t n_faces, const cs_lnum_t face_ids[])
Set face extrusion information by zone.
Definition: cs_mesh_extrude.c:2314
Definition: cs_mesh_extrude.h:54
Definition: cs_mesh_extrude.h:78
cs_lnum_t * vertex_ids
Definition: cs_mesh_extrude.h:83
cs_lnum_t n_vertices
Definition: cs_mesh_extrude.h:81
cs_coord_3_t * coord_shift
Definition: cs_mesh_extrude.h:85
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";
&n_selected_faces,
selected_faces);
false,
n_layers,
thickness,
expansion_factor,
n_selected_faces,
selected_faces);
for(int i=0; i<n_selected_elts; i++)
selected_elts[i] = n_prev_cells + i;
"solid",
n_selected_elts,
selected_elts);
}
void cs_mesh_group_cells_add(cs_mesh_t *mesh, const char *name, cs_lnum_t n_selected_cells, const cs_lnum_t selected_cell_id[])
Add selected cells to a given group.
Definition: cs_mesh_group.c:1102
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.
{
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};
for (int z_id = 0; z_id < n_zones; z_id++) {
zone_layers[z_id],
zone_thickness[z_id],
zone_expansion[z_id],
n_faces,
face_list);
}
}
void cs_mesh_boundary_layer_insert(cs_mesh_t *m, cs_mesh_extrude_vectors_t *e, cs_real_t min_volume_factor, bool interior_gc, cs_lnum_t n_fixed_vertices, const cs_lnum_t *fixed_vertex_ids)
Insert mesh boundary layers.
Definition: cs_mesh_boundary_layer.c:410
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:
{
const char criteria[] = "box[0.5, 0.5, 0, 1, 1, 0.05]";
&n_selected_elts,
selected_elts);
"source_region",
n_selected_elts,
selected_elts);
}
void cs_selector_get_cell_list(const char *criteria, cs_lnum_t *n_cells, cs_lnum_t cell_list[])
Fill a list of cells verifying a given selection criteria.
Definition: cs_selector.c:369
The user can also add groups for boundary faces:
{
const char criteria[] = "box[0.5, 0.5, 0, 1, 1, 0.05]";
&n_selected_elts,
selected_elts);
"source_region",
n_selected_elts,
selected_elts);
}
void cs_mesh_group_b_faces_add(cs_mesh_t *mesh, const char *name, cs_lnum_t n_selected_faces, const cs_lnum_t selected_face_id[])
Add selected boundary faces to a given group.
Definition: cs_mesh_group.c:1156
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]";
&n_selected_cells,
selected_cells);
true,
n_selected_cells,
selected_cells);
}
void cs_mesh_refine_simple_selected(cs_mesh_t *m, bool conforming, cs_lnum_t n_cells, const cs_lnum_t cells[])
Refine selected mesh cells.
Definition: cs_mesh_refine.c:4786
The following code shows an example of mesh refinement based on the intersection with a CAD defined in an STL file.
{
"cad_file.stl");
3,
2);
}
cs_stl_mesh_t * cs_stl_mesh_add(const char *name)
Add a new STL mesh structure.
Definition: cs_stl.c:1048
void cs_stl_file_read(cs_stl_mesh_t *stl_mesh, const char *path)
Read a binary STL file and store its content in a STL mesh structure.
Definition: cs_stl.c:1157
void cs_stl_refine(cs_stl_mesh_t *stl_mesh, int n_ref, int n_add_layer)
Refine the mesh following a given STL mesh.
Definition: cs_stl.c:1969
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.
{
}
void cs_preprocessor_data_add_file(const char *file_name, size_t n_group_renames, const char **group_rename, const double transf_matrix[3][4])
Definition: cs_preprocessor_data.c:2150
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.
{
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.}};
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.
{
const char criteria[] = "box[-250, -250, 0, 250, 250, 100]";
&n_selected_elts,
selected_elts);
char *flag;
flag[i] = 0;
}
for (
cs_lnum_t i = 0; i < n_selected_elts; i++) {
flag[selected_elts[i]] = 1;
}
}
#define CS_MESH_MODIFIED_BALANCE
Definition: cs_mesh.h:65
void cs_mesh_remove_cells(cs_mesh_t *m, char flag[], const char *group_name)
Remove flagged cells.
Definition: cs_mesh_remove.c:109
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.
{
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;
int postprocess = 0;
postprocess);
void cs_mesh_warping_set_defaults(double max_warp_angle, int postprocess)
Definition: cs_mesh_warping.c:678
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;
int *vtx_is_fixed = NULL;
feature_angle,
vtx_is_fixed);
void cs_mesh_smoother_fix_by_feature(cs_mesh_t *mesh, cs_real_t feature_angle, int vtx_is_fixed[])
Set fixed vertices flag based on feature angle criterion.
Definition: cs_mesh_smoother.c:797
void cs_mesh_smoother_unwarp(cs_mesh_t *mesh, const int vtx_is_fixed[])
Unwarping smoother.
Definition: cs_mesh_smoother.c:886
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.
unsigned *bad_cell_flag = mesh_quantities->bad_cell_flag;
double *
volume = mesh_quantities->cell_vol;
for (cell_id = 0; cell_id < n_cells_wghosts; cell_id++)
bad_vol_cells[cell_id] = 0;
n_cells_tot = n_cells;
for (cell_id = 0; cell_id < n_cells; cell_id++) {
bad_vol_cells[cell_id] = 1;
}
}
ibad = 0;
iwarning = 0;
for (cell_id = 0; cell_id < n_cells; cell_id++) {
ibad++;
iwarning++;
}
}
}
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)
(" 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");
int bft_printf(const char *const format,...)
Replacement for printf() with modifiable behavior.
Definition: bft_printf.c:140
int cs_glob_rank_id
Definition: cs_defs.c:174
unsigned long cs_gnum_t
global mesh entity number
Definition: cs_defs.h:298
#define CS_BAD_CELL_USER
Definition: cs_mesh_bad_cells.h:55
static void cs_parall_counter(cs_gnum_t cpt[], const int n)
Sum values of a counter on all default communicator processes.
Definition: cs_parall.h:93
double precision, dimension(:), pointer volume
volume of each cell
Definition: mesh.f90:152
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;
int verbosity = 1;
int visualization = 1;
float fraction = 0.10, plane = 25.;
fraction,
plane,
verbosity,
visualization);
int cs_join_add(const char *sel_criteria, float fraction, float plane, int verbosity, int visualization)
Definition: cs_join.c:1591
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).
1.00,
0.10,
1,
1,
500,
100,
30,
25,
5.0,
2.0);
void cs_join_set_advanced_param(int join_num, double mtf, double pmf, int tcm, int icm, int max_break, int max_sub_faces, int tml, int tmb, double tmr, double tmr_distrib)
Set advanced parameters for the joining algorithm.
Definition: cs_join.c:1687
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;
int visualization = 1;
float fraction = 0.10, plane = 25.;
const double translation[3] = {1.0, 0.0, 0.0};
fraction,
plane,
verbosity,
visualization,
translation);
}
int cs_join_perio_add_translation(const char *sel_criteria, double fraction, double plane, int verbosity, int visualization, const double trans[3])
Definition: cs_join_perio.c:342
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;
int visualization = 1;
float fraction = 0.10, plane = 25.;
double axis[3] = {1.0, 0, 0};
double invariant[3] = {0, 0, 0};
fraction = 0.2;
verbosity = 2;
fraction,
plane,
verbosity,
visualization,
axis,
invariant);
}
int cs_join_perio_add_rotation(const char *sel_criteria, double fraction, double plane, int verbosity, int visualization, double theta, const double axis[3], const double invariant[3])
Definition: cs_join_perio.c:394
double precision, dimension(:,:,:), allocatable theta
Definition: atimbr.f90:123
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.
{
int join_num;
int verbosity = 1;
int visualization = 1;
float fraction = 0.10, plane = 25.;
double matrix[3][4] = {{1., 0., 0., 0.5},
{0., 1., 0., 0.},
{0., 0., 1., 0.}};
fraction,
plane,
verbosity,
visualization,
}
int cs_join_perio_add_mixed(const char *sel_criteria, double fraction, double plane, int verbosity, int visualization, double matrix[3][4])
Definition: cs_join_perio.c:454
void matrix(const int *iconvp, const int *idiffp, const int *ndircp, const int *isym, const cs_real_t *thetap, const int *imucpp, const cs_real_t coefbp[], const cs_real_t cofbfp[], const cs_real_t rovsdt[], const cs_real_t i_massflux[], const cs_real_t b_massflux[], const cs_real_t i_visc[], const cs_real_t b_visc[], const cs_real_t xcpp[], cs_real_t da[], cs_real_t xa[])
Definition: cs_matrix_building.c:111
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.
{
double mtf = 1.00;
double pmf = 0.10;
int tcm = 1;
int max_break = 500;
int max_sub_face = 100;
int tml = 30;
int tmb = 25;
double tmr = 5.0;
double tmr_distrib = 2.0;
max_break, max_sub_face,
tml, tmb, tmr, tmr_distrib);
}
integer, save icm
the intersection computation mode. If its value is:
Definition: ppincl.f90:332
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:
{
const char criteria[] = "plane[0, -1, 0, 0.5, epsilon = 0.0001]"
" or plane[-1, 0, 0, 0.5, epsilon = 0.0001]";
&n_selected_faces,
selected_faces);
n_selected_faces,
selected_faces);
}
void cs_mesh_boundary_insert(cs_mesh_t *mesh, cs_lnum_t n_faces, cs_lnum_t face_id[])
Insert boundary into the mesh.
Definition: cs_mesh_boundary.c:1753
void cs_selector_get_i_face_list(const char *criteria, cs_lnum_t *n_i_faces, cs_lnum_t i_face_list[])
Fill a list of interior faces verifying a given selection criteria.
Definition: cs_selector.c:291
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.
{
const char criteria[] = "box[0.5, 0.5, 0, 1, 1, 0.05]";
&n_selected_cells,
selected_cells);
"zone_interface",
n_selected_cells,
selected_cells);
}
void cs_mesh_boundary_insert_separating_cells(cs_mesh_t *mesh, const char *group_name, cs_lnum_t n_cells, const cs_lnum_t cell_id[])
Insert a boundary into the mesh between a given set of cells and the the others.
Definition: cs_mesh_boundary.c:1818