ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/REPOS_ERICCA/magic/lib/geometrie/src/mg_arete.cpp
Revision: 672
Committed: Mon Jun 1 15:37:21 2015 UTC (9 years, 11 months ago) by francois
File size: 17473 byte(s)
Log Message:
bug sur la sonde pour les aretes periodiques+ bug de nom de fichier sur la comparaison de valeur entre deux champs d'une solution

File Contents

# User Rev Content
1 francois 283 //------------------------------------------------------------
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     // mg_arete.cpp
16     //
17     //------------------------------------------------------------
18     //------------------------------------------------------------
19     // COPYRIGHT 2000
20     // Version du 02/03/2006 � 11H22
21     //------------------------------------------------------------
22     //------------------------------------------------------------
23    
24    
25     #include "gestionversion.h"
26     #include "mg_arete.h"
27     #include "vct_arete.h"
28     //#include "message.h"
29     //#include "affiche.h"
30 francois 375 #include "mg_definition.h"
31 francois 283 #include <math.h>
32     #include "ot_mathematique.h"
33     #include "constantegeo.h"
34 francois 604 #include "fem_solution.h"
35     #include "fem_maillage.h"
36 francois 283 #define PI 3.1415926535897932384626433832795
37    
38     MG_ARETE::MG_ARETE(std::string idori,unsigned long num,MG_COSOMMET* mgcosom1,MG_COSOMMET* mgcosom2,MG_COURBE* crb,int sens):MG_ELEMENT_TOPOLOGIQUE(num,idori),cosommet1(mgcosom1),cosommet2(mgcosom2),courbe(crb),orientation(sens),vect(NULL)
39     {
40     }
41    
42     MG_ARETE::MG_ARETE(std::string idori,MG_COSOMMET* mgcosom1,MG_COSOMMET* mgcosom2,MG_COURBE* crb,int sens):MG_ELEMENT_TOPOLOGIQUE(idori),cosommet1(mgcosom1),cosommet2(mgcosom2),courbe(crb),orientation(sens),vect(NULL)
43     {
44     }
45    
46     MG_ARETE::MG_ARETE(MG_ARETE& mdd):MG_ELEMENT_TOPOLOGIQUE(mdd),cosommet1(mdd.cosommet1),cosommet2(mdd.cosommet2),courbe(mdd.courbe),orientation(mdd.orientation),vect(NULL)
47     {
48     }
49    
50     MG_ARETE::MG_ARETE(std::string idori,unsigned long num,class MG_COURBE* crb,int sens):MG_ELEMENT_TOPOLOGIQUE(num,idori),courbe(crb),orientation(sens),vect(NULL)
51     {
52     }
53     MG_ARETE::MG_ARETE(std::string idori,class MG_COURBE* crb,int sens):MG_ELEMENT_TOPOLOGIQUE(idori),courbe(crb),orientation(sens),vect(NULL)
54     {
55     }
56     void MG_ARETE::get_topologie_sousjacente(TPL_MAP_ENTITE<MG_ELEMENT_TOPOLOGIQUE*> *lst)
57     {
58     MG_SOMMET* som1=cosommet1->get_sommet();
59     MG_SOMMET* som2=cosommet2->get_sommet();
60     lst->ajouter(som1);
61     lst->ajouter(som2);
62     som1->get_topologie_sousjacente(lst);
63     som2->get_topologie_sousjacente(lst);
64     }
65    
66     MG_ARETE::~MG_ARETE()
67     {
68     // if (lst_coarete.size()!=0) afficheur << WARCOARETE << this->get_id()<< enderr;
69     segment.vide();
70     if (vect!=NULL) delete vect;
71     }
72    
73     void MG_ARETE::changer_cosommet1(class MG_COSOMMET* cosom)
74     {
75     cosommet1=cosom;
76     }
77     void MG_ARETE::changer_cosommet2(class MG_COSOMMET* cosom)
78     {
79     cosommet2=cosom;
80     }
81     MG_COSOMMET* MG_ARETE::get_cosommet1(void)
82     {
83     return cosommet1;
84     }
85     MG_COSOMMET* MG_ARETE::get_cosommet2(void)
86     {
87     return cosommet2;
88     }
89     MG_COURBE* MG_ARETE::get_courbe(void)
90     {
91     return courbe;
92     }
93     int MG_ARETE::get_orientation(void)
94     {
95     return orientation;
96     }
97    
98    
99    
100    
101     void MG_ARETE::ajouter_mg_coarete(class MG_COARETE* coarete)
102     {
103     lst_coarete.insert(lst_coarete.end(),coarete);
104     }
105    
106     int MG_ARETE::get_nb_mg_coarete(void)
107     {
108     return lst_coarete.size();
109     }
110    
111 francois 576 bool MG_ARETE::est_une_arete_element(void)
112     {
113     return false;
114     }
115 francois 283
116     void MG_ARETE::evaluer(double t,double *xyz)
117     {
118     if (orientation!=MEME_SENS) t=courbe->get_tmin()+courbe->get_tmax()-t;
119     courbe->evaluer(t,xyz);
120     }
121    
122     void MG_ARETE::deriver(double t,double *xyz)
123     {
124     if (orientation!=MEME_SENS) t=courbe->get_tmin()+courbe->get_tmax()-t;
125     courbe->deriver(t,xyz);
126     if (orientation!=MEME_SENS)
127     {
128     xyz[0]=-xyz[0];
129     xyz[1]=-xyz[1];
130     xyz[2]=-xyz[2];
131     }
132     }
133    
134     void MG_ARETE::deriver_seconde(double t,double *ddxyz,double* dxyz,double* xyz)
135     {
136     if (orientation!=MEME_SENS) t=courbe->get_tmin()+courbe->get_tmax()-t;
137     courbe->deriver_seconde(t,ddxyz,dxyz,xyz);
138     if (orientation!=MEME_SENS)
139     {
140     dxyz[0]=-dxyz[0];
141     dxyz[1]=-dxyz[1];
142     dxyz[2]=-dxyz[2];
143     }
144     }
145    
146     void MG_ARETE::inverser(double& t,double *xyz,double precision)
147     {
148     courbe->inverser(t,xyz,precision);
149     if (orientation!=MEME_SENS) t=courbe->get_tmin()+courbe->get_tmax()-t;
150     }
151    
152     double MG_ARETE::get_tmin(void)
153     {
154     return cosommet1->get_t();
155     }
156    
157     double MG_ARETE::get_tmax(void)
158     {
159     if (!courbe->est_periodique()) return cosommet2->get_t();
160     double tmin=cosommet1->get_t();
161     double tmax=cosommet2->get_t();
162     if (tmax<tmin+1e-6) tmax=tmax+courbe->get_periode();
163     return tmax;
164     }
165    
166     double MG_ARETE::get_longueur(double t1,double t2,double precis)
167     {
168     if (orientation!=MEME_SENS)
169     {
170     t1=courbe->get_tmin()+courbe->get_tmax()-t1;
171     t2=courbe->get_tmin()+courbe->get_tmax()-t2;
172     return (-courbe->get_longueur(t1,t2,precis));
173     }
174     return courbe->get_longueur(t1,t2,precis);
175    
176     }
177    
178     double MG_ARETE::get_M(double t)
179     {
180     if (orientation!=MEME_SENS) t=courbe->get_tmin()+courbe->get_tmax()-t;
181     return courbe->get_M(t);
182     }
183    
184    
185    
186     void MG_ARETE::supprimer_mg_coarete(class MG_COARETE* coarete)
187     {
188     std::vector<MG_COARETE*>::iterator i;
189     for (i=lst_coarete.begin();i!=lst_coarete.end();i++)
190     {
191     if ((*i)==coarete)
192     {
193     lst_coarete.erase(i);
194     return;
195     }
196     }
197     }
198    
199    
200    
201     MG_COARETE* MG_ARETE::get_mg_coarete(int num)
202     {
203     return lst_coarete[num];
204     }
205    
206    
207 gervaislavoie 332 TPL_LISTE_ENTITE<class MG_SEGMENT*>* MG_ARETE::get_lien_segment(void)
208 francois 283 {
209     return &segment;
210     }
211    
212     int MG_ARETE::get_dimension(void)
213     {
214     return 1;
215     }
216    
217    
218     VCT& MG_ARETE::get_vectorisation(void)
219     {
220     if (vect==NULL) vect=new VCT_ARETE(this);
221     return *vect;
222     }
223    
224     void MG_ARETE::enregistrer(std::ostream& o)
225     {
226     int nb=get_nb_ccf();
227     o << "%" << get_id() << "=ARETE("<< get_idoriginal() << ",$" << courbe->get_id() << ",$"<<cosommet1->get_id() << ",$" <<cosommet2->get_id() << "," << orientation << "," << nb;
228     if (nb!=0)
229     {
230     o << ",(";
231     for (int i=0;i<nb;i++)
232     {
233     char nom[3];
234     get_type_ccf(i,nom);
235     o << "(" << nom << "," << get_valeur_ccf(i) << ")";
236     if (i!=nb-1) o << "," ;
237     }
238     o << ")";
239     }
240     o << ");" << std::endl;
241     }
242    
243    
244     void MG_ARETE::get_intersection(double t1, double t2,double* point_iners)
245     {
246     OT_VECTEUR_3D P1_INTERSECTION;
247    
248     double U,W;
249    
250     double xyz2[3];
251     double xyz1[3];
252     double drv1[3];
253     double drv2[3];
254    
255     evaluer(t1,xyz1);
256     deriver(t1,drv1);
257    
258     evaluer(t2,xyz2);
259     deriver(t2,drv2);
260    
261     OT_VECTEUR_3D a(xyz1) ;
262     OT_VECTEUR_3D b(drv1) ;
263     OT_VECTEUR_3D c(xyz2) ;
264     OT_VECTEUR_3D d(drv2) ;
265     OT_VECTEUR_3D f=c-a;
266     d=-d;
267     double det1=b[0]*d[1]-b[1]*d[0];
268     double det2=b[0]*d[2]-b[2]*d[0];
269     double det3=b[1]*d[2]-b[2]*d[1];
270    
271     bool c_fait=false;
272    
273     if (fabs(det1)>=1e-10)
274     {
275     double u,w;
276     u=(d[1]*f[0]-d[0]*f[1])/det1;
277     w=(b[0]*f[1]-b[1]*f[0])/det1;
278     double verif_x=a[2] +u*b[2];
279     double verif_y=c[2]-w*d[2];
280     if (fabs(verif_x-verif_y)<0.000001)
281     {
282     U=u ;
283     W=w ;
284     c_fait=true;
285     }
286    
287     }
288    
289     if (!c_fait&&fabs(det2)>=1e-10)
290     {
291     double u,w;
292     u=(d[2]*f[0]-d[0]*f[2])/det2;
293     w=(b[0]*f[2]-b[2]*f[0])/det2;
294     double verif_x=a[1]+u*b[1];
295     double verif_y=c[1]-w*d[1];
296     if (fabs(verif_x-verif_y)<0.000001)
297     {
298     U=u ;
299     W=w ;
300     c_fait=true;
301     }
302     }
303    
304     if (!c_fait&&fabs(det3)>=1e-10)
305     {
306     double u,w;
307     u=(d[2]*f[1]-d[1]*f[2])/det3;
308     w=(b[1]*f[2]-b[2]*f[1])/det3;
309     double verif_x=a[0] +u*b[0];
310     double verif_y=c[0]-w*d[0];
311     if (fabs(verif_x-verif_y)<0.000001)
312     {
313     U=u ;
314     W=w ;
315     c_fait=true;
316     }
317     }
318    
319     P1_INTERSECTION=a+U*b;
320    
321     point_iners[0]=P1_INTERSECTION[0];
322     point_iners[1]=P1_INTERSECTION[1];
323     point_iners[2]=P1_INTERSECTION[2];
324    
325     }
326    
327    
328    
329     void MG_ARETE:: get_param_NURBS(int& indx_premier_ptctr,TPL_LISTE_ENTITE<double> &param)
330     {
331    
332     TPL_LISTE_ENTITE<double> param1;
333     int type=courbe->get_type_geometrique(param1);
334    
335     double tcs1=get_tmin();
336     double tcs2=get_tmax();
337     double angle=tcs1-tcs2;
338    
339     if ((type==MGCo_CIRCLE||type==MGCo_ELLIPSE)&&(tcs1!=tcs2)&&fabs(angle)!=2*PI)
340     {
341    
342     double xyz[3];
343     param.ajouter(1);
344     param.ajouter(4);
345     param.ajouter(0);
346     param.ajouter(7);
347     param.ajouter(0);
348    
349     double arc_lengh= get_longueur(tcs1,tcs2);
350     double rayon=arc_lengh/fabs(angle);
351    
352     double dif=fabs(fabs(angle)-PI);
353     if (fabs(fabs(angle)-PI)<=1e-6||fabs(angle)-PI<0.) //ouverture inf ou egale a PI
354     {
355    
356    
357     double request_lengh=arc_lengh/2;
358    
359     param.ajouter(0);
360     param.ajouter(0);
361     param.ajouter(0);
362     param.ajouter(0.5);
363     param.ajouter(0.5);
364     param.ajouter(1);
365     param.ajouter(1);
366     param.ajouter(1);
367    
368    
369     double xyz1[3];
370     double xyz1_inters[3];
371     double xyz2_inters[3];
372     double xyz2[3];
373     double xyz3[3];
374     double t_mid=get_t(tcs1,tcs2,request_lengh);
375    
376     evaluer(tcs1,xyz1);
377    
378     param.ajouter(xyz1[0]);
379     param.ajouter(xyz1[1]);
380     param.ajouter(xyz1[2]);
381     param.ajouter(1);
382    
383     get_intersection(tcs1, t_mid,xyz1_inters);
384    
385     param.ajouter(xyz1_inters[0]);
386     param.ajouter(xyz1_inters[1]);
387     param.ajouter(xyz1_inters[2]);
388     param.ajouter(0.5);
389    
390     get_intersection(t_mid, tcs2,xyz2_inters);
391    
392     param.ajouter(xyz2_inters[0]);
393     param.ajouter(xyz2_inters[1]);
394     param.ajouter(xyz2_inters[2]);
395     param.ajouter(0.5);
396    
397     evaluer(tcs2,xyz3);
398    
399     param.ajouter(xyz3[0]);
400     param.ajouter(xyz3[1]);
401     param.ajouter(xyz3[2]);
402     param.ajouter(1);
403    
404     indx_premier_ptctr=13;
405    
406     }
407    
408     if (fabs(fabs(angle)-PI)>1e-6&&fabs(angle)-PI>0.)
409     {
410    
411     double xyz[3];
412    
413 francois 290 double lengh_at_mid=arc_lengh;
414     double request_lengh=lengh_at_mid/2.;
415     double tmid=get_t(tcs1,tcs2,request_lengh);
416     double tmid1=get_t(tcs1,tmid,request_lengh/2.);
417     double tmid2=get_t(tmid,tcs2,request_lengh/2.);
418 francois 283
419     param.ajouter(0);
420     param.ajouter(0);
421     param.ajouter(0);
422     param.ajouter(0.25);
423     param.ajouter(0.25);
424     param.ajouter(0.5);
425     param.ajouter(0.5);
426     param.ajouter(0.75);
427     param.ajouter(0.75);
428     param.ajouter(1);
429     param.ajouter(1);
430     param.ajouter(1);
431    
432    
433    
434     evaluer(tcs1,xyz);
435    
436     param.ajouter(xyz[0]);
437     param.ajouter(xyz[1]);
438     param.ajouter(xyz[2]);
439     param.ajouter(1);
440    
441 francois 290
442     get_intersection(tcs1,tmid1,xyz);
443 francois 283
444     param.ajouter(xyz[0]);
445     param.ajouter(xyz[1]);
446     param.ajouter(xyz[2]);
447     param.ajouter(0.5);
448    
449 francois 290 get_intersection(tmid1, tmid,xyz);
450 francois 283
451     param.ajouter(xyz[0]);
452     param.ajouter(xyz[1]);
453     param.ajouter(xyz[2]);
454     param.ajouter(0.5);
455    
456 francois 290 evaluer(tmid,xyz);
457 francois 283
458     param.ajouter(xyz[0]);
459     param.ajouter(xyz[1]);
460     param.ajouter(xyz[2]);
461     param.ajouter(1);
462    
463 francois 290
464     get_intersection(tmid,tmid2,xyz);
465 francois 283
466     param.ajouter(xyz[0]);
467     param.ajouter(xyz[1]);
468     param.ajouter(xyz[2]);
469     param.ajouter(0.5);
470    
471    
472 francois 290 get_intersection(tmid2,tcs2,xyz);
473 francois 283
474     param.ajouter(xyz[0]);
475     param.ajouter(xyz[1]);
476     param.ajouter(xyz[2]);
477     param.ajouter(0.5);
478    
479     evaluer(tcs2,xyz);
480    
481     param.ajouter(xyz[0]);
482     param.ajouter(xyz[1]);
483     param.ajouter(xyz[2]);
484     param.ajouter(1);
485    
486     indx_premier_ptctr=17;
487     }
488    
489     }
490    
491    
492     if ((type==MGCo_ELLIPSE||type==MGCo_CIRCLE)&&OPERATEUR::egal(fabs(angle),2*PI,0.000001))
493     courbe->get_param_NURBS(indx_premier_ptctr,param);
494    
495     if (type==MGCo_BSPLINE)
496     courbe->get_param_NURBS(indx_premier_ptctr,param);
497    
498     if (type==MGCo_LINE)
499     {
500    
501     param.ajouter(1);
502     param.ajouter(2);
503     param.ajouter(0);
504     param.ajouter(2);
505     param.ajouter(0);
506    
507    
508     param.ajouter(0);
509     param.ajouter(0);
510     param.ajouter(1);
511     param.ajouter(1);
512    
513     double xyz1[3];
514     double xyz2[3];
515    
516     double tcs1=get_tmin();
517     double tcs2=get_tmax();
518    
519     evaluer(tcs1,xyz1);
520    
521     param.ajouter(xyz1[0]);
522     param.ajouter(xyz1[1]);
523     param.ajouter(xyz1[2]);
524     param.ajouter(1);
525    
526     evaluer(tcs2,xyz2);
527    
528     param.ajouter(xyz2[0]);
529     param.ajouter(xyz2[1]);
530     param.ajouter(xyz2[2]);
531     param.ajouter(1);
532    
533     indx_premier_ptctr=9;
534     }
535    
536     }
537    
538     struct integrale {
539     double ti;
540     double li;
541     };
542    
543     double MG_ARETE:: get_t(double t1,double t2,double longueur_voulue,double dt)
544     {
545     double pt[3];
546     double t,ti,tii;
547    
548     //double longueur_voulue,longueur=0;
549     double longueur=0;
550     std::vector<integrale> tab_integrale;
551     //longueur_voulue= alpha *get_longueur(t1,t2);
552    
553     integrale pas;
554     ti=t1;
555     do
556     {
557     tii=ti+(t2-t1)*dt;
558    
559     t=0.7886751345*ti+0.2113248654*tii;
560     deriver(t,pt);
561    
562     double facteur=pow(pt[0],2)+pow(pt[1],2)+pow(pt[2],2);
563     longueur=longueur+0.5*(tii-ti)*sqrt(facteur);
564     t=0.7886751345*tii+0.2113248654*ti;
565     deriver(t,pt);
566     facteur=pow(pt[0],2)+pow(pt[1],2)+pow(pt[2],2);
567     longueur=longueur+0.5*(tii-ti)*sqrt(facteur);
568     integrale pas;
569     pas.ti=ti;
570     pas.li=longueur;
571     tab_integrale.insert(tab_integrale.end(),pas);
572     ti=tii;
573     }
574     while (fabs(tab_integrale.back().li)<fabs(longueur_voulue));
575    
576     //double t_trouve=0.7886751345*tab_integrale.back().ti+0.2113248654*tab_integrale.back().ti+(t2-t1)*dt;
577     double tn_1=tab_integrale.back().ti;
578     double ln_1=tab_integrale[tab_integrale.size()-2].li;
579     double tn=ti;
580     double ln=tab_integrale.back().li;
581    
582    
583     double t_trouve=tn_1+(tn-tn_1)/(ln-ln_1)*(longueur_voulue-ln_1);
584     return t_trouve ;
585     }
586    
587 francois 604 void MG_ARETE::recupere_resultat(class FEM_SOLUTION* sol,int numchamps,char* fichier)
588     {
589     FEM_MAILLAGE *mai=sol->get_maillage();
590     MG_SOMMET* som1=cosommet1->get_sommet();
591     MG_SOMMET* som2=cosommet2->get_sommet();
592     sol->active_solution(numchamps);
593    
594     FEM_NOEUD* noeudcourant=NULL;
595     for (int i=0;i<som1->get_lien_fem_maillage()->get_nb();i++)
596     {
597     FEM_NOEUD* no=(FEM_NOEUD*)som1->get_lien_fem_maillage()->get(i);
598     if (mai->get_fem_noeudid(no->get_id())==no) noeudcourant=no;
599     }
600     FEM_NOEUD* noeudfin=NULL;
601     for (int i=0;i<som2->get_lien_fem_maillage()->get_nb();i++)
602     {
603     FEM_NOEUD* no=(FEM_NOEUD*)som2->get_lien_fem_maillage()->get(i);
604     if (mai->get_fem_noeudid(no->get_id())==no) noeudfin=no;
605     }
606     int n1=0;
607     int n2;
608     if (mai->get_degre()==1) n2=1;
609     if (mai->get_degre()==2) n2=2;
610     double *xyz=noeudcourant->get_coord();
611     double t;
612     inverser(t,xyz);
613     std::vector<double> resultat;
614     double tprec=t;
615     resultat.push_back(t);
616     resultat.push_back(0.);
617     resultat.push_back(noeudcourant->get_solution());
618     while (noeudcourant!=noeudfin)
619     {
620     FEM_ELEMENT1 *segcourant=NULL;
621     for (int i=0;i<get_lien_fem_maillage()->get_nb();i++)
622     {
623     FEM_ELEMENT1* seg=(FEM_ELEMENT1*)get_lien_fem_maillage()->get(i);
624     if (mai->get_fem_element1id(seg->get_id())==seg)
625     {
626     if ((seg->get_fem_noeud(n1)==noeudcourant)||(seg->get_fem_noeud(n2)==noeudcourant))
627     if (seg!=segcourant)
628     segcourant=seg;
629     }
630     }
631     FEM_NOEUD* autrenoeud=NULL;
632     if (segcourant->get_fem_noeud(n1)==noeudcourant) autrenoeud=segcourant->get_fem_noeud(n2);
633     if (segcourant->get_fem_noeud(n2)==noeudcourant) autrenoeud=segcourant->get_fem_noeud(n1);
634     FEM_NOEUD *noprec=noeudcourant;
635     for (int i=1;i<segcourant->get_nb_fem_noeud();i++)
636     {
637     FEM_NOEUD* no=segcourant->get_fem_noeud(i);
638     double *xyz=no->get_coord();
639     double t2;
640     inverser(t2,xyz);
641 francois 672 if (get_courbe()->est_periodique())
642     if (t2<tprec) t2=t2+get_courbe()->get_periode();
643 francois 604 double s=resultat[resultat.size()-2];
644     int pas=32;
645     for (int j=0;j<pas;j++)
646     {
647     double ti=tprec+1.0*j*(t2-tprec)/pas;
648     double tii=tprec+1.0*(j+1)*(t2-tprec)/pas;
649     double tgauss1=0.7886751345*tii+0.2113248654*ti;
650     double dxyz[3];
651     deriver(tgauss1,dxyz);
652     s=s+0.5*(tii-ti)*sqrt(dxyz[0]*dxyz[0]+dxyz[1]*dxyz[1]+dxyz[2]*dxyz[2]);
653     double tgauss2=0.7886751345*ti+0.2113248654*tii;
654     deriver(tgauss2,dxyz);
655     s=s+0.5*(tii-ti)*sqrt(dxyz[0]*dxyz[0]+dxyz[1]*dxyz[1]+dxyz[2]*dxyz[2]);
656     }
657     resultat.push_back(t2);
658     resultat.push_back(s);
659     resultat.push_back(no->get_solution());
660     noprec=no;
661     tprec=t2;
662     }
663     noeudcourant=autrenoeud;
664     }
665     FILE *in=fopen(fichier,"wt");
666     char mess[500];
667     sprintf(mess,"%s_%s",(char*)sol->get_nom().c_str(),(char*)sol->get_legende(numchamps).c_str());
668     fprintf(in,"t;s;%s;\n",mess);
669     for (int i=0;i<resultat.size();i=i+3)
670     fprintf(in,"%lf;%lf;%lf;\n",resultat[i],resultat[i+1],resultat[i+2]);
671     fclose(in);
672     }