-// MEFISTO : library to compute 2D triangulation from segmented boundaries
+// MEFISTO2: a library to compute 2D triangulation from segmented boundaries
//
-// Copyright (C) 2003 Laboratoire J.-L. Lions UPMC Paris
//
-// This library is free software; you can redistribute it and/or
-// modify it under the terms of the GNU Lesser General Public
-// License as published by the Free Software Foundation; either
-// version 2.1 of the License.
+// Copyright (C) 2006-2014 CEA/DEN, EDF R&D, OPEN CASCADE
//
-// This library is distributed in the hope that it will be useful,
-// but WITHOUT ANY WARRANTY; without even the implied warranty of
-// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
-// Lesser General Public License for more details.
+// This library is free software; you can redistribute it and/or
+// modify it under the terms of the GNU Lesser General Public
+// License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
//
-// You should have received a copy of the GNU Lesser General Public
-// License along with this library; if not, write to the Free Software
-// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
+// This library is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+// Lesser General Public License for more details.
//
-// See http://www.ann.jussieu.fr/~perronne or email Perronnet@ann.jussieu.fr
+// You should have received a copy of the GNU Lesser General Public
+// License along with this library; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
//
+// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
//
-// File : aptrte.cxx
+// File : aptrte.cxx le C++ de l'appel du trianguleur plan
// Module : SMESH
-// Author: Alain PERRONNET
+// Author : Alain PERRONNET
+// Date : 13 novembre 2006
-using namespace std;
#include "Rn.h"
#include "aptrte.h"
#include "utilities.h"
+using namespace std;
+
extern "C"
{
R aretemaxface_;
- R areteideale_( R3 xyz, R3 direction )
+ MEFISTO2D_EXPORT
+ R
+ #ifdef WIN32
+ #ifdef F2C_BUILD
+ #else
+ __stdcall
+ #endif
+ #endif
+ areteideale()//( R3 xyz, R3 direction )
{
return aretemaxface_;
}
static double cpunew, cpuold=0;
-void tempscpu_( double & tempsec )
+void
+#ifdef WIN32
+#ifdef F2C_BUILD
+#else
+ __stdcall
+#endif
+#endif
+tempscpu_( double & tempsec )
//Retourne le temps CPU utilise en secondes
{
tempsec = ( (double) clock() ) / CLOCKS_PER_SEC;
}
-void deltacpu_( R & dtcpu )
+void
+#ifdef WIN32
+#ifdef F2C_BUILD
+#else
+ __stdcall
+#endif
+#endif
+deltacpu_( R & dtcpu )
//Retourne le temps CPU utilise en secondes depuis le precedent appel
{
tempscpu_( cpunew );
}
-void aptrte( Z nutysu, R aretmx,
- Z nblf, Z * nudslf, R2 * uvslf,
- Z nbpti, R2 *uvpti,
- Z & nbst, R2 * & uvst, Z & nbt, Z * & nust,
- Z & ierr )
+void aptrte( Z nutysu, R aretmx,
+ Z nblf, Z * nudslf, R2 * uvslf,
+ Z nbpti, R2 * uvpti,
+ Z & nbst, R2 * & uvst,
+ Z & nbt, Z * & nust,
+ Z & ierr )
//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
// but : appel de la triangulation par un arbre-4 recouvrant
// ----- de triangles equilateraux
// ierr : 0 si pas d'erreur
// > 0 sinon
//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
-// auteur : Alain Perronnet Analyse Numerique Paris UPMC decembre 2001
+// auteur : Alain Perronnet Laboratoire J.-L. LIONS Paris UPMC mars 2006
//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
{
+ Z nbsttria=4; //Attention: 4 sommets stockes par triangle
+ //no st1, st2, st3, 0 (non quadrangle)
+
R d, tcpu=0;
- R3 direction=R3(0,0,0); //direction pour areteideale() inactive ici!
+// R3 direction=R3(0,0,0); //direction pour areteideale() inactive ici!
Z nbarfr=nudslf[nblf]; //nombre total d'aretes des lignes fermees
Z mxtrou = Max( 1024, nblf ); //nombre maximal de trous dans la surface
R3 *mnpxyd=NULL;
Z *mnsoar=NULL, mosoar=7, mxsoar, n1soar; //le hachage des aretes
Z *mnartr=NULL, moartr=3, mxartr, n1artr; //le no des 3 aretes des triangles
- Z *mntree=NULL, motree=9, mxtree; //L'arbre 4 de TE et nombre d'entiers par TE
+ Z *mntree=NULL, motree=9, mxtree; //L'arbre 4 de TE et nombre d'entiers par TE
Z *mnqueu=NULL, mxqueu;
Z *mn1arcf=NULL;
Z *mnarcf=NULL, mxarcf;
Z *mnarcf2=NULL;
Z *mnarcf3=NULL;
Z *mntrsu=NULL;
- Z *mndalf=NULL;
Z *mnslig=NULL;
Z *mnarst=NULL;
Z *mnlftr=NULL;
R3 comxmi[2]; //coordonnees UV Min et Maximales
R aremin, aremax; //longueur minimale et maximale des aretes
+ R airemx; //aire maximale souhaitee d'un triangle
R quamoy, quamin;
Z noar0, noar, na;
Z i, l, n, ns, ns0, ns1, ns2, nosotr[3], nt;
Z mxsomm, nbsomm, nbarpi, nbarli, ndtri0, mn;
Z moins1=-1;
-
- aretemaxface_ = aretmx;
+ Z nuds = 0;
// initialisation du temps cpu
deltacpu_( d );
// quelques reservations de tableaux pour faire les calculs
// ========================================================
- // le tableau pointeur sur la premiere arete de chaque ligne fermee
- if( mndalf!=NULL ) delete [] mndalf;
- mndalf = new Z[1+nblf];
- if( mndalf==NULL ) goto ERREUR;
- mndalf[0]=0;
-
// declaration du tableau des coordonnees des sommets de la frontiere
// puis des sommets internes ajoutes
// majoration empirique du nombre de sommets de la triangulation
- i = 4*nbarfr/10;
+ i = 4*nbarfr/10;
mxsomm = Max( 20000, 64*nbpti+i*i );
- MESSAGE( "APTRTE: Depart de la triangulation avec " );
- MESSAGE( "nutysu=" << nutysu << " aretmx=" << aretmx << " mxsomm=" << mxsomm );
+ MESSAGE( "APTRTE: Debut de la triangulation plane avec " );
+ MESSAGE( "nutysu=" << nutysu << " aretmx=" << aretmx
+ << " mxsomm=" << mxsomm );
+ MESSAGE( nbarfr << " sommets sur la frontiere et " << nbpti << " points internes");
NEWDEPART:
//mnpxyd( 3, mxsomm ) les coordonnees UV des sommets et la taille d'arete aux sommets
mnsoar = new Z[mosoar*mxsoar];
if( mnsoar==NULL ) goto ERREUR;
//initialiser le tableau mnsoar pour le hachage des aretes
- insoar_( mxsomm, mosoar, mxsoar, n1soar, mnsoar );
+ insoar( mxsomm, mosoar, mxsoar, n1soar, mnsoar );
// mnarst( mxsomm ) numero mnsoar d'une arete pour chacun des sommets
if( mnarst!=NULL ) delete [] mnarst;
mnarst = new Z[1+mxsomm];
if( mnarst==NULL ) goto ERREUR;
n = 1+mxsomm;
- azeroi_( n, mnarst );
+ azeroi( n, mnarst );
// mnslig( mxsomm ) no de sommet dans sa ligne pour chaque sommet frontalier
// ou no du point si interne forc'e par l'utilisateur
if( mnslig!=NULL ) delete [] mnslig;
mnslig = new Z[mxsomm];
if( mnslig==NULL ) goto ERREUR;
- azeroi_( mxsomm, mnslig );
+ azeroi( mxsomm, mnslig );
// initialisation des aretes frontalieres de la triangulation future
// renumerotation des sommets des aretes des lignes pour la triangulation
ns0 = nudslf[n-1];
mnpxyd[ns0].x = uvslf[ns0].x;
mnpxyd[ns0].y = uvslf[ns0].y;
- mnpxyd[ns0].z = areteideale_( mnpxyd[ns0], direction );
-// cout << "Sommet " << ns0 << ": " << mnpxyd[ns0].x
-// << " " << mnpxyd[ns0].y << " longueur arete=" << mnpxyd[ns0].z << endl;
+ mnpxyd[ns0].z = areteideale();//( mnpxyd[ns0], direction );
+// MESSAGE("Sommet " << ns0 << ": " << mnpxyd[ns0].x
+// << " " << mnpxyd[ns0].y << " longueur arete=" << mnpxyd[ns0].z);
//carre de la longueur de l'arete 1 de la ligne fermee n
d = pow( uvslf[ns0+1].x - uvslf[ns0].x, 2 )
//le numero n de la ligne du sommet et son numero ns1 dans la ligne
mnslig[ns0-1] = 1000000 * n + ns1-nudslf[n-1];
- fasoar_( ns1, ns2, moins1, moins1, n,
- mosoar, mxsoar, n1soar, mnsoar, mnarst,
- noar0, ierr );
+ fasoar( ns1, ns2, moins1, moins1, n,
+ mosoar, mxsoar, n1soar, mnsoar, mnarst,
+ noar0, ierr );
//pas de test sur ierr car pas de saturation possible a ce niveau
//le pointeur dans le hachage sur la premiere arete de la ligne fermee n
- mndalf[n] = noar0;
+ //mndalf[n] = noar0;
//la nouvelle arete est la suivante de l'arete definie juste avant
if( noar > 0 )
{
ns1 = ns2; //le numero de l'arete et le numero du premier sommet de l'arete
if( i < nbarli )
- //nbs+1 est le 2-eme sommet de l'arete i de la ligne fermee n
- ns2 = ns1+1;
+ //nbs+1 est le 2-eme sommet de l'arete i de la ligne fermee n
+ ns2 = ns1+1;
else
- //le 2-eme sommet de la derniere arete est le premier sommet de la ligne
- ns2 = ns0;
+ //le 2-eme sommet de la derniere arete est le premier sommet de la ligne
+ ns2 = ns0;
//l'arete precedente est dotee de sa suivante:celle cree ensuite
//les 2 coordonnees du sommet ns2 de la ligne
ns = ns1 - 1;
+//debut ajout 5/10/2006 ................................................
+ nuds = Max( nuds, ns ); //le numero du dernier sommet traite
+//fin ajout 5/10/2006 ................................................
mnpxyd[ns].x = uvslf[ns].x;
mnpxyd[ns].y = uvslf[ns].y;
- mnpxyd[ns].z = areteideale_( mnpxyd[ns], direction );
-// cout << "Sommet " << ns << ": " << mnpxyd[ns].x
-// << " " << mnpxyd[ns].y << " longueur arete=" << mnpxyd[ns].z << endl;
+ mnpxyd[ns].z = areteideale();//( mnpxyd[ns], direction );
+// MESSAGE("Sommet " << ns << ": " << mnpxyd[ns].x
+// << " " << mnpxyd[ns].y << " longueur arete=" << mnpxyd[ns].z);
//carre de la longueur de l'arete
d = pow( uvslf[ns2-1].x - uvslf[ns1-1].x, 2)
aremin = Min( aremin, d );
aremax = Max( aremax, d );
+//debut ajout du 5/10/2006 .............................................
+ //la longueur de l'arete ns1-ns2
+ d = sqrt( d );
+ //longueur arete = Min ( aretmx, aretes incidentes )
+ mnpxyd[ns ].z = Min( mnpxyd[ns ].z, d );
+ mnpxyd[ns2-1].z = Min( mnpxyd[ns2-1].z, d );
+//fin ajout du 5/10/2006 ...............................................
+
//le numero n de la ligne du sommet et son numero ns1 dans la ligne
mnslig[ns] = 1000000 * n + ns1-nudslf[n-1];
//ajout de l'arete dans la liste
- fasoar_( ns1, ns2, moins1, moins1, n,
- mosoar, mxsoar, n1soar, mnsoar,
- mnarst, noar, ierr );
+ fasoar( ns1, ns2, moins1, moins1, n,
+ mosoar, mxsoar, n1soar, mnsoar,
+ mnarst, noar, ierr );
//pas de test sur ierr car pas de saturation possible a ce niveau
//chainage des aretes frontalieres en position 6 du tableau mnsoar
aremin = sqrt( aremin ); //longueur minimale d'une arete des lignes fermees
aremax = sqrt( aremax ); //longueur maximale d'une arete
- aretmx = Min( aretmx, aremax ); //pour homogeneiser
- cout << "nutysu=" << nutysu << " aretmx=" << aretmx
- << " arete min=" << aremin << " arete max=" << aremax << endl;
+//debut ajout 9/11/2006 ................................................
+ // devenu un commentaire aretmx = Min( aretmx, aremax ); //pour homogeneiser
+
+ // protection contre une arete max desiree trop grande ou trop petite
+ if( aretmx > aremax*2.05 ) aretmx = aremax;
+
+ // protection contre une arete max desiree trop petite
+ if( (aremax-aremin) > (aremin+aremax)*0.05 && aretmx < aremin*0.5 )
+ aretmx =(aremin+aremax*2)/3.0;
+
+ if( aretmx < aremin && aremin > 0 )
+ aretmx = aremin;
+
+ //sauvegarde pour la fonction areteideale_
+ aretemaxface_ = aretmx;
+
+ //aire maximale souhaitee des triangles
+ airemx = aretmx * aretmx * sqrt(3.0) / 2.0; //Aire triangle equilateral
+
+ for(i=0; i<=nuds; i++ )
+ mnpxyd[i].z = Min( mnpxyd[i].z, aretmx );
+ //MESSAGE("Numero du dernier sommet frontalier=" << nuds+1);
+//fin ajout 9/11/2006 .................................................
+
+
+ MESSAGE("Sur le bord: arete min=" << aremin << " arete max=" << aremax );
+ MESSAGE("Triangulation: arete mx=" << aretmx
+ << " triangle aire mx=" << airemx );
//chainage des aretes frontalieres : la derniere arete frontaliere
mnsoar[ mosoar * noar - mosoar + 5 ] = 0;
//les 2 coordonnees du point i de sommet nbs
mnpxyd[ns1].x = uvpti[i].x;
mnpxyd[ns1].y = uvpti[i].y;
- mnpxyd[ns1].z = areteideale_( mnpxyd[ns1], direction );
+ mnpxyd[ns1].z = areteideale();//( mnpxyd[ns1], direction );
//le numero i du point interne
mnslig[ns1] = i+1;
ns1++;
if( mntree==NULL ) goto ERREUR;
//initialisation du tableau letree et ajout dans letree des sommets 1 a nbsomm
- teajte_( mxsomm, nbsomm, mnpxyd, comxmi, aretmx, mxtree, mntree, ierr );
+ comxmi[0].x = comxmi[1].x = uvslf[0].x;
+ comxmi[0].y = comxmi[1].y = uvslf[0].y;
+ teajte( mxsomm, nbsomm, mnpxyd, comxmi, aretmx, mxtree, mntree, ierr );
comxmi[0].z=0;
comxmi[1].z=0;
mnqueu = new Z[mxqueu];
if( mnqueu==NULL) goto ERREUR;
- tehote_( nutysu, nbarpi, mxsomm, nbsomm, mnpxyd,
- comxmi, aretmx,
- mntree, mxqueu, mnqueu,
- ierr );
+ tehote( nutysu, nbarpi, mxsomm, nbsomm, mnpxyd,
+ comxmi, aretmx,
+ mntree, mxqueu, mnqueu,
+ ierr );
deltacpu_( d );
tcpu += d;
- cout << "Temps de l'adaptation et l'homogeneisation de l'arbre-4 des TE="
- << d << " secondes" << endl;
+ MESSAGE("Temps de l'adaptation et l'homogeneisation de l'arbre-4 des TE="
+ << d << " secondes");
if( ierr != 0 )
{
//destruction du tableau auxiliaire et de l'arbre
// trianguler les triangles equilateraux feuilles a partir de leurs 3 sommets
// et des points de la frontiere, des points internes imposes interieurs
// ==========================================================================
- tetrte_( comxmi, aretmx, nbarpi, mxsomm, mnpxyd,
- mxqueu, mnqueu, mntree, mosoar, mxsoar, n1soar, mnsoar,
- moartr, mxartr, n1artr, mnartr, mnarst,
- ierr );
+ tetrte( comxmi, aretmx, nbarpi, mxsomm, mnpxyd,
+ mxqueu, mnqueu, mntree, mosoar, mxsoar, n1soar, mnsoar,
+ moartr, mxartr, n1artr, mnartr, mnarst,
+ ierr );
// destruction de la queue et de l'arbre devenus inutiles
delete [] mnqueu; mnqueu=NULL;
//qualites de la triangulation actuelle
qualitetrte( mnpxyd, mosoar, mxsoar, mnsoar, moartr, mxartr, mnartr,
- nbt, quamoy, quamin );
+ nbt, quamoy, quamin );
// boucle sur les aretes internes (non sur une ligne de la frontiere)
// avec echange des 2 diagonales afin de rendre la triangulation delaunay
// ======================================================================
// formation du chainage 6 des aretes internes a echanger eventuellement
- aisoar_( mosoar, mxsoar, mnsoar, na );
- tedela_( mnpxyd, mnarst,
- mosoar, mxsoar, n1soar, mnsoar, na,
- moartr, mxartr, n1artr, mnartr, n );
+ aisoar( mosoar, mxsoar, mnsoar, na );
+ tedela( mnpxyd, mnarst,
+ mosoar, mxsoar, n1soar, mnsoar, na,
+ moartr, mxartr, n1artr, mnartr, n );
MESSAGE( "Nombre d'echanges des diagonales de 2 triangles=" << n );
deltacpu_( d );
tcpu += d;
- cout << "Temps de la triangulation Delaunay par echange des diagonales="
- << d << " secondes" << endl;
+ MESSAGE("Temps de la triangulation Delaunay par echange des diagonales="
+ << d << " secondes");
//qualites de la triangulation actuelle
qualitetrte( mnpxyd, mosoar, mxsoar, mnsoar, moartr, mxartr, mnartr,
- nbt, quamoy, quamin );
+ nbt, quamoy, quamin );
// detection des aretes frontalieres initiales perdues
// triangulation frontale pour les restaurer
mnarcf2 = new Z[mxarcf];
if( mnarcf2 == NULL ) goto ERREUR;
- terefr_( nbarpi, mnpxyd,
- mosoar, mxsoar, n1soar, mnsoar,
- moartr, n1artr, mnartr, mnarst,
- mxarcf, mn1arcf, mnarcf, mnarcf1, mnarcf2,
- n, ierr );
+ terefr( nbarpi, mnpxyd,
+ mosoar, mxsoar, n1soar, mnsoar,
+ moartr, mxartr, n1artr, mnartr, mnarst,
+ mxarcf, mn1arcf, mnarcf, mnarcf1, mnarcf2,
+ n, ierr );
- MESSAGE( "Restauration de " << n << " aretes perdues de la frontiere" );
+ MESSAGE( "Restauration de " << n << " aretes perdues de la frontiere ierr=" << ierr );
deltacpu_( d );
tcpu += d;
- cout << "Temps de la recuperation des aretes perdues de la frontiere="
- << d << " secondes" << endl;
+ MESSAGE("Temps de la recuperation des aretes perdues de la frontiere="
+ << d << " secondes");
if( ierr != 0 ) goto ERREUR;
//qualites de la triangulation actuelle
qualitetrte( mnpxyd, mosoar, mxsoar, mnsoar, moartr, mxartr, mnartr,
- nbt, quamoy, quamin );
+ nbt, quamoy, quamin );
// fin de la triangulation avec respect des aretes initiales frontalieres
for (n=0; n<nblf; n++) //numero de la ligne fermee de 1 a nblf
mnlftr[n] = n+1;
- tesuex_( nblf, mnlftr,
- ndtri0, nbsomm, mnpxyd, mnslig,
- mosoar, mxsoar, mnsoar,
- moartr, mxartr, n1artr, mnartr, mnarst,
- nbt, mntrsu, ierr );
+ tesuex( nblf, mnlftr,
+ ndtri0, nbsomm, mnpxyd, mnslig,
+ mosoar, mxsoar, mnsoar,
+ moartr, mxartr, n1artr, mnartr, mnarst,
+ nbt, mntrsu, ierr );
delete [] mnlftr; mnlftr=NULL;
delete [] mntrsu; mntrsu=NULL;
deltacpu_( d );
tcpu += d;
- MESSAGE( "Temps de la suppression des triangles externes=" << d );
+ MESSAGE( "Temps de la suppression des triangles externes=" << d << "ierr=" << ierr );
if( ierr != 0 ) goto ERREUR;
//qualites de la triangulation actuelle
qualitetrte( mnpxyd, mosoar, mxsoar, mnsoar, moartr, mxartr, mnartr,
- nbt, quamoy, quamin );
+ nbt, quamoy, quamin );
// amelioration de la qualite de la triangulation par
// barycentrage des sommets internes a la triangulation
// mise en delaunay de la triangulation
// =====================================================
mnarcf3 = new Z[mxarcf];
- if( mnarcf3 == NULL ) goto ERREUR;
-
- teamqt_( nutysu,
- mnarst, mosoar, mxsoar, n1soar, mnsoar,
- moartr, mxartr, n1artr, mnartr,
- mxarcf, mnarcf2, mnarcf3,
- mn1arcf, mnarcf, mnarcf1,
- comxmi, nbarpi, nbsomm, mxsomm, mnpxyd, mnslig,
- ierr );
+ if( mnarcf3 == NULL )
+ {
+ MESSAGE ( "aptrte: MC saturee mnarcf3=" << mnarcf3 );
+ goto ERREUR;
+ }
+ teamqt( nutysu, aretmx, airemx,
+ mnarst, mosoar, mxsoar, n1soar, mnsoar,
+ moartr, mxartr, n1artr, mnartr,
+ mxarcf, mnarcf2, mnarcf3,
+ mn1arcf, mnarcf, mnarcf1,
+ nbarpi, nbsomm, mxsomm, mnpxyd, mnslig,
+ ierr );
+ if( mnarcf3 != NULL ) {delete [] mnarcf3; mnarcf3=NULL;}
if( mn1arcf != NULL ) {delete [] mn1arcf; mn1arcf=NULL;}
if( mnarcf != NULL ) {delete [] mnarcf; mnarcf =NULL;}
if( mnarcf1 != NULL ) {delete [] mnarcf1; mnarcf1=NULL;}
if( mnarcf2 != NULL ) {delete [] mnarcf2; mnarcf2=NULL;}
- if( mnarcf3 != NULL ) {delete [] mnarcf3; mnarcf3=NULL;}
deltacpu_( d );
tcpu += d;
MESSAGE( "Temps de l'amelioration de la qualite de la triangulation=" << d );
- if( ierr != 0 ) goto ERREUR;
+ if( ierr == -13 ) ierr=0; //6/10/2006 arret de l'amelioration apres boucle infinie dans caetoi
+ if( ierr != 0 ) goto ERREUR;
//qualites de la triangulation finale
qualitetrte( mnpxyd, mosoar, mxsoar, mnsoar, moartr, mxartr, mnartr,
- nbt, quamoy, quamin );
+ nbt, quamoy, quamin );
// renumerotation des sommets internes: mnarst(i)=numero final du sommet
// ===================================
if( mnartr[nt*moartr-moartr] != 0 )
{
//le numero des 3 sommets du triangle nt
- nusotr_( nt, mosoar, mnsoar, moartr, mnartr, nosotr );
+ nusotr( nt, mosoar, mnsoar, moartr, mnartr, nosotr );
//les 3 sommets du triangle sont actifs
mnarst[ nosotr[0] ] = 1;
mnarst[ nosotr[1] ] = 1;
n = mnslig[i];
if( n > 0 )
{
- if( n >= 1000000 )
- {
- //sommet d'une ligne
- //retour aux coordonnees initiales dans uvslf
- l = n / 1000000;
- n = n - 1000000 * l + nudslf[l-1] - 1;
- uvst[nbst].x = uvslf[n].x;
- uvst[nbst].y = uvslf[n].y;
- }
- else
- {
- //point utilisateur n interne impose
- //retour aux coordonnees initiales dans uvpti
- uvst[nbst].x = uvpti[n-1].x;
- uvst[nbst].y = uvpti[n-1].y;
- }
+ if( n >= 1000000 )
+ {
+ //sommet d'une ligne
+ //retour aux coordonnees initiales dans uvslf
+ l = n / 1000000;
+ n = n - 1000000 * l + nudslf[l-1] - 1;
+ uvst[nbst].x = uvslf[n].x;
+ uvst[nbst].y = uvslf[n].y;
+ }
+ else
+ {
+ //point utilisateur n interne impose
+ //retour aux coordonnees initiales dans uvpti
+ uvst[nbst].x = uvpti[n-1].x;
+ uvst[nbst].y = uvpti[n-1].y;
+ }
}
}
}
// -----------------------------------------------------
// boucle sur les triangles occupes (internes et externes)
if( nust != NULL ) delete [] nust;
- nust = new Z[4*nbt];
+ nust = new Z[nbsttria*nbt];
if( nust == NULL ) goto ERREUR;
nbt = 0;
for (i=1; i<=mxartr; i++)
if( mnartr[i*moartr-moartr] != 0 )
{
//le triangle i est interne => nosotr numero de ses 3 sommets
- nusotr_( i, mosoar, mnsoar, moartr, mnartr, nosotr );
+ nusotr( i, mosoar, mnsoar, moartr, mnartr, nosotr );
nust[nbt++] = mnarst[ nosotr[0] ];
nust[nbt++] = mnarst[ nosotr[1] ];
nust[nbt++] = mnarst[ nosotr[2] ];
nust[nbt++] = 0;
}
}
- nbt /= 4; //le nombre final de triangles de la surface
- cout << "Nombre de sommets=" << nbst
- << " Nombre de triangles=" << nbt << endl;
-
+ nbt /= nbsttria; //le nombre final de triangles de la surface
+ MESSAGE( "APTRTE: Fin de la triangulation plane avec "<<nbst<<" sommets et "
+ << nbt << " triangles" );
deltacpu_( d );
tcpu += d;
- MESSAGE( "Temps total de la triangulation=" << tcpu << " secondes" );
+ MESSAGE( "APTRTE: Temps total de la triangulation plane=" << tcpu << " secondes" );
// destruction des tableaux auxiliaires
// ------------------------------------
if( mnslig != NULL ) delete [] mnslig;
if( mnsoar != NULL ) delete [] mnsoar;
if( mnpxyd != NULL ) delete [] mnpxyd;
- if( mndalf != NULL ) delete [] mndalf;
if( mntree != NULL ) delete [] mntree;
if( mnqueu != NULL ) delete [] mnqueu;
if( mntrsu != NULL ) delete [] mntrsu;
}
else
{
- MESSAGE( "Triangulation non realisee " << ierr );
+ MESSAGE( "APTRTE: Triangulation NON REALISEE avec erreur=" << ierr );
if( ierr == 0 ) ierr=1;
goto NETTOYAGE;
}
}
-
-
-void qualitetrte( R3 *mnpxyd,
- Z & mosoar, Z & mxsoar, Z *mnsoar,
- Z & moartr, Z & mxartr, Z *mnartr,
- Z & nbtria, R & quamoy, R & quamin )
+void
+#ifdef WIN32
+#ifdef F2C_BUILD
+#else
+ __stdcall
+#endif
+#endif
+ qualitetrte( R3 *mnpxyd,
+ Z & mosoar, Z & mxsoar, Z *mnsoar,
+ Z & moartr, Z & mxartr, Z *mnartr,
+ Z & nbtria, R & quamoy, R & quamin )
// +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
// but : calculer la qualite moyenne et minimale de la triangulation
// ----- actuelle definie par les tableaux mnsoar et mnartr
// +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
{
R d, aire, qualite;
- Z nosotr[3], mn, nbtrianeg, nt;
+ Z nosotr[3], mn, nbtrianeg, nt, ntqmin;
aire = 0;
quamoy = 0;
quamin = 2.0;
nbtria = 0;
nbtrianeg = 0;
+ ntqmin = 0;
mn = -moartr;
for ( nt=1; nt<=mxartr; nt++ )
nbtria++;
//le numero des 3 sommets du triangle nt
- nusotr_( nt, mosoar, mnsoar, moartr, mnartr, nosotr );
+ nusotr( nt, mosoar, mnsoar, moartr, mnartr, nosotr );
//la qualite du triangle ns1 ns2 ns3
- qutr2d_( mnpxyd[nosotr[0]-1], mnpxyd[nosotr[1]-1], mnpxyd[nosotr[2]-1],
- qualite );
+ qutr2d( mnpxyd[nosotr[0]-1], mnpxyd[nosotr[1]-1], mnpxyd[nosotr[2]-1],
+ qualite );
//la qualite moyenne
quamoy += qualite;
//la qualite minimale
- quamin = Min( quamin, qualite );
+ if( qualite < quamin )
+ {
+ quamin = qualite;
+ ntqmin = nt;
+ }
//aire signee du triangle nt
- d = surtd2_( mnpxyd[nosotr[0]-1], mnpxyd[nosotr[1]-1], mnpxyd[nosotr[2]-1] );
+ d = surtd2( mnpxyd[nosotr[0]-1], mnpxyd[nosotr[1]-1], mnpxyd[nosotr[2]-1] );
if( d<0 )
{
- //un triangle d'aire negative de plus
- nbtrianeg++;
- cout << "ATTENTION: le triangle " << nt << " de sommets:"
- << nosotr[0] << " " << nosotr[1] << " " << nosotr[2]
- << " a une aire " << d <<"<=0" << endl;
+ //un triangle d'aire negative de plus
+ nbtrianeg++;
+ MESSAGE("ATTENTION: le triangle " << nt << " de sommets:"
+ << nosotr[0] << " " << nosotr[1] << " " << nosotr[2]
+ << " a une aire " << d <<"<=0");
}
//aire des triangles actuels
//les affichages
quamoy /= nbtria;
- cout << "Qualite moyenne=" << quamoy
+ MESSAGE("Qualite moyenne=" << quamoy
<< " Qualite minimale=" << quamin
- << " des " << nbtria << " triangles de surface totale="
- << aire << endl;
+ << " des " << nbtria << " triangles de surface plane totale="
+ << aire);
+
+ if( quamin<0.3 )
+ {
+ //le numero des 3 sommets du triangle ntqmin de qualite minimale
+ nusotr(ntqmin, mosoar, mnsoar, moartr, mnartr, nosotr );
+ MESSAGE("Triangle de qualite minimale "<<quamin<<" de sommets:"
+ <<nosotr[0]<<" "<<nosotr[1]<<" "<<nosotr[2]<<" ");
+ for (int i=0;i<3;i++)
+ MESSAGE("Sommet "<<nosotr[i]<<": x="<< mnpxyd[nosotr[i]-1].x
+ <<" y="<< mnpxyd[nosotr[i]-1].y);
+ }
if( nbtrianeg>0 )
- MESSAGE( "ATTENTION: nombre de triangles d'aire negative=" << nbtrianeg );
+ MESSAGE( "ATTENTION: "<< nbtrianeg << " TRIANGLES d'AIRE NEGATIVE" );
+
+ MESSAGE(" ");
return;
}