ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/REPOS_ERICCA/acismesh/o3d_move.cpp
Revision: 1
Committed: Mon Jun 11 22:53:07 2007 UTC (17 years, 11 months ago)
File size: 6372 byte(s)
Log Message:

File Contents

# User Rev Content
1 1 /*****************************************************************
2    
3     m3d_move.cpp Type:Func
4    
5     bouge de point 3D
6    
7     Date de creation : 6-2-1998 11 :22 :51
8     Derniere version : 6-2-1998 11 :22 :51
9    
10     Vincent FRANCOIS
11    
12     *****************************************************************/
13    
14    
15    
16    
17    
18     /**************************/
19     /* include */
20     #include <stdio.h>
21     #include <math.h>
22     #include "const.h"
23     #include "memoire.h"
24     #include "struct.h"
25     #include "struct3d.h"
26     #include "prototype.h"
27    
28     /**************************/
29     /* variables globales */
30     extern struct environnement env;
31     extern struct s_mesh *mesh;
32    
33    
34    
35     /**************************/
36     /* programme principal */
37    
38    
39    
40     void o3d_move(int nump,float *coord_res,float *crit)
41     {
42     int n1,n2,n3 ;
43     struct s_noeud *no,*no1,*no2,*no3;
44     struct s_tetra *tet;
45     int i, j ;
46     float tab_crit[1000], vec_opt[3], cmin ;
47     float tab_opt[3000] ;
48     float gamma, vn[3], xi, yi, zi, vab[3], vac[3], vbc[3], hauteur ;
49     float perimetre, coeff, alpha ;
50     float delta[3] ;
51     float vec_save[3], bmin, bmax, eps, crit_save ;
52     float calpha, calpha_moins_eps, vresu[3] ;
53    
54    
55     no=ADRESSE(nump,noeud,mesh->);
56     if (no->type!=BODY) return;
57     /* sauvegarde des coordonnees du noeud */
58     vec_save[0] = no->x;
59     vec_save[1] = no->y;
60     vec_save[2] = no->z;
61    
62     /* tetraedres connectes au noeud */
63     for (i=0;i<no->nb_tetra;i++)
64     {
65     tet=no->tetra[i];
66     if (tet->etat!=ACTIF) continue;
67     if (tet->n1==no->num)
68     {
69     n1=tet->n2;
70     n2=tet->n4;
71     n3=tet->n3;
72     }
73     if (tet->n2==no->num)
74     {
75     n1=tet->n1;
76     n2=tet->n3;
77     n3=tet->n4;
78     }
79     if (tet->n3==no->num)
80     {
81     n1=tet->n1;
82     n2=tet->n4;
83     n3=tet->n2;
84     }
85     if (tet->n4==no->num)
86     {
87     n1=tet->n1;
88     n2=tet->n2;
89     n3=tet->n3;
90     }
91     no1=ADRESSE(n1,noeud,mesh->);
92     no2=ADRESSE(n2,noeud,mesh->);
93     no3=ADRESSE(n3,noeud,mesh->);
94     /* determination du point optimal de la face nb_ele */
95     xi = (no1->x+no2->x+no3->x)/3. ;
96     yi = (no1->y+no2->y+no3->y)/3. ;
97     zi = (no1->z+no2->z+no3->z)/3. ;
98    
99     /* calcul de la normale a la face */
100     VEC(vab,no1,no2);
101     NORME(vab) ;
102     VEC(vac,no1,no3);
103     NORME(vac) ;
104     VEC(vbc,no2,no3);
105     NORME(vbc) ;
106    
107     /* calcul du critere 2 D */
108     perimetre = vab[3] + vbc[3] + vbc[3] ;
109     hauteur = (perimetre/3.) * 0.8 ;
110     PVEC(vn,vab,vac);
111     NORME(vn) ;
112     /* calcul du vecteur AP et orientation de la face */
113     /* point optimal */
114     tab_opt[3*i] = xi + hauteur * vn[0] ;
115     tab_opt[3*i+1] = yi + hauteur * vn[1] ;
116     tab_opt[3*i+2] = zi + hauteur * vn[2] ;
117     }
118     /* calcul de la qualite de chaque tetraedre */
119     cmin = 1. ;
120     for (i=0;i<no->nb_tetra;i++)
121     {
122     tet=no->tetra[i];
123     /* calcul du critere */
124     if (tet->etat==ACTIF)
125     {
126     tab_crit[i] = m3d_cal_qual(tet->n1,tet->n2,tet->n3,tet->n4) ;
127     if (tab_crit[i] < cmin) cmin = tab_crit[i] ;
128     }
129     else tab_crit[i]=0.;
130     }
131     /* calcul des coordonnees du point ideal */
132     gamma = 0. ;
133     for (i=0;i<no->nb_tetra;i++)
134     if (tab_crit[i]!=0.) gamma = gamma + 1./(tab_crit[i]*tab_crit[i]) ;/* somme des carres des inverses */
135     /* calcul du point optimal pour la boule */
136     gamma = 1./gamma ;
137     vec_opt[0] = 0. ;
138     vec_opt[1] = 0. ;
139     vec_opt[2] = 0. ;
140    
141     for (i=0;i<no->nb_tetra;i++)
142     {
143     if (tab_crit[i]!=0.)
144     {
145     coeff = gamma/(tab_crit[i]*tab_crit[i]) ;
146     vec_opt[0] = vec_opt[0] + coeff * tab_opt[3*i] ;
147     vec_opt[1] = vec_opt[1] + coeff * tab_opt[3*i+1] ;
148     vec_opt[2] = vec_opt[2] + coeff * tab_opt[3*i+2] ;
149     }
150     }
151    
152     delta[0] = vec_opt[0] - no->x ;
153     delta[1] = vec_opt[1] - no->y ;
154     delta[2] = vec_opt[2] - no->z ;
155     /* sauvegarde du critere mini */
156     crit_save = cmin ;
157    
158     /* initialisation des bornes du parametre */
159     bmin = 0. ;
160     bmax = 1. ;
161    
162     vresu[0] = vec_save[0] ;
163     vresu[1] = vec_save[1] ;
164     vresu[2] = vec_save[2] ;
165    
166     for (i=0;i<5;i++)
167     {
168     alpha = 0.5 * (bmin + bmax) ;
169     no->x = vec_save[0] + alpha * delta[0] ;
170     no->y = vec_save[1] + alpha * delta[1] ;
171     no->z = vec_save[2] + alpha * delta[2] ;
172    
173     /* evaluation de la qualite pour la valeur alpha */
174     /* on evalue la qualite cmin pour alpha */
175     calpha = 1. ;
176     for (j=0;j<no->nb_tetra;j++)
177     {
178     tet=no->tetra[j];
179     if (tet->etat!=ACTIF) continue;
180     /* calcul du critere */
181     tab_crit[j] = m3d_cal_qual(tet->n1,tet->n2,tet->n3,tet->n4) ;
182     if (tab_crit[j] < calpha) calpha = tab_crit[j] ;
183     }
184     /* la valeur du critere pour alpha est superieure, il faut la stocker */
185     if (calpha > crit_save)
186     {
187     vresu[0] = no->x ;
188     vresu[1] = no->y ;
189     vresu[2] = no->z ;
190     crit_save = calpha ;
191     }
192    
193     /* on evalue la qualite c1 pour alpha - eps */
194     eps = (bmax-bmin)/50. ;
195     no->x = vec_save[0] + (alpha - eps) * delta[0] ;
196     no->y = vec_save[1] + (alpha - eps) * delta[1] ;
197     no->z = vec_save[2] + (alpha - eps) * delta[2] ;
198     calpha_moins_eps = 1. ;
199     for (j=0;j<no->nb_tetra;j++)
200     {
201     tet=no->tetra[j];
202     if (tet->etat!=ACTIF) continue;
203     /* calcul du critere */
204     tab_crit[j] = m3d_cal_qual(tet->n1,tet->n2,tet->n3,tet->n4) ;
205     if (tab_crit[j] < calpha_moins_eps) calpha_moins_eps = tab_crit[j] ;
206     }
207    
208     no->x = vec_save[0] + alpha * delta[0] ;
209     no->y = vec_save[1] + alpha * delta[1] ;
210     no->z = vec_save[2] + alpha * delta[2] ;
211    
212     if (calpha_moins_eps < calpha)
213     {
214     bmin = alpha ;
215     /* bmax est inchange */
216     }
217     else
218     {
219     /* bmin est inchange */
220     bmax = alpha ;
221     }
222     }
223    
224     coord_res[0] = vresu[0] ;
225     coord_res[1] = vresu[1] ;
226     coord_res[2] = vresu[2] ;
227    
228     no->x=vec_save[0];
229     no->y=vec_save[1];
230     no->z=vec_save[2];
231    
232     *crit=crit_save;
233    
234     }
235    
236    
237    
238    
239    
240    
241