mmgs
mmgs.h File Reference
#include "libmmgcommon.h"
#include "libmmgs.h"
Include dependency graph for mmgs.h:
This graph shows which files directly or indirectly include this file:

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_LMAX   1024
 
#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 (MMG5_pMesh mesh)
 
int MMGS_inqua (MMG5_pMesh, MMG5_pSol)
 
int MMGS_outqua (MMG5_pMesh, MMG5_pSol)
 
int MMGS_hashTria (MMG5_pMesh)
 
int curvpo (MMG5_pMesh, MMG5_pSol)
 
int MMG5_mmgs1 (MMG5_pMesh, MMG5_pSol, int *)
 
int MMGS_mmgs2 (MMG5_pMesh, MMG5_pSol, MMG5_pSol)
 
int MMGS_bdryUpdate (MMG5_pMesh mesh)
 
int boulet (MMG5_pMesh mesh, int start, int ip, int *list)
 
int boulechknm (MMG5_pMesh mesh, int start, int ip, int *list)
 
int bouletrid (MMG5_pMesh mesh, int start, int ip, int *il1, int *l1, int *il2, int *l2, int *ip0, int *ip1)
 
int MMGS_newPt (MMG5_pMesh mesh, double c[3], double n[3])
 
void MMGS_delPt (MMG5_pMesh mesh, int ip)
 
int MMGS_newElt (MMG5_pMesh mesh)
 
int MMGS_delElt (MMG5_pMesh mesh, int iel)
 
int chkedg (MMG5_pMesh, 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, int k, int i, int ip)
 
int MMGS_split1_sim (MMG5_pMesh mesh, MMG5_pSol met, int k, int i, int *vx)
 
int MMG5_split2_sim (MMG5_pMesh mesh, MMG5_pSol met, int k, int *vx)
 
int MMGS_split3_sim (MMG5_pMesh mesh, MMG5_pSol met, int k, int *vx)
 
int MMGS_split1 (MMG5_pMesh mesh, MMG5_pSol met, int k, int i, int *vx)
 
int MMGS_split2 (MMG5_pMesh mesh, MMG5_pSol met, int k, int *vx)
 
int MMGS_split3 (MMG5_pMesh mesh, MMG5_pSol met, int k, int *vx)
 
int split1b (MMG5_pMesh mesh, int k, int8_t i, int ip)
 
int chkcol (MMG5_pMesh mesh, MMG5_pSol met, int k, int8_t i, int *list, int8_t typchk)
 
int colver (MMG5_pMesh mesh, int *list, int ilist)
 
int colver3 (MMG5_pMesh mesh, int *list)
 
int colver2 (MMG5_pMesh mesh, int *ilist)
 
int swapar (MMG5_pMesh mesh, int k, int i)
 
int chkswp (MMG5_pMesh mesh, MMG5_pSol met, int k, int i, int8_t typchk)
 
int swpedg (MMG5_pMesh mesh, MMG5_pSol met, int *list, int ilist, int8_t typchk)
 
int8_t typelt (MMG5_pPoint p[3], int8_t *ia)
 
int litswp (MMG5_pMesh mesh, int k, int8_t i, double kal)
 
int litcol (MMG5_pMesh mesh, int k, int8_t i, double kal)
 
int MMG5_mmgsChkmsh (MMG5_pMesh, int, 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, int k, int8_t i, double s, double mr[6])
 
int MMG5_intridmet (MMG5_pMesh, MMG5_pSol, int, int, double, double *, double *)
 
int setref (MMG5_pMesh, int, 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)
 
void MMGS_keep_only1Subdomain (MMG5_pMesh mesh, int nsd)
 
int MMGS_indElt (MMG5_pMesh mesh, int kel)
 
int MMGS_indPt (MMG5_pMesh mesh, int kp)
 
void MMG5_Init_parameters (MMG5_pMesh mesh)
 
double caleltsig_ani (MMG5_pMesh mesh, MMG5_pSol met, int iel)
 
double caleltsig_iso (MMG5_pMesh mesh, MMG5_pSol met, int iel)
 
int MMGS_defsiz_iso (MMG5_pMesh mesh, MMG5_pSol met)
 
int MMGS_defsiz_ani (MMG5_pMesh mesh, MMG5_pSol met)
 
void MMG5_defaultValues (MMG5_pMesh)
 
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, int k, int8_t i, int ip, double s)
 
int intmet_ani (MMG5_pMesh mesh, MMG5_pSol met, int k, int8_t i, int ip, double s)
 
int MMGS_intmet33_ani (MMG5_pMesh, MMG5_pSol, int, int8_t, int, double)
 
int MMGS_paramDisp (MMG5_pMesh mesh, int it1, int it2, double l1old, double l2old, int8_t isrid1, int8_t isrid2, int ip0, int ip1, int ip2, double step, double o[3], int8_t *isrid)
 
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, int *list, int ilist)
 
int movintpt_iso (MMG5_pMesh mesh, MMG5_pSol met, int *list, int ilist)
 
int movridpt_ani (MMG5_pMesh mesh, MMG5_pSol met, int *list, int ilist)
 
int movintpt_ani (MMG5_pMesh mesh, MMG5_pSol met, int *list, int ilist)
 
int MMGS_surfballRotation (MMG5_pMesh, MMG5_pPoint, int *, int, double r[3][3], double *)
 
int MMGS_prilen (MMG5_pMesh mesh, MMG5_pSol met, int)
 
int MMGS_set_metricAtPointsOnReqEdges (MMG5_pMesh, MMG5_pSol, int8_t)
 
static void MMGS_Set_commonFunc (void)
 

Variables

double(* MMG5_calelt )(MMG5_pMesh mesh, MMG5_pSol met, MMG5_pTria ptt)
 
int(* MMGS_defsiz )(MMG5_pMesh mesh, MMG5_pSol met)
 
int(* MMGS_gradsiz )(MMG5_pMesh mesh, MMG5_pSol met)
 
int(* MMGS_gradsizreq )(MMG5_pMesh mesh, MMG5_pSol met)
 
int(* intmet )(MMG5_pMesh mesh, MMG5_pSol met, int k, int8_t i, int ip, double s)
 
int(* movridpt )(MMG5_pMesh mesh, MMG5_pSol met, int *list, int ilist)
 
int(* movintpt )(MMG5_pMesh mesh, MMG5_pSol met, int *list, int ilist)
 

Macro Definition Documentation

◆ MMGS_ALPHAD

#define MMGS_ALPHAD   3.464101615137755 /* 6.0 / sqrt(3.0) */

◆ MMGS_BADKAL

#define MMGS_BADKAL   2.e-2

◆ MMGS_LLONG

#define MMGS_LLONG   2.0

◆ MMGS_LMAX

#define MMGS_LMAX   1024

◆ MMGS_LOPTL

#define MMGS_LOPTL   1.4

◆ MMGS_LOPTS

#define MMGS_LOPTS   0.71

◆ MMGS_LSHRT

#define MMGS_LSHRT   0.3

◆ MMGS_NPMAX

#define MMGS_NPMAX   500000

◆ MMGS_NTMAX

#define MMGS_NTMAX   1000000

◆ MMGS_NULKAL

#define MMGS_NULKAL   1.e-4

◆ MMGS_POINT_REALLOC

#define MMGS_POINT_REALLOC (   mesh,
  sol,
  ip,
  wantedGap,
  law,
  o,
  tag 
)
Value:
do \
{ \
int klink; \
assert ( mesh && mesh->point ); \
MMG5_TAB_RECALLOC(mesh,mesh->point,mesh->npmax,wantedGap,MMG5_Point, \
"larger point table",law); \
\
mesh->npnil = mesh->np+1; \
for (klink=mesh->npnil; klink<mesh->npmax-1; klink++) \
mesh->point[klink].tmp = klink+1; \
\
/* solution */ \
if ( sol->m ) { \
MMG5_ADD_MEM(mesh,(sol->size*(mesh->npmax-sol->npmax))*sizeof(double), \
"larger solution",law); \
MMG5_SAFE_REALLOC(sol->m,sol->size*(sol->npmax+1), \
sol->size*(mesh->npmax+1),double, \
"larger solution",law); \
} \
sol->npmax = mesh->npmax; \
\
/* We try again to add the point */ \
ip = MMGS_newPt(mesh,o,tag); \
if ( !ip ) {law;} \
}while(0)
MMG5_pMesh * mesh
Definition: API_functionsf_s.c:63
int MMGS_newPt(MMG5_pMesh mesh, double c[3], double n[3])
Definition: zaldy_s.c:39
MMG5_pPoint point
Definition: libmmgtypes.h:589
int npnil
Definition: libmmgtypes.h:569
int np
Definition: libmmgtypes.h:559
int npmax
Definition: libmmgtypes.h:559
Structure to store points of a MMG mesh.
Definition: libmmgtypes.h:241

Reallocation of point table and sol table and creation of point ip with coordinates o and tag tag

◆ MMGS_RETURN_AND_FREE

#define MMGS_RETURN_AND_FREE (   mesh,
  met,
  ls,
  val 
)
Value:
do \
{ \
MMG5_ARG_end) ) { \
return MMG5_LOWFAILURE; \
} \
return val; \
}while(0)
int MMGS_Free_all(const int starter,...)
Definition: API_functions_s.c:1577
#define MMG5_ARG_ppMesh
Definition: libmmgtypes.h:85
#define MMG5_ARG_end
Definition: libmmgtypes.h:162
#define MMG5_ARG_ppLs
Definition: libmmgtypes.h:95
#define MMG5_LOWFAILURE
Definition: libmmgtypes.h:51
#define MMG5_ARG_start
Definition: libmmgtypes.h:76
#define MMG5_ARG_ppMet
Definition: libmmgtypes.h:105

Free allocated pointers of mesh and sol structure and return value val

◆ MMGS_TRIA_REALLOC

#define MMGS_TRIA_REALLOC (   mesh,
  jel,
  wantedGap,
  law 
)
Value:
do \
{ \
int klink,oldSiz; \
\
oldSiz = mesh->ntmax; \
MMG5_TAB_RECALLOC(mesh,mesh->tria,mesh->ntmax,wantedGap,MMG5_Tria, \
"larger tria table",law); \
\
mesh->nenil = mesh->nt+1; \
for (klink=mesh->nenil; klink<mesh->ntmax-1; klink++) \
mesh->tria[klink].v[2] = klink+1; \
if ( mesh->adja ) { \
/* adja table */ \
MMG5_ADD_MEM(mesh,3*(mesh->ntmax-oldSiz)*sizeof(int), \
"larger adja table",law); \
MMG5_SAFE_RECALLOC(mesh->adja,3*oldSiz+5,3*mesh->ntmax+5,int \
,"larger adja table",law); \
} \
\
/* We try again to add the point */ \
jel = MMGS_newElt(mesh); \
if ( !jel ) {law;} \
}while(0)
if(!ier) exit(EXIT_FAILURE)
int MMGS_newElt(MMG5_pMesh mesh)
Definition: zaldy_s.c:71
int ntmax
Definition: libmmgtypes.h:559
int nenil
Definition: libmmgtypes.h:570
int nt
Definition: libmmgtypes.h:559
MMG5_pTria tria
Definition: libmmgtypes.h:595
int * adja
Definition: libmmgtypes.h:572
Definition: libmmgtypes.h:301

Reallocation of tria table and creation of tria jel

◆ MMGS_XPMAX

#define MMGS_XPMAX   500000

◆ MS_SIN

#define MS_SIN (   tag)    ((tag & MG_CRN) || (tag & MG_REQ) || (tag & MG_NOM))

Function Documentation

◆ boulechknm()

int boulechknm ( MMG5_pMesh  mesh,
int  start,
int  ip,
int *  list 
)
Parameters
meshpointer toward the mesh structure.
startindex of tetra to start to compute the ball.
ipindex of point in tetra start for which we want to compute the ball.
listpointer toward 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.

Here is the caller graph for this function:

◆ boulet()

int boulet ( MMG5_pMesh  mesh,
int  start,
int  ip,
int *  list 
)
Parameters
meshpointer toward the mesh structure.
startindex of triangle to start.
ipindex of point for wich we compute the ball.
listpointer toward the computed ball of ip.
Returns
the size of the computed ball or 0 if fail.

Find all triangles sharing ip, $list[0] =$ start do not stop when crossing ridge.

Here is the caller graph for this function:

◆ bouletrid()

int bouletrid ( MMG5_pMesh  mesh,
int  start,
int  ip,
int *  il1,
int *  l1,
int *  il2,
int *  l2,
int *  ip0,
int *  ip1 
)
Parameters
meshpointer toward the mesh structure.
startindex of the starting triangle.
ipindex of the looked ridge point.
il1pointer toward the first ball size.
l1pointer toward the first computed ball (associated to n1's side).
il2pointer toward the second ball size.
l2pointer toward the second computed ball (associated to n2's side).
ip0index of the first extremity of the ridge.
ip1index of the second extremity of the ridge.
Returns
0 if fail, 1 otherwise.

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.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ caleltsig_ani()

double caleltsig_ani ( MMG5_pMesh  mesh,
MMG5_pSol  met,
int  iel 
)
Parameters
meshpointer toward the mesh structure.
metpointer toward the metric structure.
ielelement index
Returns
0 if fail, -1 if orientation is reversed with regards to orientation of vertices, the computed quality otherwise.

Quality function identic to caltri_ani but puts a sign according to deviation to normal to vertices.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ caleltsig_iso()

double caleltsig_iso ( MMG5_pMesh  mesh,
MMG5_pSol  met,
int  iel 
)
Here is the caller graph for this function:

◆ chkcol()

int chkcol ( MMG5_pMesh  mesh,
MMG5_pSol  met,
int  k,
int8_t  i,
int *  list,
int8_t  typchk 
)
Parameters
meshpointer toward the mesh
metpointer toward the metric
kindex of the element in wich we collapse
iindex of the edge to collapse
listpointer toward the ball of point
typchktype of check to perform
Returns
0 if we can't move of if we fail, 1 if success

check if geometry preserved by collapsing edge i

Here is the call graph for this function:
Here is the caller graph for this function:

◆ chkedg()

int chkedg ( MMG5_pMesh  mesh,
int  iel 
)
Here is the call graph for this function:
Here is the caller graph for this function:

◆ chkmet()

int chkmet ( MMG5_pMesh  ,
MMG5_pSol   
)

◆ chknor()

int chknor ( MMG5_pMesh  mesh)
Parameters
meshpointer toward the mesh
Returns
1 if success, 0 if fail.

Check normal vectors consistency

Warning
unused
Here is the call graph for this function:
Here is the caller graph for this function:

◆ chkswp()

int chkswp ( MMG5_pMesh  mesh,
MMG5_pSol  met,
int  k,
int  i,
int8_t  typchk 
)
Here is the call graph for this function:
Here is the caller graph for this function:

◆ colver()

int colver ( MMG5_pMesh  mesh,
int *  list,
int  ilist 
)
Here is the call graph for this function:
Here is the caller graph for this function:

◆ colver2()

int colver2 ( MMG5_pMesh  mesh,
int *  ilist 
)
Here is the call graph for this function:
Here is the caller graph for this function:

◆ colver3()

int colver3 ( MMG5_pMesh  mesh,
int *  list 
)
Parameters
meshpointer toward the mesh structure.
listpointer toward the ball of the point to collapse.
Returns
1 if success, 0 if fail.

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.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ curvpo()

int curvpo ( MMG5_pMesh  ,
MMG5_pSol   
)

◆ delref()

int delref ( MMG5_pMesh  mesh)

◆ intmet_ani()

int intmet_ani ( MMG5_pMesh  mesh,
MMG5_pSol  met,
int  k,
int8_t  i,
int  ip,
double  s 
)
Parameters
meshpointer toward the mesh structure.
metpointer toward the metric structure.
kelement index.
ilocal index of edge in k.
ipglobal index of the new point in which we want to compute the metric.
sinterpolation parameter (between 0 and 1).
Returns
0 if fail, 1 otherwise.

Interpolation of anisotropic sizemap at parameter s along edge i of elt k for special storage of ridges metrics (after defsiz call).

Here is the call graph for this function:
Here is the caller graph for this function:

◆ intmet_iso()

int intmet_iso ( MMG5_pMesh  mesh,
MMG5_pSol  met,
int  k,
int8_t  i,
int  ip,
double  s 
)
Parameters
meshpointer toward the mesh structure.
metpointer toward the metric structure.
ktriangle in which we interpole the metrics.
iedge along which we interpole the metrics.
ipindex of point in which we compute the interpolated metric.
sparameter at which we compute the interpolation.
Returns
1 if success, 0 otherwise.

Linear interpolation of sizemap along edge i of tria k.

Here is the caller graph for this function:

◆ intregmet()

int intregmet ( MMG5_pMesh  mesh,
MMG5_pSol  met,
int  k,
int8_t  i,
double  s,
double  mr[6] 
)
Parameters
meshpointer toward the mesh structure.
metpointer toward the metric structure.
kelement index.
ilocal index of edge in k.
sinterpolation parameter.
mrcomputed metric.
Returns
call to MMG5_interpreg_ani (thus, 0 if fail, 1 otherwise).

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.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ litcol()

int litcol ( MMG5_pMesh  mesh,
int  k,
int8_t  i,
double  kal 
)
Here is the call graph for this function:

◆ litswp()

int litswp ( MMG5_pMesh  mesh,
int  k,
int8_t  i,
double  kal 
)
Here is the call graph for this function:

◆ MMG5_defaultValues()

void MMG5_defaultValues ( MMG5_pMesh  )

◆ MMG5_Init_parameters()

void MMG5_Init_parameters ( MMG5_pMesh  mesh)
Parameters
meshpointer toward the mesh structure.

Initialization of the input parameters.

MMG3D_IPARAM_lag is used by mmg3d only but need to be negative in the scaleMesh function

Here is the call graph for this function:
Here is the caller graph for this function:

◆ MMG5_intridmet()

int MMG5_intridmet ( MMG5_pMesh  ,
MMG5_pSol  ,
int  ,
int  ,
double  ,
double *  ,
double *   
)

◆ MMG5_memSize()

size_t MMG5_memSize ( void  )
Returns
the available memory size of the computer.

Compute the available memory size of the computer.

Here is the caller graph for this function:

◆ MMG5_mmgs1()

int MMG5_mmgs1 ( MMG5_pMesh  mesh,
MMG5_pSol  met,
int *  permNodGlob 
)
Parameters
meshpointer toward the mesh structure.
metpointer toward the metric structure.
permNodGlobif provided, strore the global permutation of nodes.
Returns
0 if failed, 1 if success.

Main adaptation routine.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ MMG5_mmgsBezierCP()

int MMG5_mmgsBezierCP ( MMG5_pMesh  mesh,
MMG5_Tria pt,
MMG5_pBezier  pb,
int8_t  ori 
)
Parameters
meshpointer toward the mesh structure.
ptpointer toward the triangle structure.
pbpointer toward the computed Bezier structure.
oritriangle orientation (unused but here for compatibility with the MMG5_bezierCP interface).
Returns
1.

Compute Bezier control points on triangle pt (cf. Vlachos)

Todo:
merge with the MMG5_mmg3dBeizerCP function and remove the pointer toward this functions.
Here is the call graph for this function:
Here is the caller graph for this function:

◆ MMG5_mmgsChkmsh()

int MMG5_mmgsChkmsh ( MMG5_pMesh  mesh,
int  severe,
int  base 
)
Parameters
meshpointer toward the mesh structure.
severelevel of performed check
baseunused argument.
Returns
0 if fail, 1 if success.

Check the mesh validity

Here is the call graph for this function:
Here is the caller graph for this function:

◆ MMG5_split2_sim()

int MMG5_split2_sim ( MMG5_pMesh  mesh,
MMG5_pSol  met,
int  k,
int *  vx 
)
Parameters
meshpointer toward the mesh structure.
metpointer toward the metric structure.
kindex of element to split.
vx$vx[i]$ is the index of the point to add on the edge i.
Returns
0 if split leads to invalid element, else 1.

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).

Here is the call graph for this function:
Here is the caller graph for this function:

◆ MMGS_analys()

int MMGS_analys ( MMG5_pMesh  mesh)
Here is the call graph for this function:
Here is the caller graph for this function:

◆ MMGS_assignEdge()

int MMGS_assignEdge ( MMG5_pMesh  mesh)
Parameters
meshpointer toward the mesh structure.
Returns
1 if success, 0.

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.

Remarks
this function handle all the provided edges.
Here is the call graph for this function:
Here is the caller graph for this function:

◆ MMGS_bdryUpdate()

int MMGS_bdryUpdate ( MMG5_pMesh  mesh)
Parameters
meshpointer toward the mesh structure.
Returns
1 if success, 0.

Copy the edge tags stored in triangles in the other triangles sharing the edge.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ MMGS_bezierInt()

int MMGS_bezierInt ( MMG5_pBezier  ,
double *  ,
double *  ,
double *  ,
double *   
)

◆ MMGS_defsiz_ani()

int MMGS_defsiz_ani ( MMG5_pMesh  mesh,
MMG5_pSol  met 
)
Parameters
meshpointer toward the mesh structure.
metpointer toward the metric stucture.
Returns
0 if fail, 1 otherwise.

Define size at points by intersecting the surfacic metric and the physical metric.

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

Here is the call graph for this function:
Here is the caller graph for this function:

◆ MMGS_defsiz_iso()

int MMGS_defsiz_iso ( MMG5_pMesh  mesh,
MMG5_pSol  met 
)
Parameters
meshpointer toward the mesh
metpointer toward the metric
Returns
1 if success, 0 if fail

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

Here is the call graph for this function:
Here is the caller graph for this function:

◆ MMGS_delElt()

int MMGS_delElt ( MMG5_pMesh  mesh,
int  iel 
)
Parameters
meshpointer toward the mesh
ielindex of the element to delete
Returns
1 if success, 0 if fail

Delete the element iel

Here is the caller graph for this function:

◆ MMGS_delPt()

void MMGS_delPt ( MMG5_pMesh  mesh,
int  ip 
)
Here is the caller graph for this function:

◆ MMGS_Free_all_var()

int MMGS_Free_all_var ( va_list  argptr)
Parameters
argptrlist 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 toward 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 toward 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).

Returns
0 if fail, 1 if success

Internal function for deallocations before return (taking a va_list as argument).

Remarks
we pass the structures by reference in order to have argument compatibility between the library call from a Fortran code and a C code.
Here is the call graph for this function:
Here is the caller graph for this function:

◆ MMGS_Free_names_var()

int MMGS_Free_names_var ( va_list  argptr)
Parameters
argptrlist 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 toward 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 toward 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).

Returns
0 if fail, 1 if success

Internal function for name deallocations before return (taking a va_list as argument).

Remarks
we pass the structures by reference in order to have argument compatibility between the library call from a Fortran code and a C code.
Here is the call graph for this function:
Here is the caller graph for this function:

◆ MMGS_Free_structures_var()

int MMGS_Free_structures_var ( va_list  argptr)
Parameters
argptrlist 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 toward 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 toward 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).

Returns
0 if fail, 1 if success

Internal function for structures deallocations before return (taking a va_list as argument).

Remarks
we pass the structures by reference in order to have argument compatibility between the library call from a Fortran code and a C code.
Here is the call graph for this function:
Here is the caller graph for this function:

◆ MMGS_gradsiz_ani()

int MMGS_gradsiz_ani ( MMG5_pMesh  mesh,
MMG5_pSol  met 
)
Parameters
meshpointer toward the mesh structure.
metpointer toward the metric structure.
Returns
1

Enforces mesh gradation by truncating metric field.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ MMGS_gradsizreq_ani()

int MMGS_gradsizreq_ani ( MMG5_pMesh  mesh,
MMG5_pSol  met 
)

◆ MMGS_hashTria()

int MMGS_hashTria ( MMG5_pMesh  mesh)
Parameters
meshpointer toward the mesh structure.
Returns
1 if success, 0 if fail.

Create adjacency table.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ MMGS_indElt()

int MMGS_indElt ( MMG5_pMesh  mesh,
int  kel 
)

find the element number in packed numerotation

Here is the caller graph for this function:

◆ MMGS_indPt()

int MMGS_indPt ( MMG5_pMesh  mesh,
int  kp 
)

find the point number in packed numerotation

Here is the caller graph for this function:

◆ MMGS_Init_mesh_var()

int MMGS_Init_mesh_var ( va_list  argptr)
Parameters
argptrlist 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 toward 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 toward a MMG5_pSol structure (that will contain the ouput metric (and the input one, if provided) and identified by the MMG5_ARG_ppMet keyword).

Returns
0 if fail, 1 if success

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).

Here is the call graph for this function:
Here is the caller graph for this function:

◆ MMGS_inqua()

int MMGS_inqua ( MMG5_pMesh  mesh,
MMG5_pSol  met 
)
Parameters
meshpointer toward the mesh structure.
metpointer toward the metric structure.
Returns
0 if the worst element has a nul quality, 1 otherwise.

Print histogram of mesh qualities for classical storage of ridges metrics (so before the the MMG5_defsiz function call).

Here is the call graph for this function:
Here is the caller graph for this function:

◆ MMGS_intmet33_ani()

int MMGS_intmet33_ani ( MMG5_pMesh  mesh,
MMG5_pSol  met,
int  k,
int8_t  i,
int  ip,
double  s 
)
Parameters
meshpointer toward the mesh structure.
metpointer toward the metric structure.
kelement index.
ilocal index of edge in k.
ipglobal index of the new point in which we want to compute the metric.
sinterpolation parameter (between 0 and 1).
Returns
0 if fail, 1 otherwise.

Interpolation of anisotropic sizemap at parameter s along edge i of elt k for classic storage of ridges metrics (before defsiz call).

Here is the call graph for this function:
Here is the caller graph for this function:

◆ MMGS_keep_only1Subdomain()

void MMGS_keep_only1Subdomain ( MMG5_pMesh  mesh,
int  nsd 
)
Parameters
meshpointer toward the mesh structure.
nsdindex of subdomain to keep.

Keep only subdomain of index nsd and remove other subdomains.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ MMGS_memOption()

int MMGS_memOption ( MMG5_pMesh  mesh)
Parameters
meshpointer toward the mesh structure
Returns
0 if fail, 1 otherwise

memory repartition for the -m option

Here is the call graph for this function:
Here is the caller graph for this function:

◆ MMGS_mmgs2()

int MMGS_mmgs2 ( MMG5_pMesh  mesh,
MMG5_pSol  sol,
MMG5_pSol  met 
)
Parameters
meshpointer toward the mesh structure.
solpointer toward the level-set
metpointer toward a metric (optionnal)
Returns
0 if fail, 1 otherwise.

Create implicit surface in mesh.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ MMGS_moveTowardPoint()

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] 
)
Parameters
meshpointer toward the mesh
p0point to move.
pneighbouring point toward which we try to move.
lloldinit length of edge p0-p
lam0first bezier basis function (order 2)
lam1second bezier basis function (order 2)
lam2third bezier basis function (order 2)
nn1normal at point p0 after relocation
nn2normal at point p0 after relocation
totangent along edge at point p0 after relocation
Returns
1 if success, 0 if fail

Update normals and tangent at ref or ridge point p0 after relocation at coordinates o.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ MMGS_newElt()

int MMGS_newElt ( MMG5_pMesh  mesh)
Here is the caller graph for this function:

◆ MMGS_newPt()

int MMGS_newPt ( MMG5_pMesh  mesh,
double  c[3],
double  n[3] 
)
Here is the caller graph for this function:

◆ MMGS_outqua()

int MMGS_outqua ( MMG5_pMesh  mesh,
MMG5_pSol  met 
)
Parameters
meshpointer toward the mesh structure.
metpointer toward the metric structure.
Returns
0 if the worst element has a nul quality, 1 otherwise.

Print histogram of mesh qualities for special storage of ridges metrics (after the MMG5_defsiz function call).

Here is the call graph for this function:
Here is the caller graph for this function:

◆ MMGS_paramDisp()

int MMGS_paramDisp ( MMG5_pMesh  mesh,
int  it1,
int  it2,
double  l1old,
double  l2old,
int8_t  isrid1,
int8_t  isrid2,
int  ip0,
int  ip1,
int  ip2,
double  step,
double  o[3],
int8_t *  isrid 
)
Parameters
meshpointer toward the mesh
it1triangle to which belongs the first edge
it2triangle to which belongs the second edge
l1oldlength of the first edge
l2oldlength of the second edge
isrid11 if the first edge is a ridge
isrid21 if the second edge is a ridge
ip0edge point that we want to move
ip1edge point connected by the ref/ridge edge to p0
ip2edge point connected by the ref/ridge edge to p0
stepdisplacement factor along the ref/ridge edge
ocoordinates of point after relocation
isrid1 if point is moved toward a ridge.
Returns
1 if success, 0 otherwise.

Infer arc length of displacement along ref or ridge edge, parameterized over edges.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ MMGS_prilen()

int MMGS_prilen ( MMG5_pMesh  mesh,
MMG5_pSol  met,
int  metRidTyp 
)
Parameters
meshpointer toward the mesh structure.
metpointer toward the metric structure.
metRidTypType of storage of ridges metrics: 0 for classic storage, 1 for special storage.
Returns
0 if fail, 1 otherwise.

Compute sizes of edges of the mesh, and displays histo.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ MMGS_Set_commonFunc()

static void MMGS_Set_commonFunc ( void  )
inlinestatic

Set common pointer functions between mmgs and mmg3d to the matching mmgs functions.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ MMGS_set_metricAtPointsOnReqEdges()

int MMGS_set_metricAtPointsOnReqEdges ( MMG5_pMesh  mesh,
MMG5_pSol  met,
int8_t  ismet 
)
Parameters
meshpointer toward the mesh
metpointer toward the metric
ismet1 if user provided metric
Returns
0 if fail, 1 otherwise

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.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ MMGS_setMeshSize_alloc()

int MMGS_setMeshSize_alloc ( MMG5_pMesh  mesh)
Parameters
meshpointer toward the mesh structure.
Returns
0 if failed, 1 otherwise.

Allocation of the array fields of the mesh.

Here is the caller graph for this function:

◆ MMGS_simbulgept()

int MMGS_simbulgept ( MMG5_pMesh  mesh,
MMG5_pSol  met,
int  k,
int  i,
int  ip 
)
Parameters
meshpointer toward the mesh structure.
metpointer toward the metric structure.
kindex of the starting triangle.
ilocal index of the edge to split in k.
ipindex of the point that we try to create.
Returns
0 if final position is invalid, 1 if all checks are ok.

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).

Remarks
Don't work for non-manifold edge.
Here is the caller graph for this function:

◆ MMGS_split1()

int MMGS_split1 ( MMG5_pMesh  mesh,
MMG5_pSol  met,
int  k,
int  i,
int *  vx 
)
Parameters
meshpointer toward the mesh structure.
metpointer toward the metric structure.
kindex of element to split.
iindex of edge to split.
vx$vx[i]$ is the index of the point to add on the edge i.
Returns
1 if success, 0 if fail.

Split element k along edge i.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ MMGS_split1_sim()

int MMGS_split1_sim ( MMG5_pMesh  mesh,
MMG5_pSol  met,
int  k,
int  i,
int *  vx 
)
Parameters
meshpointer toward the mesh structure.
metpointer toward the metric structure.
kindex of element to split.
iindex of edge to split.
vx$vx[i]$ is the index of the point to add on the edge i.
kindex of element to split.
Returns
0 if split leads to invalid element, else 1.

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).

Here is the call graph for this function:
Here is the caller graph for this function:

◆ MMGS_split2()

int MMGS_split2 ( MMG5_pMesh  mesh,
MMG5_pSol  met,
int  k,
int *  vx 
)
Parameters
meshpointer toward the mesh structure.
metpointer toward the metric structure.
kindex of element to split.
vx$vx[i]$ is the index of the point to add on the edge i.
Returns
1 if success, 0 if fail.

Split element k along the 2 edges i1 and i2.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ MMGS_split3()

int MMGS_split3 ( MMG5_pMesh  mesh,
MMG5_pSol  met,
int  k,
int *  vx 
)
Parameters
meshpointer toward the mesh structure.
metpointer toward the metric structure.
kindex of element to split.
vx$vx[i]$ is the index of the point to add on the edge i.
Returns
1 if success, 0 if fail.

Split element k along the 3 edges

Here is the call graph for this function:
Here is the caller graph for this function:

◆ MMGS_split3_sim()

int MMGS_split3_sim ( MMG5_pMesh  mesh,
MMG5_pSol  met,
int  k,
int *  vx 
)
Parameters
meshpointer toward the mesh structure.
metpointer toward the metric structure.
kindex of element to split.
vx$vx[i]$ is the index of the point to add on the edge i.
Returns
0 if split leads to invalid element, else 1.

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).

Here is the call graph for this function:
Here is the caller graph for this function:

◆ MMGS_surfballRotation()

int MMGS_surfballRotation ( MMG5_pMesh  mesh,
MMG5_pPoint  p0,
int *  list,
int  ilist,
double  r[3][3],
double *  lispoi 
)
Parameters
meshpointer toward the mesh structure.
p0starting point
listball of p0
ilistnumber of tria in the ball of p0
rrotation that send the normal at p0 onto the z vector
lipointrotated ball of point p0
Returns
1 if success, 0 otherwise.

Compute the rotation matrix that sends the tangent plane at p0 onto z=0 and apply this rotation to the ball of p0.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ MMGS_zaldy()

int MMGS_zaldy ( MMG5_pMesh  mesh)
Parameters
meshpointer toward the mesh
Returns
1 if success, 0 if fail

allocate main structure

Here is the call graph for this function:
Here is the caller graph for this function:

◆ movintpt_ani()

int movintpt_ani ( MMG5_pMesh  mesh,
MMG5_pSol  met,
int *  list,
int  ilist 
)
Parameters
meshpointer toward the mesh structure.
metpointer toward the metric structure.
listball of point.
ilistsize of the point ball.
Returns
0 if fail, 1 otherwise.

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

Here is the call graph for this function:
Here is the caller graph for this function:

◆ movintpt_iso()

int movintpt_iso ( MMG5_pMesh  mesh,
MMG5_pSol  met,
int *  list,
int  ilist 
)
Here is the call graph for this function:
Here is the caller graph for this function:

◆ movridpt_ani()

int movridpt_ani ( MMG5_pMesh  mesh,
MMG5_pSol  met,
int *  list,
int  ilist 
)
Here is the call graph for this function:
Here is the caller graph for this function:

◆ movridpt_iso()

int movridpt_iso ( MMG5_pMesh  mesh,
MMG5_pSol  met,
int *  list,
int  ilist 
)
Here is the call graph for this function:
Here is the caller graph for this function:

◆ paratmet()

int paratmet ( double  c0[3],
double  n0[3],
double  m[6],
double  c1[3],
double  n1[3],
double  mt[6] 
)

◆ setref()

int setref ( MMG5_pMesh  mesh,
int  start,
int  ref,
int  putreq 
)
Parameters
meshpointer toward the mesh
startindex of the tetra from which we start
refreference to set
putreq1 if boundary edges must be set to required
Returns
1 if success, 0 if fail

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.

Here is the call graph for this function:

◆ split1b()

int split1b ( MMG5_pMesh  mesh,
int  k,
int8_t  i,
int  ip 
)
Parameters
meshpointer toward the mesh structure.
kindex of element to split.
iindex of edge to split.
ipindex of the new point.
Returns
0 if lack of memory, 1 otherwise.

Split element k along edge i, inserting point ip and updating the adjacency relations.

Remarks
do not call this function in non-manifold case
Here is the call graph for this function:
Here is the caller graph for this function:

◆ swapar()

int swapar ( MMG5_pMesh  mesh,
int  k,
int  i 
)
Parameters
meshpoiner toward the mesh structure.
kelt index.
iindex of the elt edge to swap.
Returns
1
Warning
the quality of the resulting triangles is not checked here... It must be checked outside to prevent the creation of empty elts.
Here is the caller graph for this function:

◆ swpedg()

int swpedg ( MMG5_pMesh  mesh,
MMG5_pSol  met,
int *  list,
int  ilist,
int8_t  typchk 
)

attempt to swap any edge below quality value list goes from 0 to ilist-1.

Warning
not used
Here is the call graph for this function:

◆ typelt()

int8_t typelt ( MMG5_pPoint  p[3],
int8_t *  ia 
)

Variable Documentation

◆ intmet

int(* intmet) (MMG5_pMesh mesh, MMG5_pSol met, int k, int8_t i, int ip, double s) ( MMG5_pMesh  mesh,
MMG5_pSol  met,
int  k,
int8_t  i,
int  ip,
double  s 
)
extern

◆ MMG5_calelt

double(* MMG5_calelt) (MMG5_pMesh mesh, MMG5_pSol met, MMG5_pTria ptt) ( MMG5_pMesh  mesh,
MMG5_pSol  met,
MMG5_pTria  ptt 
)
extern

◆ MMGS_defsiz

int(* MMGS_defsiz) (MMG5_pMesh mesh, MMG5_pSol met) ( MMG5_pMesh  mesh,
MMG5_pSol  met 
)
extern

◆ MMGS_gradsiz

int(* MMGS_gradsiz) (MMG5_pMesh mesh, MMG5_pSol met) ( MMG5_pMesh  mesh,
MMG5_pSol  met 
)
extern

◆ MMGS_gradsizreq

int(* MMGS_gradsizreq) (MMG5_pMesh mesh, MMG5_pSol met) ( MMG5_pMesh  mesh,
MMG5_pSol  met 
)
extern

◆ movintpt

int(* movintpt) (MMG5_pMesh mesh, MMG5_pSol met, int *list, int ilist) ( MMG5_pMesh  mesh,
MMG5_pSol  met,
int *  list,
int  ilist 
)
extern

◆ movridpt

int(* movridpt) (MMG5_pMesh mesh, MMG5_pSol met, int *list, int ilist) ( MMG5_pMesh  mesh,
MMG5_pSol  met,
int *  list,
int  ilist 
)
extern