mmgs
mmgcommon.h File Reference
#include <assert.h>
#include <stdlib.h>
#include <stdio.h>
#include <limits.h>
#include <string.h>
#include <signal.h>
#include <ctype.h>
#include <float.h>
#include <math.h>
#include <complex.h>
#include "eigenv.h"
#include "libmmgcommon.h"
Include dependency graph for mmgcommon.h:
This graph shows which files directly or indirectly include this file:

Go to the source code of this file.

Data Structures

struct  _MMG5_Bezier
 
struct  _MMG5_hedge
 Used to hash edges (memory economy compared to MMG5_hgeom). More...
 
struct  _MMG5_Hash
 Identic as MMG5_HGeom but use _MMG5_hedge to store edges instead of MMG5_hgeom (memory economy). More...
 
struct  _MMG5_iNode_s
 
struct  _MMG5_dNode_s
 

Macros

#define POSIX
 
#define GNU
 
#define MG_VER   "5.3.11"
 
#define MG_REL   "Apr. 10, 2017"
 
#define MG_CPY   "Copyright (c) IMB-LJLL, 2004-"
 
#define MG_STR   "&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&"
 
#define MG_SMSGN(a, b)   (((double)(a)*(double)(b) > (0.0)) ? (1) : (0))
 
#define _MMG5_BOXSIZE   500
 
#define _MMG5_MEMMAX   800
 
#define MG_PLUS   2
 
#define MG_MINUS   3
 
#define _MMG5_ANGEDG   0.707106781186548 /*0.573576436351046 */
 
#define _MMG5_ANGLIM   -0.999999
 
#define _MMG5_ATHIRD   0.333333333333333
 
#define _MMG5_EPSD   1.e-30
 
#define _MMG5_EPSD2   1.0e-200
 
#define _MMG5_EPS   1.e-06
 
#define _MMG5_EPSOK   1.e-15
 
#define _MMG5_NULKAL   1.e-30
 
#define _MMG5_SQR32   0.866025403784439
 
#define A64TH   0.015625
 
#define A16TH   0.0625
 
#define A32TH   0.03125
 
#define _MMG5_MEMMIN   38
 
#define MG_MAX(a, b)   (((a) > (b)) ? (a) : (b))
 
#define MG_MIN(a, b)   (((a) < (b)) ? (a) : (b))
 
#define MG_NOTAG   (0)
 
#define MG_REF   (1 << 0)
 
#define MG_GEO   (1 << 1)
 
#define MG_REQ   (1 << 2)
 
#define MG_NOM   (1 << 3)
 
#define MG_BDY   (1 << 4)
 
#define MG_CRN   (1 << 5)
 
#define MG_NOSURF   (1 << 6)
 
#define MG_OPNBDY   (1 << 7)
 
#define MG_PARBDY   (1 << 13)
 
#define MG_NUL   (1 << 14)
 
#define MG_Vert   (1 << 0 )
 
#define MG_Tria   (1 << 1 )
 
#define MG_Tetra   (1 << 2 )
 
#define MG_VOK(ppt)   (ppt && ((ppt)->tag < MG_NUL))
 
#define MG_EOK(pt)   (pt && ((pt)->v[0] > 0))
 
#define MG_EDG(tag)   ((tag & MG_GEO) || (tag & MG_REF))
 
#define MG_SIN(tag)   ((tag & MG_CRN) || (tag & MG_REQ))
 
#define MG_SET(flag, bit)   ((flag) |= (1 << (bit)))
 
#define MG_CLR(flag, bit)   ((flag) &= ~(1 << (bit)))
 
#define MG_GET(flag, bit)   ((flag) & (1 << (bit)))
 
#define _MMG5_KA   7
 
#define _MMG5_KB   11
 
#define _LIBMMG5_RETURN(mesh, met, val)
 
#define _MMG5_CHK_MEM(mesh, size, string, law)
 
#define _MMG5_DEL_MEM(mesh, ptr, size)
 
#define _MMG5_ADD_MEM(mesh, size, message, law)
 
#define _MMG5_SAFE_FREE(ptr)
 
#define _MMG5_SAFE_CALLOC(ptr, size, type, retval)
 
#define _MMG5_SAFE_MALLOC(ptr, size, type, retval)
 
#define _MMG5_SAFE_REALLOC(ptr, size, type, message, retval)
 
#define _MMG5_SAFE_RECALLOC(ptr, prevSize, newSize, type, message, retval)
 
#define _MMG5_TAB_RECALLOC(mesh, ptr, initSize, wantedGap, type, message, law, retval)
 
#define _MMG5_INCREASE_MEM_MESSAGE()
 
#define _MMG5_SAFELL2LCAST(longlongval)   (((longlongval) > (LONG_MAX)) ? 0 : ((long)(longlongval)))
 
#define _MMG5_SAFELL2ICAST(longlongval)   (((longlongval) > (INT_MAX)) ? 0 : ((int)(longlongval)))
 
#define FORTRAN_NAME(nu, nl, pl, pc)
 Adds function definitions. More...
 
#define FORTRAN_VARIADIC(nu, nl, pl, body)
 Adds function definitions. More...
 

Typedefs

typedef _MMG5_Bezier_MMG5_pBezier
 
typedef struct _MMG5_iNode_s _MMG5_iNode
 
typedef struct _MMG5_dNode_s _MMG5_dNode
 

Functions

static void _MMG5_excfun (int sigid)
 
double _MMG5_det3pt1vec (double c0[3], double c1[3], double c2[3], double v[3])
 
double _MMG5_det4pt (double c0[3], double c1[3], double c2[3], double c3[3])
 
int _MMG5_devangle (double *n1, double *n2, double crit)
 
double _MMG5_orvol (MMG5_pPoint point, int *v)
 
int _MMG5_Add_inode (MMG5_pMesh mesh, _MMG5_iNode **liLi, int val)
 
int _MMG5_Alloc_inode (MMG5_pMesh mesh, _MMG5_iNode **node)
 
int _MMG5_Add_dnode (MMG5_pMesh mesh, _MMG5_dNode **liLi, int, double)
 
int _MMG5_Alloc_dnode (MMG5_pMesh mesh, _MMG5_dNode **node)
 
void _MMG5_bezierEdge (MMG5_pMesh, int, int, double *, double *, char, double *)
 
int _MMG5_buildridmet (MMG5_pMesh, MMG5_pSol, int, double, double, double, double *)
 
int _MMG5_buildridmetfic (MMG5_pMesh, double *, double *, double, double, double, double *)
 
int _MMG5_buildridmetnor (MMG5_pMesh, MMG5_pSol, int, double *, double *)
 
int _MMG5_paratmet (double c0[3], double n0[3], double m[6], double c1[3], double n1[3], double mt[6])
 
int _MMG5_rmtr (double r[3][3], double m[6], double mr[6])
 
int _MMG5_boundingBox (MMG5_pMesh mesh)
 
int _MMG5_boulec (MMG5_pMesh, int *, int, int i, double *tt)
 
int _MMG5_boulen (MMG5_pMesh, int *, int, int i, double *nn)
 
int _MMG5_bouler (MMG5_pMesh, int *, int, int i, int *, int *, int *, int *, int)
 
double _MMG5_caltri33_ani (MMG5_pMesh mesh, MMG5_pSol met, MMG5_pTria pt)
 
double _MMG5_caltri_ani (MMG5_pMesh mesh, MMG5_pSol met, MMG5_pTria ptt)
 
double _MMG5_caltri_iso (MMG5_pMesh mesh, MMG5_pSol met, MMG5_pTria ptt)
 
void _MMG5_defUninitSize (MMG5_pMesh mesh, MMG5_pSol met, char ismet)
 
void _MMG5_displayHisto (MMG5_pMesh, int, double *, int, int, double, int, int, double, int, double *, int *, char)
 
int _MMG5_minQualCheck (int iel, double minqual, double alpha)
 
int _MMG5_elementWeight (MMG5_pMesh, MMG5_pSol, MMG5_pTria, MMG5_pPoint, _MMG5_Bezier *, double r[3][3], double gv[2])
 
void _MMG5_fillDefmetregSys (int, MMG5_pPoint, int, _MMG5_Bezier, double r[3][3], double *, double *, double *, double *)
 
void _MMG5_Free_ilinkedList (MMG5_pMesh mesh, _MMG5_iNode *liLi)
 
void _MMG5_Free_dlinkedList (MMG5_pMesh mesh, _MMG5_dNode *liLi)
 
int _MMG5_grad2metSurf (MMG5_pMesh mesh, MMG5_pSol met, MMG5_pTria pt, int i)
 
int _MMG5_hashEdge (MMG5_pMesh mesh, _MMG5_Hash *hash, int a, int b, int k)
 
int _MMG5_hashUpdate (_MMG5_Hash *hash, int a, int b, int k)
 
int _MMG5_hashGet (_MMG5_Hash *hash, int a, int b)
 
int _MMG5_hashNew (MMG5_pMesh mesh, _MMG5_Hash *hash, int hsiz, int hmax)
 
int _MMG5_intmetsavedir (MMG5_pMesh mesh, double *m, double *n, double *mr)
 
int _MMG5_intridmet (MMG5_pMesh, MMG5_pSol, int, int, double, double *, double *)
 
int _MMG5_mmgIntmet33_ani (double *, double *, double *, double)
 
int _MMG5_mmgIntextmet (MMG5_pMesh, MMG5_pSol, int, double *, double *)
 
long long _MMG5_memSize (void)
 
void _MMG5_mmgDefaultValues (MMG5_pMesh mesh)
 
int _MMG5_mmgHashTria (MMG5_pMesh mesh, int *adja, _MMG5_Hash *, int chkISO)
 
void _MMG5_mmgInit_parameters (MMG5_pMesh mesh)
 
void _MMG5_mmgUsage (char *prog)
 
int _MMG5_nonUnitNorPts (MMG5_pMesh, int, int, int, double *)
 
double _MMG5_nonorsurf (MMG5_pMesh mesh, MMG5_pTria pt)
 
int _MMG5_norpts (MMG5_pMesh, int, int, int, double *)
 
int _MMG5_nortri (MMG5_pMesh mesh, MMG5_pTria pt, double *n)
 
void _MMG5_printTria (MMG5_pMesh mesh, char *fileName)
 
int _MMG5_rotmatrix (double n[3], double r[3][3])
 
int _MMG5_invmat (double *m, double *mi)
 
int _MMG5_invmatg (double m[9], double mi[9])
 
double _MMG5_ridSizeInNormalDir (MMG5_pMesh, int, double *, _MMG5_pBezier, double, double)
 
double _MMG5_ridSizeInTangentDir (MMG5_pMesh, MMG5_pPoint, int, int *, double, double)
 
int _MMG5_scaleMesh (MMG5_pMesh mesh, MMG5_pSol met)
 
int _MMG5_scotchCall (MMG5_pMesh mesh, MMG5_pSol sol)
 
int _MMG5_solveDefmetregSys (MMG5_pMesh, double r[3][3], double *, double *, double *, double *, double, double, double)
 
int _MMG5_solveDefmetrefSys (MMG5_pMesh, MMG5_pPoint, int *, double r[3][3], double *, double *, double *, double *, double, double, double)
 
double _MMG5_surftri_ani (MMG5_pMesh mesh, MMG5_pSol met, MMG5_pTria ptt)
 
double _MMG5_surftri33_ani (MMG5_pMesh, MMG5_pTria, double *, double *, double *)
 
double _MMG5_surftri_iso (MMG5_pMesh mesh, MMG5_pSol met, MMG5_pTria ptt)
 
int _MMG5_sys33sym (double a[6], double b[3], double r[3])
 
int _MMG5_unscaleMesh (MMG5_pMesh mesh, MMG5_pSol met)
 
int _MMG5_interpreg_ani (MMG5_pMesh, MMG5_pSol, MMG5_pTria, char, double, double *mr)
 
int _MMG5_interp_iso (double *ma, double *mb, double *mp, double t)
 
int _MMG5_intersecmet22 (MMG5_pMesh mesh, double *m, double *n, double *mr)
 
int _MMG5_countLocalParamAtTri (MMG5_pMesh, _MMG5_iNode **)
 
int _MMG5_writeLocalParamAtTri (MMG5_pMesh, _MMG5_iNode *, FILE *)
 
double MMG2_quickarea (double a[2], double b[2], double c[2])
 
int MMG5_loadMshMesh_part1 (MMG5_pMesh mesh, const char *filename, FILE **inm, long *posNodes, long *posElts, long **posNodeData, int *bin, int *iswp, int *nelts, int *nsols)
 
int MMG5_loadMshMesh_part2 (MMG5_pMesh mesh, MMG5_pSol *sol, FILE **inm, const long posNodes, const long posElts, const long *posNodeData, const int bin, const int iswp, const int nelts)
 
int MMG5_saveMshMesh (MMG5_pMesh, MMG5_pSol *, const char *, const int)
 
int MMG5_loadSolHeader (const char *, int, FILE **, int *, int *, int *, int *, int *, int *, int **, long *)
 
int MMG5_chkMetricType (MMG5_pMesh mesh, int *type, FILE *inm)
 
void MMG5_readFloatSol3D (MMG5_pSol, FILE *, int, int, int)
 
void MMG5_readDoubleSol3D (MMG5_pSol, FILE *, int, int, int)
 
int MMG5_saveSolHeader (MMG5_pMesh, const char *, FILE **, int, int *, int, int, int, int *, int *)
 
void MMG5_writeDoubleSol3D (MMG5_pMesh, MMG5_pSol, FILE *, int, int, int)
 
void MMG5_printMetStats (MMG5_pMesh mesh, MMG5_pSol met)
 
void MMG5_printSolStats (MMG5_pMesh mesh, MMG5_pSol *sol)
 
void MMG5_chooseOutputFormat (MMG5_pMesh mesh, int *msh)
 
void _MMG5_Set_commonFunc ()
 

Variables

static const unsigned char _MMG5_inxt2 [6] = {1,2,0,1,2}
 
static const unsigned char _MMG5_iprv2 [3] = {2,0,1}
 
int(* _MMG5_chkmsh )(MMG5_pMesh, int, int)
 
int(* _MMG5_bezierCP )(MMG5_pMesh, MMG5_Tria *, _MMG5_pBezier, char)
 
double(* _MMG5_lenSurfEdg )(MMG5_pMesh mesh, MMG5_pSol sol, int, int, char)
 
int(* _MMG5_indElt )(MMG5_pMesh mesh, int kel)
 
int(* _MMG5_indPt )(MMG5_pMesh mesh, int kp)
 

Macro Definition Documentation

◆ _LIBMMG5_RETURN

#define _LIBMMG5_RETURN (   mesh,
  met,
  val 
)
Value:
do \
{ \
signal(SIGABRT,SIG_DFL); \
signal(SIGFPE,SIG_DFL); \
signal(SIGILL,SIG_DFL); \
signal(SIGSEGV,SIG_DFL); \
signal(SIGTERM,SIG_DFL); \
signal(SIGINT,SIG_DFL); \
mesh->npi = mesh->np; \
mesh->nti = mesh->nt; \
mesh->nai = mesh->na; \
mesh->nei = mesh->ne; \
met->npi = met->np; \
return(val); \
}while(0)
double val
Definition: mmgcommon.h:428
int nt
Definition: libmmgtypes.h:480
int np
Definition: libmmgtypes.h:480
int na
Definition: libmmgtypes.h:480
int ne
Definition: libmmgtypes.h:480
MMG5_pMesh * mesh
Definition: API_functionsf_s.c:63

Reset the customized signals and set the internal counters of points, edges, tria and tetra to the suitable value (needed by users to recover their mesh using the API)

◆ _MMG5_ADD_MEM

#define _MMG5_ADD_MEM (   mesh,
  size,
  message,
  law 
)
Value:
do \
{ \
(mesh)->memCur += (long long)(size); \
_MMG5_CHK_MEM(mesh,size,message,law); \
}while(0)
! int size
Definition: libmmgtypesf.h:583
! long long memCur
Definition: libmmgtypesf.h:525
MMG5_pMesh * mesh
Definition: API_functionsf_s.c:63

Increment memory counter memCur and check if we don't overflow the maximum authorizied memory memMax.

◆ _MMG5_ANGEDG

#define _MMG5_ANGEDG   0.707106781186548 /*0.573576436351046 */

◆ _MMG5_ANGLIM

#define _MMG5_ANGLIM   -0.999999

◆ _MMG5_ATHIRD

#define _MMG5_ATHIRD   0.333333333333333

◆ _MMG5_BOXSIZE

#define _MMG5_BOXSIZE   500

size of box for renumbering with scotch.

◆ _MMG5_CHK_MEM

#define _MMG5_CHK_MEM (   mesh,
  size,
  string,
  law 
)
Value:
do \
{ \
if ( ((mesh)->memCur) > ((mesh)->memMax) || \
((mesh)->memCur < 0 )) { \
fprintf(stderr," ## Error:"); \
fprintf(stderr," unable to allocate %s.\n",string); \
fprintf(stderr," ## Check the mesh size or "); \
fprintf(stderr,"increase maximal authorized memory with the -m option.\n"); \
(mesh)->memCur -= (long long)(size); \
law; \
} \
}while(0)
! int size
Definition: libmmgtypesf.h:583
! long long memCur
Definition: libmmgtypesf.h:525
! long long memMax
Definition: libmmgtypesf.h:524
MMG5_pMesh * mesh
Definition: API_functionsf_s.c:63

Check if used memory overflow maximal authorized memory. Execute the command law if lack of memory.

◆ _MMG5_DEL_MEM

#define _MMG5_DEL_MEM (   mesh,
  ptr,
  size 
)
Value:
do \
{ \
(mesh)->memCur -= (long long)(size); \
free(ptr); \
ptr = NULL; \
}while(0)
! int size
Definition: libmmgtypesf.h:583
! long long memCur
Definition: libmmgtypesf.h:525
MMG5_pMesh * mesh
Definition: API_functionsf_s.c:63

Free pointer ptr of mesh structure and compute the new used memory. size is the size of the pointer

◆ _MMG5_EPS

#define _MMG5_EPS   1.e-06

◆ _MMG5_EPSD

#define _MMG5_EPSD   1.e-30

◆ _MMG5_EPSD2

#define _MMG5_EPSD2   1.0e-200

◆ _MMG5_EPSOK

#define _MMG5_EPSOK   1.e-15

◆ _MMG5_INCREASE_MEM_MESSAGE

#define _MMG5_INCREASE_MEM_MESSAGE ( )
Value:
do \
{ \
printf(" ## Check the mesh size or increase maximal"); \
printf(" authorized memory with the -m option.\n"); \
} while(0)

Error message when lack of memory

◆ _MMG5_KA

#define _MMG5_KA   7

Key for hash tables.

◆ _MMG5_KB

#define _MMG5_KB   11

Key for hash tables.

◆ _MMG5_MEMMAX

#define _MMG5_MEMMAX   800

Maximal memory used if available memory compitation fail.

◆ _MMG5_MEMMIN

#define _MMG5_MEMMIN   38

minimal memory needed to store the mesh/sol names

◆ _MMG5_NULKAL

#define _MMG5_NULKAL   1.e-30

◆ _MMG5_SAFE_CALLOC

#define _MMG5_SAFE_CALLOC (   ptr,
  size,
  type,
  retval 
)
Value:
do \
{ \
ptr = (type *)calloc((size),sizeof(type)); \
if ( !ptr ) { \
perror(" ## Memory problem: calloc"); \
return retval; \
} \
}while(0)
! int size
Definition: libmmgtypesf.h:583
! int type
Definition: libmmgtypesf.h:529
MMG5_pMesh char int int * retval
Definition: API_functionsf_s.c:584

Safe allocation with calloc

◆ _MMG5_SAFE_FREE

#define _MMG5_SAFE_FREE (   ptr)
Value:
do \
{ \
free(ptr); \
ptr = NULL; \
}while(0)

Safe deallocation

◆ _MMG5_SAFE_MALLOC

#define _MMG5_SAFE_MALLOC (   ptr,
  size,
  type,
  retval 
)
Value:
do \
{ \
ptr = (type *)malloc((size)*sizeof(type)); \
if ( !ptr ) { \
perror(" ## Memory problem: malloc"); \
return retval; \
} \
}while(0)
! int size
Definition: libmmgtypesf.h:583
! int type
Definition: libmmgtypesf.h:529
MMG5_pMesh char int int * retval
Definition: API_functionsf_s.c:584

Safe allocation with malloc

◆ _MMG5_SAFE_REALLOC

#define _MMG5_SAFE_REALLOC (   ptr,
  size,
  type,
  message,
  retval 
)
Value:
do \
{ \
type* tmp; \
tmp = (type *)realloc((ptr),(size)*sizeof(type)); \
if ( !tmp ) { \
_MMG5_SAFE_FREE(ptr); \
perror(" ## Memory problem: realloc"); \
return retval; \
} \
\
(ptr) = tmp; \
}while(0)
! int size
Definition: libmmgtypesf.h:583
! int type
Definition: libmmgtypesf.h:529
MMG5_pMesh char int int * retval
Definition: API_functionsf_s.c:584
! int tmp
Definition: libmmgtypesf.h:248

Safe reallocation

◆ _MMG5_SAFE_RECALLOC

#define _MMG5_SAFE_RECALLOC (   ptr,
  prevSize,
  newSize,
  type,
  message,
  retval 
)
Value:
do \
{ \
type* tmp; \
\
tmp = (type *)realloc((ptr),(newSize)*sizeof(type)); \
if ( !tmp ) { \
_MMG5_SAFE_FREE(ptr); \
perror(" ## Memory problem: realloc"); \
return retval; \
} \
\
(ptr) = tmp; \
if ( newSize > prevSize ) { \
memset(&((ptr)[prevSize]),0,((newSize)-(prevSize))*sizeof(type)); \
} \
}while(0)
! int type
Definition: libmmgtypesf.h:529
if(!ier) exit(EXIT_FAILURE)
MMG5_pMesh char int int * retval
Definition: API_functionsf_s.c:584
! int tmp
Definition: libmmgtypesf.h:248

safe reallocation with memset at 0 for the new values of tab

◆ _MMG5_SAFELL2ICAST

#define _MMG5_SAFELL2ICAST (   longlongval)    (((longlongval) > (INT_MAX)) ? 0 : ((int)(longlongval)))

◆ _MMG5_SAFELL2LCAST

#define _MMG5_SAFELL2LCAST (   longlongval)    (((longlongval) > (LONG_MAX)) ? 0 : ((long)(longlongval)))

◆ _MMG5_SQR32

#define _MMG5_SQR32   0.866025403784439

◆ _MMG5_TAB_RECALLOC

#define _MMG5_TAB_RECALLOC (   mesh,
  ptr,
  initSize,
  wantedGap,
  type,
  message,
  law,
  retval 
)
Value:
do \
{ \
int gap; \
(long long) (wantedGap*initSize*sizeof(type)) ) { \
gap = (int)((mesh->memMax-mesh->memCur)/sizeof(type)); \
if(gap<1) { \
fprintf(stderr," ## Error:"); \
fprintf(stderr," unable to allocate %s.\n",message); \
fprintf(stderr," ## Check the mesh size or "); \
fprintf(stderr,"increase maximal authorized memory with the -m option.\n"); \
law; \
} \
} \
else \
gap = (int)(wantedGap*initSize); \
_MMG5_ADD_MEM(mesh,gap*sizeof(type),message,law); \
_MMG5_SAFE_RECALLOC((ptr),initSize+1,initSize+gap+1,type,message,retval); \
initSize = initSize+gap; \
}while(0);
long long memCur
Definition: libmmgtypes.h:475
! int type
Definition: libmmgtypesf.h:529
if(!ier) exit(EXIT_FAILURE)
! double gap
Definition: libmmgtypesf.h:526
MMG5_pMesh char int int * retval
Definition: API_functionsf_s.c:584
long long memMax
Definition: libmmgtypes.h:474
MMG5_pMesh * mesh
Definition: API_functionsf_s.c:63
#define _MMG5_ADD_MEM(mesh, size, message, law)
Definition: mmgcommon.h:185

Reallocation of ptr of type type at size (initSize+wantedGap*initSize) if possible or at maximum available size if not. Execute the command law if reallocation failed. Memset to 0 for the new values of table.

◆ A16TH

#define A16TH   0.0625

◆ A32TH

#define A32TH   0.03125

◆ A64TH

#define A64TH   0.015625

◆ FORTRAN_NAME

#define FORTRAN_NAME (   nu,
  nl,
  pl,
  pc 
)
Value:
void nu pl; \
void nl pl \
{ nu pc; } \
void nl##_ pl \
{ nu pc; } \
void nl##__ pl \
{ nu pc; } \
void nu pl

Adds function definitions.

Parameters
nufunction name in upper case.
nlfunction name in lower case.
pltype of arguments.
pcname of arguments.
Note
Macro coming from Scotch library.

Adds function definitions with upcase, underscore and double underscore to match any fortran compiler.

◆ FORTRAN_VARIADIC

#define FORTRAN_VARIADIC (   nu,
  nl,
  pl,
  body 
)
Value:
void nu pl \
{ body } \
void nl pl \
{ body } \
void nl##_ pl \
{ body } \
void nl##__ pl \
{ body } \

Adds function definitions.

Parameters
nufunction name in upper case.
nlfunction name in lower case.
pltype of arguments.
bodybody of the function.

Adds function definitions with upcase, underscore and double underscore to match any fortran compiler.

◆ GNU

#define GNU

◆ MG_BDY

#define MG_BDY   (1 << 4)

16 boundary entity

◆ MG_CLR

#define MG_CLR (   flag,
  bit 
)    ((flag) &= ~(1 << (bit)))

bit number bit is set to 0

◆ MG_CPY

#define MG_CPY   "Copyright (c) IMB-LJLL, 2004-"

◆ MG_CRN

#define MG_CRN   (1 << 5)

32 corner

◆ MG_EDG

#define MG_EDG (   tag)    ((tag & MG_GEO) || (tag & MG_REF))

Edge or Ridge

◆ MG_EOK

#define MG_EOK (   pt)    (pt && ((pt)->v[0] > 0))

Element OK

◆ MG_GEO

#define MG_GEO   (1 << 1)

2 geometric ridge

◆ MG_GET

#define MG_GET (   flag,
  bit 
)    ((flag) & (1 << (bit)))

return bit number bit value

◆ MG_MAX

#define MG_MAX (   a,
  b 
)    (((a) > (b)) ? (a) : (b))

◆ MG_MIN

#define MG_MIN (   a,
  b 
)    (((a) < (b)) ? (a) : (b))

◆ MG_MINUS

#define MG_MINUS   3

◆ MG_NOM

#define MG_NOM   (1 << 3)

8 non manifold

◆ MG_NOSURF

#define MG_NOSURF   (1 << 6)

64 freezed boundary

◆ MG_NOTAG

#define MG_NOTAG   (0)

◆ MG_NUL

#define MG_NUL   (1 << 14)

16384 vertex removed

◆ MG_OPNBDY

#define MG_OPNBDY   (1 << 7)

128 open boundary

◆ MG_PARBDY

#define MG_PARBDY   (1 << 13)

8192 parallel boundary

◆ MG_PLUS

#define MG_PLUS   2

◆ MG_REF

#define MG_REF   (1 << 0)

1 edge reference

◆ MG_REL

#define MG_REL   "Apr. 10, 2017"

◆ MG_REQ

#define MG_REQ   (1 << 2)

4 required entity

◆ MG_SET

#define MG_SET (   flag,
  bit 
)    ((flag) |= (1 << (bit)))

bit number bit is set to 1

◆ MG_SIN

#define MG_SIN (   tag)    ((tag & MG_CRN) || (tag & MG_REQ))

Corner or Required

◆ MG_SMSGN

#define MG_SMSGN (   a,
  b 
)    (((double)(a)*(double)(b) > (0.0)) ? (1) : (0))

Check if a and b have the same sign

◆ MG_STR

#define MG_STR   "&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&"

◆ MG_Tetra

#define MG_Tetra   (1 << 2 )

4 local parameter applied over tetrahedron

◆ MG_Tria

#define MG_Tria   (1 << 1 )

2 local parameter applied over triangle

◆ MG_VER

#define MG_VER   "5.3.11"

◆ MG_Vert

#define MG_Vert   (1 << 0 )

1 local parameter applied over vertex

◆ MG_VOK

#define MG_VOK (   ppt)    (ppt && ((ppt)->tag < MG_NUL))

Vertex OK

◆ POSIX

#define POSIX

Typedef Documentation

◆ _MMG5_dNode

typedef struct _MMG5_dNode_s _MMG5_dNode

◆ _MMG5_iNode

typedef struct _MMG5_iNode_s _MMG5_iNode

◆ _MMG5_pBezier

Function Documentation

◆ _MMG5_Add_dnode()

int _MMG5_Add_dnode ( MMG5_pMesh  mesh,
_MMG5_dNode **  liLi,
int  k,
double  val 
)
inline
Parameters
meshpointer toward the mesh structure (for count of used memory).
liLipointer toward the address of the root of the linked list.
kinteger value to add to the linked list.
valreal value to add to the linked list.
Returns
1 if the node is inserted, 0 if the node is not inserted, -1 if fail.

Add a node with integer value k and real value val to a sorted linked list with unique entries.

Remarks
as the linked list had unique entries, we don't insert a node if it exists.
Here is the call graph for this function:

◆ _MMG5_Add_inode()

int _MMG5_Add_inode ( MMG5_pMesh  mesh,
_MMG5_iNode **  liLi,
int  val 
)
inline
Parameters
meshpointer toward the mesh structure (for count of used memory).
liLipointer toward the address of the root of the linked list.
valvalue to add to the linked list.
Returns
1 if the node is inserted, 0 if the node is not inserted, -1 if fail.

Add a node with value val to a sorted linked list with unique entries.

Remarks
as the linked list had unique entries, we don't insert a node if it exists.
Here is the call graph for this function:
Here is the caller graph for this function:

◆ _MMG5_Alloc_dnode()

int _MMG5_Alloc_dnode ( MMG5_pMesh  mesh,
_MMG5_dNode **  node 
)
inline
Parameters
meshpointer toward the mesh structure (for count of used memory).
nodepointer toward a _MMG5_dNode (cell for linked list)
Returns
1 if we can alloc the node node, 0 otherwise.

Node allocation.

Here is the caller graph for this function:

◆ _MMG5_Alloc_inode()

int _MMG5_Alloc_inode ( MMG5_pMesh  mesh,
_MMG5_iNode **  node 
)
inline
Parameters
meshpointer toward the mesh structure (for count of used memory).
nodepointer toward a _MMG5_iNode (cell for linked list)
Returns
1 if we can alloc the node node, 0 otherwise.

Node allocation.

Here is the caller graph for this function:

◆ _MMG5_bezierEdge()

void _MMG5_bezierEdge ( MMG5_pMesh  ,
int  ,
int  ,
double *  ,
double *  ,
char  ,
double *   
)

◆ _MMG5_boulec()

int _MMG5_boulec ( MMG5_pMesh  mesh,
int *  adjt,
int  start,
int  ip,
double *  tt 
)
Parameters
meshpointer toward the mesh structure.
adjtpointer toward the table of triangle adjacency.
startindex of triangle where we start to work.
ipindex of vertex where the tangent is computed.
ttpointer toward the computed tangent.
Returns
0 if fail, 1 otherwise.

Compute the tangent to the curve at point ip.

Here is the caller graph for this function:

◆ _MMG5_boulen()

int _MMG5_boulen ( MMG5_pMesh  mesh,
int *  adjt,
int  start,
int  ip,
double *  nn 
)
Parameters
meshpointer toward the mesh structure.
adjtpointer toward the table of triangle adjacency.
startindex of triangle where we start to work.
ipindex of vertex where the normal is computed.
nnpointer toward the computed tangent.
Returns
0 if fail, 1 otherwise.

Compute average normal of triangles sharing P without crossing ridge.

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

◆ _MMG5_bouler()

int _MMG5_bouler ( MMG5_pMesh  mesh,
int *  adjt,
int  start,
int  ip,
int *  list,
int *  listref,
int *  ng,
int *  nr,
int  lmax 
)
Parameters
meshpointer toward the mesh structure.
adjtpointer toward the table of triangle adjacency.
startindex of triangle where we start to work.
ipindex of vertex on which we work.
listpointer toward the computed list of GEO vertices incident to ip.
listrefpointer toward the corresponding edge references
ngpointer toward the number of ridges.
nrpointer toward the number of reference edges.
lmaxmaxmum size for the ball of the point ip.
Returns
The number of edges incident to the vertex ip.

Store edges and count the number of ridges and reference edges incident to the vertex ip.

Here is the caller graph for this function:

◆ _MMG5_boundingBox()

int _MMG5_boundingBox ( MMG5_pMesh  mesh)
Parameters
meshpointer toward the mesh structure.
Returns
1 if success, 0 if fail (computed bounding box too small).

Compute the mesh bounding box and fill the min, max and delta fields of the _MMG5_info structure.

Here is the caller graph for this function:

◆ _MMG5_buildridmet()

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

◆ _MMG5_buildridmetfic()

int _MMG5_buildridmetfic ( MMG5_pMesh  ,
double *  ,
double *  ,
double  ,
double  ,
double  ,
double *   
)

◆ _MMG5_buildridmetnor()

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

◆ _MMG5_caltri33_ani()

double _MMG5_caltri33_ani ( MMG5_pMesh  mesh,
MMG5_pSol  met,
MMG5_pTria  pt 
)
Parameters
meshpointer toward the mesh structure.
metpointer toward the meric structure.
ptpointer toward the triangle structure.
Returns
The computed quality.

Compute the quality of the surface triangle ptt with respect to an anisotropic metric and a classic storage of the ridges metrics.

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

◆ _MMG5_caltri_ani()

double _MMG5_caltri_ani ( MMG5_pMesh  mesh,
MMG5_pSol  met,
MMG5_pTria  ptt 
)
inline
Parameters
meshpointer toward the mesh structure.
metpointer toward the meric structure.
pttpointer toward the triangle structure.
Returns
The computed quality.

Compute the quality of the surface triangle ptt with respect to an anisotropic metric.

Warning
The quality is computed as if the triangle is a "straight" triangle.
Here is the call graph for this function:
Here is the caller graph for this function:

◆ _MMG5_caltri_iso()

double _MMG5_caltri_iso ( MMG5_pMesh  mesh,
MMG5_pSol  met,
MMG5_pTria  ptt 
)
inline
Parameters
meshpointer toward the mesh structure.
metpointer toward the meric structure.
pttpointer toward the triangle structure.
Returns
The computed quality.

Compute the quality of the surface triangle ptt with respect to an isotropic metric.

Here is the caller graph for this function:

◆ _MMG5_countLocalParamAtTri()

int _MMG5_countLocalParamAtTri ( MMG5_pMesh  mesh,
_MMG5_iNode **  bdryRefs 
)
inline
Parameters
meshpointer toward the mesh structure.
bdryRefspointer toward the list of the boundary references.
Returns
npar, the number of local parameters at triangles if success, 0 otherwise.

Count the local default values at triangles and fill the list of the boundary references.

Count the number of different boundary references and list it

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

◆ _MMG5_defUninitSize()

void _MMG5_defUninitSize ( MMG5_pMesh  mesh,
MMG5_pSol  met,
char  ismet 
)
Parameters
meshpointer toward the mesh structure.
metpointer toward the metric structure.
ismethas the user provided a metric?

Search for points with unintialized metric and define anisotropic size at this points.

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

◆ _MMG5_det3pt1vec()

double _MMG5_det3pt1vec ( double  c0[3],
double  c1[3],
double  c2[3],
double  v[3] 
)
inline

Compute 3 * 3 determinant : det(c1-c0,c2-c0,v)

Here is the caller graph for this function:

◆ _MMG5_det4pt()

double _MMG5_det4pt ( double  c0[3],
double  c1[3],
double  c2[3],
double  c3[3] 
)
inline

Compute 3 * 3 determinant : det(c1-c0,c2-c0,c3-c0)

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

◆ _MMG5_devangle()

int _MMG5_devangle ( double *  n1,
double *  n2,
double  crit 
)
Parameters
n1first normal
n2second normal
Returns
1 if success, 0 if fail

Check if the angle between n1 and n2 is larger than the ridge criterion. If yes, return 1, 0 otherwise (ridge creation).

◆ _MMG5_displayHisto()

void _MMG5_displayHisto ( MMG5_pMesh  mesh,
int  ned,
double *  avlen,
int  amin,
int  bmin,
double  lmin,
int  amax,
int  bmax,
double  lmax,
int  nullEdge,
double *  bd,
int *  hl,
char  shift 
)
Parameters
meshpointer toward the mesh structure.
nededges number.
avlenpointer toward the average edges lengths.
aminindex of first extremity of the smallest edge.
bminindex of second extremity of the smallest edge.
lminsmallest edge length.
amaxindex of first extremity of the largest edge.
bmaxindex of second extremity of the largest edge.
lmaxlargest edge length.
nullEdgenumber of edges for which we are unable to compute the length
bdpointer toward the table of the quality span.
hlpointer toward the table that store the number of edges for eac
shiftvalue to shift the target lenght interval span of quality

Display histogram of edge length.

Here is the caller graph for this function:

◆ _MMG5_elementWeight()

int _MMG5_elementWeight ( MMG5_pMesh  mesh,
MMG5_pSol  met,
MMG5_pTria  pt,
MMG5_pPoint  p0,
_MMG5_Bezier pb,
double  r[3][3],
double  gv[2] 
)
Parameters
meshpointer toward the mesh structure.
metpointer toward the metric structure.
ptpointer toward the tria on which we integrate.
p0pointer toward the point that we want to move.
pbbezier patch of the triangle.
rrotation matrix that sends the normal at point p0 to e_z.
gvcentre of mass that we want to update using the computed element weight.
Returns
0 if fail, 1 otherwise.

Compute integral of sqrt(T^J(xi) M(P(xi)) J(xi)) * P(xi) over the triangle.

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

◆ _MMG5_excfun()

static void _MMG5_excfun ( int  sigid)
inlinestatic

Inlined functions for libraries and executables

Parameters
sigidsignal number.

Signal handling: specify error messages depending from catched signal.

Here is the caller graph for this function:

◆ _MMG5_fillDefmetregSys()

void _MMG5_fillDefmetregSys ( int  ,
MMG5_pPoint  ,
int  ,
_MMG5_Bezier  ,
double  r[3][3],
double *  ,
double *  ,
double *  ,
double *   
)

◆ _MMG5_Free_dlinkedList()

void _MMG5_Free_dlinkedList ( MMG5_pMesh  mesh,
_MMG5_dNode liLi 
)
inline
Parameters
meshpointer toward the mesh structure (for count of used memory).
liLipointer toward the root of the linked list.

Free the memory used by the linked list whose root is liLi.

◆ _MMG5_Free_ilinkedList()

void _MMG5_Free_ilinkedList ( MMG5_pMesh  mesh,
_MMG5_iNode liLi 
)
inline
Parameters
meshpointer toward the mesh structure (for count of used memory).
liLipointer toward the root of the linked list.

Free the memory used by the linked list whose root is liLi.

Here is the caller graph for this function:

◆ _MMG5_grad2metSurf()

int _MMG5_grad2metSurf ( MMG5_pMesh  mesh,
MMG5_pSol  met,
MMG5_pTria  pt,
int  i 
)
Parameters
meshpointer toward the mesh.
metpointer toward the metric structure.
ptpointer toward a triangle.
iedge index in triangle pt.
Returns
-1 if no gradation is needed, else index of graded point.

Enforces gradation of metric in one extremity of edge i in tria pt with respect to the other, along the direction of the associated support curve first, then along the normal direction.

Warning
The gradation along the direction normal to the surface is made in an "isotropic way".
Here is the call graph for this function:
Here is the caller graph for this function:

◆ _MMG5_hashEdge()

int _MMG5_hashEdge ( MMG5_pMesh  mesh,
_MMG5_Hash hash,
int  a,
int  b,
int  k 
)
Parameters
meshpointer toward the mesh structure.
hashpointer toward the hash table of edges.
aindex of the first extremity of the edge.
bindex of the second extremity of the edge.
kindex of point along the edge.
Returns
1 if success, 0 if fail.

Add edge $[a;b]$ to the hash table.

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

◆ _MMG5_hashGet()

int _MMG5_hashGet ( _MMG5_Hash hash,
int  a,
int  b 
)
Parameters
hashpointer toward the hash table of edges.
aindex of the first extremity of the edge.
bindex of the second extremity of the edge.
Returns
the index of point stored along $[a;b]$.

Find the index of point stored along $[a;b]$.

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

◆ _MMG5_hashNew()

int _MMG5_hashNew ( MMG5_pMesh  mesh,
_MMG5_Hash hash,
int  hsiz,
int  hmax 
)
Parameters
meshpointer toward the mesh structure.
hashpointer toward the hash table of edges.
hsizinitial size of hash table.
hmaxmaximal size of hash table.
Returns
1 if success, 0 if fail.

Hash edges or faces.

Here is the caller graph for this function:

◆ _MMG5_hashUpdate()

int _MMG5_hashUpdate ( _MMG5_Hash hash,
int  a,
int  b,
int  k 
)
Parameters
meshpointer toward the mesh structure.
hashpointer toward the hash table of edges.
aindex of the first extremity of the edge.
bindex of the second extremity of the edge.
knew index of point along the edge.
Returns
1 if success, 0 if fail (edge is not found).

Update the index of the point stored along the edge $[a;b]$

Here is the caller graph for this function:

◆ _MMG5_interp_iso()

int _MMG5_interp_iso ( double *  ma,
double *  mb,
double *  mp,
double  t 
)
Parameters
mapointer on a metric
mbpointer on a metric
mppointer on the computed interpolated metric
tinterpolation parameter (comprise between 0 and 1)

Linear interpolation of isotropic sizemap along an edge

◆ _MMG5_interpreg_ani()

int _MMG5_interpreg_ani ( MMG5_pMesh  ,
MMG5_pSol  ,
MMG5_pTria  ,
char  ,
double  ,
double *  mr 
)

◆ _MMG5_intersecmet22()

int _MMG5_intersecmet22 ( MMG5_pMesh  mesh,
double *  m,
double *  n,
double *  mr 
)
Parameters
meshpointer toward the mesh structure.
mpointer toward a $(2x2)$ metric.
npointer toward a $(2x2)$ metric.
mrcomputed $(2x2)$ metric.
Returns
0 if fail, 1 otherwise.

Compute the intersected (2 x 2) metric from metrics m and n : take simultaneous reduction, and proceed to truncation in sizes.

Here is the caller graph for this function:

◆ _MMG5_intmetsavedir()

int _MMG5_intmetsavedir ( MMG5_pMesh  mesh,
double *  m,
double *  n,
double *  mr 
)
Parameters
meshpointer toward the mesh structure.
mpointer toward the first metric to intersect.
npointer toward the second metric to intersect.
mrpointer toward the computed intersected metric.
Returns
1.

Compute the intersected (2 x 2) metric between metrics m and n, PRESERVING the directions of m. Result is stored in mr.

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_invmat()

int _MMG5_invmat ( double *  m,
double *  mi 
)
Parameters
mpointer toward a 3x3 symetric matrix
mipointer toward the computed 3x3 matrix.

Invert m (3x3 symetric matrix) and store the result on mi

◆ _MMG5_invmatg()

int _MMG5_invmatg ( double  m[9],
double  mi[9] 
)
Parameters
minitial matrix.
miinverted matrix.

Invert 3x3 non-symmetric matrix.

Here is the caller graph for this function:

◆ _MMG5_memSize()

long long _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_minQualCheck()

int _MMG5_minQualCheck ( int  iel,
double  minqual,
double  alpha 
)
Parameters
ielindex of the worst tetra of the mesh
minqualquality of the worst tetra of the mesh (normalized by alpha)
alphanormalisation parameter for the quality
Returns
1 if success, 0 if fail (the quality is lower than _MMG5_NULKAL).

Print warning or error messages depending on the quality of the worst tetra of the mesh.

Here is the caller graph for this function:

◆ _MMG5_mmgDefaultValues()

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

Print the default parameters values.

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

◆ _MMG5_mmgHashTria()

int _MMG5_mmgHashTria ( MMG5_pMesh  mesh,
int *  adjt,
_MMG5_Hash hash,
int  chkISO 
)
Parameters
meshpointer toward the mesh structure.
adjtpointer toward the adjacency table of the surfacic mesh.
hashpointer toward the edge hash table.
chkISOflag to say if we check ISO references (so if we come from mmg3d).
Returns
1 if success, 0 otherwise.

Create surface adjacency

Remarks
the ph->s field computation is useless in mmgs.
Here is the call graph for this function:
Here is the caller graph for this function:

◆ _MMG5_mmgInit_parameters()

void _MMG5_mmgInit_parameters ( MMG5_pMesh  mesh)

◆ _MMG5_mmgIntextmet()

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

◆ _MMG5_mmgIntmet33_ani()

int _MMG5_mmgIntmet33_ani ( double *  m,
double *  n,
double *  mr,
double  s 
)
Parameters
minput metric.
ninput metric.
mrcomputed output metric.
sparameter coordinate for the interpolation of metrics m and n.
Returns
0 if fail, 1 otherwise.

Compute the interpolated $(3 x 3)$ metric from metrics m and n, at parameter s : $ mr = (1-s)*m +s*n $, both metrics being expressed in the simultaneous reduction basis: linear interpolation of sizes.

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

◆ _MMG5_mmgUsage()

void _MMG5_mmgUsage ( char *  prog)
Parameters
*progpointer toward the program name.

Print help for common options of mmg3d and mmgs.

Here is the caller graph for this function:

◆ _MMG5_nonorsurf()

double _MMG5_nonorsurf ( MMG5_pMesh  mesh,
MMG5_pTria  pt 
)
inline
Parameters
meshpointer toward the mesh stucture.
pttriangle for which we compute the surface.
Returns
the computed surface

Compute non-oriented surface area of a triangle.

Here is the call graph for this function:

◆ _MMG5_nonUnitNorPts()

int _MMG5_nonUnitNorPts ( MMG5_pMesh  mesh,
int  ip1,
int  ip2,
int  ip3,
double *  n 
)
inline
Parameters
meshpointer toward the mesh stucture.
ip1first point of face.
ip2second point of face.
ip3third point of face.
npointer to store the computed normal.
Returns
1

Compute non-normalized face normal given three points on the surface.

Here is the caller graph for this function:

◆ _MMG5_norpts()

int _MMG5_norpts ( MMG5_pMesh  mesh,
int  ip1,
int  ip2,
int  ip3,
double *  n 
)
inline
Parameters
meshpointer toward the mesh stucture.
ip1first point of face.
ip2second point of face.
ip3third point of face.
npointer to store the computed normal.
Returns
1 if success, 0 otherwise.

Compute normalized face normal given three points on the surface.

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

◆ _MMG5_nortri()

int _MMG5_nortri ( MMG5_pMesh  mesh,
MMG5_pTria  pt,
double *  n 
)
inline
Parameters
meshpointer toward the mesh stucture.
ptpointer toward the triangle structure.
npointer to store the computed normal.
Returns
1

Compute triangle normal.

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

◆ _MMG5_orvol()

double _MMG5_orvol ( MMG5_pPoint  point,
int *  v 
)
inline
Parameters
pointPointer toward the points array
vpointer toward the point indices
Returns
the oriented volume of tetra

Compute oriented volume of a tetrahedron

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

◆ _MMG5_paratmet()

int _MMG5_paratmet ( double  c0[3],
double  n0[3],
double  m[6],
double  c1[3],
double  n1[3],
double  mt[6] 
)
Parameters
c0table of the coordinates of the starting point.
n0normal at the starting point.
mmetric to be transported.
c1table of the coordinates of the ending point.
n1normal at the ending point.
mtcomputed metric.
Returns
0 if fail, 1 otherwise.

Parallel transport of a metric tensor field, attached to point c0, with normal n0, to point c1, with normal n1.

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

◆ _MMG5_printTria()

void _MMG5_printTria ( MMG5_pMesh  mesh,
char *  fileName 
)
Parameters
meshpointer toward the mesh structure.
fileNamepointer toward the file name.

Debug function (not use in clean code): write mesh->tria structure in file.

◆ _MMG5_ridSizeInNormalDir()

double _MMG5_ridSizeInNormalDir ( MMG5_pMesh  ,
int  ,
double *  ,
_MMG5_pBezier  ,
double  ,
double   
)

◆ _MMG5_ridSizeInTangentDir()

double _MMG5_ridSizeInTangentDir ( MMG5_pMesh  mesh,
MMG5_pPoint  p0,
int  idp,
int *  iprid,
double  isqhmin,
double  isqhmax 
)
Parameters
meshpointer toward the mesh structure.
p0pointer toward the point at which we define the metric.
idpglobal index of the point at which we define the metric.
ipridpointer toward the two extremities of the ridge.
isqhminminimum edge size.
isqhmaxmaximum edge size.
Returns
the computed ridge size in the tangent direction.

Compute the specific size of a ridge in the direction of the tangent of the ridge.

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

◆ _MMG5_rmtr()

int _MMG5_rmtr ( double  r[3][3],
double  m[6],
double  mr[6] 
)
inline
Here is the caller graph for this function:

◆ _MMG5_rotmatrix()

int _MMG5_rotmatrix ( double  n[3],
double  r[3][3] 
)
inline
Parameters
npointer toward the vector that we want to send on the third vector of canonical basis.
rcomputed rotation matrix.

Compute rotation matrix that sends vector n to the third vector of canonical basis.

Here is the caller graph for this function:

◆ _MMG5_scaleMesh()

int _MMG5_scaleMesh ( MMG5_pMesh  mesh,
MMG5_pSol  met 
)
Parameters
meshpointer toward the mesh structure.
metpointer toward the metric or solution structure.
Returns
1 if success, 0 if fail (computed bounding box too small or one af the anisotropic input metric is not valid).

Scale the mesh and the size informations between 0 and 1. Compute a default value for the hmin/hmax parameters if needed.

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

◆ _MMG5_scotchCall()

int _MMG5_scotchCall ( MMG5_pMesh  mesh,
MMG5_pSol  met 
)
Parameters
meshpointer toward the mesh structure.
metpointer toward the solution structure.
Returns
0 if _MMG5_renumbering fail (non conformal mesh), 1 otherwise (renumerotation success of renumerotation fail but the mesh is still conformal).

Call scotch renumbering.

Here is the caller graph for this function:

◆ _MMG5_Set_commonFunc()

void _MMG5_Set_commonFunc ( )

◆ _MMG5_solveDefmetrefSys()

int _MMG5_solveDefmetrefSys ( MMG5_pMesh  ,
MMG5_pPoint  ,
int *  ,
double  r[3][3],
double *  ,
double *  ,
double *  ,
double *  ,
double  ,
double  ,
double   
)

◆ _MMG5_solveDefmetregSys()

int _MMG5_solveDefmetregSys ( MMG5_pMesh  ,
double  r[3][3],
double *  ,
double *  ,
double *  ,
double *  ,
double  ,
double  ,
double   
)

◆ _MMG5_surftri33_ani()

double _MMG5_surftri33_ani ( MMG5_pMesh  ,
MMG5_pTria  ,
double *  ,
double *  ,
double *   
)

◆ _MMG5_surftri_ani()

double _MMG5_surftri_ani ( MMG5_pMesh  mesh,
MMG5_pSol  met,
MMG5_pTria  ptt 
)
Parameters
meshpointer toward the mesh structure.
metpointer toward the metric structure.
pttpointer toward the triangle structure.
Returns
The double of the triangle area.

Compute the double of the area of the surface triangle ptt with respect to the anisotropic metric met (for special storage of ridges metrics).

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

◆ _MMG5_surftri_iso()

double _MMG5_surftri_iso ( MMG5_pMesh  mesh,
MMG5_pSol  met,
MMG5_pTria  ptt 
)
Parameters
meshpointer toward the mesh structure.
metpointer toward the meric structure.
pttpointer toward the triangle structure.
Returns
The computed area.

Compute the area of the surface triangle ptt with respect to the isotropic metric met.

◆ _MMG5_sys33sym()

int _MMG5_sys33sym ( double  a[6],
double  b[3],
double  r[3] 
)
inline
Parameters
amatrix to invert.
blast member.
rvector of unknowns.
Returns
0 if fail, 1 otherwise.

Solve $ 3\times 3$ symmetric system $ A . r = b $.

Here is the caller graph for this function:

◆ _MMG5_unscaleMesh()

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

Unscale the mesh and the size informations to their initial sizes.

Here is the caller graph for this function:

◆ _MMG5_writeLocalParamAtTri()

int _MMG5_writeLocalParamAtTri ( MMG5_pMesh  mesh,
_MMG5_iNode bdryRefs,
FILE *  out 
)
inline
Parameters
meshpointer toward the mesh structure.
bdryRefspointer toward the list of the boundary references.
nparnumber of local param at triangles.
outpointer toward the file in which to write.
Returns
1 if success, 0 otherwise.

Write the local default values at triangles in the parameter file.

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

◆ MMG2_quickarea()

double MMG2_quickarea ( double  a[2],
double  b[2],
double  c[2] 
)
Parameters
apoint coordinates
bpoint coor
cpoint coor

Compute tria area.

Here is the caller graph for this function:

◆ MMG5_chkMetricType()

int MMG5_chkMetricType ( MMG5_pMesh  mesh,
int *  type,
FILE *  inm 
)
Parameters
meshpointer toward the mesh structure.
typetype of the metric
inmmetric file
Returns
1 if success, -1 if fail

Check if the type of the metric is compatible with the remeshing mode. If not, deallocate the type array and close the metric file.

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

◆ MMG5_chooseOutputFormat()

void MMG5_chooseOutputFormat ( MMG5_pMesh  mesh,
int *  msh 
)
Parameters
meshpointer toward the mesh structure.
meshpointer toward the msh value.

Update the msh value if we detect that the user want to force output at Gmsh or Medit format.

Here is the caller graph for this function:

◆ MMG5_loadMshMesh_part1()

int MMG5_loadMshMesh_part1 ( MMG5_pMesh  mesh,
const char *  filename,
FILE **  inm,
long *  posNodes,
long *  posElts,
long **  posNodeData,
int *  bin,
int *  iswp,
int *  nelts,
int *  nsols 
)
Parameters
meshpointer toward the mesh
filenamepointer toward the name of file
inmpointer toward the file pointer
posNodespointer toward the position of nodes data in file
posEltspointer toward the position of elts data in file
posNodeDatapointer toward the list of the positions of data in file
bin1 if binary format
neltsnumber of elements in file
nsolnumber of data in file
Returns
1 if success, 0 if file is not found, -1 if fail.

Begin to read mesh at MSH file format. Read the mesh size informations.

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

◆ MMG5_loadMshMesh_part2()

int MMG5_loadMshMesh_part2 ( MMG5_pMesh  mesh,
MMG5_pSol sol,
FILE **  inm,
const long  posNodes,
const long  posElts,
const long *  posNodeData,
const int  bin,
const int  iswp,
const int  nelts 
)
Parameters
meshpointer toward the mesh
solpointer toward the solutions array
inmpointer toward the file pointer
posNodesposition of nodes data in file
posEltsposition of elts data in file
posNodeDataposition of solution data in file
bin1 if binary format
neltsnumber of elements in file
Returns
1 if success, 0 if fail.

End to read mesh and solution array at MSH file format after the mesh/solution array alloc.

Second step: read the nodes and elements
Read the solution at nodes

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

◆ MMG5_loadSolHeader()

int MMG5_loadSolHeader ( const char *  filename,
int  meshDim,
FILE **  inm,
int *  ver,
int *  bin,
int *  iswp,
int *  np,
int *  dim,
int *  nsols,
int **  type,
long *  posnp 
)
Parameters
filenamename of file.
meshDimmesh dimenson.
inmallocatable pointer toward the FILE structure
verfile version (1=simple precision, 2=double)
bin1 if the file is a binary
iswp1 or 0 depending on the endianness (binary only)
npnumber of solutions of each type
dimsolution dimension
nsolsnumber of solutions of different types in the file
typetype of solutions
posnppointer toward the position of the point list in the file
Returns
-1 data invalid or we fail, 0 no file, 1 ok.

Open the "filename" solution file and read the file header.

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

◆ MMG5_printMetStats()

void MMG5_printMetStats ( MMG5_pMesh  mesh,
MMG5_pSol  met 
)
Parameters
meshpointer toward the mesh structure.
metpointer toward the metric structure.

print metric statistics

Here is the caller graph for this function:

◆ MMG5_printSolStats()

void MMG5_printSolStats ( MMG5_pMesh  mesh,
MMG5_pSol sol 
)
Parameters
meshpointer toward the mesh structure.
solpointer toward the solutions array.

print solutions statistics

Here is the caller graph for this function:

◆ MMG5_readDoubleSol3D()

void MMG5_readDoubleSol3D ( MMG5_pSol  sol,
FILE *  inm,
int  bin,
int  iswp,
int  pos 
)
Parameters
solpointer toward an allocatable sol structure.
inmpointer toward the solution file
bin1 if binary file
iswpEndianess
indexof the readed solution

Read the solution value for vertex of index pos in double precision.

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

◆ MMG5_readFloatSol3D()

void MMG5_readFloatSol3D ( MMG5_pSol  sol,
FILE *  inm,
int  bin,
int  iswp,
int  pos 
)
Parameters
solpointer toward an allocatable sol structure.
inmpointer toward the solution file
bin1 if binary file
iswpEndianess
indexof the readed solution

Read the solution value for vertex of index pos in floating precision.

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

◆ MMG5_saveMshMesh()

int MMG5_saveMshMesh ( MMG5_pMesh  mesh,
MMG5_pSol sol,
const char *  filename,
int  metricData 
)
Parameters
meshpointer toward the mesh structure.
solpointer toward an array of solutions.
filenamename of file.
metricData1 if the data saved is a metric (if only 1 data)
Returns
0 if failed, 1 otherwise.

Write mesh and a list of solutions at MSH file format (.msh extension). Write binary file for .mshb extension.and ASCII for .msh one.

First step: Count the number of elements of each type


Second step: save the elements at following format: "idx type tagNumber tag0 tag1... v0_elt v1_elt..."

Write solution
Save the solution at following format: "idx sol"

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

◆ MMG5_saveSolHeader()

int MMG5_saveSolHeader ( MMG5_pMesh  mesh,
const char *  filename,
FILE **  inm,
int  ver,
int *  bin,
int  np,
int  dim,
int  nsols,
int *  type,
int *  size 
)
Parameters
meshpointer toward the mesh structure.
filenamename of file.
inmallocatable pointer toward the FILE structure.
verfile version (1=simple precision, 2=double).
bin1 if the file is a binary.
npnumber of solutions of each type.
dimsolution dimension.
nsolsnumber of solutions of different types in the file.
typetype of solutions.
sizesize of solutions.
Returns
0 if unable to open the file, 1 if success.

Open the "filename" solution file and read the file header.

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

◆ MMG5_writeDoubleSol3D()

void MMG5_writeDoubleSol3D ( MMG5_pMesh  mesh,
MMG5_pSol  sol,
FILE *  inm,
int  bin,
int  pos,
int  metricData 
)
Parameters
meshpointer toward the mesh structure
solpointer toward an allocatable sol structure.
inmpointer toward the solution file
bin1 if binary file
posof the writted solution
metricData1 if the data saved is a metric (if only 1 data)

Write the solution value for vertex of index pos in double precision.

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

Variable Documentation

◆ _MMG5_bezierCP

int(* _MMG5_bezierCP) (MMG5_pMesh,MMG5_Tria *, _MMG5_pBezier,char)

◆ _MMG5_chkmsh

int(* _MMG5_chkmsh) (MMG5_pMesh, int, int)

◆ _MMG5_indElt

int(* _MMG5_indElt) (MMG5_pMesh mesh, int kel)

◆ _MMG5_indPt

int(* _MMG5_indPt) (MMG5_pMesh mesh, int kp)

◆ _MMG5_inxt2

const unsigned char _MMG5_inxt2[6] = {1,2,0,1,2}
static

next vertex of triangle: {1,2,0}

◆ _MMG5_iprv2

const unsigned char _MMG5_iprv2[3] = {2,0,1}
static

previous vertex of triangle: {2,0,1}

◆ _MMG5_lenSurfEdg

double(* _MMG5_lenSurfEdg) (MMG5_pMesh mesh, MMG5_pSol sol,int,int, char)