ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/REPOS_ERICCA/magic/lib/mailleur/src/mailleur3d_optimisation.cpp
Revision: 339
Committed: Wed May 30 18:20:47 2012 UTC (12 years, 11 months ago) by francois
File size: 23739 byte(s)
Log Message:
Correction bug mailleur bloc + correction bug inversion avec open cascade + preparation pour element  XFEM

File Contents

# User Rev Content
1 francois 283 //------------------------------------------------------------
2     //------------------------------------------------------------
3     // MAGiC
4     // Jean Christophe Cuilli�e et Vincent FRANCOIS
5     // D�artement de G�ie M�anique - UQTR
6     //------------------------------------------------------------
7     // Le projet MAGIC est un projet de recherche du d�artement
8     // de g�ie m�anique de l'Universit�du Qu�ec �
9     // Trois Rivi�es
10     // Les librairies ne peuvent �re utilis�s sans l'accord
11     // des auteurs (contact : francois@uqtr.ca)
12     //------------------------------------------------------------
13     //------------------------------------------------------------
14     //
15     // mailleur3d.cpp
16     //
17     //------------------------------------------------------------
18     //------------------------------------------------------------
19     // COPYRIGHT 2000
20     // Version du 02/03/2006 �11H23
21     //------------------------------------------------------------
22     //------------------------------------------------------------
23    
24    
25     #include "gestionversion.h"
26 francois 287 #include "mailleur3d_optimisation.h"
27 francois 283 #include "ot_mathematique.h"
28     #include "m3d_tetra.h"
29 francois 287 #include "m3d_triangle.h"
30 francois 283 #include "m3d_noeud.h"
31 francois 287 #include "mg_maillage.h"
32 francois 288 #include "const.h"
33 francois 283 //#include "mg_gestionnaire.h"
34     #include <math.h>
35    
36     //---------------------------------------------------------------------------
37    
38     //#pragma package(smart_init)
39    
40 francois 287
41     MAILLEUR3D_OPTIMISATION::MAILLEUR3D_OPTIMISATION(MG_MAILLAGE* mgmai,int niv):mg_maillage(mgmai),niveau_optimisation(niv),affichageactif(0)
42 francois 283 {
43 francois 287 o3d_data();
44     o3d_data2();
45     }
46    
47     MAILLEUR3D_OPTIMISATION::~MAILLEUR3D_OPTIMISATION()
48     {
49     }
50    
51    
52    
53     double MAILLEUR3D_OPTIMISATION::get_volume(MG_TETRA* tet)
54     {
55 francois 283 double *xyz1=tet->get_noeud1()->get_coord();
56     double *xyz2=tet->get_noeud2()->get_coord();
57     double *xyz3=tet->get_noeud3()->get_coord();
58     double *xyz4=tet->get_noeud4()->get_coord();
59     OT_VECTEUR_3D ab(xyz1,xyz2);
60     OT_VECTEUR_3D ac(xyz1,xyz3);
61     OT_VECTEUR_3D ad(xyz1,xyz4);
62     double vol=(ab&ac)*ad;
63     vol=vol/6.;
64     return vol;
65     }
66    
67 francois 287 double MAILLEUR3D_OPTIMISATION::get_volume(double *xyz1,double *xyz2,double *xyz3,double *xyz4)
68 francois 283 {
69     OT_VECTEUR_3D ab(xyz1,xyz2);
70     OT_VECTEUR_3D ac(xyz1,xyz3);
71     OT_VECTEUR_3D ad(xyz1,xyz4);
72     double vol=(ab&ac)*ad;
73     vol=vol/6.;
74     return vol;
75     }
76    
77 francois 288 void MAILLEUR3D_OPTIMISATION::optimise(MG_VOLUME* mgvol)
78     {
79     char mess[255];
80     sprintf(mess," Niveau d'optimisation : %d",niveau_optimisation);
81     if (affichageactif==1) affiche(mess);
82     int nbaoptimiser;
83     int nbaoptimiserapres=mg_maillage->get_nb_mg_tetra();
84     do
85     {
86     nbaoptimiser=nbaoptimiserapres;
87     optimise(mgvol,nbaoptimiserapres);
88     if (nbaoptimiserapres==0)
89     {
90     nbaoptimiser=0;
91     sprintf(mess," Tout est optimise");
92     if (affichageactif==1) affiche(mess);
93     }
94    
95     }
96     while (nbaoptimiserapres!=nbaoptimiser);
97     }
98    
99    
100    
101 francois 287 void MAILLEUR3D_OPTIMISATION::optimise(MG_VOLUME* mgvol,int& nbmauvais)
102 francois 283 {
103     static int passe=0;
104     passe++;
105     //preparation
106 francois 339 //LISTE_MG_TETRA::iterator it;
107     TPL_SET<MG_ELEMENT_MAILLAGE*>::ITERATEUR it;
108 francois 283 int nbtet=0;
109 francois 339 for (MG_TETRA* tet=(MG_TETRA*)mgvol->get_lien_maillage()->get_premier(it);tet;tet=(MG_TETRA*)mgvol->get_lien_maillage()->get_suivant(it))
110     //for (MG_TETRA* tet=mg_maillage->get_premier_tetra(it);tet;tet=mg_maillage->get_suivant_tetra(it))
111 francois 283 {
112 francois 339 if (mg_maillage->get_mg_tetraid(tet->get_id())==NULL) continue;
113     M3D_TETRA* mtet=(M3D_TETRA*)tet;
114 francois 283 double qual=OPERATEUR::qualite_tetra(mtet->get_noeud1()->get_coord(),mtet->get_noeud2()->get_coord(),mtet->get_noeud3()->get_coord(),mtet->get_noeud4()->get_coord());
115     mtet->change_qualite(qual);
116     if (qual<0.1*niveau_optimisation) nbtet++;
117     ajouter_ordre_tetra(mtet,0);
118     }
119 francois 287
120 francois 283 //on commence
121 francois 339 if (nbtet==0) return;
122 francois 283 for (int phase=0;phase<2;phase++)
123     {
124     char mess[255];
125     if (phase==0) sprintf(mess," Phase %d : %d a optimiser",2*passe-1+phase,nbtet);
126     if (phase==1) sprintf(mess," Phase %d",2*passe-1+phase);
127     if (affichageactif==1) affiche(mess);
128     while (lst_tetra[phase].size()>0)
129     {
130     ORDRE_TETRA::iterator i=lst_tetra[phase].begin();
131     M3D_TETRA* tet=(*i).second;
132     supprimer_ordre_tetra(tet);
133     double crit_avant=tet->get_qualite();
134 francois 339 if (crit_avant<0.1*niveau_optimisation)
135 francois 283 {
136     double crit[4]={-1.,-1.,-1.,-1.};
137     M3D_NOEUD* no[4];
138     no[0]=(M3D_NOEUD*)tet->get_noeud1();
139     no[1]=(M3D_NOEUD*)tet->get_noeud2();
140     no[2]=(M3D_NOEUD*)tet->get_noeud3();
141     no[3]=(M3D_NOEUD*)tet->get_noeud4();
142     double x[4],y[4],z[4];
143     bouge_point(mgvol,no[0],crit[0],x[0],y[0],z[0]);
144     bouge_point(mgvol,no[1],crit[1],x[1],y[1],z[1]);
145     bouge_point(mgvol,no[2],crit[2],x[2],y[2],z[2]);
146     bouge_point(mgvol,no[3],crit[3],x[3],y[3],z[3]);
147     std::vector<double> vecteur_crit(crit,crit+4);
148     std::vector<double>::iterator it1 = max_element(vecteur_crit.begin(), vecteur_crit.end());
149     double crit_opt=*it1;
150     if (crit_opt>crit_avant)
151     {
152     int num=it1-vecteur_crit.begin();
153     no[num]->change_x(x[num]);
154     no[num]->change_y(y[num]);
155     no[num]->change_z(z[num]);
156     int nb=no[num]->get_lien_tetra()->get_nb();
157     for (int j=0;j<nb;j++)
158     {
159     M3D_TETRA* mtet=(M3D_TETRA*)no[num]->get_lien_tetra()->get(j);
160     double qual=OPERATEUR::qualite_tetra(mtet->get_noeud1()->get_coord(),mtet->get_noeud2()->get_coord(),mtet->get_noeud3()->get_coord(),mtet->get_noeud4()->get_coord());
161     mtet->change_qualite(qual);
162     }
163     crit_avant=crit_opt;
164     }
165     }
166     int coquille_non_optimise=1;
167 francois 339 if (crit_avant<0.1*niveau_optimisation)
168 francois 283 {
169     COQUILLE coque[6];
170     double crit[6]={-1,-1,-1,-1,-1,-1};
171     remaille_coquille(tet->get_noeud1(),tet->get_noeud2(),crit[0],coque[0]);
172     remaille_coquille(tet->get_noeud1(),tet->get_noeud3(),crit[1],coque[1]);
173     remaille_coquille(tet->get_noeud1(),tet->get_noeud4(),crit[2],coque[2]);
174     remaille_coquille(tet->get_noeud2(),tet->get_noeud3(),crit[3],coque[3]);
175     remaille_coquille(tet->get_noeud2(),tet->get_noeud3(),crit[4],coque[4]);
176     remaille_coquille(tet->get_noeud3(),tet->get_noeud4(),crit[5],coque[5]);
177     std::vector<double> vecteur_crit(crit,crit+6);
178     std::vector<double>::iterator it1 = max_element(vecteur_crit.begin(), vecteur_crit.end());
179     double crit_opt=*it1;
180     if (crit_opt>crit_avant)
181     {
182     coquille_non_optimise=0;
183     int num=it1-vecteur_crit.begin();
184     for (int i=0;i<coque[num].taille;i++)
185     {
186     M3D_TETRA *tettmp=(M3D_TETRA*)coque[num].tet[i];
187     if (tet!=tettmp) supprimer_ordre_tetra(tettmp);
188     mg_maillage->supprimer_mg_tetraid(coque[num].tet[i]->get_id());
189     }
190     for (int i=0;i<2*coque[num].taille-4;i++)
191     {
192     MG_NOEUD* noeud1=coque[num].new_tetra[4*i];
193     MG_NOEUD* noeud2=coque[num].new_tetra[4*i+1];
194     MG_NOEUD* noeud3=coque[num].new_tetra[4*i+2];
195     MG_NOEUD* noeud4=coque[num].new_tetra[4*i+3];
196     MG_TRIANGLE* triangle1=mg_maillage->get_mg_triangle(noeud1->get_id(),noeud3->get_id(),noeud2->get_id());
197     MG_TRIANGLE* triangle2=mg_maillage->get_mg_triangle(noeud1->get_id(),noeud2->get_id(),noeud4->get_id());
198     MG_TRIANGLE* triangle3=mg_maillage->get_mg_triangle(noeud2->get_id(),noeud3->get_id(),noeud4->get_id());
199     MG_TRIANGLE* triangle4=mg_maillage->get_mg_triangle(noeud1->get_id(),noeud4->get_id(),noeud3->get_id());
200     if (triangle1==NULL) triangle1=insere_triangle(mgvol,noeud1,noeud3,noeud2,MAILLEUR_AUTO);
201     if (triangle2==NULL) triangle2=insere_triangle(mgvol,noeud1,noeud2,noeud4,MAILLEUR_AUTO);
202     if (triangle3==NULL) triangle3=insere_triangle(mgvol,noeud2,noeud3,noeud4,MAILLEUR_AUTO);
203     if (triangle4==NULL) triangle4=insere_triangle(mgvol,noeud1,noeud4,noeud3,MAILLEUR_AUTO);
204     M3D_TETRA* mtet=new M3D_TETRA(mgvol,noeud1,noeud2,noeud3,noeud4,triangle1,triangle2,triangle3,triangle4,MAILLEUR_AUTO);
205     mg_maillage->ajouter_mg_tetra(mtet);
206     double qual=OPERATEUR::qualite_tetra(mtet->get_noeud1()->get_coord(),mtet->get_noeud2()->get_coord(),mtet->get_noeud3()->get_coord(),mtet->get_noeud4()->get_coord());
207     mtet->change_qualite(qual);
208     if (phase==0)
209     ajouter_ordre_tetra(mtet,1);
210     }
211    
212    
213     }
214    
215     }
216     if ((phase==0) && (coquille_non_optimise==1))
217     ajouter_ordre_tetra(tet,1);
218     }
219    
220     }
221     int nb=0;
222 francois 339 //for (MG_TETRA* tet=mg_maillage->get_premier_tetra(it);tet;tet=mg_maillage->get_suivant_tetra(it))
223     for (MG_TETRA* tet=(MG_TETRA*)mgvol->get_lien_maillage()->get_premier(it);tet;tet=(MG_TETRA*)mgvol->get_lien_maillage()->get_suivant(it))
224 francois 283 {
225 francois 339 if (mg_maillage->get_mg_tetraid(tet->get_id())==NULL) continue;
226 francois 283 M3D_TETRA* mtet=(M3D_TETRA*)tet;
227     if (mtet->get_qualite()<0.1*niveau_optimisation) nb++;
228     }
229    
230     nbmauvais=nb;
231    
232     }
233    
234    
235    
236    
237    
238 francois 287 int MAILLEUR3D_OPTIMISATION::bouge_point(MG_VOLUME* mgvol,M3D_NOEUD* noeud,double& crit,double& x,double& y, double& z)
239 francois 283 {
240 francois 288 if (mgvol!=NULL)
241     if (noeud->get_lien_topologie()!=mgvol) return 0;
242     if (mgvol==NULL)
243     if (noeud->get_etat()==INACTIF) return 0;
244 francois 283 if (noeud->get_origine()!=MAILLEUR_AUTO) return 0;
245    
246     double xopt=0.;
247     double yopt=0.;
248     double zopt=0.;
249     double qual_dep=1.;
250     int nb_tet=noeud->get_lien_tetra()->get_nb();
251     double tab_coord[3000];
252     double gamma=0;
253     for (int i=0;i<nb_tet;i++)
254     {
255     M3D_TETRA* tet=(M3D_TETRA*)noeud->get_lien_tetra()->get(i);
256     qual_dep=std::min(qual_dep,tet->get_qualite());
257     MG_NOEUD* no1=tet->get_noeud1();
258     MG_NOEUD* no2=tet->get_noeud2();
259     MG_NOEUD* no3=tet->get_noeud3();
260     MG_NOEUD* no4=tet->get_noeud4();
261     MG_NOEUD *noeud1,*noeud2,*noeud3;
262     if (no1==noeud)
263     {
264     noeud1=no2;
265     noeud2=no4;
266     noeud3=no3;
267     }
268     if (no2==noeud)
269     {
270     noeud1=no1;
271     noeud2=no3;
272     noeud3=no4;
273     }
274     if (no3==noeud)
275     {
276     noeud1=no1;
277     noeud2=no4;
278     noeud3=no2;
279     }
280     if (no4==noeud)
281     {
282     noeud1=no1;
283     noeud2=no2;
284     noeud3=no3;
285     }
286     double xyzi[3];
287     double *xyz1=noeud1->get_coord();
288     double *xyz2=noeud2->get_coord();
289     double *xyz3=noeud3->get_coord();
290     xyzi[0]=(xyz1[0]+xyz2[0]+xyz3[0])/3.;
291     xyzi[1]=(xyz1[1]+xyz2[1]+xyz3[1])/3.;
292     xyzi[2]=(xyz1[2]+xyz2[2]+xyz3[2])/3.;
293     OT_VECTEUR_3D v12(xyz1,xyz2);
294     OT_VECTEUR_3D v13(xyz1,xyz3);
295     OT_VECTEUR_3D v23(xyz2,xyz3);
296     OT_VECTEUR_3D normal=v12&v13;
297     double perimetre=v12.get_longueur()+v13.get_longueur()+v23.get_longueur();
298     double hauteur = (perimetre/3.) * 0.8 ;
299     normal.norme();
300     tab_coord[3*i]=xyzi[0]+normal.get_x()*hauteur;
301     tab_coord[3*i+1]=xyzi[1]+normal.get_y()*hauteur;
302     tab_coord[3*i+2]=xyzi[2]+normal.get_z()*hauteur;
303     if (tet->get_qualite()> 1e-10) gamma=gamma+1./tet->get_qualite()/tet->get_qualite();
304     }
305     gamma=1./gamma;
306     for (int i=0;i<nb_tet;i++)
307     {
308     M3D_TETRA* tet=(M3D_TETRA*)noeud->get_lien_tetra()->get(i);
309     if (tet->get_qualite()> 1e-10)
310     {
311     double alpha=gamma/tet->get_qualite()/tet->get_qualite();
312     xopt=xopt+alpha*tab_coord[3*i];
313     yopt=yopt+alpha*tab_coord[3*i+1];
314     zopt=zopt+alpha*tab_coord[3*i+2];
315     }
316    
317     }
318     double delta[3]={xopt-noeud->get_x(),yopt-noeud->get_y(),zopt-noeud->get_z()};
319     double bmin=0.,bmax=1.;
320     double vieuxx=noeud->get_x();
321     double vieuxy=noeud->get_y();
322     double vieuxz=noeud->get_z();
323     double qualini=0.;
324     for (int iteration=0;iteration<5;iteration++)
325     {
326     double alpha=(bmin+bmax)*0.5;
327     noeud->change_x(vieuxx+alpha*delta[0]);
328     noeud->change_y(vieuxy+alpha*delta[1]);
329     noeud->change_z(vieuxz+alpha*delta[2]);
330     double qualcour=1.;
331     for (int i=0;i<nb_tet;i++)
332     {
333     M3D_TETRA* mtet=(M3D_TETRA*)noeud->get_lien_tetra()->get(i);
334     double qual=OPERATEUR::qualite_tetra(mtet->get_noeud1()->get_coord(),mtet->get_noeud2()->get_coord(),mtet->get_noeud3()->get_coord(),mtet->get_noeud4()->get_coord());
335     qualcour=std::min(qualcour,qual);
336     }
337     double alpha_eps=(bmin+bmax)*0.5-(bmax-bmin)/50.;
338     noeud->change_x(vieuxx+alpha_eps*delta[0]);
339     noeud->change_y(vieuxy+alpha_eps*delta[1]);
340     noeud->change_z(vieuxz+alpha_eps*delta[2]);
341     double qualcour_eps=1.;
342     for (int i=0;i<nb_tet;i++)
343     {
344     M3D_TETRA* mtet=(M3D_TETRA*)noeud->get_lien_tetra()->get(i);
345     double qual=OPERATEUR::qualite_tetra(mtet->get_noeud1()->get_coord(),mtet->get_noeud2()->get_coord(),mtet->get_noeud3()->get_coord(),mtet->get_noeud4()->get_coord());
346     qualcour_eps=std::min(qualcour_eps,qual);
347     }
348     if (qualcour>qualcour_eps) bmin =alpha;
349     else bmax=alpha;
350     qualini=std::max(qualini,qualcour);
351     }
352     noeud->change_x(vieuxx);
353     noeud->change_y(vieuxy);
354     noeud->change_z(vieuxz);
355     if (qualini>qual_dep)
356     {
357     x=vieuxx+(bmin+bmax)*0.5*delta[0];
358     y=vieuxy+(bmin+bmax)*0.5*delta[1];
359     z=vieuxz+(bmin+bmax)*0.5*delta[2];
360     crit=qualini;
361     return 1;
362     }
363     return 0;
364     }
365    
366    
367    
368 francois 287 void MAILLEUR3D_OPTIMISATION::remaille_coquille(MG_NOEUD* noeud1,MG_NOEUD* noeud2, double& crit, COQUILLE& coque)
369 francois 283 {
370     //recherche du segment a supprimer
371     MG_NOEUD* noeud;
372     if (noeud1->get_id()>noeud2->get_id()) noeud=noeud2;
373     else noeud=noeud1;
374     int nb=noeud->get_lien_petit_segment()->get_nb();
375     for (int i=0;i<nb;i++)
376     {
377     MG_SEGMENT* seg=noeud->get_lien_petit_segment()->get(i);
378     if ( ((seg->get_noeud1()==noeud1) && (seg->get_noeud2()==noeud2)) || ((seg->get_noeud1()==noeud2) && (seg->get_noeud2()==noeud1)))
379 francois 288 {
380     if (mg_maillage->get_mg_geometrie()==NULL)
381     {
382     if (seg->get_dimension_topo_null()==2)
383     {
384     coque.taille=0;
385     crit=0.;
386     return;
387     }
388     }
389     else
390     if ((seg->get_lien_topologie()->get_dimension()!=3) || (seg->get_origine()!=MAILLEUR_AUTO))
391 francois 283 {
392     coque.taille=0;
393     crit=0.;
394     return;
395     }
396 francois 288 }
397 francois 283 }
398     //recherche des tetra qui s enroule autour du segment
399     int nb1=noeud1->get_lien_tetra()->get_nb();
400     int nb2=noeud2->get_lien_tetra()->get_nb();
401     MG_TETRA* coq[100];
402     int nb_coq=0;
403     for (int i=0;i<nb1;i++)
404     for (int j=0;j<nb2;j++)
405     if (noeud1->get_lien_tetra()->get(i)==noeud2->get_lien_tetra()->get(j))
406     {
407     coq[nb_coq]=noeud1->get_lien_tetra()->get(i);
408     nb_coq++;
409     }
410     if ((nb_coq<4) || (nb_coq>10))
411     {
412     coque.taille=0;
413     crit=0.;
414     return;
415     }
416     //recherche des noeuds qui compose la coquille
417     MG_NOEUD* tet_noeud[50];
418     int nb_tet_noeud=0;
419     coque.taille=nb_coq;
420     for (int i=0;i<nb_coq;i++)
421     {
422     coque.tet[i]=coq[i];
423     if ((coq[i]->get_noeud1()!=noeud1) && (coq[i]->get_noeud1()!=noeud2)) tet_noeud[nb_tet_noeud++]=coq[i]->get_noeud1();
424     if ((coq[i]->get_noeud2()!=noeud1) && (coq[i]->get_noeud2()!=noeud2)) tet_noeud[nb_tet_noeud++]=coq[i]->get_noeud2();
425     if ((coq[i]->get_noeud3()!=noeud1) && (coq[i]->get_noeud3()!=noeud2)) tet_noeud[nb_tet_noeud++]=coq[i]->get_noeud3();
426     if ((coq[i]->get_noeud4()!=noeud1) && (coq[i]->get_noeud4()!=noeud2)) tet_noeud[nb_tet_noeud++]=coq[i]->get_noeud4();
427     coque.volume=coque.volume+fabs(get_volume(coq[i]));
428     }
429     //ordonnancement des noeuds qui compose la coquille en polygone
430     MG_NOEUD* polygone[20];
431     int nb_poly=2;
432     polygone[0]=tet_noeud[0];
433     tet_noeud[0]=NULL;
434     polygone[1]=tet_noeud[1];
435     tet_noeud[1]=NULL;
436     while (nb_poly!=nb_coq+1)
437     {
438     for (int j=0;j<nb_coq;j++)
439     {
440     if (tet_noeud[2*j]==polygone[nb_poly-1])
441     {
442     polygone[nb_poly++]=tet_noeud[2*j+1];
443     tet_noeud[2*j]=NULL;
444     tet_noeud[2*j+1]=NULL;
445     }
446     if (tet_noeud[2*j+1]==polygone[nb_poly-1])
447     {
448     polygone[nb_poly++]=tet_noeud[2*j];
449     tet_noeud[2*j]=NULL;
450     tet_noeud[2*j+1]=NULL;
451     }
452    
453     }
454     }
455     //etude du positionnement du polygone
456     double vol1=get_volume(polygone[0]->get_coord(),polygone[1]->get_coord(),polygone[nb_poly-1]->get_coord(),noeud1->get_coord());
457     double vol2=get_volume(polygone[0]->get_coord(),polygone[nb_poly-1]->get_coord(),polygone[1]->get_coord(),noeud2->get_coord());
458     if (vol1*vol2<0.)
459     {
460     crit=0.;
461     coque.taille=0;
462     return;
463     }
464     MG_NOEUD *noeuda=noeud1;
465     MG_NOEUD *noeudb=noeud2;
466     if (vol1<0.)
467     {
468     noeuda=noeud2;
469     noeudb=noeud1;
470     }
471     //Examen des solutions
472     int nb_solution,nb_triangle;
473     if (nb_coq==4) {
474     nb_solution=2;
475     nb_triangle=4;
476     }
477     if (nb_coq==5) {
478     nb_solution=5;
479     nb_triangle=10;
480     }
481     if (nb_coq==6) {
482     nb_solution=14;
483     nb_triangle=20;
484     }
485     if (nb_coq==7) {
486     nb_solution=42;
487     nb_triangle=35;
488     }
489     if (nb_coq==8) {
490     nb_solution=132;
491     nb_triangle=56;
492     }
493     if (nb_coq==9) {
494     nb_solution=429;
495     nb_triangle=84;
496     }
497     if (nb_coq==10) {
498     nb_solution=1430;
499     nb_triangle=120;
500     }
501     double crit_triangle[120];
502     for (int i=0;i<nb_triangle;i++)
503     {
504     double crit1=OPERATEUR::qualite_tetra(polygone[tab_face[nb_coq-4][i][0]]->get_coord(),polygone[tab_face[nb_coq-4][i][1]]->get_coord(),polygone[tab_face[nb_coq-4][i][2]]->get_coord(),noeuda->get_coord());
505     double crit2=OPERATEUR::qualite_tetra(polygone[tab_face[nb_coq-4][i][0]]->get_coord(),polygone[tab_face[nb_coq-4][i][2]]->get_coord(),polygone[tab_face[nb_coq-4][i][1]]->get_coord(),noeudb->get_coord());
506     crit_triangle[i]=std::min(crit1,crit2);
507     }
508     /* examen de chaque solution */
509     double crit_opt=0.;
510     int numero_solution= -1;
511     double crit_solution[1430];
512     for (int i=0;i<nb_solution;i++)
513     {
514     double volume=0.;
515     for (int j=0;j<nb_coq-2;j++)
516     {
517     MG_NOEUD* noa=polygone[tab_face[nb_coq-4][tab_solution[nb_coq-4][i][j]][0]];
518     MG_NOEUD* nob=polygone[tab_face[nb_coq-4][tab_solution[nb_coq-4][i][j]][1]];
519     MG_NOEUD* noc=polygone[tab_face[nb_coq-4][tab_solution[nb_coq-4][i][j]][2]];
520     MG_NOEUD* nod=noeuda;
521     volume=volume+fabs(get_volume(noa->get_coord(),nob->get_coord(),noc->get_coord(),nod->get_coord()));
522     nod=noeudb;
523     volume=volume+fabs(get_volume(noa->get_coord(),nob->get_coord(),noc->get_coord(),nod->get_coord()));
524     }
525     double eps=0.0018*pow(volume,0.666666666);
526     if (OPERATEUR::egal(volume,coque.volume,eps))
527     {
528     crit_solution[i]=1.;
529     for (int j=0;j<nb_coq-2;j++)
530     crit_solution[i]=std::min(crit_solution[i],crit_triangle[tab_solution[nb_coq-4][i][j]]);
531     }
532     else crit_solution[i]=0.;
533     if (crit_opt<crit_solution[i])
534     {
535     crit_opt=crit_solution[i];
536     numero_solution=i;
537     }
538    
539     }
540     if (numero_solution==(-1))
541     {
542     crit=0.;
543     coque.taille=0;
544     return;
545     }
546     crit=crit_opt;
547     for (int j=0;j<nb_coq-2;j++)
548     {
549     MG_NOEUD* noa=polygone[tab_face[nb_coq-4][tab_solution[nb_coq-4][numero_solution][j]][0]];
550     MG_NOEUD* nob=polygone[tab_face[nb_coq-4][tab_solution[nb_coq-4][numero_solution][j]][1]];
551     MG_NOEUD* noc=polygone[tab_face[nb_coq-4][tab_solution[nb_coq-4][numero_solution][j]][2]];
552     coque.new_tetra[8*j]=noa;
553     coque.new_tetra[8*j+1]=nob;
554     coque.new_tetra[8*j+2]=noc;
555     coque.new_tetra[8*j+3]=noeuda;
556     coque.new_tetra[8*j+4]=noa;
557     coque.new_tetra[8*j+5]=noc;
558     coque.new_tetra[8*j+6]=nob;
559     coque.new_tetra[8*j+7]=noeudb;
560     }
561     }
562    
563    
564 francois 287 void MAILLEUR3D_OPTIMISATION::ajouter_ordre_tetra(M3D_TETRA* tet,int num)
565 francois 283 {
566 francois 339 double val;
567     if (num==0) val=tet->get_qualite();
568     if (num==1) val=1./tet->get_qualite();
569     std::pair<double,M3D_TETRA*> tmp(val,tet);
570 francois 283 ORDRE_TETRA::iterator p=lst_tetra[num].insert(tmp);
571     lst_tetraid[num][tet->get_id()]=p;
572     }
573    
574 francois 287 void MAILLEUR3D_OPTIMISATION::supprimer_ordre_tetra(M3D_TETRA* tet)
575 francois 283 {
576     int num=0;
577     ORDRE_TETRA_PARID::iterator it=lst_tetraid[num].find(tet->get_id());
578     if (it==lst_tetraid[num].end())
579     {
580     num=1;
581     it=lst_tetraid[num].find(tet->get_id());
582     }
583     if (it==lst_tetraid[num].end()) return;
584     lst_tetra[num].erase(it->second);
585     lst_tetraid[num].erase(it);
586     }
587    
588    
589 francois 287 void MAILLEUR3D_OPTIMISATION::change_niveau_optimisation(int num)
590 francois 283 {
591     niveau_optimisation=num;
592     }
593    
594 francois 287 int MAILLEUR3D_OPTIMISATION::get_niveau_optimisation(void)
595 francois 283 {
596     return niveau_optimisation;
597     }
598    
599 francois 287 void MAILLEUR3D_OPTIMISATION::active_affichage(void (*fonc)(char*))
600     {
601     affiche=fonc;
602     affichageactif=1;
603     }
604    
605    
606    
607    
608     class MG_TRIANGLE* MAILLEUR3D_OPTIMISATION::insere_triangle(class MG_VOLUME* mgvol,MG_NOEUD* noeud1,MG_NOEUD* noeud2,MG_NOEUD* noeud3,int origine)
609     {
610     MG_TRIANGLE* triangle=mg_maillage->get_mg_triangle(noeud1->get_id(),noeud3->get_id(),noeud2->get_id());
611     if (triangle!=NULL) return NULL;
612     // insertion du triangle
613     MG_SEGMENT* segment1=mg_maillage->get_mg_segment(noeud1->get_id(),noeud2->get_id());
614     MG_SEGMENT* segment2=mg_maillage->get_mg_segment(noeud1->get_id(),noeud3->get_id());
615     MG_SEGMENT* segment3=mg_maillage->get_mg_segment(noeud2->get_id(),noeud3->get_id());
616     if (segment1==NULL) segment1=cree_segment(mgvol,noeud1,noeud2,origine);
617     if (segment2==NULL) segment2=cree_segment(mgvol,noeud1,noeud3,origine);
618     if (segment3==NULL) segment3=cree_segment(mgvol,noeud2,noeud3,origine);
619     MG_TRIANGLE* tri=cree_triangle(mgvol,noeud1,noeud2,noeud3,segment1,segment2,segment3,origine);
620     return tri;
621     }
622    
623    
624    
625    
626    
627    
628     MG_TRIANGLE* MAILLEUR3D_OPTIMISATION::cree_triangle(class MG_VOLUME* mgvol,MG_NOEUD* noeud1,MG_NOEUD* noeud2,MG_NOEUD* noeud3,MG_SEGMENT* segment1,MG_SEGMENT* segment2,MG_SEGMENT* segment3,int origine)
629     {
630     M3D_TRIANGLE* tri=new M3D_TRIANGLE(mgvol,noeud1,noeud2,noeud3,segment1,segment2,segment3,origine);
631     mg_maillage->ajouter_mg_triangle(tri);
632     return tri;
633     }
634    
635    
636    
637     MG_SEGMENT* MAILLEUR3D_OPTIMISATION::cree_segment(class MG_VOLUME* mgvol,MG_NOEUD* noeud1,MG_NOEUD* noeud2,int origine)
638     {
639     MG_SEGMENT* seg=mg_maillage->ajouter_mg_segment(mgvol,noeud1,noeud2,origine);
640 francois 288 seg->change_dimension_topo_null(3);
641 francois 287 return seg;
642     }