Generator: Mesh Generation
Version: 2.4 (29/03/2017)
Author: Onera

Preamble
This module can be used to generate meshes from geometries.
A mesh
can be stored as an array (as defined in the Converter documentation)
or in a zone node of a CGNS/python tree (pyTree).
This module is part of Cassiopee, a free open-source
pre- and post-processor for CFD simulations.
When using the array interface, import the Generator module:
import Generator as G
Then, a is an array and A is a list of arrays.
When using the pyTree interface, import the module:
import Generator.PyTree as G
Then, a is a zone node and A is a list of zone nodes or a pyTree.
Basic grid generation
G.cart: create a structured Cartesian mesh with ni x nj x nk
points starting from point (xo,yo,zo) and of step (hi,hj,hk):
a = G.cart((xo,yo,zo), (hi,hj,hk), (ni,nj,nk))
Example of use: Cartesian mesh (array),
Cartesian mesh (pyTree).
G.cartHexa: create an unstructured hexahedral mesh defined from a Cartesian grid of ni x nj x nk points
starting from point (xo,yo,zo)
and of step (hi,hj,hk). Type of elements are 'QUAD' for 2D arrays and 'HEXA' for 3D arrays:
a = G.cartHexa((xo,yo,zo), (hi,hj,hk), (ni,nj,nk))
Example of use: hexahedral Cartesian mesh (array),
hexahedral Cartesian mesh (pyTree).
G.cartTetra: create an unstructured tetrahedral mesh defined from a Cartesian grid of ni x nj x nk points starting from point (xo,yo,zo)
and of step (hi,hj,hk). Type of elements are 'TRI' for 2D arrays and 'TETRA' for 3D arrays:
a = G.cartTetra((xo,yo,zo), (hi,hj,hk), (ni,nj,nk))
Example of use: tetrahedral Cartesian mesh (array),
tetrahedral Cartesian mesh (pyTree).
G.cartPenta: create an unstructured prismatic mesh defined from a regular Cartesian mesh.
The initial Cartesian mesh is defined by ni x nj x nk points starting from point (xo,yo,zo)
and of step (hi,hj,hk). Type of elements is 'PENTA':
a = G.cartPenta((xo,yo,zo), (hi,hj,hk), (ni,nj,nk))
Example of use: prismatic mesh (array),
prismatic mesh (pyTree).
G.cartPyra: create an unstructured pyramidal mesh defined from a regular Cartesian mesh.
The initial Cartesian mesh is defined by ni x nj x nk points starting from point (xo,yo,zo)
and of step (hi,hj,hk). Type of elements is 'PYRA':
a = G.cartPyra((xo,yo,zo), (hi,hj,hk), (ni,nj,nk))
Example of use: pyramidal mesh(array),
pyramidal mesh (pyTree).
G.cartNGon: create a NGON mesh defined from a regular Cartesian mesh.
The initial Cartesian mesh is defined by ni x nj x nk points starting from point (xo,yo,zo)
and of step (hi,hj,hk). Type of elements is 'NGON':
a = G.cartNGon((xo,yo,zo), (hi,hj,hk), (ni,nj,nk))
Example of use: NGON Cartesian mesh(array),
NGON Cartesian mesh (pyTree).
G.cylinder: create a regular cylindrical grid (or a portion of cylinder between tetas and tetae)
with ni x nj x nk points, of center-bottom point (xo,yo,zo), of
inner radius R1, outer radius R2 and height H. For a direct mesh,
use tetae < tetas:
a = G.cylinder((xo,yo,zo), R1, R2, tetas, tetae, H, (ni,nj,nk))
Example of use: regular cylindrical mesh (array),
regular cylindrical mesh (pyTree).
G.cylinder2: create an irregular cylindrical grid (or a portion of cylinder between
tetas and tetae) with ni x nj x nk points, of center-bottom point
(xo,yo,zo), of inner radius R1, outer radius R2, height H and with distributions in
r, teta, z. Distributions are arrays defining 1D meshes (x and i varying)
giving a distribution in [0,1]. Their number of points gives ni, nj, nk:
a = G.cylinder2((xo,yo,zo), R1, R2, tetas, tetae, H, arrayR, arrayTeta, arrayZ)
Example of use: irregular cylindrical mesh (array),
irregular cylindrical mesh (pyTree).
G.cylinder3: create an irregular cylindrical grid (or a portion of cylinder between
tetas and tetae) from a xz plane mesh defined by a and a
teta distribution defined by arrayTeta:
b = G.cylinder3(a, tetas, tetae, arrayTeta)
Example of use: irregular cylindrical mesh (array),
irregular cylindrical mesh (pyTree).
General purpose grid generators
G.delaunay: create a 2D Delaunay type mesh from an array. The array can be a 2D structured array, or an unstructured array of type 'NODE', 'TRI' or 'QUAD'.
Tol is a geometric tolerance. Points nearer than tol are merged. If keepBB is set to 1, the bounding box is kept in the final triangulation:
b = G.delaunay(a, tol=1.e-10, keepBB=0)
Example of use: Delaunay mesh from a 2D Cartesian grid (array),
Delaunay mesh from a 2D Cartesian grid (pyTree).
G.constrainedDelaunay: create a constrained Delaunay triangulation of the convex hull of a contour c.
Contour must be a BAR-array and must be in the plane (x,y).
Tol is a geometric tolerance.
Points nearer than tol are merged.
If keepBB is set to 1, the bounding box is kept in the final triangulation:
b = G.constrainedDelaunay(c, tol=1.e-10, keepBB=0)
Example of use: constrained Delaunay mesh from a polyline (array),
constrained Delaunay mesh from a polyline (pyTree).
G.checkDelaunay: check if the Delaunay triangulation defined in tri is inside the contour c:
c2 = G.checkDelaunay(c, tri)
Example of use: check the Delaunay triangulation inside the contour(array),
check the Delaunay triangulation inside the contour(pyTree).
G.T3mesher2D: create a 2D Delaunay mesh given a BAR defined in a. If
triangulateOnly=1 then only points of a are triangulated,
if triangulateOnly=0, then interior points are inserted:
b = G.T3mesher2D(a, triangulateOnly=0)
Example of use: create a Delaunay mesh in a circle (array),
create a Delaunay mesh in a circle (pyTree).
G.TetraMesher: create a 3D tetra mesh given a TRI surface defined in a. If the TRI surface has external normals, tetras are filled inside
the surface.
If algo=0, netgen is used, if algo=1, tetgen is used:
b = G.tetraMesher(a, algo=1)
Example of use: create a tetra mesh in a cube (array),
create a tetra mesh in a cube (pyTree).
G.TFI: generate a mesh by transfinite interpolation (TFI). Generated mesh can be 2D or 3D structured, or unstructured TRI or PENTA mesh.
Warning: the boundaries can be in a different order from the examples below, except for the PENTA TFI meshes.
2D structured mesh is built from imin, imax, jmin, jmax boundaries.
3D structured mesh is built from imin, imax, jmin, jmax, kmin, kmax boundaries.
Dimensions must be equal for each pair (imin,imax), (jmin,jmax)...
TRI mesh is built from imin, jmin, diag boundaries. Each boundary is a structured array with the same dimension.
PENTA mesh is built from Tmin, Tmax triangles boundary and imin, imax, diag boundaries. Tmin, Tmax must be structured triangles of dimension nxn.
imin, jmin, diag must be structured n*p arrays:
2D-struct: a = G.TFI([imin, imax, jmin, jmax])
3D-struct: a = G.TFI([imin, imax, jmin, jmax, kmin, kmax])
TRI: a = G.TFI([imin, jmin, diag])
PENTA: a = G.TFI([Tmin, Tmax, imin, imax, diag])
Example of use: TFI example (array),
TFI example (pyTree).
G.TFITri: generate three meshes by transfinite interpolation around three given curves a1, a2, a3. N3-N2+N1 must be odd:
A = G.TFITri(a1, a2, a3)
Example of use: TFITri example (array),
TFITri example (pyTree).
G.TFIO: generate five meshes by transfinite interpolation around one given curves a. The number of points of a must be odd:
A = G.TFIO(a)
Example of use: TFIO example (array),
TFIO example (pyTree).
G.TFIHalfO: generate four meshes by transfinite interpolation around two given curves a1 and a2 forming a half-O. N1 and N2 must be odd:
A = G.TFIHalfO(a)
Example of use: TFIHalfO example (array),
TFIHalfO example (pyTree).
G.TFIMono: generate one mesh by transfinite interpolation around two given curves a1 and a2 forming a half-O. N1-N2 must be even:
A = G.TFIMono(a)
Example of use: TFIMono example (array),
TFIMono example (pyTree).
G.hyper2D: generate an hyperbolic mesh (2D) of "C" or "O" type from a from
a line defined by line a and from a distribution defined by distrib.
The resulting mesh is nearly orthogonal:
b = G.hyper2D(line, distrib, "C")
Example of use: hyperbolic mesh generation (array),
hyperbolic mesh generation (pyTree).
G.PolyLine.polyLineMesher: generate a 2D mesh around a 2D polyline:
B = G.PolyLine.polyLineMesher(a, h, hf, density)
where a is the input polyline (BAR-array), h is the height of the mesh,
hf is the height of the first cell and density is the number of points
per unity of length.
In the 'array' version, it returns a list where B[0] is the list
of generated meshes, B[1] is the list of wall boundaries, B[2] is
the list of overlap boundaries, B[3] is h, B[4] is density
(eventually modified by the mesher).
In the pyTree version, it returns a list [zones,hs,densities], where zones is a list of zones of a CGNS python tree, containing the blocks,
wall boundaries, match and overlap boundaries; hs is the list of heights (modified if necessary),
and densities the list of densities (also modified if necessary).
Example of use: generates meshes around a polyline (array),
generates meshes around a polyline (pyTree).
G.PolyC1.polyC1Mesher: generate a 2D mesh around a 2D polyC1 curve:
B = G.PolyC1.polyC1Mesher(A, h, hf, density, splitCrit=10.)
where A is a list of i-arrays each representing a
C1 curve. All i-arrays put together must represent a polyC1 curve.
SplitCrit is a curvature radius triggering split.
Other arguments are similar to polyLineMesher. The function return
is also similar to polyLineMesher.
Example of use: generates meshes around a polyC1 curve (array),
generates meshes around a polyC1 curve (pyTree).
G.pointedHat: create a structured mesh from a curve defined by a i-array and a point. For the pyTree version: if a contains a solution, it is not taken into
account in b:
b = G.pointedHat(a, (x,y,z))
Example of use: pointed hat (array),
pointed hat (pyTree).
G.stitchedHat: create a stitched mesh from a curve defined by a i-array. The surface
is stitched in the middle. Tol is the accuracy of the search, tol2 is
a merging tolerance and offx, offy, off z an optional offset.
For the pyTree version: if a contains a solution, it is not taken into account in b:
b = G.stitchedHat(a, (offx,offy,offz), tol=1.e-6, tol2=1.e-5)
Example of use: stitched hat (array).
G.surfaceWalk: surface extrusion starting from a curve,
resulting into a surface mesh.
dj is the distribution of points in the extrusion direction starting from c,
niter the number of smoothing iterations.
check=1 means that the extrusion stops at the layer before cells intersect
alphaRef is the deviation angle wrt 180 degrees enabling to stop the extrusion before it crosses a sharp edge on the surface.
toldist is a tolerance below which points are considered matching.
Constraints can be set as 1D zones:
walk = G.surfaceWalk(surfaces, c, dj, constraints=[], niter=0, alphaRef=180., check=0, toldist=1.e-6)
Example of use: surface walk (array),
surface walk (pyTree).
G.collarMesh*: create a collar mesh at junction(s) between two surfaces s1 and s2 in
union or difference assembly, using a distribution along the surface dj
and a distribution in the normal direction to the wall dk. niterj and
niterk are the number of smoothing iterations for j and k directions.
ext is the extension of the collar mesh for difference assembly.
type is the assembly type, and can be 'union' or 'difference'.
alphaRef is the deviation angle wrt 180 degrees above which the walk is stopped.
contour is the starting contour to create the collar grids, constraints1 and constraints2 are 1D zones
defining the curves the collar grid must follow on surfaces s1 and s2 respectively.
toldist is the matching point tolerance. Parameter 'topology' can be 'overset' or 'extruded', only useful in case of difference.
Topology set to 'overset' results in two overlapping collar grids, whereas it results in a collar grid extruded from the
surface grid in the other case:
A = G.collarMesh(s1, s2, dj,dk, niterj=100, niterk=100, ext=5, alphaRef=30., type='union', contour=[], constraints1=[], constraints2=[], toldist=1.e-10, topology='overset')
Example of use: collar mesh generation (array),
collar mesh generation (pyTree).
Cartesian grid generators
G.gencartmb: simple Cartesian generator.
Create a set of Cartesian grids (B) around a list of body grids (A).
Those grids are patched with a ratio of 2. The user controls the number of levels, and the number of points for each level of grid.
h is the spatial step on the finest level.
Dfar is the maximal distance to the body.
nlvl is a list that provides the number of points per level (nlvl[0]: finest grid), except for the finest level:
B = G.gencartmb(A, h, Dfar, nlvl)
Example of use: multiblock Cartesian grid set (array),
multiblock Cartesian grid set (pyTree).
G.octree: create a QUAD quadtree mesh in 2D or an HEXA octree mesh in 3D starting
from a list of bodies and snears. Each parameter snear is the required spatial
step of the octree near the corresponding body; dfar is the extension
of the octree mesh in all the directions; balancing=1 means that the octree
is balanced, i.e. adjacent elements are at worst twice as big/small;
levelMax is the maximum number of levels required.
If ratio=2, then a classical octree mesh is built. If ratio=3, a 27-tree mesh is built,
in which case the spacing ratio is 3 (and not 2) between two adjacent elements:
b = G.octree(surfs, snears, dfar=5., balancing=0, levelMax=1000, ratio=2)
Example of use: generates a quadtree starting from a circle (array),
generates a quad tree starting from a circle (pyTree).
G.octree2Struct: convert an octree or a quadtree mesh into a set of Cartesian grids.
Parameter ext is the extension of Cartesian grids in all the directions;
vmin can be an integer defining the number of points in each Cartesian grid, or a list of integers, defining the
number of points per refinement level. In that case, the first element of the list of vmin defines the
finest level. Specifying all the levels is not mandatory.
If optimized=1, the ext value is reduced by -1 at overlap borders for the coarsest grid for minimum overlapping.
If merged=1, Cartesian grids are merged in order to reduce the number of created grids.
If AMR=1, a set of AMR zones are generated.
Parameter sizeMax can be used when merging is applied: in that case, the number of points per grid does not exceed sizeMax.
Warning: to obtain multigrid blocks, optimized must be set to 0:
b = G.octree2Struct(octree, vmin=15, ext=0, optimized=1, merged=1, AMR=0, sizeMax=1000000)
Example of use: builds a Cartesian mesh from an octree (array),
builds a Cartesian mesh from an octree (pyTree).
G.adaptOctree: adapt an unstructured octree with respect
to an indicator field located at element centers.
If 'indicator' is strictly positive for an element, then the element
must be refined as many times as required by the indicator number.
If 'indicator' is strictly negative, the element is coarsened if
possible as many times as required by the indicator number.
If 'indicator' is 0., the element remains unchanged.
balancing=1 means that the octree is balanced after adaptation.
If ratio=2, then a classical octree mesh is built. If ratio=3, a 27-tree
mesh is built, in which case the spacing ratio is 3 (and not 2) between
two adjacent elements. For array interface indicator is an array,
for pyTree version, indicator is the name of field
stored as a solution located at centers:
res = G.adaptOctree(octree, indicator, balancing=1, ratio=2)
Example of use: adapt an octree (array),
adapt an octree (pyTree).
G.expandLayer:
Expand the layer of given level for an octree unstructured mesh. If corners=1, expand also in corners directions:
res = G.expandLayer(octree, level=0, corners=0, balancing=0)
Example of use: expand the layer of finest level of an octree (array),
expand the layer of finest level of an octree (pyTree).
Operations on meshes
G.close: close a mesh defined by array a. Points that are distant of tol maximum to one another are merged:
b = G.close(a, tol=1.e-12) .or. B = G.close(A, tol=1.e-12)
Example of use: mesh closing (array),
mesh closing (pyTree).
G.selectInsideElts: select elements of a TRI-array, whose centers are inside the given
list of curves, defined by BAR-arrays:
b = G.selectInsideElts(a, curves)
Example of use: select inside elements (array).
G.map: map a distribution on a curve or on a structured surface:
b = G.map(a, distrib)
Map a i-array distribution in a direction (dir=1,2,3) in a surface or volume mesh:
b = G.map(a, distrib, dir)
Example of use: distribution mapping (array),
distribution mapping (pyTree).
G.mapSplit: split a i-array and map a distribution on the splitted i-array.
SplitCrit is the curvature radius triggering split:
b = G.mapSplit(a, distrib, splitCrit=100)
Example of use: splitting and distribution mapping(array),
splitting and distribution mapping (pyTree).
G.refine: refine a structured array. The original distribution is kept
but the number of points is multiplied by power. Dir is the direction
of refinement (1, 2, 3). If dir=0, refine in all directions:
b = G.refine(a, power, dir)
Example of use: refine an array (array),
refine a pyTree (pyTree).
G.mapCurvature: map a structured array following the curvature. N is the final
number of points. Dir is the direction
of remeshing (1, 2, 3):
b = G.mapCurvature(a, N, power, dir)
Example of use: map an array following curvature (array),
map a pyTree (pyTree).
G.densify: densify a i-array or a BAR-array with a new discretization step h.
Discretization points from the original array are kept:
b = G.densify(a, h)
Example of use: densify a 1D-array (array),
densify a 1D-pyTree (pyTree).
G.grow: grow a surface array of one layer. Vector is the node displacement.
For the array version, vector is defined by an array. For the PyTree
version, vector = ['v1','v2','v3'] where variables 'v1', 'v2', 'v3'
are defined as solutions in a, located at nodes:
b = G.grow(a, vector)
Example of use: grow a mesh of one layer (array),
grow a mesh of one layer (pyTree).
G.stack: stack two structured meshes (with the same nixnj) into a single mesh:
c = G.stack(a, b)
Example of use: stack two layers (array),
stack two layers (pyTree).
G.addNormalLayers: normal extrusion from a surface mesh. d is a 1D distribution providing the height of each layer.
If check=1, the extrusion stops before negative volume cells are created.
Niter specifies the number of iterations for normals smoothing:
b = G.addNormalLayers(a, d, check=0, niter=0)
Example of use: add normal layers to a surface,
normal extrusion (pyTree).
G.TTM: smooth a mesh using elliptic generator:
b = G.TTM(a, niter=100)
Example of use: TTM example (array),
TTM example (pyTree).
G.snapFront: snap a mesh to a surface S.
A front must be defined in a by a cellN field.
Points of this front are snapped to the surface. If optimized=0,
the exterior front cellN=1 is snapped, else if optimized=1 optimized front cellN=1 is snapped, else if optimized=2, front cellN=0 is snapped:
b = G.snapFront(a, S, optimized=1)
Example of use: snap front on a circular surface (array),
snap front on a circular surface (pyTree).
G.snapSharpEdges: snap a mesh to a surface S,
constrained by sharp edges and corners. if step != None,
sharp edges are refined with this step.
Sharp Edges are calculated depending on angle:
b = G.snapSharpEdges(a, S, step=None, angle=30.)
Example of use: snap a mesh on a polyline (array),
snap a mesh on a polyline (pyTree).
Operations on surface meshes
G.fittingPlaster: fit a surface structured patch to a curve a. BumpFactor controls
the curvature of the patch:
b = G.fittingPlaster(a, bumpFactor=0.)
Example of use: create a patch from a curve (array),
create a patch from a curve (pyTree).
G.gapfixer: fill a gap defined by a BAR contour a
drawn on a surface c. You can force the generated mesh
to pass through HardPoints (NODES). If refine=0, no inside
points are added:
b = G.gapfixer(a, c, hardPoints=None, refine=1)
Example of use: fill the gap in a circle (array),
fix a gap in a circle (pyTree).
G.gapsmanager: fill multiple gaps in a set of surface components A. Also, eliminate
overlap regions between components if any. Normals for all patches must be pointed outwards.
Set mode=0 for nodal mesh, 1 for center mesh, and 2 otherwise.
Set coplanar=1 if all components are lying on a same plane:
B = G.gapsmanager(A, mode=0, coplanar=0)
Example of use: fill multiple gaps (array),
fill multiple gaps (pyTree).
Information on generated meshes
G.check: check regularity, orthogonality for a mesh defined by an array:
G.check(a)
Example of use: mesh check (array),
mesh check (pyTree).
G.barycenter: return the barycenter of a, with optional weight:
bary = G.barycenter(a, weight=None) .or. bary = G.barycenter(A, weight=None)
Example of use: barycenter of array (array),
barycenter of zone (pyTree).
G.bbox: return the bounding box [xmin, ymin, zmin, xmax, ymax, zmax] of a or A:
bb = G.bbox(a) .or. bb = G.bbox(A)
Example of use: bounding box (array),
bounding box (pyTree).
G.bboxOfCells: return the bounding box of each cell of a. The bounding box field is
located at centers of cells:
bb = G.bboxOfCells(a) .or. BB = G.bboxOfCells(A)
Example of use: bounding box of all cells (array),
bounding box of all cells (pyTree).
G.BB: return the bounding box of a as an array or a zone. If method is 'AABB', then it computes
the Axis-Aligned Bounding-Box, if method is 'OBB' then it computes the Oriented Bounding-Box. The argument
weighting may be 0, and the OBB is computed using a Cloud-Point approach, or 1, and it is computed using a
Surface-Weighting approach. If weighting=1, then the provided array must be a surface composed of triangles:
b = G.BB(a, method='AABB', weighting=0) .or. B = G.BB(A, method='AABB', weighting=0)
Example of use: bounding box of array (array),
bounding box of zone (pyTree).
G.CEBBIntersection: test the Cartesian Elements Bounding Box (CEBB) intersection between a1
and a2. Tolerance is a float given by tol. Return 0 if no intersection, 1 otherwise:
intersect = G.CEBBIntersection(a1, a2, tol=1.e-10)
Example of use: CEBB intersection (array),
CEBB intersection (pyTree).
G.bboxIntersection: test if a1 and a2 intersects. Three options
are available: method='AABB' (intersection between two Axis-Aligned Bounding Boxes, by default); method='OBB' (intersection between
two Oriented Bounding Boxes, the most general case); method='AABBOBB' (intersection between an AABB -a1- and an OBB -a2-).;
If a1 and a2 are directly the corresponding bounding boxes, the user may switch isBB=True in order to avoid recalculating them.
Return 0 if no intersection, 1 otherwise:
intersect = G.bboxIntersection(a1, a2, tol=1.e-6, isBB=False, method='AABB')
Example of use: bounding box intersection (array),
bounding box intersection (pyTree).
G.checkPointInCEBB: test if a given point is in the CEBB of a:
inside = G.checkPointInCEBB(a, (x,y,z))
Example of use: point in CEBB check (array),
Point in CEBB check (pyTree).
G.getVolumeMap: return the volume field of an array. Volume is located at centers of cells:
b = G.getVolumeMap(a) .or. B = G.getVolumeMap(A)
Example of use: volume field (array),
volume field (pyTree).
G.getNormalMap: return the surface normals field of a surface array. It is located at
centers of cells:
b = G.getNormalMap(a) .or. B = G.getNormalMap(A)
Example of use: surface normals field (array),
surface normals field (pyTree).
G.getSmoothNormalMap: return the smoothed surface normals field of a surface array, located at nodes.
niter is the number of smoothing operations, and eps is a smoothing weight :
b = G.getSmoothNormalMap(a, niter=2,eps=0.4) .or. B = G.getNormalMap(A, niter=2,eps=0.4)
Example of use: smoothed surface normals field (array),
smoothed surface normals field (pyTree).
G.getOrthogonalityMap: return the orthogonality map of an array.
The orthogonality map corresponds to the maximum deviation of
all dihedral angles of an element.
The orthogonality map is expressed in degree and located at centers:
b = G.getOrthogonalityMap(a) .or. B = G.getOrthogonalityMap(A)
Example of use: orthogonality map (array),
orthogonality map (pyTree).
G.getRegularityMap: return the regularity map of an array.
The regularity map corresponds to the maximum deviation of the volume ratio
of an element and all its neigbouring cells.
The regularity map is located at centers:
b = G.getRegularityMap(a) .or. B = G.getRegularityMap(A)
Example of use: regularity map (array),
regularity map (pyTree).
G.getTriQualityMap: return the quality map of a TRI array.
The triangle quality is a value between 0. (degenerated triangle) and 1. (equilateral triangle).
The quality map is located at centers:
b = G.getTriQualityMap(a) .or. B = G.getTriQualityMap(A)
Example of use: quality map of a TRI array (array),
quality map of TRI zone (pyTree).
G.getCellPlanarity: return a measure of cell planarity for each cell. It is located at
centers of cells:
b = G.getCellPlanarity(a) .or. B = G.getCellPlanarity(A)
Example of use: cell planarity (array),
cell planarity (pyTree).
G.getCircumCircleMap: return the map of circum circle radius of any cell of a 'TRI' array:
b = G.getCircumCircleMap(a) .or. B = G.getCircumCircleMap(A)
Example of use: circum circle radius map (array),
circum circle radius map (pyTree).
G.getInCircleMap: return the map of inscribed circle radius of any cell of a 'TRI' array:
b = G.getInCircleMap(a) .or. B = G.getInCircleMap(A)
Example of use: inscribed circle radius map (array),
inscribed circle radius map (pyTree).
G.getEdgeRatio: return the ratio between the longest and the smallest edges of a cell:
r = G.getEdgeRatio(a) .or. B = G.getEdgeRatio(A)
Example of use: ratio between edges of a mesh (array),
ratio between edges of a mesh (pyTree).
G.getMaxEdgeLength: return the length of
the longer edge of each cell:
r = G.getMaxLength(a) .or. B = G.getMaxLength(A)
Example of use: max length between edges of a mesh (array),
max length between edges of a mesh (pyTree).
Operations on distributions
Distributions are Cartesian meshes that can be used to be
mapped on profiles to make curvilinear meshes for instance.
Distributions in 2D (x,y) distribution
represent the length and height of each cell.
Distributions in 3D (x,y,z) represents the lengths in the three topological
directions of each cell. Distributions can be modified by the enforce
functions.
The three following functions are also available in Y and Z directions: replace suffix 'X' by 'Y' or 'Z' in enforceX, enforcePlusX or enforceMoinsX functions.
G.enforceX: enforce a region around a line x=x0. The size of the cell around
the line is enforcedh. "supp" points are suppressed from the starting
distribution on the left and right side. "add" points are added on
the left and add points are added on the right.
Add exactely add points:
b = G.enforceX(a, x0, enforcedh, (supp,add))
Adjust add in order to have a monotonic distribution:
b = G.enforceX(a, x0, enforcedh, supp, add)
Example of use: point addition on both sides of line x (array),
point addition in 2 sides of line x (pyTree).
Example of use: monotonic enforce (array),
monotonic enforce (pyTree).
G.enforceMoinsX: same as before but with a one sided distribution (left).
This can be usefull to create a boundary layer distribution in an Euler mesh.
b = G.enforceMoinsX(a, enforcedh, (supp,add)) .or. b = G.enforceMoinsX(a, enforcedh, supp, add)
Example of use: point addition on the left side of line x (array),
point addition on the left side of line x (pyTree).
G.enforcePlusX: same as before but with a one sided distribution (right):
b = G.enforcePlusX(a, enforcedh, (supp,add)) .or. b = G.enforcePlusX(a, enforcedh, supp, add)
Example of use: point addition on the right side of line x (array),
point addition on the right side of line x (pyTree).
G.enforceLine: enforce a curvilinear line defined by the array line in a distribution defined by the array a:
b = G.enforceLine(a, line, enforcedh, (supp,add))
Example of use: line enforce (array),
line enforce (pyTree).
G.enforcePoint: enforce a point in the distribution. The index of enforced point is returned:
ind = G.enforcePoint(a, x0)
Example of use: point addition (array),
point addition (pyTree).
G.enforceCurvature: enforce the curvature of an i-curve in a distribution defined by a.
Power reflects the power of stretching:
b = G.enforceCurvature(a, curve, power=0.5)
Example of use: enforce curvature in a distribution (array),
enforce curvature in a distribution (pyTree).
G.addPointInDistribution: add a point in a distribution at index ind:
b = G.addPointInDistribution(a, ind)
Example of use: point addition in a distribution (array),
point addition in a distribution (pyTree).
More general examples of use
Return to main userguide