Skip to content

Increased diagnostic output and added a few lines of documentation #238

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Draft
wants to merge 10 commits into
base: develop
Choose a base branch
from
Draft
Show file tree
Hide file tree
Changes from 1 commit
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
33 changes: 19 additions & 14 deletions src/mmg3d/boulep_3d.c
Original file line number Diff line number Diff line change
Expand Up @@ -1723,15 +1723,16 @@ int MMG3D_coquilFaceFirstLoop(MMG5_pMesh mesh,MMG5_int start,MMG5_int na,MMG5_in
pradj = (*adj);
pri = i;

/* travel through new tetra */
// travel through new tetra.
// ier=face tag if boundary, 0 if not, -1 if error
// Will update *adj, *piv, iface.
ier = MMG5_coquilTravel(mesh,na,nb,adj,piv,&iface,&i);

/* fill the shell */
list[(*ilist)] = 6*(int64_t)pradj +pri;
(*ilist)++;

/* overflow */
if ( (*ilist) > MMG3D_LMAX-2 ) {
if ( (*ilist) > MMG3D_LMAX-2 ) { // overflow
if ( !mmgErr0 ) {
fprintf(stderr,"\n ## Warning: %s: problem in remesh process."
" Coquil of edge %" MMG5_PRId "-%" MMG5_PRId " contains too many elts.\n",
Expand All @@ -1743,16 +1744,11 @@ int MMG3D_coquilFaceFirstLoop(MMG5_pMesh mesh,MMG5_int start,MMG5_int na,MMG5_in
return -1;
}

if ( ier<0 ) return -1;
else if ( !ier ) continue;
if ( ier<0 ) return -1; // eror
else if ( !ier ) continue; // not a boundary

if ( !(*it2) ) {
*it2 = 4*pradj+iface;
(*nbdy)++;
}
else {
(*nbdy)++;
}
if ( !(*it2) ) *it2 = 4*pradj+iface;
(*nbdy)++;

} while ( (*adj) && ((*adj) != start) );

Expand Down Expand Up @@ -1831,7 +1827,11 @@ void MMG3D_coquilFaceSecondLoopInit(MMG5_pMesh mesh,MMG5_int piv,int8_t *iface,
*
* Find all tets sharing edge \a ia of tetra \a start, and stores boundary faces
* when met. \f$ it1 \f$ and \f$ it2 = 6*iel + iface\f$, \a iel = index of
* tetra, \a iface = index of face in tetra.
* tetra, \a iface = index of face in tetra. This function can print a warning
* or error message when it finds that the edge has more than one boundary
* face. This is an error condition if the edge is supposed to be a manifold
* edge. Indeed this function is supposed not to be called for non-manifold
* edges, i.e. edges where multiple boundaries join.
*
* \warning Don't work if \a ia has only one boundary face in its shell.
*/
Expand All @@ -1845,6 +1845,10 @@ int MMG5_coquilface(MMG5_pMesh mesh,MMG5_int start,int8_t iface,int ia,int64_t *

pt = &mesh->tetra[start];

/* MMG5_coquilface is called only on edges marked as manifold, check this */
assert ( pt->xt );
assert ( !(mesh->xtetra[pt->xt].tag[ia] & MG_NOM) );

na = pt->v[ MMG5_iare[ia][0] ];
nb = pt->v[ MMG5_iare[ia][1] ];

Expand Down Expand Up @@ -1883,7 +1887,8 @@ int MMG5_coquilface(MMG5_pMesh mesh,MMG5_int start,int8_t iface,int ia,int64_t *
printf(" ## Warning: %s: you have %d boundary triangles in the closed shell"
" of a manifold edge.\n",__func__,nbdy);
printf(" Problem may occur during remesh process.\n");
mmgWarn0 = 1;
MMG5_show_tet_location(mesh, pt, start);
if(0) mmgWarn0 = 1; // disabled, I want to see how many there are!
Comment on lines +1893 to +1894
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think that it is not obvious so I prefer to explain this a little: the number of times the message is displayed does not reflect the number of edges for which an error is detected.

The function that may raise this warning is called only:

  1. when attempting to swap a boundary edge;
  2. when attempting to split such an edge.

In consequence:

  • It is possible to have an erroneous edge without seeing this warning;
  • a single erroneous edge may raise lot of warnings because:
    • this edge may be seen by the swapping and/or splitting operator;
    • the edge will probably be seen during multiple remeshing iterations;
    • at each iteration and for each operator we will see the edge from all the elements to which it belongs (and as the edge is not modified, we will try to do something and fail each time we will see it);

As the number of warning that is printed is not really meaningful, and as lot of users doesn't need to have the message more than once, I propose to keep this proposition of higher level of verbosity only in debug mode (-d command line option of Mmg):

if ( mesh->info.ddebug ) {
  // print indices and coordinates of the vertices of one tet.
  MMG5_show_tet_location(mesh, pt, start);
} else {
  // disable warnings for this error as the message has already been printed once
  mmgWarn0 = 1;  
} 


/* MMG5_coquilface is called only on edges marked as manifold, check this */
assert ( pt->xt );
Expand Down
2 changes: 2 additions & 0 deletions src/mmg3d/libmmg3d.c
Original file line number Diff line number Diff line change
Expand Up @@ -972,6 +972,8 @@ int MMG3D_packMesh(MMG5_pMesh mesh,MMG5_pSol sol,MMG5_pSol met) {
return 1;
}

// Do all the work in a remeshing run (apart from input etc)
//
int MMG3D_mmg3dlib(MMG5_pMesh mesh,MMG5_pSol met) {
MMG5_pSol sol=NULL; // unused
mytime ctim[TIMEMAX];
Expand Down
6 changes: 6 additions & 0 deletions src/mmg3d/libmmg3d.h
Original file line number Diff line number Diff line change
Expand Up @@ -3057,4 +3057,10 @@ LIBMMG3D_EXPORT int MMG3D_loadVtuMesh_and_allData(MMG5_pMesh mesh,MMG5_pSol *sol
}
#endif


/* -------------------------------------- Mark's hacks --------------------------------------------- */

void MMG5_show_tet_location(MMG5_pMesh mesh, MMG5_pTetra pt, int iel);


#endif
22 changes: 22 additions & 0 deletions src/mmg3d/mmg3d.c
Original file line number Diff line number Diff line change
Expand Up @@ -21,6 +21,7 @@
** =============================================================================
*/


/**
* \file mmg3d/mmg3d.c
* \brief Main file for MMG3D executable: perform 3d mesh adaptation.
Expand Down Expand Up @@ -538,3 +539,24 @@ int main(int argc,char *argv[]) {
/* free mem */
MMG5_RETURN_AND_FREE(mesh,met,ls,disp,ier);
}



/* -------------------------------------- Mark's hacks --------------------------------------------- */

// Print the indices and coordinates of the vertices of one tet. This is to
// inform the user about the location of a problem when the tet index is not
// helpful, for example because it does not correpond to the input or output
// mesh.
//
void MMG5_show_tet_location(MMG5_pMesh mesh, MMG5_pTetra pt, int iel)
{
fprintf(stderr, " ## tet index %d\n", iel);
for(int j=0; j<4; j++){
double U[3], *S = mesh->point[pt->v[j]].c; // unscaled and scaled coords
for(int i=0; i<3; i++) U[i] = S[i]*mesh->info.delta + mesh->info.min[i];
fprintf(stderr, " ## vertex %d at (%f,%f,%f)\n", pt->v[j], U[0], U[1], U[2]);
}
}


Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thanks, it can be useful.

Just for your information, we also have debugging tools to print the packed indices of tetra and points (which means the indices they will have in the output mesh):

  • the MMG3D_indElt(mesh,k) function gives you the index of tetra k inside the output mesh;
  • the MMG3D_indPt(mesh,k) function gives you the index of point k inside the output mesh;

Note that it is useful and trustable only if you save the mesh immediately after calling these functions (because tetra and points indices will be modified if we continue the remeshing so, even if the targetted point/tetra are not touched, the mesh packing will be different).

83 changes: 61 additions & 22 deletions src/mmg3d/mmg3d1.c
Original file line number Diff line number Diff line change
Expand Up @@ -625,16 +625,24 @@ MMG5_int MMG5_swpmsh(MMG5_pMesh mesh,MMG5_pSol met,MMG3D_pPROctree PROctree, int

ret = MMG5_coquilface(mesh,k,i,ia,list,&it1,&it2,0);
ilist = ret / 2;
if ( ret < 0 ) return -1;
if ( ret < 0 ){
MMG5_show_tet_location(mesh, pt, k);
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Please, keep this in debug mode only (-d command line option, mesh->info.ddebug variable inside the code) :

if ( mesh->info.ddebug ) {
   MMG5_show_tet_location(mesh, pt, k);
}

return -1;
}
/* CAUTION: trigger collapse with 2 elements */
if ( ilist <= 1 ) continue;
ier = MMG5_chkswpbdy(mesh,met,list,ilist,it1,it2,typchk);
if ( ier < 0 )
if ( ier < 0 ){
MMG5_show_tet_location(mesh, pt, k);
Comment on lines +635 to +636
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Please, keep this in debug mode only (-d command line option, mesh->info.ddebug variable inside the code) :

if ( mesh->info.ddebug ) {
   MMG5_show_tet_location(mesh, pt, k);
}

return -1;
}
else if ( ier ) {
ier = MMG5_swpbdy(mesh,met,list,ret,it1,PROctree,typchk);
if ( ier > 0 ) ns++;
else if ( ier < 0 ) return -1;
else if ( ier < 0 ){
MMG5_show_tet_location(mesh, pt, k);
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Please, keep this in debug mode only (-d command line option, mesh->info.ddebug variable inside the code) :

if ( mesh->info.ddebug ) {
   MMG5_show_tet_location(mesh, pt, k);
}

return -1;
}
break;
}
}
Expand Down Expand Up @@ -694,11 +702,17 @@ MMG5_int MMG5_swptet(MMG5_pMesh mesh,MMG5_pSol met,double crit,double declic,
}

nconf = MMG5_chkswpgen(mesh,met,k,i,&ilist,list,crit,typchk);
if ( nconf<0 ) return -1;
if ( nconf<0 ){
MMG5_show_tet_location(mesh, pt, k);
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Please, keep this in debug mode only (-d command line option, mesh->info.ddebug variable inside the code) :

if ( mesh->info.ddebug ) {
   MMG5_show_tet_location(mesh, pt, k);
}

return -1;
}
else if ( nconf ) {
ier = MMG5_swpgen(mesh,met,nconf,ilist,list,PROctree,typchk);
if ( ier > 0 ) ns++;
else if ( ier < 0 ) return -1;
else if ( ier < 0 ){
MMG5_show_tet_location(mesh, pt, k);
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Please, keep this in debug mode only (-d command line option, mesh->info.ddebug variable inside the code) :

if ( mesh->info.ddebug ) {
   MMG5_show_tet_location(mesh, pt, k);
}

return -1;
}
break;
}
}
Expand Down Expand Up @@ -807,34 +821,41 @@ MMG5_int MMG5_movtet(MMG5_pMesh mesh,MMG5_pSol met, MMG3D_pPROctree PROctree,
if( !ier ) continue;
else if ( ier>0 )
ier = MMG5_movbdynompt(mesh,met,PROctree,listv,ilistv,lists,ilists,improveVolSurf);
else
else{
MMG5_show_tet_location(mesh, pt, k);
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Please, keep this in debug mode only (-d command line option, mesh->info.ddebug variable inside the code) :

if ( mesh->info.ddebug ) {
   MMG5_show_tet_location(mesh, pt, k);
}

return -1;
}
}
}
else if ( ppt->tag & MG_GEO ) {
ier=MMG5_boulesurfvolp(mesh,k,i0,i,listv,&ilistv,lists,&ilists,0);
if ( !ier ) continue;
else if ( ier>0 )
ier = MMG5_movbdyridpt(mesh,met,PROctree,listv,ilistv,lists,ilists,improveVolSurf);
else
else{
MMG5_show_tet_location(mesh, pt, k);
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Please, keep this in debug mode only (-d command line option, mesh->info.ddebug variable inside the code) :

if ( mesh->info.ddebug ) {
   MMG5_show_tet_location(mesh, pt, k);
}

return -1;
}
}
else if ( ppt->tag & MG_REF ) {
ier=MMG5_boulesurfvolp(mesh,k,i0,i,listv,&ilistv,lists,&ilists,0);
if ( !ier )
continue;
else if ( ier>0 )
ier = MMG5_movbdyrefpt(mesh,met,PROctree,listv,ilistv,lists,ilists,improveVolSurf);
else
else{
MMG5_show_tet_location(mesh, pt, k);
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Please, keep this in debug mode only (-d command line option, mesh->info.ddebug variable inside the code) :

if ( mesh->info.ddebug ) {
   MMG5_show_tet_location(mesh, pt, k);
}

return -1;
}
}
else {
ier=MMG5_boulesurfvolp(mesh,k,i0,i,listv,&ilistv,lists,&ilists,0);
if ( !ier )
continue;
else if ( ier<0 )
else if ( ier<0 ){
MMG5_show_tet_location(mesh, pt, k);
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Please, keep this in debug mode only (-d command line option, mesh->info.ddebug variable inside the code) :

if ( mesh->info.ddebug ) {
   MMG5_show_tet_location(mesh, pt, k);
}

return -1;

}
n = &(mesh->xpoint[ppt->xp].n1[0]);

/* If the orientation of the tetra face is
Expand All @@ -846,7 +867,10 @@ MMG5_int MMG5_movtet(MMG5_pMesh mesh,MMG5_pSol met, MMG3D_pPROctree PROctree,
}
ier = MMG5_movbdyregpt(mesh,met,PROctree,listv,ilistv,
lists,ilists,improveSurf,improveVolSurf);
if (ier < 0 ) return -1;
if (ier < 0 ){
MMG5_show_tet_location(mesh, pt, k);
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Please, keep this in debug mode only (-d command line option, mesh->info.ddebug variable inside the code) :

if ( mesh->info.ddebug ) {
   MMG5_show_tet_location(mesh, pt, k);
}

return -1;
}
else if ( ier ) ns++;
}
}
Expand Down Expand Up @@ -964,14 +988,18 @@ static int MMG5_coltet(MMG5_pMesh mesh,MMG5_pSol met,int8_t typchk) {
else {
if ( mesh->adja[4*(k-1)+1+i] ) continue;
if (MMG5_boulesurfvolpNom(mesh,k,ip,i,
list,&ilist,lists,&ilists,&refmin,&refplus,p0->tag & MG_NOM) < 0 )
list,&ilist,lists,&ilists,&refmin,&refplus,p0->tag & MG_NOM) < 0 ){
MMG5_show_tet_location(mesh, pt, k);
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Please, keep this in debug mode only (-d command line option, mesh->info.ddebug variable inside the code) :

if ( mesh->info.ddebug ) {
   MMG5_show_tet_location(mesh, pt, k);
}

return -1;
}
}
}
else {
if (MMG5_boulesurfvolp(mesh,k,ip,i,
list,&ilist,lists,&ilists,p0->tag & MG_NOM) < 0 )
list,&ilist,lists,&ilists,p0->tag & MG_NOM) < 0 ){
MMG5_show_tet_location(mesh, pt, k);
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Please, keep this in debug mode only (-d command line option, mesh->info.ddebug variable inside the code) :

if ( mesh->info.ddebug ) {
   MMG5_show_tet_location(mesh, pt, k);
}

return -1;
}
}
}
else {
Expand Down Expand Up @@ -1095,14 +1123,18 @@ static int MMG5_coltet(MMG5_pMesh mesh,MMG5_pSol met,int8_t typchk) {
else {
if ( mesh->adja[4*(k-1)+1+i] ) continue;
if (MMG5_boulesurfvolpNom(mesh,k,ip,i,
list,&ilist,lists,&ilists,&refmin,&refplus,p0->tag & MG_NOM) < 0 )
list,&ilist,lists,&ilists,&refmin,&refplus,p0->tag & MG_NOM) < 0 ){
MMG5_show_tet_location(mesh, pt, k);
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Please, keep this in debug mode only (-d command line option, mesh->info.ddebug variable inside the code) :

if ( mesh->info.ddebug ) {
   MMG5_show_tet_location(mesh, pt, k);
}

return -1;
}
}
}
else {
if (MMG5_boulesurfvolp(mesh,k,ip,i,
list,&ilist,lists,&ilists,p0->tag & MG_NOM) < 0 )
list,&ilist,lists,&ilists,p0->tag & MG_NOM) < 0 ){
MMG5_show_tet_location(mesh, pt, k);
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Please, keep this in debug mode only (-d command line option, mesh->info.ddebug variable inside the code) :

if ( mesh->info.ddebug ) {
   MMG5_show_tet_location(mesh, pt, k);
}

return -1;
}
}
}
else {
Expand Down Expand Up @@ -1134,13 +1166,19 @@ static int MMG5_coltet(MMG5_pMesh mesh,MMG5_pSol met,int8_t typchk) {

if ( ilist > 0 ) {
ier = MMG5_colver(mesh,met,list,ilist,iq,typchk);
if ( ier < 0 ) return -1;
if ( ier < 0 ){
MMG5_show_tet_location(mesh, pt, k);
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Please, keep this in debug mode only (-d command line option, mesh->info.ddebug variable inside the code) :

if ( mesh->info.ddebug ) {
   MMG5_show_tet_location(mesh, pt, k);
}

return -1;
}
else if ( ier ) {
MMG3D_delPt(mesh,ier);
break;
}
}
else if (ilist < 0 ) return -1;
else if (ilist < 0 ){
MMG5_show_tet_location(mesh, pt, k);
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Please, keep this in debug mode only (-d command line option, mesh->info.ddebug variable inside the code) :

if ( mesh->info.ddebug ) {
   MMG5_show_tet_location(mesh, pt, k);
}

return -1;
}
}
if ( ier ) {
p1->flag = base;
Expand Down Expand Up @@ -2593,6 +2631,7 @@ MMG3D_anatets_iso(MMG5_pMesh mesh,MMG5_pSol met,int8_t typchk) {
fprintf(stderr,"\n ## Warning: %s: surfacic pattern: unable to find"
" a valid split for at least 1 point. Point(s) deletion.\n",
__func__ );
MMG5_show_tet_location(mesh, pt, k);
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Please, keep this in debug mode only (-d command line option, mesh->info.ddebug variable inside the code) :

if ( mesh->info.ddebug ) {
   MMG5_show_tet_location(mesh, pt, k);
}

mmgWarn2 = 1;
}

Expand Down Expand Up @@ -3161,7 +3200,7 @@ int MMG5_anatet(MMG5_pMesh mesh,MMG5_pSol met,int8_t typchk, int patternMode) {
ier = MMG3D_anatets(mesh,met,typchk);

if ( ier < 0 ) {
fprintf(stderr,"\n ## Unable to complete surface mesh. Exit program.\n");
fprintf(stderr,"\n ## Unable to complete surface mesh in MMG5_anatet. Exit program.\n");
return 0;
}
ns += ier;
Expand All @@ -3172,7 +3211,7 @@ int MMG5_anatet(MMG5_pMesh mesh,MMG5_pSol met,int8_t typchk, int patternMode) {
/* analyze internal tetras */
ier = MMG5_anatetv(mesh,met,typchk);
if ( ier < 0 ) {
fprintf(stderr,"\n ## Unable to complete volume mesh. Exit program.\n");
fprintf(stderr,"\n ## Unable to complete volume mesh in MMG5_anatet. Exit program.\n");
return 0;
}
ns += ier;
Expand All @@ -3188,7 +3227,7 @@ int MMG5_anatet(MMG5_pMesh mesh,MMG5_pSol met,int8_t typchk, int patternMode) {
if ( !mesh->info.noinsert ) {
nc = MMG5_coltet(mesh,met,typchk);
if ( nc < 0 ) {
fprintf(stderr,"\n ## Unable to collapse mesh. Exiting.\n");
fprintf(stderr,"\n ## Unable to collapse mesh in MMG5_anatet: nc=%d. Exiting.\n", nc);
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

if we pass here, nc is equal to -1 no ?

return 0;
}
}
Expand All @@ -3198,14 +3237,14 @@ int MMG5_anatet(MMG5_pMesh mesh,MMG5_pSol met,int8_t typchk, int patternMode) {
if ( !mesh->info.noswap ) {
ier = MMG5_swpmsh(mesh,met,NULL,typchk);
if ( ier < 0 ) {
fprintf(stderr,"\n ## Unable to improve mesh. Exiting.\n");
fprintf(stderr,"\n ## Unable to improve mesh in MMG5_anatet. Exiting.\n");
return 0;
}
nf += ier;

ier = MMG5_swptet(mesh,met,MMG3D_LSWAPIMPROVE,MMG3D_SWAP06,NULL,typchk,mesh->mark-2);
if ( ier < 0 ) {
fprintf(stderr,"\n ## Unable to improve mesh. Exiting.\n");
fprintf(stderr,"\n ## Unable to improve mesh in MMG5_anatet. Exiting.\n");
return 0;
}
nf += ier;
Expand Down
6 changes: 4 additions & 2 deletions src/mmg3d/quality_3d.c
Original file line number Diff line number Diff line change
Expand Up @@ -91,7 +91,9 @@ int MMG3D_tetraQual(MMG5_pMesh mesh, MMG5_pSol met,int8_t metRidTyp) {

/* Here the quality is not normalized by alpha, thus we need to
* normalized it */
return MMG5_minQualCheck(iel,minqual,MMG3D_ALPHAD);
int r = MMG5_minQualCheck(iel,minqual,MMG3D_ALPHAD);
if(r==0) MMG5_show_tet_location(mesh, pt, iel);
return r;
}

/**
Expand Down Expand Up @@ -476,7 +478,7 @@ int MMG3D_displayQualHisto(MMG5_int ne,double max,double avg,double min,MMG5_int
fprintf(stdout," (LES)");
fprintf(stdout," %" MMG5_PRId "\n",ne);

fprintf(stdout," BEST %8.6f AVRG. %8.6f WRST. %8.6f (%" MMG5_PRId ")\n",
fprintf(stdout," BEST %8.6f AVRG. %8.6f WRST. %.6g (%" MMG5_PRId ")\n",
max,avg / ne,min,iel);

return ( MMG3D_displayQualHisto_internal(ne,max,avg,min,iel,good,med,his,
Expand Down