Mmg
Simplicial remeshers (mesh adaptation, isovalue discretization, lagrangian movement)
|
#include "libmmgcommon_private.h"
Go to the source code of this file.
Macros | |
#define | MMGS_ALPHAD 3.464101615137755 /* 6.0 / sqrt(3.0) */ |
#define | MMGS_LOPTL 1.4 |
#define | MMGS_LOPTS 0.71 |
#define | MMGS_LLONG 2.0 |
#define | MMGS_LSHRT 0.3 |
#define | MMGS_BADKAL 2.e-2 |
#define | MMGS_NULKAL 1.e-4 |
#define | MMGS_NPMAX 500000 |
#define | MMGS_NTMAX 1000000 |
#define | MMGS_XPMAX 500000 |
#define | MS_SIN(tag) ((tag & MG_CRN) || (tag & MG_REQ) || (tag & MG_NOM)) |
#define | MMGS_RETURN_AND_FREE(mesh, met, ls, val) |
#define | MMGS_POINT_REALLOC(mesh, sol, ip, wantedGap, law, o, tag) |
#define | MMGS_TRIA_REALLOC(mesh, jel, wantedGap, law) |
Functions | |
int | MMGS_Init_mesh_var (va_list argptr) |
int | MMGS_Free_all_var (va_list argptr) |
int | MMGS_Free_structures_var (va_list argptr) |
int | MMGS_Free_names_var (va_list argptr) |
int | MMGS_zaldy (MMG5_pMesh mesh) |
int | MMGS_assignEdge (MMG5_pMesh mesh) |
int | MMGS_analys_for_norver (MMG5_pMesh mesh) |
int | MMGS_analys (MMG5_pMesh mesh) |
int | MMGS_regver (MMG5_pMesh mesh) |
int | MMGS_inqua (MMG5_pMesh, MMG5_pSol) |
int | MMGS_outqua (MMG5_pMesh, MMG5_pSol) |
int | MMGS_hashTria (MMG5_pMesh) |
int | MMGS_setadj (MMG5_pMesh mesh) |
int | curvpo (MMG5_pMesh, MMG5_pSol) |
int | MMG5_mmgs1 (MMG5_pMesh, MMG5_pSol, MMG5_int *) |
int | MMGS_mmgs2 (MMG5_pMesh, MMG5_pSol, MMG5_pSol) |
int | MMGS_bdryUpdate (MMG5_pMesh mesh) |
int | boulechknm (MMG5_pMesh mesh, MMG5_int start, int ip, MMG5_int *list) |
int | bouletrid (MMG5_pMesh mesh, MMG5_int start, MMG5_int ip, int *il1, MMG5_int *l1, int *il2, MMG5_int *l2, MMG5_int *ip0, MMG5_int *ip1) |
MMG5_int | MMGS_newPt (MMG5_pMesh mesh, double c[3], double n[3]) |
void | MMGS_delPt (MMG5_pMesh mesh, MMG5_int ip) |
MMG5_int | MMGS_newElt (MMG5_pMesh mesh) |
int | MMGS_delElt (MMG5_pMesh mesh, MMG5_int iel) |
int | chkedg (MMG5_pMesh, MMG5_int) |
int | MMG5_mmgsBezierCP (MMG5_pMesh, MMG5_Tria *, MMG5_pBezier, int8_t ori) |
int | MMGS_bezierInt (MMG5_pBezier, double *, double *, double *, double *) |
int | MMGS_simbulgept (MMG5_pMesh mesh, MMG5_pSol met, MMG5_int k, int i, MMG5_int ip) |
int | MMGS_split1_sim (MMG5_pMesh mesh, MMG5_pSol met, MMG5_int k, int i, MMG5_int *vx) |
int | MMG5_split2_sim (MMG5_pMesh mesh, MMG5_pSol met, MMG5_int k, MMG5_int *vx) |
int | MMGS_split3_sim (MMG5_pMesh mesh, MMG5_pSol met, MMG5_int k, MMG5_int *vx) |
int | MMGS_split1 (MMG5_pMesh mesh, MMG5_pSol met, MMG5_int k, int i, MMG5_int *vx) |
int | MMGS_split2 (MMG5_pMesh mesh, MMG5_pSol met, MMG5_int k, MMG5_int *vx) |
int | MMGS_split3 (MMG5_pMesh mesh, MMG5_pSol met, MMG5_int k, MMG5_int *vx) |
int | split1b (MMG5_pMesh mesh, MMG5_int k, int8_t i, MMG5_int ip) |
int | chkcol (MMG5_pMesh mesh, MMG5_pSol met, MMG5_int k, int8_t i, MMG5_int *list, int8_t typchk, double(*MMGS_lenEdg)(MMG5_pMesh, MMG5_pSol, MMG5_int, MMG5_int, int8_t), double(*MMGS_caltri)(MMG5_pMesh, MMG5_pSol, MMG5_pTria)) |
int | colver (MMG5_pMesh mesh, MMG5_int *list, int ilist) |
int | colver3 (MMG5_pMesh mesh, MMG5_int *list) |
int | colver2 (MMG5_pMesh mesh, MMG5_int *ilist) |
int | swapar (MMG5_pMesh mesh, MMG5_int k, int i) |
int | chkswp (MMG5_pMesh mesh, MMG5_pSol met, MMG5_int k, int i, int8_t typchk, double(*MMGS_lenEdg)(MMG5_pMesh, MMG5_pSol, MMG5_int, MMG5_int, int8_t), double(*MMGS_caltri)(MMG5_pMesh, MMG5_pSol, MMG5_pTria)) |
int | swpedg (MMG5_pMesh mesh, MMG5_pSol met, MMG5_int *list, int ilist, int8_t typchk) |
int8_t | typelt (MMG5_pPoint p[3], int8_t *ia) |
int | litswp (MMG5_pMesh mesh, MMG5_int k, int8_t i, double kal) |
int | litcol (MMG5_pMesh mesh, MMG5_int k, int8_t i, double kal) |
int | MMG5_mmgsChkmsh (MMG5_pMesh, int, MMG5_int) |
int | paratmet (double c0[3], double n0[3], double m[6], double c1[3], double n1[3], double mt[6]) |
int | intregmet (MMG5_pMesh mesh, MMG5_pSol met, MMG5_int k, int8_t i, double s, double mr[6]) |
int | MMG5_intridmet (MMG5_pMesh, MMG5_pSol, MMG5_int, MMG5_int, double, double *, double *) |
int | setref (MMG5_pMesh, MMG5_int, MMG5_int, int) |
int | delref (MMG5_pMesh) |
int | chkmet (MMG5_pMesh, MMG5_pSol) |
int | chknor (MMG5_pMesh) |
size_t | MMG5_memSize (void) |
int | MMGS_memOption (MMG5_pMesh mesh) |
int | MMGS_setMeshSize_alloc (MMG5_pMesh mesh) |
int | MMG5_mmgsRenumbering (int, MMG5_pMesh, MMG5_pSol, MMG5_pSol, MMG5_int *) |
void | MMGS_keep_only1Subdomain (MMG5_pMesh mesh, MMG5_int nsd) |
MMG5_int | MMGS_indElt (MMG5_pMesh mesh, MMG5_int kel) |
MMG5_int | MMGS_indPt (MMG5_pMesh mesh, MMG5_int kp) |
double | caleltsig_ani (MMG5_pMesh mesh, MMG5_pSol met, MMG5_int iel) |
double | caleltsig_iso (MMG5_pMesh mesh, MMG5_pSol met, MMG5_int iel) |
int | MMGS_defsiz_iso (MMG5_pMesh mesh, MMG5_pSol met) |
int | MMGS_defsiz_ani (MMG5_pMesh mesh, MMG5_pSol met) |
int | MMGS_doSol_iso (MMG5_pMesh, MMG5_pSol) |
int | MMGS_doSol_ani (MMG5_pMesh, MMG5_pSol) |
int | MMGS_gradsiz_ani (MMG5_pMesh mesh, MMG5_pSol met) |
int | MMGS_gradsizreq_ani (MMG5_pMesh mesh, MMG5_pSol met) |
int | intmet_iso (MMG5_pMesh mesh, MMG5_pSol met, MMG5_int k, int8_t i, MMG5_int ip, double s) |
int | intmet_ani (MMG5_pMesh mesh, MMG5_pSol met, MMG5_int k, int8_t i, MMG5_int ip, double s) |
int | MMGS_intmet33_ani (MMG5_pMesh, MMG5_pSol, MMG5_int, int8_t, MMG5_int, double) |
int | MMGS_paramDisp (MMG5_pMesh, MMG5_int, int8_t, MMG5_int, MMG5_int, double, double[3]) |
int | MMGS_moveTowardPoint (MMG5_pMesh mesh, MMG5_pPoint p0, MMG5_pPoint p, double llold, double lam0, double lam1, double lam2, double nn1[3], double nn2[3], double to[3]) |
int | movridpt_iso (MMG5_pMesh mesh, MMG5_pSol met, MMG5_int *list, int ilist) |
int | movintpt_iso (MMG5_pMesh mesh, MMG5_pSol met, MMG5_int *list, int ilist) |
int | movridpt_ani (MMG5_pMesh mesh, MMG5_pSol met, MMG5_int *list, int ilist) |
int | movintpt_ani (MMG5_pMesh mesh, MMG5_pSol met, MMG5_int *list, int ilist) |
int | MMGS_surfballRotation (MMG5_pMesh, MMG5_pPoint, MMG5_int *, int, double r[3][3], double *, double n[3]) |
int | MMGS_prilen (MMG5_pMesh mesh, MMG5_pSol met, int) |
int | MMGS_set_metricAtPointsOnReqEdges (MMG5_pMesh, MMG5_pSol, int8_t) |
#define MMGS_ALPHAD 3.464101615137755 /* 6.0 / sqrt(3.0) */ |
Definition at line 34 of file libmmgs_private.h.
#define MMGS_BADKAL 2.e-2 |
Definition at line 41 of file libmmgs_private.h.
#define MMGS_LLONG 2.0 |
Definition at line 38 of file libmmgs_private.h.
#define MMGS_LOPTL 1.4 |
Definition at line 36 of file libmmgs_private.h.
#define MMGS_LOPTS 0.71 |
Definition at line 37 of file libmmgs_private.h.
#define MMGS_LSHRT 0.3 |
Definition at line 39 of file libmmgs_private.h.
#define MMGS_NPMAX 500000 |
Definition at line 44 of file libmmgs_private.h.
#define MMGS_NTMAX 1000000 |
Definition at line 45 of file libmmgs_private.h.
#define MMGS_NULKAL 1.e-4 |
Definition at line 42 of file libmmgs_private.h.
Reallocation of point table and sol table and creation of point ip with coordinates o and tag tag
Definition at line 66 of file libmmgs_private.h.
#define MMGS_RETURN_AND_FREE | ( | mesh, | |
met, | |||
ls, | |||
val | |||
) |
Free allocated pointers of mesh and sol structure and return value val
Definition at line 53 of file libmmgs_private.h.
#define MMGS_TRIA_REALLOC | ( | mesh, | |
jel, | |||
wantedGap, | |||
law | |||
) |
Reallocation of tria table and creation of tria jel
Definition at line 94 of file libmmgs_private.h.
#define MMGS_XPMAX 500000 |
Definition at line 46 of file libmmgs_private.h.
Definition at line 49 of file libmmgs_private.h.
int boulechknm | ( | MMG5_pMesh | mesh, |
MMG5_int | start, | ||
int | ip, | ||
MMG5_int * | list | ||
) |
mesh | pointer to the mesh structure. |
start | index of tetra to start to compute the ball. |
ip | index of point in tetra start for which we want to compute the ball. |
list | pointer to the computed ball of point. |
Find all triangles sharing ip, \(list[0] = start\) . Do not stop when crossing ridge. Check whether resulting configuration is manifold.
Definition at line 51 of file boulep_s.c.
int bouletrid | ( | MMG5_pMesh | mesh, |
MMG5_int | start, | ||
MMG5_int | ip, | ||
int * | il1, | ||
MMG5_int * | l1, | ||
int * | il2, | ||
MMG5_int * | l2, | ||
MMG5_int * | ip0, | ||
MMG5_int * | ip1 | ||
) |
mesh | pointer to the mesh structure. |
start | index of the starting triangle. |
ip | index of the looked ridge point. |
il1 | pointer to the first ball size. |
l1 | pointer to the first computed ball (associated to n1's side). |
il2 | pointer to the second ball size. |
l2 | pointer to the second computed ball (associated to n2's side). |
global | ip0 index of the first extremity of the ridge. |
global | ip1 index of the second extremity of the ridge. |
Computation of the two balls of a ridge point: the list l1 is associated to normal n1's side. ip0 and ip1 are the indices of the 2 ending point of the ridge. Both lists are returned enumerated in direct order.
Definition at line 201 of file boulep_s.c.
double caleltsig_ani | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met, | ||
MMG5_int | iel | ||
) |
mesh | pointer to the mesh structure. |
met | pointer to the metric structure. |
iel | element index |
Quality function identic to caltri_ani but puts a sign according to deviation to normal to vertices.
Definition at line 54 of file quality_s.c.
double caleltsig_iso | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met, | ||
MMG5_int | iel | ||
) |
int chkcol | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met, | ||
MMG5_int | k, | ||
int8_t | i, | ||
MMG5_int * | list, | ||
int8_t | typchk, | ||
double(*)(MMG5_pMesh, MMG5_pSol, MMG5_int, MMG5_int, int8_t) | MMGS_lenEdg, | ||
double(*)(MMG5_pMesh, MMG5_pSol, MMG5_pTria) | MMGS_caltri | ||
) |
mesh | pointer to the mesh |
met | pointer to the metric |
k | index of the element in wich we collapse |
i | index of the edge to collapse |
list | pointer to the ball of point |
typchk | type of check to perform |
MMGS_lenEdg | pointer to the suitable fct to compute edge lengths depending on presence of input metric, metric type (iso/aniso) and typchk value (i.e. stage of adaptation) |
MMGS_caltri | pointer to the suitable fct to compute tria quality depending on presence of input metric, metric type (iso/aniso) and typchk value (i.e. stage of adaptation) |
check if geometry preserved by collapsing edge i
Definition at line 59 of file colver_s.c.
int chkedg | ( | MMG5_pMesh | mesh, |
MMG5_int | iel | ||
) |
int chkmet | ( | MMG5_pMesh | , |
MMG5_pSol | |||
) |
int chknor | ( | MMG5_pMesh | mesh | ) |
mesh | pointer to the mesh |
Check normal vectors consistency
Definition at line 242 of file chkmsh_s.c.
int chkswp | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met, | ||
MMG5_int | k, | ||
int | i, | ||
int8_t | typchk, | ||
double(*)(MMG5_pMesh, MMG5_pSol, MMG5_int, MMG5_int, int8_t) | MMGS_lenEdg, | ||
double(*)(MMG5_pMesh, MMG5_pSol, MMG5_pTria) | MMGS_caltri | ||
) |
mesh | pointer to the mesh |
met | pointer to the metric |
k | index of the element in wich we perform the edge swap |
i | index of the edge to swap |
typchk | type of check to perform |
MMGS_lenEdg | pointer to the suitable fct to compute edge lengths depending on presence of input metric, metric type (iso/aniso) and typchk value (i.e. stage of adaptation) |
MMGS_caltri | pointer to the suitable fct to compute tria quality depending on presence of input metric, metric type (iso/aniso) and typchk value (i.e. stage of adaptation) |
Check whether edge i of triangle k should be swapped for geometric approximation purposes
Definition at line 58 of file swapar_s.c.
int colver | ( | MMG5_pMesh | mesh, |
MMG5_int * | list, | ||
int | ilist | ||
) |
Definition at line 280 of file colver_s.c.
int colver2 | ( | MMG5_pMesh | mesh, |
MMG5_int * | ilist | ||
) |
Definition at line 433 of file colver_s.c.
int colver3 | ( | MMG5_pMesh | mesh, |
MMG5_int * | list | ||
) |
mesh | pointer to the mesh structure. |
list | pointer to the ball of the point to collapse. |
Collapse edge \(list[0]\%3\) in tet \(list[0]/3\) ( \( ip->i1\) ) for a ball of the collapsed point of size 3: the collapsed point is removed.
Definition at line 370 of file colver_s.c.
int curvpo | ( | MMG5_pMesh | , |
MMG5_pSol | |||
) |
int delref | ( | MMG5_pMesh | mesh | ) |
Definition at line 39 of file gentools_s.c.
int intmet_ani | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met, | ||
MMG5_int | k, | ||
int8_t | i, | ||
MMG5_int | ip, | ||
double | s | ||
) |
mesh | pointer to the mesh structure. |
met | pointer to the metric structure. |
k | element index. |
i | local index of edge in k. |
ip | global index of the new point in which we want to compute the metric. |
s | interpolation parameter (between 0 and 1). |
Interpolation of anisotropic sizemap at parameter s along edge i of elt k for special storage of ridges metrics (after defsiz call).
Definition at line 104 of file intmet_s.c.
int intmet_iso | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met, | ||
MMG5_int | k, | ||
int8_t | i, | ||
MMG5_int | ip, | ||
double | s | ||
) |
mesh | pointer to the mesh structure. |
met | pointer to the metric structure. |
k | triangle in which we interpole the metrics. |
i | edge along which we interpole the metrics. |
ip | index of point in which we compute the interpolated metric. |
s | parameter at which we compute the interpolation. |
Linear interpolation of sizemap along edge i of tria k.
Definition at line 77 of file intmet_s.c.
int intregmet | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met, | ||
MMG5_int | k, | ||
int8_t | i, | ||
double | s, | ||
double | mr[6] | ||
) |
mesh | pointer to the mesh structure. |
met | pointer to the metric structure. |
k | element index. |
i | local index of edge in k. |
s | interpolation parameter. |
mr | computed metric. |
Metric interpolation on edge i in elt it at parameter \( 0 <= s0 <= 1 \) from p1 result is stored in mr. edge \( p_1-p_2 \) must not be a ridge.
Definition at line 57 of file intmet_s.c.
int litcol | ( | MMG5_pMesh | mesh, |
MMG5_int | k, | ||
int8_t | i, | ||
double | kal | ||
) |
int litswp | ( | MMG5_pMesh | mesh, |
MMG5_int | k, | ||
int8_t | i, | ||
double | kal | ||
) |
int MMG5_intridmet | ( | MMG5_pMesh | , |
MMG5_pSol | , | ||
MMG5_int | , | ||
MMG5_int | , | ||
double | , | ||
double * | , | ||
double * | |||
) |
size_t MMG5_memSize | ( | void | ) |
int MMG5_mmgs1 | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met, | ||
MMG5_int * | permNodGlob | ||
) |
mesh | pointer to the mesh structure. |
met | pointer to the metric structure. |
permNodGlob | if provided, strore the global permutation of nodes. |
Main adaptation routine.
Definition at line 1462 of file mmgs1.c.
int MMG5_mmgsBezierCP | ( | MMG5_pMesh | mesh, |
MMG5_Tria * | pt, | ||
MMG5_pBezier | pb, | ||
int8_t | ori | ||
) |
mesh | pointer to the mesh structure. |
pt | pointer to the triangle structure. |
pb | pointer to the computed Bezier structure. |
ori | triangle orientation (unused but here for compatibility with the MMG5_bezierCP interface). |
Compute Bezier control points on triangle pt (cf. [4]).
Definition at line 54 of file bezier_s.c.
int MMG5_mmgsChkmsh | ( | MMG5_pMesh | mesh, |
int | severe, | ||
MMG5_int | base | ||
) |
mesh | pointer to the mesh structure. |
severe | level of performed check |
base | unused argument. |
Check the mesh validity
Definition at line 48 of file chkmsh_s.c.
int MMG5_mmgsRenumbering | ( | int | boxVertNbr, |
MMG5_pMesh | mesh, | ||
MMG5_pSol | sol, | ||
MMG5_pSol | fields, | ||
MMG5_int * | permNodGlob | ||
) |
boxVertNbr | number of vertices by box. |
mesh | pointer to the mesh structure. |
sol | pointer to the solution structure |
fields | pointer to an array of solution fields |
permNodGlob | array to store the global permutation of nodes (if provided) |
Modifies the node indicies to prevent from cache missing.
Definition at line 80 of file librnbg_s.c.
int MMG5_split2_sim | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met, | ||
MMG5_int | k, | ||
MMG5_int * | vx | ||
) |
mesh | pointer to the mesh structure. |
met | pointer to the metric structure. |
k | index of element to split. |
vx | \(vx[i]\) is the index of the point to add on the edge i. |
Simulate the splitting of element k along the 2 edges i1 and i2. Check that the new triangles are not empty (otherwise we can create a 0 surface triangle).
Definition at line 383 of file split_s.c.
int MMGS_analys | ( | MMG5_pMesh | mesh | ) |
Definition at line 1094 of file analys_s.c.
int MMGS_analys_for_norver | ( | MMG5_pMesh | mesh | ) |
mesh | pointer to the mesh structure. |
Preprocessing stage: mesh analysis.
Definition at line 1042 of file analys_s.c.
int MMGS_assignEdge | ( | MMG5_pMesh | mesh | ) |
mesh | pointer to the mesh structure. |
Copy the properties (ref and tag) of the declared edges to the triangles, where they are assigned to the individual corners of the triangle. First a hash is created for rapid lookup of the edges. Then in a loop over all edges of all triangles, the hash is probed for each edge, and if it exists its properties are copied. Thus, declared edges that do not occur in any triangle will be silently ignored.
Definition at line 113 of file hash_s.c.
int MMGS_bdryUpdate | ( | MMG5_pMesh | mesh | ) |
int MMGS_bezierInt | ( | MMG5_pBezier | , |
double * | , | ||
double * | , | ||
double * | , | ||
double * | |||
) |
int MMGS_defsiz_ani | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met | ||
) |
mesh | pointer to the mesh structure. |
met | pointer to the metric stucture. |
Define size at points by intersecting the surfacic metric and the physical metric.
The output metric is:
Step 1: Set metric at points belonging to a required edge: compute the metric as the mean of the length of the required eges passing through the point
Step 2: Travel all the points (via triangles) in the mesh and set metric tensor
search for unintialized metric
Remark: as non manifold points are marked as CRN and REQ, we first try to call defmetsin that fails (because MMG5_boulet don't work for non-manifod points), then we pass here and we set hmax on non-manifold points
Definition at line 735 of file anisosiz_s.c.
int MMGS_defsiz_iso | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met | ||
) |
mesh | pointer to the mesh |
met | pointer to the metric |
Define isotropic size map at all vertices of the mesh, associated with geometric approx ; by convention, p0->h stores desired length at point p0
Step 1: Set metric at points belonging to a required edge: compute the metric as the mean of the length of the required eges passing through the point
Step 2: size at non required internal points
Step 3: Minimum size feature imposed by the hausdorff distance
Definition at line 142 of file isosiz_s.c.
int MMGS_delElt | ( | MMG5_pMesh | mesh, |
MMG5_int | iel | ||
) |
void MMGS_delPt | ( | MMG5_pMesh | mesh, |
MMG5_int | ip | ||
) |
int MMGS_doSol_ani | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met | ||
) |
mesh | pointer to the mesh |
met | pointer to the metric |
Compute the unit metric tensor at mesh vertices (from edges passing through the vertices).
Increment base marker to detect points already processed
Step 1: treat non-manifold points: travel triangle edges and add edge contribution to extremities (we have to do that because ball of non-manifold points can't be computed)
Corner point (no normal at vertex): if the corner defines an angle we can compute the 3D unit tensor (non singular), in the other cases (corner along a flat surface, for example because it is at the intersection of 3 specific edges or because it is provided by the user) we will project the edges onto the tangent plane and compute the 2D unit tensor (as for required and regular points)
Required point (no normal at vertex): projection of edges onto the tangent plane and computation of the 2D unit tensor
Ridge point (2 normals): normally we can compute the 3D unit tensor. If the angle is too flat, the computation fails and we try to compute the 2D unit tensor on the tangent plane of the point.
Regular or reference point: projection of edges onto the tangent plane and computation of the 2D unit metric tensor
Step 3: check computed metric
Step 4: computation of hmin/hmax (if needed) and metric truncature
Definition at line 1334 of file libmmgs_tools.c.
int MMGS_doSol_iso | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met | ||
) |
mesh | pointer to the mesh |
met | pointer to the metric |
Compute isotropic size map according to the mean of the length of the edges passing through a point.
Definition at line 859 of file libmmgs_tools.c.
int MMGS_Free_all_var | ( | va_list | argptr | ) |
argptr | list of the mmg structures that must be deallocated. Each structure must follow one of the MMG5_ARG preprocessor variable that allow to identify it. |
argptr contains at least a pointer to a MMG5_pMesh structure (that will contain the mesh and identified by the MMG5_ARG_ppMesh keyword).
To call the MMGS_mmgslib function, you must also provide a pointer to a MMG5_pSol structure (that will contain the ouput metric (and the input one, if provided) and identified by the MMG5_ARG_ppMet keyword).
To call the MMGS_mmgsls function, you must also provide a pointer toward a MMG5_pSol structure (that will contain the level-set function and identified by the MMG5_ARG_ppLs keyword).
Internal function for deallocations before return (taking a va_list as argument).
Definition at line 225 of file variadic_s.c.
int MMGS_Free_names_var | ( | va_list | argptr | ) |
argptr | list of the mmg structures for whose we want to deallocate the name. Each structure must follow one of the MMG5_ARG preprocessor variable that allow to identify it. |
argptr contains at least a pointer to a MMG5_pMesh structure (that will contain the mesh and identified by the MMG5_ARG_ppMesh keyword).
To call the MMGS_mmgslib function, you must also provide a pointer to a MMG5_pSol structure (that will contain the ouput metric (and the input one, if provided) and identified by the MMG5_ARG_ppMet keyword).
To call the MMGS_mmgsls function, you must also provide a pointer toward a MMG5_pSol structure (that will contain the level-set function and identified by the MMG5_ARG_ppLs keyword).
Internal function for name deallocations before return (taking a va_list as argument).
Definition at line 426 of file variadic_s.c.
int MMGS_Free_structures_var | ( | va_list | argptr | ) |
argptr | list of the mmg structures that must be deallocated. Each structure must follow one of the MMG5_ARG* preprocessor variable that allow to identify it. |
argptr contains at least a pointer to a MMG5_pMesh structure (that will contain the mesh and identified by the MMG5_ARG_ppMesh keyword).
To call the MMGS_mmgslib function, you must also provide a pointer to a MMG5_pSol structure (that will contain the ouput metric (and the input one, if provided) and identified by the MMG5_ARG_ppMet keyword).
To call the MMGS_mmgsls function, you must also provide a pointer toward a MMG5_pSol structure (that will contain the level-set function and identified by the MMG5_ARG_ppLs keyword).
Internal function for structures deallocations before return (taking a va_list as argument).
Definition at line 324 of file variadic_s.c.
int MMGS_gradsiz_ani | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met | ||
) |
mesh | pointer to the mesh structure. |
met | pointer to the metric structure. |
Enforces mesh gradation by truncating metric field.
Definition at line 835 of file anisosiz_s.c.
int MMGS_gradsizreq_ani | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met | ||
) |
int MMGS_hashTria | ( | MMG5_pMesh | mesh | ) |
MMG5_int MMGS_indElt | ( | MMG5_pMesh | mesh, |
MMG5_int | kel | ||
) |
find the element number in packed numerotation
Definition at line 123 of file gentools_s.c.
MMG5_int MMGS_indPt | ( | MMG5_pMesh | mesh, |
MMG5_int | kp | ||
) |
find the point number in packed numerotation
Definition at line 139 of file gentools_s.c.
int MMGS_Init_mesh_var | ( | va_list | argptr | ) |
argptr | list of the mmg structures that must be initialized. Each structure must follow one of the MMG5_ARG* preprocessor variable that allow to identify it. |
argptr contains at least a pointer to a MMG5_pMesh structure (that will contain the mesh and identified by the MMG5_ARG_ppMesh keyword).
To call the MMGS_mmgslib function, you must also provide a pointer to a MMG5_pSol structure (that will contain the ouput metric (and the input one, if provided) and identified by the MMG5_ARG_ppMet keyword).
To call the MMGS_mmgsls function, you must also provide a pointer toward a MMG5_pSol structure (that will contain the level-set function and identified by the MMG5_ARG_ppLs keyword).
Internal function for structure allocations (taking a va_list argument).
Definition at line 149 of file variadic_s.c.
int MMGS_inqua | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met | ||
) |
mesh | pointer to the mesh structure. |
met | pointer to the metric structure. |
Print histogram of mesh qualities for classical storage of ridges metrics (so before the the MMG5_defsiz function call).
Definition at line 399 of file quality_s.c.
int MMGS_intmet33_ani | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met, | ||
MMG5_int | k, | ||
int8_t | i, | ||
MMG5_int | ip, | ||
double | s | ||
) |
mesh | pointer to the mesh structure. |
met | pointer to the metric structure. |
k | element index. |
i | local index of edge in k. |
ip | global index of the new point in which we want to compute the metric. |
s | interpolation parameter (between 0 and 1). |
Interpolation of anisotropic sizemap at parameter s along edge i of elt k for classic storage of ridges metrics (before defsiz call).
Definition at line 144 of file intmet_s.c.
void MMGS_keep_only1Subdomain | ( | MMG5_pMesh | mesh, |
MMG5_int | nsd | ||
) |
mesh | pointer to the mesh structure. |
nsd | index of subdomain to keep. |
Keep only subdomain of index nsd and remove other subdomains.
Definition at line 161 of file gentools_s.c.
int MMGS_memOption | ( | MMG5_pMesh | mesh | ) |
int MMGS_mmgs2 | ( | MMG5_pMesh | mesh, |
MMG5_pSol | sol, | ||
MMG5_pSol | met | ||
) |
int MMGS_moveTowardPoint | ( | MMG5_pMesh | mesh, |
MMG5_pPoint | p0, | ||
MMG5_pPoint | p, | ||
double | llold, | ||
double | lam0, | ||
double | lam1, | ||
double | lam2, | ||
double | nn1[3], | ||
double | nn2[3], | ||
double | to[3] | ||
) |
mesh | pointer to the mesh |
p0 | point to move. |
p | neighbouring point toward which we try to move. |
llold | init length of edge p0-p |
lam0 | first bezier basis function (order 2) |
lam1 | second bezier basis function (order 2) |
lam2 | third bezier basis function (order 2) |
nn1 | normal at point p0 after relocation |
nn2 | normal at point p0 after relocation |
to | tangent along edge at point p0 after relocation |
Update normals and tangent at ref or ridge point p0 after relocation at coordinates o.
Definition at line 555 of file movpt_s.c.
MMG5_int MMGS_newElt | ( | MMG5_pMesh | mesh | ) |
MMG5_int MMGS_newPt | ( | MMG5_pMesh | mesh, |
double | c[3], | ||
double | n[3] | ||
) |
int MMGS_outqua | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met | ||
) |
mesh | pointer to the mesh structure. |
met | pointer to the metric structure. |
Print histogram of mesh qualities for special storage of ridges metrics (after the MMG5_defsiz function call).
Definition at line 467 of file quality_s.c.
int MMGS_paramDisp | ( | MMG5_pMesh | mesh, |
MMG5_int | it, | ||
int8_t | isrid, | ||
MMG5_int | ip0, | ||
MMG5_int | ip, | ||
double | step, | ||
double | o[3] | ||
) |
mesh | pointer to the mesh |
it | triangle to which belongs the edge along which we move |
isrid | 1 if the edge is a ridge |
ip0 | edge point that we want to move |
ip | edge point connected by the ref/ridge edge to p0 |
step | displacement factor along the ref/ridge edge |
o | coordinates of point after relocation |
Infer arc length of displacement along ref or ridge edge, parameterized over edges.
Definition at line 375 of file movpt_s.c.
int MMGS_prilen | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met, | ||
int | metRidTyp | ||
) |
mesh | pointer to the mesh structure. |
met | pointer to the metric structure. |
metRidTyp | Type of storage of ridges metrics: 0 for classic storage, 1 for special storage. |
Compute sizes of edges of the mesh, and displays histo.
Definition at line 283 of file quality_s.c.
int MMGS_regver | ( | MMG5_pMesh | mesh | ) |
mesh | pointer to a MMG5 mesh structure. |
Regularization procedure for vertices coordinates, dual Laplacian for a surface mesh.
Definition at line 817 of file analys_s.c.
int MMGS_set_metricAtPointsOnReqEdges | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met, | ||
int8_t | ismet | ||
) |
mesh | pointer to the mesh |
met | pointer to the metric |
ismet | 1 if user provided metric |
Compute the metric at points on trequired adges as the mean of the lengths of the required eges to which belongs the point. The processeed points are marked with flag 3.
Definition at line 90 of file isosiz_s.c.
int MMGS_setadj | ( | MMG5_pMesh | mesh | ) |
mesh | pointer to the mesh |
topology: set adjacent, detect Moebius, flip faces, count connected comp.
Definition at line 47 of file analys_s.c.
int MMGS_setMeshSize_alloc | ( | MMG5_pMesh | mesh | ) |
int MMGS_simbulgept | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met, | ||
MMG5_int | k, | ||
int | i, | ||
MMG5_int | ip | ||
) |
mesh | pointer to the mesh structure. |
met | pointer to the metric structure. |
k | index of the starting triangle. |
i | local index of the edge to split in k. |
ip | index of the point that we try to create. |
Simulate the creation of the point ip, to be inserted at an edge. Check that the new triangles are not empty (otherwise we can create a 0 surface triangle).
Definition at line 162 of file split_s.c.
int MMGS_split1 | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met, | ||
MMG5_int | k, | ||
int | i, | ||
MMG5_int * | vx | ||
) |
mesh | pointer to the mesh structure. |
met | pointer to the metric structure. |
k | index of element to split. |
i | index of edge to split. |
vx | \(vx[i]\) is the index of the point to add on the edge i. |
Split element k along edge i.
Definition at line 109 of file split_s.c.
int MMGS_split1_sim | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met, | ||
MMG5_int | k, | ||
int | i, | ||
MMG5_int * | vx | ||
) |
mesh | pointer to the mesh structure. |
met | pointer to the metric structure. |
k | index of element to split. |
i | index of edge to split. |
vx | \(vx[i]\) is the index of the point to add on the edge i. |
k | index of element to split. |
Simulate the splitting of element k along edge i. Check that the new triangles are not empty (otherwise we can create a 0 surface triangle).
Definition at line 52 of file split_s.c.
int MMGS_split2 | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met, | ||
MMG5_int | k, | ||
MMG5_int * | vx | ||
) |
mesh | pointer to the mesh structure. |
met | pointer to the metric structure. |
k | index of element to split. |
vx | \(vx[i]\) is the index of the point to add on the edge i. |
Split element k along the 2 edges i1 and i2.
Definition at line 459 of file split_s.c.
int MMGS_split3 | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met, | ||
MMG5_int | k, | ||
MMG5_int * | vx | ||
) |
mesh | pointer to the mesh structure. |
met | pointer to the metric structure. |
k | index of element to split. |
vx | \(vx[i]\) is the index of the point to add on the edge i. |
Split element k along the 3 edges
Definition at line 620 of file split_s.c.
int MMGS_split3_sim | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met, | ||
MMG5_int | k, | ||
MMG5_int * | vx | ||
) |
mesh | pointer to the mesh structure. |
met | pointer to the metric structure. |
k | index of element to split. |
vx | \(vx[i]\) is the index of the point to add on the edge i. |
Simulate the splitting of element k along the 3 edges. Check that the new triangles are not empty (otherwise we can create a 0 surface triangle).
Definition at line 532 of file split_s.c.
int MMGS_surfballRotation | ( | MMG5_pMesh | mesh, |
MMG5_pPoint | p0, | ||
MMG5_int * | list, | ||
int | ilist, | ||
double | r[3][3], | ||
double * | lispoi, | ||
double | n[3] | ||
) |
mesh | pointer to the mesh structure. |
p0 | starting point |
list | ball of p0 |
ilist | number of tria in the ball of p0 |
r | rotation that send the normal at p0 onto the z vector |
lipoint | rotated ball of point p0 |
n | normal at point p0 |
Compute the rotation matrix that sends the tangent plane at p0 onto z=0 and apply this rotation to the ball of p0.
Definition at line 529 of file anisosiz_s.c.
int MMGS_zaldy | ( | MMG5_pMesh | mesh | ) |
int movintpt_ani | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met, | ||
MMG5_int * | list, | ||
int | ilist | ||
) |
mesh | pointer to the mesh structure. |
met | pointer to the metric structure. |
list | ball of point. |
ilist | size of the point ball. |
Compute movement of an internal point whose ball is passed.
Step 2 : Compute gradient towards optimal position = centre of mass of the ball, projected to tangent plane
Step 3 : locate new point in the ball, and compute its barycentric coordinates
Step 4 : come back to original problem, and compute patch in triangle iel
Definition at line 49 of file anisomovpt_s.c.
int movintpt_iso | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met, | ||
MMG5_int * | list, | ||
int | ilist | ||
) |
mesh | pointer to the mesh structure. |
met | pointer to the metric structure. |
list | pointer to the ball of the point. |
ilist | size of the ball. |
Move internal point whose volumic is passed.
Definition at line 52 of file movpt_s.c.
int movridpt_ani | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met, | ||
MMG5_int * | list, | ||
int | ilist | ||
) |
Definition at line 252 of file anisomovpt_s.c.
int movridpt_iso | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met, | ||
MMG5_int * | list, | ||
int | ilist | ||
) |
int paratmet | ( | double | c0[3], |
double | n0[3], | ||
double | m[6], | ||
double | c1[3], | ||
double | n1[3], | ||
double | mt[6] | ||
) |
int setref | ( | MMG5_pMesh | mesh, |
MMG5_int | start, | ||
MMG5_int | ref, | ||
int | putreq | ||
) |
mesh | pointer to the mesh |
start | index of the tetra from which we start |
ref | reference to set |
putreq | 1 if boundary edges must be set to required |
Start from triangle start, and pile up triangles by adjacency, till a GEO or REF curve is met ; pass all references of travelled faces to ref ; putreq = 1 if boundary edges met must be set to MG_REQ, 0 otherwise.
Definition at line 64 of file gentools_s.c.
int split1b | ( | MMG5_pMesh | mesh, |
MMG5_int | k, | ||
int8_t | i, | ||
MMG5_int | ip | ||
) |
mesh | pointer to the mesh structure. |
k | index of element to split. |
i | index of edge to split. |
ip | index of the new point. |
Split element k along edge i, inserting point ip and updating the adjacency relations.
Definition at line 283 of file split_s.c.
int swapar | ( | MMG5_pMesh | mesh, |
MMG5_int | k, | ||
int | i | ||
) |
mesh | poiner toward the mesh structure. |
k | elt index. |
i | index of the elt edge to swap. |
Definition at line 318 of file swapar_s.c.
int swpedg | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met, | ||
MMG5_int * | list, | ||
int | ilist, | ||
int8_t | typchk | ||
) |
int8_t typelt | ( | MMG5_pPoint | p[3], |
int8_t * | ia | ||
) |
Definition at line 525 of file quality_s.c.