mmgs
mmgs.h File Reference
#include "libmmgs.h"
#include "mmgcommon.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, val)
 
#define _MMGS_POINT_REALLOC(mesh, sol, ip, wantedGap, law, o, tag, retval)
 
#define _MMGS_TRIA_REALLOC(mesh, jel, wantedGap, law, retval)
 

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 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 _MMGS_mmgs2 (MMG5_pMesh, MMG5_pSol)
 
int boulet (MMG5_pMesh mesh, int start, int ip, int *list)
 
int boulechknm (MMG5_pMesh mesh, int start, int ip, int *list)
 
int boulep (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, char 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, char i, int ip)
 
int chkcol (MMG5_pMesh mesh, MMG5_pSol met, int k, char i, int *list, char 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, char typchk)
 
int swpedg (MMG5_pMesh mesh, MMG5_pSol met, int *list, int ilist, char typchk)
 
char typelt (MMG5_pPoint p[3], char *ia)
 
int litswp (MMG5_pMesh mesh, int k, char i, double kal)
 
int litcol (MMG5_pMesh mesh, int k, char 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, char 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)
 
long long _MMG5_memSize (void)
 
int _MMGS_memOption (MMG5_pMesh mesh)
 
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 gradsiz_iso (MMG5_pMesh mesh, MMG5_pSol met)
 
int gradsiz_ani (MMG5_pMesh mesh, MMG5_pSol met)
 
int intmet_iso (MMG5_pMesh mesh, MMG5_pSol met, int k, char i, int ip, double s)
 
int intmet_ani (MMG5_pMesh mesh, MMG5_pSol met, int k, char i, int ip, double s)
 
int _MMGS_intmet33_ani (MMG5_pMesh, MMG5_pSol, int, char, int, double)
 
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_prilen (MMG5_pMesh mesh, MMG5_pSol met, int)
 
static void _MMGS_Set_commonFunc ()
 

Variables

double(* _MMG5_calelt )(MMG5_pMesh mesh, MMG5_pSol met, MMG5_pTria ptt)
 
int(* _MMG5_defsiz )(MMG5_pMesh mesh, MMG5_pSol met)
 
int(* gradsiz )(MMG5_pMesh mesh, MMG5_pSol met)
 
int(* intmet )(MMG5_pMesh mesh, MMG5_pSol met, int k, char 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,
  retval 
)
Value:
do \
{ \
int klink; \
"larger point table",law,retval); \
\
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*(mesh->npmax+1),double, \
"larger solution",retval); \
} \
sol->npmax = mesh->npmax; \
\
/* We try again to add the point */ \
ip = _MMGS_newPt(mesh,o,tag); \
if ( !ip ) {law;} \
}while(0)
MMG5_pPoint point
Definition: libmmgtypes.h:505
#define _MMG5_TAB_RECALLOC(mesh, ptr, initSize, wantedGap, type, message, law, retval)
Definition: mmgcommon.h:253
int npmax
Definition: libmmgtypes.h:480
int npnil
Definition: libmmgtypes.h:490
! int npmax
Definition: libmmgtypesf.h:530
int np
Definition: libmmgtypes.h:480
Structure to store points of a MMG mesh.
Definition: libmmgtypes.h:205
MMG5_pMesh char int int * retval
Definition: API_functionsf_s.c:584
int _MMGS_newPt(MMG5_pMesh mesh, double c[3], double n[3])
Definition: zaldy_s.c:39
! int16_t tag
Definition: libmmgtypesf.h:252
MMG5_pMesh * mesh
Definition: API_functionsf_s.c:63

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,
  val 
)
Value:
do \
{ \
MMG5_ARG_end) ) { \
return MMG5_LOWFAILURE; \
} \
return(val); \
}while(0)
#define MMG5_ARG_start
Definition: libmmgtypes.h:73
double val
Definition: mmgcommon.h:427
#define MMG5_LOWFAILURE
Definition: libmmgtypes.h:48
#define MMG5_ARG_end
Definition: libmmgtypes.h:159
int MMGS_Free_all(const int starter,...)
Definition: API_functions_s.c:1274
MMG5_pMesh * mesh
Definition: API_functionsf_s.c:63
#define MMG5_ARG_ppMet
Definition: libmmgtypes.h:102
#define MMG5_ARG_ppMesh
Definition: libmmgtypes.h:82

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

◆ _MMGS_TRIA_REALLOC

#define _MMGS_TRIA_REALLOC (   mesh,
  jel,
  wantedGap,
  law,
  retval 
)
Value:
do \
{ \
int klink,oldSiz; \
\
oldSiz = mesh->ntmax; \
_MMG5_TAB_RECALLOC(mesh,mesh->tria,mesh->ntmax,wantedGap,MMG5_Tria, \
"larger tria table",law,retval); \
\
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*mesh->nt+5,3*mesh->ntmax+5,int \
,"larger adja table",retval); \
} \
\
/* We try again to add the point */ \
jel = _MMGS_newElt(mesh); \
if ( !jel ) {law;} \
}while(0)
int _MMGS_newElt(MMG5_pMesh mesh)
Definition: zaldy_s.c:71
int ntmax
Definition: libmmgtypes.h:480
! int ntmax
Definition: libmmgtypesf.h:530
int nt
Definition: libmmgtypes.h:480
if(!ier) exit(EXIT_FAILURE)
int nenil
Definition: libmmgtypes.h:491
MMG5_pMesh char int int * retval
Definition: API_functionsf_s.c:584
int * adja
Definition: libmmgtypes.h:493
MMG5_pTria tria
Definition: libmmgtypes.h:511
MMG5_pMesh * mesh
Definition: API_functionsf_s.c:63
Definition: libmmgtypes.h:261

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

◆ _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 (stored in the Info structure).

MMG5_IPARAM_verbose = 1















MMG*_IPARAM_iso = 0














MMG5_IPARAM_mem = -1













MMG5_IPARAM_debug = 0












MMG5_IPARAM_npar = 0











MMG5_IPARAM_noinsert = 0










MMG5_IPARAM_noswap = 0









MMG5_IPARAM_nomove = 0








MMG5_IPARAM nmat = 0







MMG5_DPARAM_angleDetection = _MMG5_ANGEDG






MMG5_DPARAM_hmin = 0.001 $\times$ bounding box size;





MMG5_DPARAM_hmax = double of the bounding box size




MMG5_DPARAM_hsiz= -1.



MMG5_DPARAM_hausd = 0.01


MMG5_DPARAM_hgrad = 1.3

MMG5_PPARAM = NULL
MMG3D_IPARAM_lag = -1 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()

long long _MMG5_memSize ( void  )
Returns
the available memory size of the computer.

Compute the available memory size of the computer.

◆ _MMG5_mmgs1()

int _MMG5_mmgs1 ( MMG5_pMesh  ,
MMG5_pSol   
)
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,
char  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_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.

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

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_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,
char  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_memOption()

int _MMGS_memOption ( MMG5_pMesh  mesh)

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 
)
Parameters
meshpointer toward the mesh structure.
solpointer toward the solution structure
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_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_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 ( )
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_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_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:

◆ assignEdge()

int 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:

◆ 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:

◆ boulep()

int boulep ( MMG5_pMesh  mesh,
int  start,
int  ip,
int *  list 
)
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 
)
inline
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 
)
inline
Here is the caller graph for this function:

◆ chkcol()

int chkcol ( MMG5_pMesh  mesh,
MMG5_pSol  met,
int  k,
char  i,
int *  list,
char  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  ,
int   
)
Here is the call graph for this function:
Here is the caller graph for this function:

◆ chkmet()

int chkmet ( MMG5_pMesh  mesh,
MMG5_pSol  met 
)
Returns
0 if fail , 1 if success

Check metric consistency.

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

◆ 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:

◆ chkswp()

int chkswp ( MMG5_pMesh  mesh,
MMG5_pSol  met,
int  k,
int  i,
char  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  )

◆ gradsiz_ani()

int 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:

◆ gradsiz_iso()

int gradsiz_iso ( MMG5_pMesh  mesh,
MMG5_pSol  met 
)

Enforces mesh gradations by truncating size map

Here is the caller graph for this function:

◆ intmet_ani()

int intmet_ani ( MMG5_pMesh  mesh,
MMG5_pSol  met,
int  k,
char  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,
char  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,
char  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,
char  i,
double  kal 
)
Here is the call graph for this function:

◆ litswp()

int litswp ( MMG5_pMesh  mesh,
int  k,
char  i,
double  kal 
)
Here is the call 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 1 : Rotation matrix that sends normal at p0 to e_z
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,
char  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,
char  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()

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

Variable Documentation

◆ _MMG5_calelt

double(* _MMG5_calelt) (MMG5_pMesh mesh, MMG5_pSol met, MMG5_pTria ptt)

◆ _MMG5_defsiz

int(* _MMG5_defsiz) (MMG5_pMesh mesh, MMG5_pSol met)

◆ gradsiz

int(* gradsiz) (MMG5_pMesh mesh, MMG5_pSol met)

◆ intmet

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

◆ movintpt

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

◆ movridpt

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