ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/REPOS_ERICCA/magic/lib/mailleur/src/mailleur3d_intersection.cpp
Revision: 469
Committed: Wed Nov 27 20:02:22 2013 UTC (11 years, 5 months ago) by francois
Original Path: magic/lib/mailleur_auto/src/mailleur3d_intersection.cpp
File size: 20183 byte(s)
Log Message:
ajout d'un generateur de carte de taille constante pour les nouvelles cartes et decoupage de la librairie mailleur en 2 pour cause d'appel cyclique

File Contents

# User Rev Content
1 5 //------------------------------------------------------------
2     //------------------------------------------------------------
3     // MAGiC
4     // Jean Christophe Cuillière et Vincent FRANCOIS
5     // Département de Génie Mécanique - UQTR
6     //------------------------------------------------------------
7     // Le projet MAGIC est un projet de recherche du département
8     // de génie mécanique de l'Université du Québec à
9     // Trois Rivières
10     // Les librairies ne peuvent être utilisées sans l'accord
11     // des auteurs (contact : francois@uqtr.ca)
12     //------------------------------------------------------------
13     //------------------------------------------------------------
14     //
15     // mailleur3d_intersection.cpp
16     //
17     //------------------------------------------------------------
18     //------------------------------------------------------------
19     // COPYRIGHT 2000
20     // Version du 02/03/2006 à 11H23
21     //------------------------------------------------------------
22     //------------------------------------------------------------
23    
24    
25     #include "gestionversion.h"
26    
27    
28     //#pragma hdrstop
29    
30     #include <math.h>
31     #include "mailleur3d.h"
32     #include "ot_mathematique.h"
33    
34    
35    
36     int MAILLEUR3D::inter_tetra_triangle(MG_NOEUD* noeud1,MG_NOEUD* noeud2,MG_NOEUD* noeud3,MG_NOEUD* noeud4,MG_NOEUD* noeud5,MG_NOEUD* noeud6,MG_NOEUD* noeud7)
37     {
38    
39     // filtre sur numero
40     int commun1;
41     int commun2;
42     int commun3;
43     if ( (noeud5==noeud1) || (noeud5==noeud2) || (noeud5==noeud3) || (noeud5==noeud4) ) commun1=1; else commun1=0;
44     if ( (noeud6==noeud1) || (noeud6==noeud2) || (noeud6==noeud3) || (noeud6==noeud4) ) commun2=1; else commun2=0;
45     if ( (noeud7==noeud1) || (noeud7==noeud2) || (noeud7==noeud3) || (noeud7==noeud4) ) commun3=1; else commun3=0;
46     int nb_commun=commun1+commun2+commun3;
47     if (nb_commun==3) return(false);
48     // filtre sur position
49     int inter;
50     if (nb_commun>0) inter=true;
51     else
52     {
53     double bxmin=std::min(noeud1->get_x(),noeud2->get_x());
54     bxmin=std::min(bxmin,noeud3->get_x());
55     bxmin=std::min(bxmin,noeud4->get_x());
56     double bxmax=std::max(noeud1->get_x(),noeud2->get_x());
57     bxmax=std::max(bxmax,noeud3->get_x());
58     bxmax=std::max(bxmax,noeud4->get_x());
59     double bymin=std::min(noeud1->get_y(),noeud2->get_y());
60     bymin=std::min(bymin,noeud3->get_y());
61     bymin=std::min(bymin,noeud4->get_y());
62     double bymax=std::max(noeud1->get_y(),noeud2->get_y());
63     bymax=std::max(bymax,noeud3->get_y());
64     bymax=std::max(bymax,noeud4->get_y());
65     double bzmin=std::min(noeud1->get_z(),noeud2->get_z());
66     bzmin=std::min(bzmin,noeud3->get_z());
67     bzmin=std::min(bzmin,noeud4->get_z());
68     double bzmax=std::max(noeud1->get_z(),noeud2->get_z());
69     bzmax=std::max(bzmax,noeud3->get_z());
70     bzmax=std::max(bzmax,noeud4->get_z());
71     BOITE_3D boite_tetra(bxmin,bymin,bzmin,bxmax,bymax,bzmax);
72     double b2xmin=std::min(noeud5->get_x(),noeud6->get_x());
73     b2xmin=std::min(b2xmin,noeud7->get_x());
74     double b2xmax=std::max(noeud5->get_x(),noeud6->get_x());
75     b2xmax=std::max(b2xmax,noeud7->get_x());
76     double b2ymin=std::min(noeud5->get_y(),noeud6->get_y());
77     b2ymin=std::min(b2ymin,noeud7->get_y());
78     double b2ymax=std::max(noeud5->get_y(),noeud6->get_y());
79     b2ymax=std::max(b2ymax,noeud7->get_y());
80     double b2zmin=std::min(noeud5->get_z(),noeud6->get_z());
81     b2zmin=std::min(b2zmin,noeud7->get_z());
82     double b2zmax=std::max(noeud5->get_z(),noeud6->get_z());
83     b2zmax=std::max(b2zmax,noeud7->get_z());
84     BOITE_3D boite_triangle(b2xmin,b2ymin,b2zmin,b2xmax,b2ymax,b2zmax);
85     inter=boite_tetra*boite_triangle;
86     }
87     if (inter==true)
88     {
89     inter=false;
90     OT_VECTEUR_3D n1n2(noeud2->get_x()-noeud1->get_x(),noeud2->get_y()-noeud1->get_y(),noeud2->get_z()-noeud1->get_z());
91     OT_VECTEUR_3D n1n3(noeud3->get_x()-noeud1->get_x(),noeud3->get_y()-noeud1->get_y(),noeud3->get_z()-noeud1->get_z());
92     OT_VECTEUR_3D n1n4(noeud4->get_x()-noeud1->get_x(),noeud4->get_y()-noeud1->get_y(),noeud4->get_z()-noeud1->get_z());
93     OT_VECTEUR_3D n2n3(noeud3->get_x()-noeud2->get_x(),noeud3->get_y()-noeud2->get_y(),noeud3->get_z()-noeud2->get_z());
94     OT_VECTEUR_3D n2n4(noeud4->get_x()-noeud2->get_x(),noeud4->get_y()-noeud2->get_y(),noeud4->get_z()-noeud2->get_z());
95     OT_VECTEUR_3D n1n5(noeud5->get_x()-noeud1->get_x(),noeud5->get_y()-noeud1->get_y(),noeud5->get_z()-noeud1->get_z());
96     OT_VECTEUR_3D n1n6(noeud6->get_x()-noeud1->get_x(),noeud6->get_y()-noeud1->get_y(),noeud6->get_z()-noeud1->get_z());
97     OT_VECTEUR_3D n1n7(noeud7->get_x()-noeud1->get_x(),noeud7->get_y()-noeud1->get_y(),noeud7->get_z()-noeud1->get_z());
98     OT_VECTEUR_3D n2n5(noeud5->get_x()-noeud2->get_x(),noeud5->get_y()-noeud2->get_y(),noeud5->get_z()-noeud2->get_z());
99     OT_VECTEUR_3D n2n6(noeud6->get_x()-noeud2->get_x(),noeud6->get_y()-noeud2->get_y(),noeud6->get_z()-noeud2->get_z());
100     OT_VECTEUR_3D n2n7(noeud7->get_x()-noeud2->get_x(),noeud7->get_y()-noeud2->get_y(),noeud7->get_z()-noeud2->get_z());
101     OT_VECTEUR_3D normal1=n1n3&n1n2;
102     OT_VECTEUR_3D normal2=n1n2&n1n4;
103     OT_VECTEUR_3D normal3=n2n3&n2n4;
104     OT_VECTEUR_3D normal4=n1n4&n1n3;
105     double eps=0.33333333333*(n1n3.diff()+n1n2.diff()+n1n5.diff());
106 francois 61 eps=18.*eps*eps*eps*1e-6;
107 5 if (commun1==0)
108     if (normal1*n1n5<eps)
109     if (normal2*n1n5<eps)
110     if (normal3*n2n5<eps)
111     if (normal4*n1n5<eps) return true;
112     eps=0.33333333333*(n1n3.diff()+n1n2.diff()+n1n6.diff());
113 francois 61 eps=18.*eps*eps*eps*1e-6;
114 5 if (commun2==0)
115     if (normal1*n1n6<eps)
116     if (normal2*n1n6<eps)
117     if (normal3*n2n6<eps)
118     if (normal4*n1n6<eps) return true;
119     eps=0.33333333333*(n1n3.diff()+n1n2.diff()+n1n7.diff());
120 francois 61 eps=18.*eps*eps*eps*1e-6;
121 5 if (commun3==0)
122     if (normal1*n1n7<eps)
123     if (normal2*n1n7<eps)
124     if (normal3*n2n7<eps)
125     if (normal4*n1n7<eps) return true;
126     if (inter_triangle_triangle(noeud1,noeud3,noeud2,noeud5,noeud6,noeud7)==true) return true;
127     if (inter_triangle_triangle(noeud1,noeud2,noeud4,noeud5,noeud6,noeud7)==true) return true;
128     if (inter_triangle_triangle(noeud2,noeud3,noeud4,noeud5,noeud6,noeud7)==true) return true;
129     if (inter_triangle_triangle(noeud1,noeud4,noeud3,noeud5,noeud6,noeud7)==true) return true;
130     }
131    
132     return inter;
133     }
134    
135    
136     int MAILLEUR3D::inter_triangle_triangle(MG_NOEUD* noeud1,MG_NOEUD* noeud2,MG_NOEUD* noeud3,MG_NOEUD* noeud4,MG_NOEUD* noeud5,MG_NOEUD* noeud6)
137     {
138     // filtre sur numero
139     int nb_commun=0;
140     if ( (noeud4==noeud1) || (noeud4==noeud2) || (noeud4==noeud3) ) nb_commun++;
141     if ( (noeud5==noeud1) || (noeud5==noeud2) || (noeud5==noeud3) ) nb_commun++;
142     if ( (noeud6==noeud1) || (noeud6==noeud2) || (noeud6==noeud3) ) nb_commun++;
143     if (nb_commun==3) return(false);
144     // filtre sur position
145     double bxmin=std::min(noeud1->get_x(),noeud2->get_x());
146     bxmin=std::min(bxmin,noeud3->get_x());
147     double bxmax=std::max(noeud1->get_x(),noeud2->get_x());
148     bxmax=std::max(bxmax,noeud3->get_x());
149     double bymin=std::min(noeud1->get_y(),noeud2->get_y());
150     bymin=std::min(bymin,noeud3->get_y());
151     double bymax=std::max(noeud1->get_y(),noeud2->get_y());
152     bymax=std::max(bymax,noeud3->get_y());
153     double bzmin=std::min(noeud1->get_z(),noeud2->get_z());
154     bzmin=std::min(bzmin,noeud3->get_z());
155     double bzmax=std::max(noeud1->get_z(),noeud2->get_z());
156     bzmax=std::max(bzmax,noeud3->get_z());
157     BOITE_3D boite_triangle1(bxmin,bymin,bzmin,bxmax,bymax,bzmax);
158     double b2xmin=std::min(noeud5->get_x(),noeud6->get_x());
159     b2xmin=std::min(b2xmin,noeud4->get_x());
160     double b2xmax=std::max(noeud5->get_x(),noeud6->get_x());
161     b2xmax=std::max(b2xmax,noeud4->get_x());
162     double b2ymin=std::min(noeud5->get_y(),noeud6->get_y());
163     b2ymin=std::min(b2ymin,noeud4->get_y());
164     double b2ymax=std::max(noeud5->get_y(),noeud6->get_y());
165     b2ymax=std::max(b2ymax,noeud4->get_y());
166     double b2zmin=std::min(noeud5->get_z(),noeud6->get_z());
167     b2zmin=std::min(b2zmin,noeud4->get_z());
168     double b2zmax=std::max(noeud5->get_z(),noeud6->get_z());
169     b2zmax=std::max(b2zmax,noeud4->get_z());
170     BOITE_3D boite_triangle2(b2xmin,b2ymin,b2zmin,b2xmax,b2ymax,b2zmax);
171     int inter=boite_triangle1*boite_triangle2;
172     if (inter==true)
173     {
174     inter=false;
175     OT_VECTEUR_3D n1n2(noeud2->get_x()-noeud1->get_x(),noeud2->get_y()-noeud1->get_y(),noeud2->get_z()-noeud1->get_z());
176     OT_VECTEUR_3D n1n3(noeud3->get_x()-noeud1->get_x(),noeud3->get_y()-noeud1->get_y(),noeud3->get_z()-noeud1->get_z());
177     OT_VECTEUR_3D n1n4(noeud4->get_x()-noeud1->get_x(),noeud4->get_y()-noeud1->get_y(),noeud4->get_z()-noeud1->get_z());
178     OT_VECTEUR_3D n1n5(noeud5->get_x()-noeud1->get_x(),noeud5->get_y()-noeud1->get_y(),noeud5->get_z()-noeud1->get_z());
179     OT_VECTEUR_3D n1n6(noeud6->get_x()-noeud1->get_x(),noeud6->get_y()-noeud1->get_y(),noeud6->get_z()-noeud1->get_z());
180     OT_VECTEUR_3D normal=n1n2&n1n3;
181     double ps1=normal*n1n4;
182     double ps2=normal*n1n5;
183     double ps3=normal*n1n6;
184     double eps1=0.333333333*(n1n2.diff()+n1n3.diff()+n1n4.diff());
185     double eps2=0.333333333*(n1n2.diff()+n1n3.diff()+n1n5.diff());
186     double eps3=0.333333333*(n1n2.diff()+n1n3.diff()+n1n6.diff());
187 francois 61 eps1=18.*eps1*eps1*eps1*1e-6;
188     eps2=18.*eps2*eps2*eps2*1e-6;
189     eps3=18.*eps3*eps3*eps3*1e-6;
190 5 if ((ps1<-eps1) && (ps2<-eps2) && (ps3<-eps3) ) return false;
191     if ((ps1>eps1) && (ps2>eps2) && (ps3>eps3) ) return false;
192     if (inter_segment_triangle(noeud1,noeud2,noeud3,noeud4,noeud5)==true) return true;
193     if (inter_segment_triangle(noeud1,noeud2,noeud3,noeud5,noeud6)==true) return true;
194     if (inter_segment_triangle(noeud1,noeud2,noeud3,noeud4,noeud6)==true) return true;
195     if (inter_segment_triangle(noeud4,noeud5,noeud6,noeud1,noeud2)==true) return true;
196     if (inter_segment_triangle(noeud4,noeud5,noeud6,noeud1,noeud3)==true) return true;
197     if (inter_segment_triangle(noeud4,noeud5,noeud6,noeud2,noeud3)==true) return true;
198     }
199     return inter;
200    
201    
202     }
203    
204     int MAILLEUR3D::inter_segment_triangle(MG_NOEUD* noeud1,MG_NOEUD* noeud2,MG_NOEUD* noeud3,MG_NOEUD* noeud4,MG_NOEUD* noeud5)
205     {
206     // filtre sur numero
207     int nb_commun=0;
208     if ( (noeud4==noeud1) || (noeud4==noeud2) || (noeud4==noeud3) ) nb_commun++;
209     if ( (noeud5==noeud1) || (noeud5==noeud2) || (noeud5==noeud3) ) nb_commun++;
210     if (nb_commun==2) return(false);
211    
212     int inter=false;
213     OT_VECTEUR_3D g1(noeud2->get_x()-noeud1->get_x(),noeud2->get_y()-noeud1->get_y(),noeud2->get_z()-noeud1->get_z());
214     OT_VECTEUR_3D g2(noeud3->get_x()-noeud1->get_x(),noeud3->get_y()-noeud1->get_y(),noeud3->get_z()-noeud1->get_z());
215     OT_VECTEUR_3D g3(noeud5->get_x()-noeud4->get_x(),noeud5->get_y()-noeud4->get_y(),noeud5->get_z()-noeud4->get_z());
216     OT_MATRICE_3D systeme(g1,g2,g3);
217     double det=systeme.get_determinant();
218     double eps=0.333333333333*(g1.diff()+g2.diff()+g3.diff());
219     eps=eps*eps*0.0018;
220 francois 35 if (OPERATEUR::egal(det,0.0,eps)==true)
221 5 {
222 francois 35 OT_VECTEUR_3D g3b(noeud5->get_x()-noeud1->get_x(),noeud5->get_y()-noeud1->get_y(),noeud5->get_z()-noeud1->get_z());
223     OT_MATRICE_3D systeme2(g1,g2,g3b);
224     double det2=systeme2.get_determinant();
225     double eps2=0.333333333333*(g1.diff()+g2.diff()+g3b.diff());
226 francois 61 eps2=18.*eps2*eps2*eps2*1e-6;
227 francois 35 if (OPERATEUR::egal(det2,0.0,eps)==true) // cas 2D
228     {
229     if (inter_segment_segment(noeud1,noeud2,noeud4,noeud5)==true) return true;
230     if (inter_segment_segment(noeud1,noeud3,noeud4,noeud5)==true) return true;
231     if (inter_segment_segment(noeud2,noeud3,noeud4,noeud5)==true) return true;
232     }
233     else return false;
234 5 }
235     else
236     {
237     if ( (noeud4==noeud1) || (noeud4==noeud2) || (noeud4==noeud3) ) return false;
238     if ( (noeud5==noeud1) || (noeud5==noeud2) || (noeud5==noeud3) ) return false;
239     double x,y,z;
240     x=1.0/det*(g2.get_y()*g3.get_z()-g2.get_z()*g3.get_y());
241     y=1.0/det*(g3.get_x()*g2.get_z()-g2.get_x()*g3.get_z());
242     z=1.0/det*(g2.get_x()*g3.get_y()-g2.get_y()*g3.get_x());
243     OT_VECTEUR_3D g1b(x,y,z);
244     x=1.0/det*(g3.get_y()*g1.get_z()-g1.get_y()*g3.get_z());
245     y=1.0/det*(g1.get_x()*g3.get_z()-g3.get_x()*g1.get_z());
246     z=1.0/det*(g3.get_x()*g1.get_y()-g1.get_x()*g3.get_y());
247     OT_VECTEUR_3D g2b(x,y,z);
248     x=1.0/det*(g1.get_y()*g2.get_z()-g1.get_z()*g2.get_y());
249     y=1.0/det*(g2.get_x()*g1.get_z()-g1.get_x()*g2.get_z());
250     z=1.0/det*(g1.get_x()*g2.get_y()-g1.get_y()*g2.get_x());
251     OT_VECTEUR_3D g3b(x,y,z);
252     OT_VECTEUR_3D n1n4(noeud4->get_x()-noeud1->get_x(),noeud4->get_y()-noeud1->get_y(),noeud4->get_z()-noeud1->get_z());
253     double alpha1=n1n4*g1b;
254     double alpha2=n1n4*g2b;
255     double alpha3=-(n1n4*g3b);
256     double alpha4=1-alpha1-alpha2;
257 francois 61 double eps=0.000001;
258 5 if ((alpha1>-eps) && (alpha1<1.+eps))
259     if ((alpha2>-eps) && (alpha2<1.+eps))
260     if ((alpha3>-eps) && (alpha3<1.+eps))
261     if ((alpha4>-eps) && (alpha4<1.+eps)) return true;
262     return false;
263     }
264     return inter;
265     }
266    
267    
268    
269     int MAILLEUR3D::noeud_est_triangle(MG_NOEUD* noeud1,MG_NOEUD* noeud2,MG_NOEUD* noeud3,double x,double y,double z)
270     {
271     OT_VECTEUR_3D n1n2(noeud2->get_x()-noeud1->get_x(),noeud2->get_y()-noeud1->get_y(),noeud2->get_z()-noeud1->get_z());
272     OT_VECTEUR_3D n1n3(noeud3->get_x()-noeud1->get_x(),noeud3->get_y()-noeud1->get_y(),noeud3->get_z()-noeud1->get_z());
273     OT_VECTEUR_3D n=n1n2&n1n3;
274     OT_VECTEUR_3D n1p(x-noeud1->get_x(),y-noeud1->get_y(),z-noeud1->get_z());
275     OT_VECTEUR_3D n2p(x-noeud2->get_x(),y-noeud2->get_y(),z-noeud2->get_z());
276     OT_VECTEUR_3D n3p(x-noeud3->get_x(),y-noeud3->get_y(),z-noeud3->get_z());
277     OT_VECTEUR_3D pv1=n1n2&n1p;
278     double ps1=pv1*n;
279     double eps=n1n2*n1n2;
280     eps=eps*eps*4*0.0001;
281     if (ps1>(-eps))
282     {
283     OT_VECTEUR_3D n2n3(noeud3->get_x()-noeud2->get_x(),noeud3->get_y()-noeud2->get_y(),noeud3->get_z()-noeud2->get_z());
284     OT_VECTEUR_3D pv2=n2n3&n2p;
285     double ps2=pv2*n;
286     if (ps2>(-eps))
287     {
288     OT_VECTEUR_3D pv3=n1n3&n3p;
289     double ps3=(-1)*(pv3*n);
290     if (ps3>(-eps)) return true;
291     }
292    
293     }
294     return false;
295     }
296    
297    
298    
299    
300    
301    
302     // PARTIE DE CODE RECUPERER DE MON PROJET DE DEA
303     // TRES ILLISIBLE MAIS FONCTIONNEL
304     // POIL AU MARRON
305    
306     #define PSCA(a,b) (a[0]*b[0]+a[1]*b[1]+a[2]*b[2])
307     #define EGAL(x,y,eps) (float)fabs((double)(x-y))<eps
308     #define DETER(a,b,c,d) (a*d-b*c)
309    
310     int MAILLEUR3D::inter_segment_segment(MG_NOEUD* noeud1,MG_NOEUD* noeud2,MG_NOEUD* noeud3,MG_NOEUD* noeud4)
311     {
312    
313     double ab[3];
314     double nm[3];
315     double am[3];
316     ab[0]=noeud2->get_x()-noeud1->get_x();
317     ab[1]=noeud2->get_y()-noeud1->get_y();
318     ab[2]=noeud2->get_z()-noeud1->get_z();
319     nm[0]=noeud3->get_x()-noeud4->get_x();
320     nm[1]=noeud3->get_y()-noeud4->get_y();
321     nm[2]=noeud3->get_z()-noeud4->get_z();
322     am[0]=noeud3->get_x()-noeud1->get_x();
323     am[1]=noeud3->get_y()-noeud1->get_y();
324     am[2]=noeud3->get_z()-noeud1->get_z();
325     int equation[4];
326     equation[0]=1; /* etat de l'equation 0 */
327     equation[1]=1;
328     equation[2]=1;
329     equation[3]=3; /* cette variable comporte le bilan du nombre d'equation */
330     double eps2=PSCA(ab,ab);
331     double eps=sqrt(eps2);
332     eps=eps*0.0001;
333     eps2=eps2*0.0001;
334     /* recherche du nombre d'equation -> inter franche ou para ou confondu */
335     if ( (EGAL(ab[0],0,eps)) && (EGAL(nm[0],0,eps)) )
336     if (EGAL(am[0],0,eps)) equation[0]=0; else return false;
337     if ( (EGAL(ab[1],0,eps)) && (EGAL(nm[1],0,eps)) )
338     if (EGAL(am[1],0,eps)) equation[1]=0; else return false;
339     if ( (EGAL(ab[2],0,eps)) && (EGAL(nm[2],0,eps)) )
340     if (EGAL(am[2],0,eps)) equation[2]=0; else return false;
341     equation[3]=equation[0]+equation[1]+equation[2];
342     if (equation[3]==3)
343     {
344     double det=DETER(ab[0],nm[0],ab[1],nm[1]);
345     if (fabs(det)>eps2)
346     {
347     det=1/det;
348     double sol1=det*DETER(am[0],nm[0],am[1],nm[1]);
349     double sol2=det*DETER(ab[0],am[0],ab[1],am[1]);
350     if ( (float)fabs((double)(sol1*ab[2]-sol2*nm[2]-am[2]))>eps2) return false;
351     return(examine_solution(sol1,sol2,1));
352     }
353     else
354     {
355     equation[0]=0;
356     equation[3]=2;
357     /* on verifie la compatibilite des deux equations dont le det est nul*/
358     double tmp;
359     if (ab[0]!=0) tmp=ab[1]*am[0]/ab[0]; else tmp=nm[1]*am[0]/nm[0];
360     if (!(EGAL(tmp,am[1],eps))) return false;
361     }
362     }
363     if (equation[3]==2)
364     {
365     /* on repere les equations qui existent */
366     int ne1;
367     int ne2;
368     if (equation[0]!=0)
369     {
370     ne1=0;
371     if (equation[1]!=0) ne2=1; else ne2=2;
372     }
373     else
374     {
375     ne1=1;
376     ne2=2;
377     }
378    
379     double det=DETER(ab[ne1],nm[ne1],ab[ne2],nm[ne2]);
380     if (fabs(det)>eps2)
381     {
382     det=1/det;
383     double sol1=det*DETER(am[ne1],nm[ne1],am[ne2],nm[ne2]);
384     double sol2=det*DETER(ab[ne1],am[ne1],ab[ne2],am[ne2]);
385     return(examine_solution(sol1,sol2,1));
386     }
387     else
388     {
389     equation[ne1]=0;
390     equation[3]=1;
391     /* on verifie la compatibilite des deux equations dont le det est nul */
392     double tmp;
393     if (ab[ne1]!=0) tmp=ab[ne2]*am[ne1]/ab[ne1]; else tmp=nm[ne2]*am[ne1]/nm[ne1];
394     if (!(EGAL(tmp,am[ne2],eps))) return false;
395     }
396    
397     }
398     if (equation[3]==1)
399     {
400     /* on repere l' equation qui existe */
401     int ne1;
402     if (equation[0]!=0) ne1=0; else
403     if (equation[1]!=0) ne1=1; else ne1=2;
404     double an[3];
405     an[0]=noeud4->get_x()-noeud1->get_x();
406     an[1]=noeud4->get_y()-noeud1->get_y();
407     an[2]=noeud4->get_z()-noeud1->get_z();
408     double tmp=1/ab[ne1];
409     double sol1=am[ne1]*tmp;
410     double sol2=an[ne1]*tmp;
411     return(examine_solution(sol1,sol2,2));
412     }
413     return 0;
414     }
415    
416    
417    
418     int MAILLEUR3D::examine_solution(double sol1,double sol2,int type)
419     {
420     double epsilon=0.0001;
421    
422     if (type==1)
423     {
424     if ( (sol1>epsilon) && ((sol1)<(1-epsilon)) && (sol2>epsilon) && ((sol2)<(1-epsilon)) ) return true;
425     if ( ( (EGAL(sol1,0,epsilon)) || (EGAL(sol1,1,epsilon))) && ( (sol2>epsilon) && ((sol2)<(1-epsilon)) ) ) return true;
426     if ( ( (EGAL(sol2,0,epsilon)) || (EGAL(sol2,1,epsilon))) && ( (sol1>epsilon) && ((sol1)<(1-epsilon)) ) ) return true;
427     if ( (sol1>epsilon) && ((sol1)<(1-epsilon)) && (sol2>(-0.1-epsilon)) && ((sol2)<(1.1-epsilon)) ) return true;
428     if ( (sol2>epsilon) && ((sol2)<(1-epsilon)) && (sol1>(-0.1-epsilon)) && ((sol1)<(1.1-epsilon)) ) return true;
429    
430     }
431     if (type==2)
432     {
433     if ( (sol1>epsilon) && ((sol1)<(1-epsilon)) ) return true;
434     if ( (sol2>epsilon) && ((sol2)<(1-epsilon)) ) return true;
435     if ( ((sol1)>(1+epsilon)) && ((-sol2)>epsilon) ) return true;
436     if ( ((sol2)>(1+epsilon)) && ((-sol1)>epsilon) ) return true;
437     }
438     return false;
439     }
440     #undef EGAL
441     #undef PSCA
442     #undef DETER
443     ///FIN DU CODE DE DEA
444     //---------------------------------------------------------------------------
445    
446     //#pragma package(smart_init)