mmg2d
|
#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"
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.8" |
#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) |
#define _LIBMMG5_RETURN | ( | mesh, | |
met, | |||
val | |||
) |
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)
Increment memory counter memCur and check if we don't overflow the maximum authorizied memory memMax.
#define _MMG5_ANGEDG 0.707106781186548 /*0.573576436351046 */ |
#define _MMG5_ANGLIM -0.999999 |
#define _MMG5_ATHIRD 0.333333333333333 |
#define _MMG5_BOXSIZE 500 |
size of box for renumbering with scotch.
Check if used memory overflow maximal authorized memory. Execute the command law if lack of memory.
Free pointer ptr of mesh structure and compute the new used memory. size is the size of the pointer
#define _MMG5_EPS 1.e-06 |
#define _MMG5_EPSD 1.e-30 |
#define _MMG5_EPSD2 1.0e-200 |
#define _MMG5_EPSOK 1.e-15 |
#define _MMG5_INCREASE_MEM_MESSAGE | ( | ) |
Error message when lack of memory
#define _MMG5_KA 7 |
Key for hash tables.
#define _MMG5_KB 11 |
Key for hash tables.
#define _MMG5_MEMMAX 800 |
Maximal memory used if available memory compitation fail.
#define _MMG5_MEMMIN 38 |
minimal memory needed to store the mesh/sol names
#define _MMG5_NULKAL 1.e-30 |
Safe allocation with calloc
#define _MMG5_SAFE_FREE | ( | ptr | ) |
Safe deallocation
Safe allocation with malloc
safe reallocation with memset at 0 for the new values of tab
#define _MMG5_SAFELL2ICAST | ( | longlongval | ) | (((longlongval) > (INT_MAX)) ? 0 : ((int)(longlongval))) |
#define _MMG5_SAFELL2LCAST | ( | longlongval | ) | (((longlongval) > (LONG_MAX)) ? 0 : ((long)(longlongval))) |
#define _MMG5_SQR32 0.866025403784439 |
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.
#define A16TH 0.0625 |
#define A32TH 0.03125 |
#define A64TH 0.015625 |
#define FORTRAN_NAME | ( | nu, | |
nl, | |||
pl, | |||
pc | |||
) |
Adds function definitions.
nu | function name in upper case. |
nl | function name in lower case. |
pl | type of arguments. |
pc | name of arguments. |
Adds function definitions with upcase, underscore and double underscore to match any fortran compiler.
#define FORTRAN_VARIADIC | ( | nu, | |
nl, | |||
pl, | |||
body | |||
) |
Adds function definitions.
nu | function name in upper case. |
nl | function name in lower case. |
pl | type of arguments. |
body | body of the function. |
Adds function definitions with upcase, underscore and double underscore to match any fortran compiler.
#define GNU |
#define MG_BDY (1 << 4) |
16 boundary entity
#define MG_CPY "Copyright (c) IMB-LJLL, 2004-" |
#define MG_CRN (1 << 5) |
32 corner
#define MG_EOK | ( | pt | ) | (pt && ((pt)->v[0] > 0)) |
Element OK
#define MG_GEO (1 << 1) |
2 geometric ridge
#define MG_MINUS 3 |
#define MG_NOM (1 << 3) |
8 non manifold
#define MG_NOSURF (1 << 6) |
64 freezed boundary
#define MG_NOTAG (0) |
#define MG_NUL (1 << 14) |
16384 vertex removed
#define MG_OPNBDY (1 << 7) |
128 open boundary
#define MG_PARBDY (1 << 13) |
8192 parallel boundary
#define MG_PLUS 2 |
#define MG_REF (1 << 0) |
1 edge reference
#define MG_REL "Apr. 10, 2017" |
#define MG_REQ (1 << 2) |
4 required entity
Check if a and b have the same sign
#define MG_STR "&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&" |
#define MG_Tetra (1 << 2 ) |
4 local parameter applied over tetrahedron
#define MG_Tria (1 << 1 ) |
2 local parameter applied over triangle
#define MG_VER "5.3.8" |
#define MG_Vert (1 << 0 ) |
1 local parameter applied over vertex
#define POSIX |
typedef struct _MMG5_dNode_s _MMG5_dNode |
typedef struct _MMG5_iNode_s _MMG5_iNode |
typedef _MMG5_Bezier* _MMG5_pBezier |
|
inline |
mesh | pointer toward the mesh structure (for count of used memory). |
liLi | pointer toward the address of the root of the linked list. |
k | integer value to add to the linked list. |
val | real value to add to the linked list. |
Add a node with integer value k and real value val to a sorted linked list with unique entries.
|
inline |
mesh | pointer toward the mesh structure (for count of used memory). |
liLi | pointer toward the address of the root of the linked list. |
val | value to add to the linked list. |
Add a node with value val to a sorted linked list with unique entries.
|
inline |
mesh | pointer toward the mesh structure (for count of used memory). |
node | pointer toward a _MMG5_dNode (cell for linked list) |
Node allocation.
|
inline |
mesh | pointer toward the mesh structure (for count of used memory). |
node | pointer toward a _MMG5_iNode (cell for linked list) |
Node allocation.
void _MMG5_bezierEdge | ( | MMG5_pMesh | , |
int | , | ||
int | , | ||
double * | , | ||
double * | , | ||
char | , | ||
double * | |||
) |
int _MMG5_boulec | ( | MMG5_pMesh | mesh, |
int * | adjt, | ||
int | start, | ||
int | ip, | ||
double * | tt | ||
) |
mesh | pointer toward the mesh structure. |
adjt | pointer toward the table of triangle adjacency. |
start | index of triangle where we start to work. |
ip | index of vertex where the tangent is computed. |
tt | pointer toward the computed tangent. |
Compute the tangent to the curve at point ip.
int _MMG5_boulen | ( | MMG5_pMesh | mesh, |
int * | adjt, | ||
int | start, | ||
int | ip, | ||
double * | nn | ||
) |
mesh | pointer toward the mesh structure. |
adjt | pointer toward the table of triangle adjacency. |
start | index of triangle where we start to work. |
ip | index of vertex where the normal is computed. |
nn | pointer toward the computed tangent. |
Compute average normal of triangles sharing P without crossing ridge.
int _MMG5_bouler | ( | MMG5_pMesh | mesh, |
int * | adjt, | ||
int | start, | ||
int | ip, | ||
int * | list, | ||
int * | listref, | ||
int * | ng, | ||
int * | nr, | ||
int | lmax | ||
) |
mesh | pointer toward the mesh structure. |
adjt | pointer toward the table of triangle adjacency. |
start | index of triangle where we start to work. |
ip | index of vertex on which we work. |
list | pointer toward the computed list of GEO vertices incident to ip. |
listref | pointer toward the corresponding edge references |
ng | pointer toward the number of ridges. |
nr | pointer toward the number of reference edges. |
lmax | maxmum size for the ball of the point ip. |
Store edges and count the number of ridges and reference edges incident to the vertex ip.
int _MMG5_boundingBox | ( | MMG5_pMesh | mesh | ) |
mesh | pointer toward the mesh structure. |
Compute the mesh bounding box and fill the min, max and delta fields of the _MMG5_info structure.
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 * | |||
) |
double _MMG5_caltri33_ani | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met, | ||
MMG5_pTria | pt | ||
) |
mesh | pointer toward the mesh structure. |
met | pointer toward the meric structure. |
pt | pointer toward the triangle structure. |
Compute the quality of the surface triangle ptt with respect to an anisotropic metric and a classic storage of the ridges metrics.
|
inline |
mesh | pointer toward the mesh structure. |
met | pointer toward the meric structure. |
ptt | pointer toward the triangle structure. |
Compute the quality of the surface triangle ptt with respect to an anisotropic metric.
|
inline |
mesh | pointer toward the mesh structure. |
met | pointer toward the meric structure. |
ptt | pointer toward the triangle structure. |
Compute the quality of the surface triangle ptt with respect to an isotropic metric.
|
inline |
mesh | pointer toward the mesh structure. |
bdryRefs | pointer toward the list of the boundary references. |
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
void _MMG5_defUninitSize | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met, | ||
char | ismet | ||
) |
mesh | pointer toward the mesh structure. |
met | pointer toward the metric structure. |
ismet | has the user provided a metric? |
Search for points with unintialized metric and define anisotropic size at this points.
|
inline |
Compute 3 * 3 determinant : det(c1-c0,c2-c0,v)
|
inline |
Compute 3 * 3 determinant : det(c1-c0,c2-c0,c3-c0)
int _MMG5_devangle | ( | double * | n1, |
double * | n2, | ||
double | crit | ||
) |
n1 | first normal |
n2 | second normal |
Check if the angle between n1 and n2 is larger than the ridge criterion. If yes, return 1, 0 otherwise (ridge creation).
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 | ||
) |
mesh | pointer toward the mesh structure. |
ned | edges number. |
avlen | pointer toward the average edges lengths. |
amin | index of first extremity of the smallest edge. |
bmin | index of second extremity of the smallest edge. |
lmin | smallest edge length. |
amax | index of first extremity of the largest edge. |
bmax | index of second extremity of the largest edge. |
lmax | largest edge length. |
nullEdge | number of edges for which we are unable to compute the length |
bd | pointer toward the table of the quality span. |
hl | pointer toward the table that store the number of edges for eac |
shift | value to shift the target lenght interval span of quality |
Display histogram of edge length.
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] | ||
) |
mesh | pointer toward the mesh structure. |
met | pointer toward the metric structure. |
pt | pointer toward the tria on which we integrate. |
p0 | pointer toward the point that we want to move. |
pb | bezier patch of the triangle. |
r | rotation matrix that sends the normal at point p0 to e_z. |
gv | centre of mass that we want to update using the computed element weight. |
Compute integral of sqrt(T^J(xi) M(P(xi)) J(xi)) * P(xi) over the triangle.
|
inlinestatic |
Inlined functions for libraries and executables
sigid | signal number. |
Signal handling: specify error messages depending from catched signal.
void _MMG5_fillDefmetregSys | ( | int | , |
MMG5_pPoint | , | ||
int | , | ||
_MMG5_Bezier | , | ||
double | r[3][3], | ||
double * | , | ||
double * | , | ||
double * | , | ||
double * | |||
) |
|
inline |
mesh | pointer toward the mesh structure (for count of used memory). |
liLi | pointer toward the root of the linked list. |
Free the memory used by the linked list whose root is liLi.
|
inline |
mesh | pointer toward the mesh structure (for count of used memory). |
liLi | pointer toward the root of the linked list. |
Free the memory used by the linked list whose root is liLi.
int _MMG5_grad2metSurf | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met, | ||
MMG5_pTria | pt, | ||
int | i | ||
) |
mesh | pointer toward the mesh. |
met | pointer toward the metric structure. |
pt | pointer toward a triangle. |
i | edge index in triangle pt. |
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.
int _MMG5_hashEdge | ( | MMG5_pMesh | mesh, |
_MMG5_Hash * | hash, | ||
int | a, | ||
int | b, | ||
int | k | ||
) |
mesh | pointer toward the mesh structure. |
hash | pointer toward the hash table of edges. |
a | index of the first extremity of the edge. |
b | index of the second extremity of the edge. |
k | index of point along the edge. |
Add edge to the hash table.
int _MMG5_hashGet | ( | _MMG5_Hash * | hash, |
int | a, | ||
int | b | ||
) |
hash | pointer toward the hash table of edges. |
a | index of the first extremity of the edge. |
b | index of the second extremity of the edge. |
Find the index of point stored along .
int _MMG5_hashNew | ( | MMG5_pMesh | mesh, |
_MMG5_Hash * | hash, | ||
int | hsiz, | ||
int | hmax | ||
) |
mesh | pointer toward the mesh structure. |
hash | pointer toward the hash table of edges. |
hsiz | initial size of hash table. |
hmax | maximal size of hash table. |
Hash edges or faces.
int _MMG5_hashUpdate | ( | _MMG5_Hash * | hash, |
int | a, | ||
int | b, | ||
int | k | ||
) |
mesh | pointer toward the mesh structure. |
hash | pointer toward the hash table of edges. |
a | index of the first extremity of the edge. |
b | index of the second extremity of the edge. |
k | new index of point along the edge. |
Update the index of the point stored along the edge
int _MMG5_interp_iso | ( | double * | ma, |
double * | mb, | ||
double * | mp, | ||
double | t | ||
) |
ma | pointer on a metric |
mb | pointer on a metric |
mp | pointer on the computed interpolated metric |
t | interpolation parameter (comprise between 0 and 1) |
Linear interpolation of isotropic sizemap along an edge
int _MMG5_interpreg_ani | ( | MMG5_pMesh | , |
MMG5_pSol | , | ||
MMG5_pTria | , | ||
char | , | ||
double | , | ||
double * | mr | ||
) |
int _MMG5_intersecmet22 | ( | MMG5_pMesh | mesh, |
double * | m, | ||
double * | n, | ||
double * | mr | ||
) |
mesh | pointer toward the mesh structure. |
m | pointer toward a ![]() |
n | pointer toward a ![]() |
mr | computed ![]() |
Compute the intersected (2 x 2) metric from metrics m and n : take simultaneous reduction, and proceed to truncation in sizes.
int _MMG5_intmetsavedir | ( | MMG5_pMesh | mesh, |
double * | m, | ||
double * | n, | ||
double * | mr | ||
) |
mesh | pointer toward the mesh structure. |
m | pointer toward the first metric to intersect. |
n | pointer toward the second metric to intersect. |
mr | pointer toward the computed intersected metric. |
Compute the intersected (2 x 2) metric between metrics m and n, PRESERVING the directions of m. Result is stored in mr.
int _MMG5_intridmet | ( | MMG5_pMesh | , |
MMG5_pSol | , | ||
int | , | ||
int | , | ||
double | , | ||
double * | , | ||
double * | |||
) |
int _MMG5_invmat | ( | double * | m, |
double * | mi | ||
) |
m | pointer toward a 3x3 symetric matrix |
mi | pointer toward the computed 3x3 matrix. |
Invert m (3x3 symetric matrix) and store the result on mi
int _MMG5_invmatg | ( | double | m[9], |
double | mi[9] | ||
) |
m | initial matrix. |
mi | inverted matrix. |
Invert 3x3 non-symmetric matrix.
long long _MMG5_memSize | ( | void | ) |
Compute the available memory size of the computer.
int _MMG5_minQualCheck | ( | int | iel, |
double | minqual, | ||
double | alpha | ||
) |
iel | index of the worst tetra of the mesh |
minqual | quality of the worst tetra of the mesh (normalized by alpha) |
alpha | normalisation parameter for the quality |
Print warning or error messages depending on the quality of the worst tetra of the mesh.
void _MMG5_mmgDefaultValues | ( | MMG5_pMesh | mesh | ) |
mesh | pointer toward the mesh structure. |
Print the default parameters values.
int _MMG5_mmgHashTria | ( | MMG5_pMesh | mesh, |
int * | adjt, | ||
_MMG5_Hash * | hash, | ||
int | chkISO | ||
) |
mesh | pointer toward the mesh structure. |
adjt | pointer toward the adjacency table of the surfacic mesh. |
hash | pointer toward the edge hash table. |
chkISO | flag to say if we check ISO references (so if we come from mmg3d). |
Create surface adjacency
void _MMG5_mmgInit_parameters | ( | MMG5_pMesh | mesh | ) |
int _MMG5_mmgIntextmet | ( | MMG5_pMesh | , |
MMG5_pSol | , | ||
int | , | ||
double * | , | ||
double * | |||
) |
int _MMG5_mmgIntmet33_ani | ( | double * | m, |
double * | n, | ||
double * | mr, | ||
double | s | ||
) |
m | input metric. |
n | input metric. |
mr | computed output metric. |
s | parameter coordinate for the interpolation of metrics m and n. |
Compute the interpolated metric from metrics m and n, at parameter s :
, both metrics being expressed in the simultaneous reduction basis: linear interpolation of sizes.
void _MMG5_mmgUsage | ( | char * | prog | ) |
*prog | pointer toward the program name. |
Print help for common options of mmg3d and mmgs.
|
inline |
mesh | pointer toward the mesh stucture. |
pt | triangle for which we compute the surface. |
Compute non-oriented surface area of a triangle.
|
inline |
mesh | pointer toward the mesh stucture. |
ip1 | first point of face. |
ip2 | second point of face. |
ip3 | third point of face. |
n | pointer to store the computed normal. |
Compute non-normalized face normal given three points on the surface.
|
inline |
mesh | pointer toward the mesh stucture. |
ip1 | first point of face. |
ip2 | second point of face. |
ip3 | third point of face. |
n | pointer to store the computed normal. |
Compute normalized face normal given three points on the surface.
|
inline |
mesh | pointer toward the mesh stucture. |
pt | pointer toward the triangle structure. |
n | pointer to store the computed normal. |
Compute triangle normal.
|
inline |
point | Pointer toward the points array |
v | pointer toward the point indices |
Compute oriented volume of a tetrahedron
int _MMG5_paratmet | ( | double | c0[3], |
double | n0[3], | ||
double | m[6], | ||
double | c1[3], | ||
double | n1[3], | ||
double | mt[6] | ||
) |
c0 | table of the coordinates of the starting point. |
n0 | normal at the starting point. |
m | metric to be transported. |
c1 | table of the coordinates of the ending point. |
n1 | normal at the ending point. |
mt | computed metric. |
Parallel transport of a metric tensor field, attached to point c0, with normal n0, to point c1, with normal n1.
void _MMG5_printTria | ( | MMG5_pMesh | mesh, |
char * | fileName | ||
) |
mesh | pointer toward the mesh structure. |
fileName | pointer toward the file name. |
Debug function (not use in clean code): write mesh->tria structure in file.
double _MMG5_ridSizeInNormalDir | ( | MMG5_pMesh | , |
int | , | ||
double * | , | ||
_MMG5_pBezier | , | ||
double | , | ||
double | |||
) |
double _MMG5_ridSizeInTangentDir | ( | MMG5_pMesh | mesh, |
MMG5_pPoint | p0, | ||
int | idp, | ||
int * | iprid, | ||
double | isqhmin, | ||
double | isqhmax | ||
) |
mesh | pointer toward the mesh structure. |
p0 | pointer toward the point at which we define the metric. |
idp | global index of the point at which we define the metric. |
iprid | pointer toward the two extremities of the ridge. |
isqhmin | minimum edge size. |
isqhmax | maximum edge size. |
Compute the specific size of a ridge in the direction of the tangent of the ridge.
|
inline |
|
inline |
n | pointer toward the vector that we want to send on the third vector of canonical basis. |
r | computed rotation matrix. |
Compute rotation matrix that sends vector n to the third vector of canonical basis.
int _MMG5_scaleMesh | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met | ||
) |
mesh | pointer toward the mesh structure. |
met | pointer toward the metric or solution structure. |
Scale the mesh and the size informations between 0 and 1. Compute a default value for the hmin/hmax parameters if needed.
int _MMG5_scotchCall | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met | ||
) |
mesh | pointer toward the mesh structure. |
met | pointer toward the solution structure. |
Call scotch renumbering.
void _MMG5_Set_commonFunc | ( | ) |
int _MMG5_solveDefmetrefSys | ( | MMG5_pMesh | , |
MMG5_pPoint | , | ||
int * | , | ||
double | r[3][3], | ||
double * | , | ||
double * | , | ||
double * | , | ||
double * | , | ||
double | , | ||
double | , | ||
double | |||
) |
int _MMG5_solveDefmetregSys | ( | MMG5_pMesh | , |
double | r[3][3], | ||
double * | , | ||
double * | , | ||
double * | , | ||
double * | , | ||
double | , | ||
double | , | ||
double | |||
) |
double _MMG5_surftri33_ani | ( | MMG5_pMesh | , |
MMG5_pTria | , | ||
double * | , | ||
double * | , | ||
double * | |||
) |
double _MMG5_surftri_ani | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met, | ||
MMG5_pTria | ptt | ||
) |
mesh | pointer toward the mesh structure. |
met | pointer toward the metric structure. |
ptt | pointer toward the triangle structure. |
Compute the double of the area of the surface triangle ptt with respect to the anisotropic metric met (for special storage of ridges metrics).
double _MMG5_surftri_iso | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met, | ||
MMG5_pTria | ptt | ||
) |
mesh | pointer toward the mesh structure. |
met | pointer toward the meric structure. |
ptt | pointer toward the triangle structure. |
Compute the area of the surface triangle ptt with respect to the isotropic metric met.
|
inline |
a | matrix to invert. |
b | last member. |
r | vector of unknowns. |
Solve symmetric system
.
int _MMG5_unscaleMesh | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met | ||
) |
mesh | pointer toward the mesh structure. |
met | pointer toward the metric or solution structure. |
Unscale the mesh and the size informations to their initial sizes.
|
inline |
mesh | pointer toward the mesh structure. |
bdryRefs | pointer toward the list of the boundary references. |
npar | number of local param at triangles. |
out | pointer toward the file in which to write. |
Write the local default values at triangles in the parameter file.
double MMG2_quickarea | ( | double | a[2], |
double | b[2], | ||
double | c[2] | ||
) |
a | point coordinates |
b | point coor |
c | point coor |
Compute tria area.
int MMG5_chkMetricType | ( | MMG5_pMesh | mesh, |
int * | type, | ||
FILE * | inm | ||
) |
mesh | pointer toward the mesh structure. |
type | type of the metric |
inm | metric file |
Check if the type of the metric is compatible with the remeshing mode. If not, deallocate the type array and close the metric file.
void MMG5_chooseOutputFormat | ( | MMG5_pMesh | mesh, |
int * | msh | ||
) |
mesh | pointer toward the mesh structure. |
mesh | pointer toward the msh value. |
Update the msh value if we detect that the user want to force output at Gmsh or Medit format.
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 | ||
) |
mesh | pointer toward the mesh |
filename | pointer toward the name of file |
inm | pointer toward the file pointer |
posNodes | pointer toward the position of nodes data in file |
posElts | pointer toward the position of elts data in file |
posNodeData | pointer toward the list of the positions of data in file |
bin | 1 if binary format |
nelts | number of elements in file |
nsol | number of data in file |
Begin to read mesh at MSH file format. Read the mesh size informations.
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 | ||
) |
mesh | pointer toward the mesh |
sol | pointer toward the solutions array |
inm | pointer toward the file pointer |
posNodes | position of nodes data in file |
posElts | position of elts data in file |
posNodeData | position of solution data in file |
bin | 1 if binary format |
nelts | number of elements in file |
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
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 | ||
) |
filename | name of file. |
meshDim | mesh dimenson. |
inm | allocatable pointer toward the FILE structure |
ver | file version (1=simple precision, 2=double) |
bin | 1 if the file is a binary |
iswp | 1 or 0 depending on the endianness (binary only) |
np | number of solutions of each type |
dim | solution dimension |
nsols | number of solutions of different types in the file |
type | type of solutions |
posnp | pointer toward the position of the point list in the file |
Open the "filename" solution file and read the file header.
void MMG5_printMetStats | ( | MMG5_pMesh | mesh, |
MMG5_pSol | met | ||
) |
mesh | pointer toward the mesh structure. |
met | pointer toward the metric structure. |
print metric statistics
void MMG5_printSolStats | ( | MMG5_pMesh | mesh, |
MMG5_pSol * | sol | ||
) |
mesh | pointer toward the mesh structure. |
sol | pointer toward the solutions array. |
print solutions statistics
void MMG5_readDoubleSol3D | ( | MMG5_pSol | sol, |
FILE * | inm, | ||
int | bin, | ||
int | iswp, | ||
int | pos | ||
) |
sol | pointer toward an allocatable sol structure. |
inm | pointer toward the solution file |
bin | 1 if binary file |
iswp | Endianess |
index | of the readed solution |
Read the solution value for vertex of index pos in double precision.
void MMG5_readFloatSol3D | ( | MMG5_pSol | sol, |
FILE * | inm, | ||
int | bin, | ||
int | iswp, | ||
int | pos | ||
) |
sol | pointer toward an allocatable sol structure. |
inm | pointer toward the solution file |
bin | 1 if binary file |
iswp | Endianess |
index | of the readed solution |
Read the solution value for vertex of index pos in floating precision.
int MMG5_saveMshMesh | ( | MMG5_pMesh | mesh, |
MMG5_pSol * | sol, | ||
const char * | filename, | ||
int | metricData | ||
) |
mesh | pointer toward the mesh structure. |
sol | pointer toward an array of solutions. |
filename | name of file. |
metricData | 1 if the data saved is a metric (if only 1 data) |
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"
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 | ||
) |
mesh | pointer toward the mesh structure. |
filename | name of file. |
inm | allocatable pointer toward the FILE structure. |
ver | file version (1=simple precision, 2=double). |
bin | 1 if the file is a binary. |
np | number of solutions of each type. |
dim | solution dimension. |
nsols | number of solutions of different types in the file. |
type | type of solutions. |
size | size of solutions. |
Open the "filename" solution file and read the file header.
void MMG5_writeDoubleSol3D | ( | MMG5_pMesh | mesh, |
MMG5_pSol | sol, | ||
FILE * | inm, | ||
int | bin, | ||
int | pos, | ||
int | metricData | ||
) |
mesh | pointer toward the mesh structure |
sol | pointer toward an allocatable sol structure. |
inm | pointer toward the solution file |
bin | 1 if binary file |
pos | of the writted solution |
metricData | 1 if the data saved is a metric (if only 1 data) |
Write the solution value for vertex of index pos in double precision.
int(* _MMG5_bezierCP) (MMG5_pMesh,MMG5_Tria *, _MMG5_pBezier,char) |
int(* _MMG5_chkmsh) (MMG5_pMesh, int, int) |
int(* _MMG5_indElt) (MMG5_pMesh mesh, int kel) |
int(* _MMG5_indPt) (MMG5_pMesh mesh, int kp) |
|
static |
next vertex of triangle: {1,2,0}
|
static |
previous vertex of triangle: {2,0,1}
double(* _MMG5_lenSurfEdg) (MMG5_pMesh mesh, MMG5_pSol sol,int,int, char) |