/[ff3d]/ff3d/solver/SurfaceMeshGenerator.cpp
ViewVC logotype

Diff of /ff3d/solver/SurfaceMeshGenerator.cpp

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 1.8 by dobj, Thu Mar 20 15:47:15 2003 UTC revision 1.9 by dobj, Tue Mar 25 18:41:01 2003 UTC
# Line 89  private: Line 89  private:
89    //! list of cut edges and intersection vertex.    //! list of cut edges and intersection vertex.
90    std::map<Edge::Pair, size_t> edgesRef;    std::map<Edge::Pair, size_t> edgesRef;
91    
92  std::vector<ReferenceCounting<EdgePairList> > triangleListes;    std::vector<ReferenceCounting<EdgePairList> > triangleListes;
93    
94    EdgePairList localTriangleList;    EdgePairList localTriangleList;
95      EdgePairList localTriangleListother;
96    
97    //! This list is used to provide an orientation to the Surfacic Elements.    //! This list is used to provide an orientation to the Surfacic Elements.
98    std::list<TinyVector<3> > InsideVertexList;    std::list<TinyVector<3> > InsideVertexList;
99    
100    //! list of cut edges and intersection vertex.    //! list of cut edges and intersection vertex.
101   //VerticesList verticesList;   //VerticesList verticesList;
102  std::vector<ReferenceCounting<VerticesList> > verticesListes;    std::vector<ReferenceCounting<VerticesList> > verticesListes;
   
103    HexaTetraAssociation hexaTetraAssociation;    HexaTetraAssociation hexaTetraAssociation;
104    
105    MotherCellList cellList;    MotherCellList cellList;
# Line 114  std::vector<ReferenceCounting<VerticesLi Line 114  std::vector<ReferenceCounting<VerticesLi
114    void __setBoundaryType(const Boundary& B,    void __setBoundaryType(const Boundary& B,
115                           POVRayReferences& vrRefList,                           POVRayReferences& vrRefList,
116                           TinyVector<6,bool>& faces);                           TinyVector<6,bool>& faces);
117      void __InTriangle(const Vertex& P,
118                TinyVector<3>& P0,
119                TinyVector<3>& V12,
120                TinyVector<3>& V02,
121                size_t& ncompt);
122      void __inout(size_t iter,
123                 const size_t nn,
124                 int& ttn0,
125                   int& ttn1,
126                   int& ttn2,
127                 size_t& re0,
128                 size_t& re1,
129                 size_t& re2,
130                 size_t& ii0,
131                 size_t& ii1,
132                 size_t& ii2);
133      
134      void __CalculPoint(const Vertex& V0, const Vertex& V1 , const Vertex& V2,
135                         const Vertex& Vo0, const Vertex& Vo1 , const Vertex& Vo2,
136                         int nn,int& ttn0,int& ttn1 ,int& ttn2,
137                         Vertex * vvn0 , Vertex* vvn1 , Vertex* vvn2);
138      void __IntersectionPoint(const Scene& scene ,const Reference::Representation& ref,                                                       const Reference::Representation& otherRef,
139                               std::map<Edge::Pair,int>& triangleedge,
140                               std::map<Edge::Pair, Vertex*>&  listofVertex);
141      
142  public:  public:
143    Internals()    Internals()
144    {    {
# Line 201  __createTriangles(const Shape& S, Line 225  __createTriangles(const Shape& S,
225      const CartesianHexahedron& currentCell      const CartesianHexahedron& currentCell
226        = static_cast<const CartesianHexahedron&>(*(*i).first);        = static_cast<const CartesianHexahedron&>(*(*i).first);
227      std::vector<TetrahedronByEdges>& T = *(((*i).second).first);      std::vector<TetrahedronByEdges>& T = *(((*i).second).first);
228          
229      for (std::vector<TetrahedronByEdges>::iterator currentT = T.begin();      for (std::vector<TetrahedronByEdges>::iterator currentT = T.begin();
230           currentT < T.end(); ++currentT) {           currentT < T.end(); ++currentT) {
231        std::vector<const Edge*> cutEdges;        std::vector<const Edge*> cutEdges;
# Line 324  __createTriangles(const Shape& S, Line 348  __createTriangles(const Shape& S,
348  }  }
349    
350    
351    void SurfaceMeshGenerator::Internals::__InTriangle(const Vertex& P,TinyVector<3>& P0,TinyVector<3>& V12,TinyVector<3>& V02,size_t& compt)
352    {
353    
354        compt=0;
355        real_t det,eps=-1e-10;
356        det=V02[0]*V12[1]-V02[1]*V12[0];
357        if(det<eps or det>-eps)
358        {
359            //std::cout<<"det "<<det<<"\n";
360            TinyMatrix<2,2> mat;
361            mat(0,0)=V02[0];
362            mat(1,0)=V02[1];
363            mat(0,1)=V12[0];
364            mat(1,1)=V12[1];
365            TinyVector<2> bary,secondmb;
366            secondmb[0]=P0[0];
367            secondmb[1]=P0[1];
368            bary=mat.invert()*secondmb;
369            real_t temp=(1.-bary[0]-bary[1]);
370            if(bary[0]>=eps && bary[1]>=eps && temp>=eps)
371            {
372                compt+=1;
373                //  std::cout<<bary[0]<<" "<<bary[1]<<" "<<1.-bary[0]-bary[1]<<" "<<compt<<"\n";
374         //std::cout<<verif[0]<<" "<<verif[1]<<" \n";
375         //std::cout<<P0[0]<<" "<<P0[2]<<" \n";
376            }
377        }
378        else
379            compt+=1;
380        {
381            det=V02[0]*V12[2]-V02[2]*V12[0];
382            if(det<eps or det>-eps)
383            {
384                TinyMatrix<2,2> mat;
385                mat(0,0)=V02[0];
386                mat(1,0)=V02[2];
387                mat(1,1)=V12[2];
388                mat(0,1)=V12[0];
389                TinyVector<2> bary,secondmb;
390                secondmb[0]=P0[0];
391                secondmb[1]=P0[2];
392                bary=mat.invert()*secondmb;
393                TinyVector<2> verif=mat*bary;
394                if(bary[0]>=eps && bary[1]>=eps && (1-bary[0]-bary[1])>=eps)
395                {
396                    compt+=1;
397                    //std::cout<<bary[0]<<" "<<bary[1]<<" "<<1.-bary[0]-bary[1]<<"\n";
398      //std::cout<<verif[0]<<" "<<verif[1]<<" \n";
399      //std::cout<<P0[0]<<" "<<P0[2]<<" \n";
400    
401    
402                }
403            }
404            else
405                compt+=1;
406        }
407        {
408            det=V02[1]*V12[2]-V02[2]*V12[1];
409            if(det<eps or det>-eps)
410            {
411                TinyMatrix<2,2> mat;
412                mat(0,0)=V02[1];
413                mat(1,0)=V02[2];
414                mat(1,1)=V12[2];
415                mat(0,1)=V12[1];
416                TinyVector<2> bary,secondmb;
417                secondmb[0]=P0[1];
418                secondmb[1]=P0[2];
419                bary=mat.invert()*secondmb;
420                TinyVector<2> verif=mat*bary;
421                if(bary[0]>=eps && bary[1]>=eps && (1-bary[0]-bary[1])>=eps)
422                {
423                    compt+=1;
424                    //std::cout<<bary[0]<<" "<<bary[1]<<" "<<1.-bary[0]-bary[1]<<" "<<compt<<"\n";
425      //std::cout<<verif[0]<<" "<<verif[1]<<" \n";
426      //std::cout<<P0[0]<<" "<<P0[1]<<" \n";
427    
428                }
429            }
430            else
431                compt+=1;
432        }
433    }
434    
435    void SurfaceMeshGenerator::Internals::__CalculPoint(const Vertex& V0, const Vertex& V1 , const Vertex& V2,
436                       const Vertex& Vo0, const Vertex& Vo1 , const Vertex& Vo2,
437                                                        int nn,int& ttn0,int& ttn1 ,int& ttn2,
438                                                        Vertex * vvn0 , Vertex * vvn1 , Vertex * vvn2)
439    {
440        TinyVector<3>  V01,V12,V02;
441        V01 = V1.Vector3() - V0.Vector3();
442        V12 = V2.Vector3() - V1.Vector3();
443        V02 = V2.Vector3() - V0.Vector3();
444        
445        //equation of plan
446        TinyVector<3>  Vo01,Vo12,Vo02;
447        Vo01=Vo1.Vector3()-Vo0.Vector3();
448        Vo12=Vo2.Vector3()-Vo1.Vector3();
449        Vo02=Vo2.Vector3()-Vo0.Vector3();
450    
451        TinyVector<3> coeffo;
452        coeffo=Vo01^Vo02;
453        real_t d2 = -(coeffo*Vo0.Vector3());
454        real_t pscalar,k;
455        size_t num;
456        {
457            Vertex P;
458            real_t norm1=Norm(V01);
459            pscalar=coeffo*V01;
460            num=0;
461            if(std::abs(pscalar)>1e-13) {
462                k=-(coeffo[0]*V0[0]+coeffo[1]*V0[1]+coeffo[2]*V0[2]+d2)/pscalar;
463                P=V0+k*V01;
464                TinyVector<3> P0=P.Vector3()-Vo2.Vector3(),P1=Vo12*(-1.),P2=Vo02*(-1.);
465                real_t norm2=Norm(P.Vector3()-V0.Vector3());
466                real_t norm3=Norm(P.Vector3()-V1.Vector3());
467                if(norm2<=norm1 and norm3<=norm1)   {
468                    __InTriangle(P,P0,P1,P2,num);
469                    if(num==3)
470                    {
471                        //std::cout<<"droite1 "<<P.Vector3()[0]<<" "<<P.Vector3()[1]<<" "<<P.Vector3()[2]<<"  \n";
472                        ttn0=1;
473                        vvn0=new Vertex(P.Vector3());
474                        //std::cout<<"P "<<nn<<" "<< vv[0][0].Vector3()[0]<<" "<<(vv[0][0]).Vector3()[1]<<" "<<(vv[0][0]).Vector3()[2]<<"\n";
475          //std::cout<<"verif "<<coeffo[0]*P.Vector3()[0]+coeffo[1]*P.Vector3()[1]+coeffo[2]*P.Vector3()[2]+d2<<"\n";
476          //std::cout<<"verif "<<coeff[0]*P.Vector3()[0]+coeff[1]*P.Vector3()[1]+coeff[2]*P.Vector3()[2]+d<<"\n";
477                    }
478                }
479            }
480        }
481        {
482            Vertex P;
483            pscalar=coeffo*V02;
484            real_t norm1=Norm(V02);
485            if(std::abs(pscalar)>1e-13) {
486                k=-(coeffo[0]*V0[0]+coeffo[1]*V0[1]+coeffo[2]*V0[2]+d2)/pscalar;
487                P=V0+k*V02;
488                TinyVector<3> P0=P.Vector3()-Vo2.Vector3(),P1=Vo12*(-1.),P2=Vo02*(-1.);
489                real_t norm2=Norm(P.Vector3()-V0.Vector3());
490                real_t norm3=Norm(P.Vector3()-V2.Vector3());
491                //std::cout<<"norm "<<norm3<<" "<<norm2<<" "<<norm1<<"\n";
492                if(norm2<=norm1 and norm3<=norm1) {
493                    __InTriangle(P,P0,P1,P2,num);
494                    if(num==3) {
495                       // std::cout<<"droite2 "<<P.Vector3()[0]<<" "<<P.Vector3()[1]<<" "<<P.Vector3()[2]<<"  \n";
496                        ttn1=1;
497                        vvn1=new Vertex(P.Vector3());
498                    }
499                }
500            }
501        }
502        {
503            Vertex P;
504            pscalar=coeffo*V12;
505            real_t norm1=Norm(V12);
506            if(std::abs(pscalar)>1e-13) {
507                k=-(coeffo*V1.Vector3() + d2)/pscalar;
508                P=V1+k*V12;
509                TinyVector<3> P0=P.Vector3()-Vo2.Vector3(),P1=Vo12*(-1.),P2=Vo02*(-1.);
510                real_t norm2=Norm(P.Vector3()-V1.Vector3());
511                real_t norm3=Norm(P.Vector3()-V2.Vector3());
512    
513                if(norm2<=norm1 and norm3<=norm1) {
514                    __InTriangle(P,P0,P1,P2,num);
515                    if(num==3) {
516    //std::cout<<"droite3 "<<P.Vector3()[0]<<" "<<P.Vector3()[1]<<" "<<P.Vector3()[2]<<"  \n";
517    
518                        ttn2=1;
519                        vvn2=new Vertex(P.Vector3());
520                    }
521                }
522            }
523        }
524    }
525    
526    void SurfaceMeshGenerator::Internals::
527    __inout(size_t iter,const size_t nn,int& ttn0,
528                   int& ttn1,
529                   int& ttn2,
530            size_t& re0,size_t& re1,size_t& re2,size_t& ii0,size_t& ii1,size_t& ii2)
531    {
532    
533        if(re0==1 and re1==1)
534        {
535            if(ttn0!=0 and iter>1)
536    std::cout<<"-- cas indetermine\n";
537        }
538        else if((re0+re1==1))
539        {
540            if(ttn0==0)
541            {
542                //std::cout<<"-- 1 0\n";
543                re0=0;
544                re1=0;
545            }
546            else
547            {
548                //std::cout<<"-- 1 in 0 out on sait traite\n";
549                ii0=re0;
550                ii1=re1;
551            }
552        }
553        else if((re0==2 and re1==0) or
554                (re0==0 and re1==2))
555        {
556            if(ttn0==0)
557            {
558                if(ttn0+ttn1+ttn2==0)
559                {
560                    re0=2;
561                    re1=2;
562                }
563                else
564                {
565    std::cout<<"-- 0 2 pas coupe\n";
566    std::cout<<ttn0<<" "<<ttn1<<" "<<ttn2<<"\n";
567    std::cout<<re0<<" "<<re1<<"\n";
568                    re0=2;
569                    re1=2;
570                }
571              /*      if(re0==2 and (ttn1==1 or ttn2==1))
572                      {
573                      re0=-nn;
574                      ii0=-nn;
575                      }
576                      else
577                      {
578                      re1=-nn;
579                      ii1=-nn;
580                      }*/
581              //  re0=0;
582              //re1=0;
583            }
584            else
585            {
586              //std::cout<<"-- 0 2 on sait traite\n";
587                ii0=re0;
588                ii1=re1;
589            }
590        }
591        else if((re0+re1==3))
592        {
593            if(ttn0==0)
594            {
595              //std::cout<<"-- 1 2 pas coupe\n";
596                re0=2;
597                re1=2;
598            }
599            else
600            {
601              //std::cout<<"-- 1 2 on sait traite\n";
602                if(re0==1)
603                {
604                    ii0=0;
605                    ii1=1;
606                    re0=0;
607                }
608                else
609                {
610                    ii0=1;
611                    ii1=0;
612                    re1=0;
613                }
614            }
615        }
616      //     else
617      //std::cout<<"--0 out 0 out rien a faire\n";
618    
619        if(re0==1 and re2==1)
620        {
621            if(ttn1!=0 and iter>1)
622    std::cout<<"--* cas indetermine\n";
623        }
624        else if((re0==1 and re2==0) or
625                (re0==0 and re2==1))
626        {
627            if(ttn1==0)
628            {
629                //std::cout<<"--* 1 O non inter\n";
630                re0=0;
631                re2=0;
632            }
633            else
634            {
635              //std::cout<<"--* 1 in 0 out on sait traite\n";
636                ii0=re0;
637                ii2=re2;
638            }
639        }
640        else if((re0==2 and re2==0) or
641                (re0==0 and re2==2))
642        {
643            if(ttn1==0)
644            {
645                if(ttn0+ttn1+ttn2==0)
646                {
647                    re0=2;
648                    re2=2;
649                }
650                else
651                {
652    
653    std::cout<<"--* 0 2 pas coupe\n";
654    std::cout<<ttn0<<" "<<ttn1<<" "<<ttn2<<"\n";
655    std::cout<<re0<<" "<<re2<<"\n";
656                    re0=2;
657                    re2=2;
658                }
659              /*   if(re0==2 and (ttn0==1 or ttn2==1))
660                   {
661                   re0=-nn;
662                   ii0=-nn;
663                   }
664                   else
665                   {
666                   re2=-nn;
667                   ii2=-nn;
668                   }*/
669    
670              //      re0=0;
671              //    re2=0;
672            }
673            else
674            {
675              //std::cout<<"--* 0 2 on sait traite\n";
676                ii0=re0;
677                ii2=re2;
678            }
679        }
680        else if((re0+re2==3))
681        {
682            if(ttn1==0)
683            {
684              //std::cout<<"--* 1 2 pas coupe\n";
685                re0=2;
686                re2=2;
687            }
688            else
689            {
690              //std::cout<<"--* 1 2 on sait traite\n";
691                if(re0==1)
692                {
693                    ii0=0;
694                    ii2=1;
695                    re0=0;
696                }
697                else
698                {
699                    ii0=1;
700                    ii2=0;
701                    re2=0;
702                }
703            }      
704        }      
705      //       else
706      //std::cout<<"--* 0 out 0 out rien a faire\n";
707    
708        if(re2==1 and re1==1)
709        {
710            if(ttn2!=0 and iter>1)
711    std::cout<<"--+ cas indetermine\n";
712        }
713        else if((re2+re1==1))
714        {
715            if(ttn2==0)
716            {
717              //std::cout<<"--+ 1 0 non inter\n";
718                re2=0;
719                re1=0;
720            }
721            else
722            {
723              //std::cout<<"--+ 1 in 0 out on sait traite\n";
724                ii1=re1;
725                ii2=re2;
726            }
727        }
728        else if((re2==2 and re1==0) or
729                (re2==0 and re1==2))
730        {
731            if(ttn2==0)
732            {
733                if(ttn0+ttn1+ttn2==0)
734                {
735                    re1=2;
736                    re2=2;
737                }
738                else
739                {
740    
741    std::cout<<"--+ 0 2 pas coupe\n";
742    std::cout<<ttn0<<" "<<ttn1<<" "<<ttn2<<"\n";
743    std::cout<<re1<<" "<<re2<<"\n";
744                    re1=2;
745                    re2=2;
746                }
747              /*   if(re1==2 and (ttn1==1 or ttn0==1))
748                   {
749                   re1=-nn;
750                   ii1=-nn;
751                   }
752                   else
753                   {
754                   re2=-nn;
755                   ii2=-nn;
756                   }*/
757              //re1=0;
758              //re2=0;
759            }
760            else
761            {
762              //std::cout<<"--+ 0 2 on sait traite\n";
763                ii1=re1;
764                ii2=re2;
765            }
766        }
767        else if((re2+re1==3))
768        {
769            if(ttn2==0)
770            {
771                //std::cout<<"--+ 1 2 pas coupe\n";
772                re2=2;
773                re1=2;
774            }
775            else
776            {
777              //std::cout<<"--+ 1 2 on sait traite\n";
778                if(re1==1)
779                {
780                    ii1=0;
781                    ii2=1;
782                    re1=0;
783                }
784                else
785                {
786                    ii1=1;
787                    ii2=0;
788                    re2=0;
789                }
790            }
791        }
792      //       else
793      //std::cout<<"--+ 0 out 0 out rien a faire\n";
794    
795        if(/*iter>1 and*/ (((re2==1 and re1==1 and re0==1) or re0+re1+re2==6)
796                           and ttn0+ttn1+ttn2==0))
797        {
798            re0=2;
799            re1=2;
800            re2=2;
801            ii0=1;
802            ii1=1;
803            ii2=1;
804        }
805        if(iter>2 and re0+re1+re2==0 and ttn0+ttn1+ttn2!=0)
806        {
807    std::cout<<"pbs de precision?? tr out cut\n";
808            ii0=0;
809            ii1=0;
810            ii2=0;
811        }
812      /*if(iter==1)
813        {
814          if(ttn0+ttn1+ttn2==0 and ((re0==0 and re1==1 and re2==1) or (re0==1 and re1==0 and re2==1)
815                                          or (re0==1 and re1==1 and re2==0) or (re0==0 and re1==0 and re2==1)
816                                          or (re0==0 and re1==1 and re2==0) or (re0==1 and re1==0 and re2==0)))
817            {
818              re0=0;
819              re1=0;
820              re2=0;
821            }
822        }*/
823    
824      if(ii0==2)
825        ii0=1;
826      if(ii1==2)
827        ii1=1;
828      if(ii2==2)
829        ii2=1;
830    
831    }
832    
833    void SurfaceMeshGenerator::Internals::
834    __IntersectionPoint(const Scene& scene ,const Reference::Representation& ref
835                        ,const Reference::Representation& otherRef,
836                        std::map<Edge::Pair,int>& triangleedge,
837                        std::map<Edge::Pair, Vertex*>&  listofVertex)
838    {
839    
840        int nn=0;
841    EdgePairList::iterator currentTriangle=(*triangleListes[0]).end();
842    EdgePairList::iterator currentTriangleother=(*triangleListes[1]).end();
843        --currentTriangle;
844        --currentTriangleother;
845        
846        for(SurfaceMeshGenerator::Internals::HexaTetraAssociation::iterator i
847            = hexaTetraAssociation.begin();
848            i != hexaTetraAssociation.end(); ++i) {
849    
850        
851    std::vector<TetrahedronByEdges>& T = *(((*i).second).first);
852    //      ffout(0) << " hexa*****************************************\n";
853    
854            for (std::vector<TetrahedronByEdges>::iterator currentT = T.begin();
855          currentT < T.end(); ++currentT) {
856    
857    std::vector<const Edge*> cutEdges;
858    std::vector<const Edge*> cutEdgesother;
859    
860                for(size_t ie=0; ie<6; ++ie) {
861                    const Edge& e = currentT->edge(ie);
862                    if ((e(0).Ref().isSet(ref)) xor (e(1).Ref().isSet(ref))) {
863                        cutEdges.push_back(&e);
864                    }
865                    if ((e(0).Ref().isSet(otherRef)) xor (e(1).Ref().isSet(otherRef))) {
866                        cutEdgesother.push_back(&e);
867                    }
868                    
869                }
870    //          ffout(0) <<cutEdges.size()<<" "<<cutEdgesother.size()<< " tetra__________*****************************************\n";
871    
872                switch (cutEdges.size()) {
873                    case 3: {
874                        localTriangleList.push_front(TinyVector<3,Edge::Pair>());
875    EdgePairList::iterator currentT = localTriangleList.begin();
876                        for(size_t i=0; i<3; ++i)
877                                    (*currentT)[i] = (*currentTriangle)[i];
878                        --currentTriangle;
879                        break;
880                       }
881                        
882                    case 4: {
883                    {                   localTriangleList.push_front(TinyVector<3,Edge::Pair>());
884    EdgePairList::iterator currentT = localTriangleList.begin();
885                        for(size_t i=0; i<3; ++i)
886                            (*currentT)[i] = (*currentTriangle)[i];
887                    }
888                        --currentTriangle;
889                        
890                        localTriangleList.push_front(TinyVector<3,Edge::Pair>());
891    EdgePairList::iterator currentT = localTriangleList.begin();
892    
893                        for(size_t i=0; i<3; ++i)
894                            (*currentT)[i] = (*currentTriangle)[i];
895                        --currentTriangle;
896                        break;
897                    }
898                    case 1:
899                    case 2: {
900                        ffout(0) << "?!\n";
901                        break;
902                    }
903                }
904    
905                switch (cutEdgesother.size()) {
906                    case 3: {
907                     localTriangleListother.push_front(TinyVector<3,Edge::Pair>());
908    EdgePairList::iterator currentT = localTriangleListother.begin();
909                        for(size_t i=0; i<3; ++i)
910                            (*currentT)[i] = (*currentTriangleother)[i];
911                        --currentTriangleother;
912                        break;
913                    }
914    
915                    case 4: {
916                    {   localTriangleListother.push_front(TinyVector<3,Edge::Pair>());
917    EdgePairList::iterator currentT = localTriangleListother.begin();
918                        for(size_t i=0; i<3; ++i)
919                            (*currentT)[i] = (*currentTriangleother)[i];
920                    }
921                        --currentTriangleother;
922    
923                        localTriangleListother.push_front(TinyVector<3,Edge::Pair>());
924    EdgePairList::iterator currentT = localTriangleListother.begin();
925    
926                        for(size_t i=0; i<3; ++i)
927                            (*currentT)[i] = (*currentTriangleother)[i];
928                        --currentTriangleother;
929                        break;
930                    }
931                    case 1:
932                    case 2: {
933                        ffout(0) << "?!\n";
934                        break;
935                    }
936                }
937            }
938    
939    
940            for(SurfaceMeshGenerator::Internals::EdgePairList::iterator currentTriangle1
941         = localTriangleList.begin();
942         currentTriangle1 != localTriangleList.end();
943         ++currentTriangle1)
944            {
945            const Edge::Pair& e0=(*currentTriangle1)[0];
946            const Edge::Pair& e1=(*currentTriangle1)[1];
947            const Edge::Pair& e2=(*currentTriangle1)[2];
948            //ffout(0) << (*e0.first).Vector3()[0]<<"\n";
949            const Vertex& V0 = *((*verticesListes[0])[e0]);
950            const Vertex& V1 = *((*verticesListes[0])[e1]);
951            const Vertex& V2 = *((*verticesListes[0])[e2]);
952            //tt[nn] = 0;
953            (triangleedge)[e0]=0;
954            (triangleedge)[e1]=0;
955            (triangleedge)[e2]=0;
956            for(SurfaceMeshGenerator::Internals::EdgePairList::iterator currentTriangle
957         = localTriangleListother.begin();
958         currentTriangle != localTriangleListother.end();
959         ++currentTriangle) {
960                    //std::cout<<"*************************\n";
961                const Edge::Pair& eo0=(*currentTriangle)[0];
962                const Edge::Pair& eo1=(*currentTriangle)[1];
963                const Edge::Pair& eo2=(*currentTriangle)[2];
964                    const Vertex& Vo0 = *((*verticesListes[1])[eo0]);
965                    const Vertex& Vo1 = *((*verticesListes[1])[eo1]);
966                    const Vertex& Vo2 = *((*verticesListes[1])[eo2]);
967    
968                    __CalculPoint(V0,V1,V2,Vo0,Vo1,Vo2,nn,
969                    (triangleedge)[e0],(triangleedge)[e1],(triangleedge)[e2]
970                    ,((listofVertex)[e0]),((listofVertex)[e1]),((listofVertex)[e2]));
971                }
972            
973            ++nn;
974            //ffout(0) << "nn "<<nn<<"\n";
975            }
976            localTriangleList.clear();
977            localTriangleListother.clear();
978        }
979    /*    nn=0;
980    
981        for(SurfaceMeshGenerator::Internals::EdgePairList::iterator currentTriangle1
982            = (*triangleListes[0]).begin();
983            currentTriangle1 != (*triangleListes[0]).end();
984            ++currentTriangle1)
985        {
986            const Edge::Pair& e0=(*currentTriangle1)[0];
987            const Edge::Pair& e1=(*currentTriangle1)[1];
988            const Edge::Pair& e2=(*currentTriangle1)[2];
989            const Vertex& V0 = *((*verticesListes[0])[e0]);
990            const Vertex& V1 = *((*verticesListes[0])[e1]);
991            const Vertex& V2 = *((*verticesListes[0])[e2]);
992    
993            //TinyVector<3>  V01,V12,V02;
994            //V01 = V1.Vector3() - V0.Vector3();
995            //V12 = V2.Vector3() - V1.Vector3();
996            //V02 = V2.Vector3() - V0.Vector3();
997            //  const TinyVector<3> coeff=V01^V02;
998     //  const real_t d=-(coeff[0]*V0[0]+coeff[1]*V0[1]+coeff[2]*V0[2]);
999     //      std::cout<<"--------------------------------------------------------------\n";
1000    
1001            tt[nn] = 0;
1002    
1003            // Pour Cecile :
1004     // Je crois que tu ne pouvais pas espérer que ça marche pour plus de 2 objets
1005     // puisque tu écris explicitement d'utiliser la liste [1]
1006    
1007            //       for(SurfaceMeshGenerator::Internals::EdgePairList::iterator currentTriangle
1008     //         = (*(*__internals).triangleListes[1]).begin();
1009     //       currentTriangle != (*(*__internals).triangleListes[1]).end();
1010     //       ++currentTriangle)
1011    
1012            for (size_t currentObjectNumber = 1; currentObjectNumber < scene.nbObjects(); ++currentObjectNumber) {
1013    
1014    SurfaceMeshGenerator::Internals::EdgePairList& otherTrianglesList
1015                = *triangleListes[currentObjectNumber];
1016    
1017    SurfaceMeshGenerator::Internals::VerticesList&
1018                otherVerticesList = (*verticesListes[currentObjectNumber]);
1019    
1020                for(SurfaceMeshGenerator::Internals::EdgePairList::iterator currentTriangle
1021             = otherTrianglesList.begin();
1022             currentTriangle != otherTrianglesList.end();
1023             ++currentTriangle)
1024                {
1025                    //std::cout<<"*************************\n";
1026                    const Vertex& Vo0 = *(otherVerticesList[(*currentTriangle)[0]]);
1027                    const Vertex& Vo1 = *(otherVerticesList[(*currentTriangle)[1]]);
1028                    const Vertex& Vo2 = *(otherVerticesList[(*currentTriangle)[2]]);
1029    
1030                    __CalculPoint(V0,V1,V2,Vo0,Vo1,Vo2,nn,tt[nn],vv[nn]);
1031                }
1032            }
1033            nn++;
1034        }*/
1035    }
1036    
1037    
1038  void  void
1039  SurfaceMeshGenerator::  SurfaceMeshGenerator::
1040  generatesAllSufaceMeshes(const Problem& problem,  generatesAllSufaceMeshes(const Problem& problem,
# Line 543  void plot2(const SurfaceMeshOfTriangles& Line 1254  void plot2(const SurfaceMeshOfTriangles&
1254            Structured3DMesh& SMesh);            Structured3DMesh& SMesh);
1255    
1256  void plotmedit(const SurfaceMeshOfTriangles& s_mesh,  void plotmedit(const SurfaceMeshOfTriangles& s_mesh,
1257            std::set<const Cell*>& toTreatHexahedra,size_t nobj);                 std::set<const Cell*>& toTreatHexahedra,
1258  void Bary(const Vertex& P,TinyVector<3>& P0,TinyVector<3>& V12,TinyVector<3>& V02,size_t& ncompt);                 size_t nobj);
1259  void inout(size_t iter,const size_t nn,TinyVector<3,int>  ttn,size_t& re0,size_t& re1,size_t& re2,size_t& ii0,size_t& ii1,size_t& ii2);  
1260    
1261    /*! Generates a surfacic mesh using only characteristic function of an object    /*! Generates a surfacic mesh using only characteristic function of an object
1262      and the Structured 3D Mesh.      and the Structured 3D Mesh.
# Line 571  generateSurfacicMesh(const Object& O, Line 1282  generateSurfacicMesh(const Object& O,
1282    const Object& object = Omega.object();    const Object& object = Omega.object();
1283    const Scene& scene = Omega.scene();    const Scene& scene = Omega.scene();
1284    
   
1285    switch (object.shapeType()) {    switch (object.shapeType()) {
1286    case Shape::union_: {    case Shape::union_: {
1287      const Union& U = static_cast<const Union&>(*(object.shape()));      const Union& U = static_cast<const Union&>(*(object.shape()));
# Line 584  generateSurfacicMesh(const Object& O, Line 1294  generateSurfacicMesh(const Object& O,
1294      const Intersection& I = static_cast<const Intersection&>(*(object.shape()));      const Intersection& I = static_cast<const Intersection&>(*(object.shape()));
1295      for(Intersection::const_iterator i = I.begin();      for(Intersection::const_iterator i = I.begin();
1296          i != I.end(); ++i) {          i != I.end(); ++i) {
1297        otherReferences[(*(*i)).reference()] = *i;          otherReferences[(*(*i)).reference()] = *i;
1298            //ffout(0) << (*(*i)).reference() << "\n";
1299      }      }
1300      break;      break;
1301    }    }
# Line 594  generateSurfacicMesh(const Object& O, Line 1305  generateSurfacicMesh(const Object& O,
1305    }    }
1306    
1307    otherReferences.erase(O.reference());    otherReferences.erase(O.reference());
1308    
1309    //  for(ObjectReferences::iterator i=otherReferences.begin();i!=otherReferences.end();++i)
1310      //    ffout(0) << (*i).first <<"\n";
1311      
1312    size_t nobject = scene.nbObjects();    size_t nobject = scene.nbObjects();
1313    
1314     (*__internals).verticesListes.resize(nobject);    (*__internals).verticesListes.resize(nobject);
1315    (*__internals).triangleListes.resize(nobject);    (*__internals).triangleListes.resize(nobject);
1316    
1317      std::cout<<"nb object "<<nobject<<"\n";
1318    
1319    for(size_t no=0 ; no<nobject ; ++no)    for(size_t no=0 ; no<nobject ; ++no)
1320    {      {
1321        (*__internals).triangleListes[no]= new SurfaceMeshGenerator::Internals::EdgePairList;        (*__internals).triangleListes[no]= new SurfaceMeshGenerator::Internals::EdgePairList;
1322        (*__internals).verticesListes[no]= new SurfaceMeshGenerator::Internals::VerticesList;        (*__internals).verticesListes[no]= new SurfaceMeshGenerator::Internals::VerticesList;
1323        }
   
 }  
1324    
1325    for (size_t i = 0; i<SMesh.nbEdges(); i++) {    for (size_t i = 0; i<SMesh.nbEdges(); i++) {
1326      Edge& e = SMesh.edge(i);      Edge& e = SMesh.edge(i);
1327  #warning THIS TEST IS NOT ENOUGH! THINK TO THE BORDER!  #warning THIS TEST IS NOT ENOUGH! THINK TO THE BORDER!
1328      if((e(0).Ref().insideDomain() xor e(1).Ref().insideDomain())      if((e(0).Ref().insideDomain() xor e(1).Ref().insideDomain())
1329         and ((e(0).Ref().isSet(ref)) xor (e(1).Ref().isSet(ref))))         and((e(0).Ref().isSet(ref)) xor (e(1).Ref().isSet(ref))))
1330        {        {
1331          (*(*__internals).verticesListes[0])[e] = new Vertex((*__internals).__splitEdge(e, ref, S));          (*(*__internals).verticesListes[0])[e] = new Vertex((*__internals).__splitEdge(e, ref, S));
1332          //(*__internals).verticesList[e] = new Vertex((*__internals).__splitEdge(e, ref, S));          //(*__internals).verticesList[e] = new Vertex((*__internals).__splitEdge(e, ref, S));
# Line 626  generateSurfacicMesh(const Object& O, Line 1340  generateSurfacicMesh(const Object& O,
1340    /// a triangulation has been generated inside the hexahedron    /// a triangulation has been generated inside the hexahedron
1341    std::set<const Cell*> meshedHexahedra;    std::set<const Cell*> meshedHexahedra;
1342    
1343  //  if ((*__internals).verticesList.size()==0) {    if ((*(*__internals).verticesListes[0]).size()==0) {
 if ((*(*__internals).verticesListes[0]).size()==0) {  
1344    
1345  fferr(2) << __FILE__ << ':' << __LINE__ << ": "      fferr(2) << __FILE__ << ':' << __LINE__ << ": "
1346               << "meshing object " << S << " a problem occured ...\n"               << "meshing object " << S << " a problem occured ...\n"
1347               << "Check that:\n"               << "Check that:\n"
1348               << "- the object is not to small compare to the background mesh\n"               << "- the object is not to small compare to the background mesh\n"
# Line 671  fferr(2) << __FILE__ << ':' << __LINE__ Line 1384  fferr(2) << __FILE__ << ':' << __LINE__
1384    
1385          for (size_t l=0; l<H.numberOfEdges(); l++) {          for (size_t l=0; l<H.numberOfEdges(); l++) {
1386            Edge& e = H.edge(l);            Edge& e = H.edge(l);
1387  //        if ((*__internals).verticesList.count(e) > 0) {            if ((*(*__internals).verticesListes[0]).count(e) > 0) {
             if ((*(*__internals).verticesListes[0]).count(e) > 0) {  
1388              intersected = true;              intersected = true;
1389              cuttedHexahedra.insert(&currentCell);              cuttedHexahedra.insert(&currentCell);
1390              break;              break;
# Line 930  fferr(2) << __FILE__ << ':' << __LINE__ Line 1642  fferr(2) << __FILE__ << ':' << __LINE__
1642                  (*(*__internals).triangleListes[0]).push_front(*t);                  (*(*__internals).triangleListes[0]).push_front(*t);
1643    
1644              }              }
1645              (*__internals).hexaTetraAssociation[&currentCell]  //          (*__internals).hexaTetraAssociation[&currentCell]
1646                  = std::pair<ReferenceCounting<std::vector<TetrahedronByEdges> >,  //              = std::pair<ReferenceCounting<std::vector<TetrahedronByEdges> >,
1647                  ReferenceCounting<std::vector<Edge> > >  //              ReferenceCounting<std::vector<Edge> > >
1648                  (pT,pAddedEdges);  //              (pT,pAddedEdges);
1649                            
1650            } else {            } else {
1651              (*__internals).hexaTetraAssociation[&currentCell]              (*__internals).hexaTetraAssociation[&currentCell]
# Line 947  fferr(2) << __FILE__ << ':' << __LINE__ Line 1659  fferr(2) << __FILE__ << ':' << __LINE__
1659          }          }
1660        }        }
1661    
1662  //(*(*__internals).triangleListes[0]).clear();  /*  (*(*__internals).triangleListes[0]).clear();
1663  //(*(*__internals).verticesListes[0]).clear();    (*(*__internals).verticesListes[0]).clear();
1664    
1665  ////  (*__internals).verticesList.clear();    (*__internals).__createTriangles(S,ref,0);
1666      int num=0;
1667    
1668      size_t nobj= (*(*__internals).triangleListes[0]).size();
1669      std::cout<<"nb tr objet1 "<<nobj<<" \n";
1670    
1671  //(*__internals).__createTriangles(S,ref,0);    for(ObjectReferences::iterator i = otherReferences.begin();
1672          i != otherReferences.end(); ++i) {
1673        num++;
1674        std::cout<<"num "<<num<<"\n";
1675        for(size_t n = 0; n<scene.nbObjects(); ++n) {
1676          Reference::Representation otherRef;
1677          otherRef.resize(scene.nbObjects());
1678        //  ffout(0) << scene.object(n).reference() << ":"<< ((*i).first) << '\n';
1679          if (scene.object(n).reference() == ((*i).first)) {
1680            otherRef[n] = true;
1681            //std::cout<<"here\n";
1682            (*__internals).__createTriangles(*(scene.object(n).shape()),
1683                                             otherRef,num);
1684          }
1685        }
1686      }*/
1687    
1688      size_t nobjother=(*(*__internals).triangleListes[1]).size();
1689      std::cout<<"nb tr objet2 "<<nobjother<<" \n";
1690      if(scene.nbObjects()>2)
1691        std::cout<<"nb tr objet3 "<<(*(*__internals).triangleListes[2]).size()<<" \n";
1692    /*
1693      // *******************************************************************************
1694    
1695    SurfaceMeshGenerator::Internals::VerticesList listofVertex;
1696    //  std::map<Edge::Pair, Vertex*>  listofVertex;
1697      std::map<Edge::Pair,int> triangleedge;
1698    
 //size_t nobj= (*(*__internals).triangleListes[0]).size();  
 //std::cout<<"nb tr objet1 "<<nobj<<" \n";  
   
 //int num=0;  
  //  for(ObjectReferences::iterator i = otherReferences.begin();  
   //     i != otherReferences.end(); ++i) {  
    //    num++;  
     // for(size_t n = 0; n<scene.nbObjects(); ++n) {  
      //  Reference::Representation otherRef;  
       // otherRef.resize(scene.nbObjects());  
   
  //      if (scene.object(n).reference() == ((*i).first)) {  
  //     otherRef[n] = true;  
  //     (*__internals).__createTriangles(*(scene.object(n).shape()),  
  //                                      otherRef,num);  
    //    }  
     // }  
   // }  
 //size_t nobjother=(*(*__internals).triangleListes[1]).size();  
 //std::cout<<"nb tr objet2 "<<nobjother<<" \n";  
   
 /*  // *******************************************************************************  
 //std::map<Edge::Pair, Vertex*>  listofVertex;  
 //std::map<Edge::Pair,int> triangleedge;  
1699  std::vector< TinyVector<3,int> > tt;  std::vector< TinyVector<3,int> > tt;
1700  std::vector< std::vector<Vertex> > vv;    std::vector< std::vector<Vertex> > vv;
1701     vv.resize((*(*__internals).triangleListes[0]).size());    vv.resize((*(*__internals).triangleListes[0]).size());
1702     tt.resize((*(*__internals).triangleListes[0]).size());    tt.resize((*(*__internals).triangleListes[0]).size());
1703     for(size_t kk=0 ;kk<nobj;kk++)    for(size_t kk=0 ;kk<nobj;kk++)
1704         vv[kk].resize(3);      vv[kk].resize(3);
1705     int nn=0;  
1706        num=0;
1707    for(SurfaceMeshGenerator::Internals::EdgePairList::iterator currentTriangle1    for(ObjectReferences::iterator i = otherReferences.begin();
1708        = (*(*__internals).triangleListes[0]).begin();        i != otherReferences.end(); ++i) {
1709        currentTriangle1 != (*(*__internals).triangleListes[0]).end();        num++;
1710        ++currentTriangle1)        for(size_t n = 0; n<scene.nbObjects(); ++n) {
1711    {  Reference::Representation otherRef;
1712        const Edge::Pair e0=(*currentTriangle1)[0];            otherRef.resize(scene.nbObjects());
1713        const Edge::Pair e1=(*currentTriangle1)[1];          //  ffout(0) << scene.object(n).reference() << ":"<< ((*i).first) << '\n';
1714        const Edge::Pair e2=(*currentTriangle1)[2];            if (scene.object(n).reference() == ((*i).first)) {
1715        const Vertex& V0 = *((*(*__internals).verticesListes[0])[e0]);                otherRef[n] = true;
1716        const Vertex& V1 = *((*(*__internals).verticesListes[0])[e1]);                (*__internals).__IntersectionPoint(scene,ref,otherRef,triangleedge,listofVertex);
       const Vertex& V2 = *((*(*__internals).verticesListes[0])[e2]);  
   
       //equation of the plan  
       TinyVector<3>  V01,V12,V02;  
       V01=V1.Vector3()-V0.Vector3();  
       V12=V2.Vector3()-V1.Vector3();  
       V02=V2.Vector3()-V0.Vector3();  
     //  const TinyVector<3> coeff=V01^V02;  
     //  const real_t d=-(coeff[0]*V0[0]+coeff[1]*V0[1]+coeff[2]*V0[2]);    
 //      std::cout<<"--------------------------------------------------------------\n";  
   
       tt[nn][0]=0;  
       tt[nn][1]=0;  
       tt[nn][2]=0;  
       for(SurfaceMeshGenerator::Internals::EdgePairList::iterator currentTriangle  
           = (*(*__internals).triangleListes[1]).begin();  
           currentTriangle != (*(*__internals).triangleListes[1]).end();  
           ++currentTriangle)  
       {  
           //std::cout<<"*************************\n";  
           const Vertex& Vo0 = *((*(*__internals).verticesListes[1])[(*currentTriangle)[0]]);  
           const Vertex& Vo1 = *((*(*__internals).verticesListes[1])[(*currentTriangle)[1]]);  
           const Vertex& Vo2 = *((*(*__internals).verticesListes[1])[(*currentTriangle)[2]]);  
   
           //equation of plan  
           TinyVector<3>  Vo01,Vo12,Vo02;  
           Vo01=Vo1.Vector3()-Vo0.Vector3();  
           Vo12=Vo2.Vector3()-Vo1.Vector3();  
           Vo02=Vo2.Vector3()-Vo0.Vector3();  
             
           TinyVector<3> coeffo;  
           coeffo=Vo01^Vo02;  
           real_t d2=-(coeffo[0]*Vo0[0]+coeffo[1]*Vo0[1]+coeffo[2]*Vo0[2]);  
 //std::cout<<"verif "<<coeffo[0]*Vo0.Vector3()[0]+coeffo[1]*Vo0.Vector3()[1]+coeffo[2]*Vo0.Vector3()[2]+d2<<"\n";  
 //std::cout<<"verif "<<coeffo[0]*Vo1.Vector3()[0]+coeffo[1]*Vo1.Vector3()[1]+coeffo[2]*Vo1.Vector3()[2]+d2<<"\n";  
 //std::cout<<"verif "<<coeffo[0]*Vo2.Vector3()[0]+coeffo[1]*Vo2.Vector3()[1]+coeffo[2]*Vo2.Vector3()[2]+d2<<"\n";  
           real_t pscalar,k;  
           size_t num;  
           {  
               Vertex P;  
               real_t norm1=Norm(V01);  
               pscalar=coeffo*V01;  
               num=0;  
               if(pscalar>1e-13 or pscalar<-1e-13)  
               {  
                   k=-(coeffo[0]*V0[0]+coeffo[1]*V0[1]+coeffo[2]*V0[2]+d2)/pscalar;  
                   P=V0+k*V01;  
                   TinyVector<3> P0=P.Vector3()-Vo2.Vector3(),P1=Vo12*(-1.),P2=Vo02*(-1.);  
                   real_t norm2=Norm(P.Vector3()-V0.Vector3());  
                   real_t norm3=Norm(P.Vector3()-V1.Vector3());  
                   if(norm2<=norm1 and norm3<=norm1)  
                   {  
                       Bary(P,P0,P1,P2,num);  
                       if(num==3)  
                       {  
 //std::cout<<"droite1 "<<P.Vector3()[0]<<" "<<P.Vector3()[1]<<" "<<P.Vector3()[2]<<"  \n";  
                           tt[nn][0]=1;  
                           vv[nn][0]=P;  
 //std::cout<<"P "<<nn<<" "<< vv[0][0].Vector3()[0]<<" "<<(vv[0][0]).Vector3()[1]<<" "<<(vv[0][0]).Vector3()[2]<<"\n";  
 //std::cout<<"verif "<<coeffo[0]*P.Vector3()[0]+coeffo[1]*P.Vector3()[1]+coeffo[2]*P.Vector3()[2]+d2<<"\n";  
 //std::cout<<"verif "<<coeff[0]*P.Vector3()[0]+coeff[1]*P.Vector3()[1]+coeff[2]*P.Vector3()[2]+d<<"\n";  
                       }  
                   }  
               }  
           }  
           {  
               Vertex P;  
               pscalar=coeffo*V02;  
               real_t norm1=Norm(V02);  
               if(pscalar>1e-13 or pscalar<-1e-13)  
               {  
                   k=-(coeffo[0]*V0[0]+coeffo[1]*V0[1]+coeffo[2]*V0[2]+d2)/pscalar;  
                   P=V0+k*V02;  
                   TinyVector<3> P0=P.Vector3()-Vo2.Vector3(),P1=Vo12*(-1.),P2=Vo02*(-1.);  
                   real_t norm2=Norm(P.Vector3()-V0.Vector3());  
                   real_t norm3=Norm(P.Vector3()-V2.Vector3());  
                   //std::cout<<"norm "<<norm3<<" "<<norm2<<" "<<norm1<<"\n";  
                   if(norm2<=norm1 and norm3<=norm1)  
                   {  
                       Bary(P,P0,P1,P2,num);  
                       if(num==3)  
                       {  
 //std::cout<<"droite2 "<<P.Vector3()[0]<<" "<<P.Vector3()[1]<<" "<<P.Vector3()[2]<<"  \n";  
                           tt[nn][1]=1;  
                           (vv[nn])[1]=P;  
                       }  
                   }  
               }  
           }  
           {  
               Vertex P;  
               pscalar=coeffo*V12;  
               real_t norm1=Norm(V12);  
               if(pscalar>1e-13 or pscalar<-1e-13)  
               {  
                   k=-(coeffo[0]*V1[0]+coeffo[1]*V1[1]+coeffo[2]*V1[2]+d2)/pscalar;  
                   P=V1+k*V12;  
                   TinyVector<3> P0=P.Vector3()-Vo2.Vector3(),P1=Vo12*(-1.),P2=Vo02*(-1.);  
                   real_t norm2=Norm(P.Vector3()-V1.Vector3());  
                   real_t norm3=Norm(P.Vector3()-V2.Vector3());  
                     
                   if(norm2<=norm1 and norm3<=norm1)  
                   {  
                       Bary(P,P0,P1,P2,num);  
                       if(num==3)  
                       {  
                           tt[nn][2]=1;  
                           (vv[nn])[2]=P;  
   
                       }  
   
                   }  
               }  
1717            }            }
1718        }        }
       nn++;  
1719    }    }
1720    //std::cout<<"Papres "<<nn<<" "<<nobj<<" "<<(vv[0][0]).Vector3()[0]<<" "<<(vv[0][0]).Vector3()[1]<<" "<<(vv[0][0]).Vector3()[2]<<"\n";    
   
1721    // create triangles    // create triangles
1722      std::vector<Vertex> vertexplot;
1723  std::vector<Vertex> vertexplot;    {
1724  {    nn=0;      int nn=0;
1725      size_t ii0,ii1,ii2;      size_t ii0,ii1,ii2;
1726      for(SurfaceMeshGenerator::Internals::EdgePairList::const_iterator currentTriangle      for(SurfaceMeshGenerator::Internals::EdgePairList::const_iterator currentTriangle
1727          = (*(*__internals).triangleListes[0]).begin();            = (*(*__internals).triangleListes[0]).begin();
1728          currentTriangle != (*(*__internals).triangleListes[0]).end();          currentTriangle != (*(*__internals).triangleListes[0]).end();
1729          ++currentTriangle)          ++currentTriangle) {
1730      {        const Edge& e0=(*currentTriangle)[0];
1731          const Edge& e0=(*currentTriangle)[0];        const Edge& e1=(*currentTriangle)[1];
1732          const Edge& e1=(*currentTriangle)[1];        const Edge& e2=(*currentTriangle)[2];
1733          const Edge& e2=(*currentTriangle)[2];        
1734  //      std::cout<<"ref "<<nn<<" "<<(*__internals).edgesRef[e0]<<" "<<(*__internals).edgesRef[e1]<<" "<<(*__internals).edgesRef[e2]<<"\n";       // if(tn0+tn1+tn2!=0) {
1735  //std::cout<<"tt "<<tt[nn][0]<<" "<<tt[nn][1]<<" "<<tt[nn][2]<<"\n";  //std::cout<<"-------------------------------------------------"<<"\n";
1736    //std::cout<<"tt "<<nn<<" "<<tn0<<" "<<tn1<<" "<<tn2<<"\n";
1737    //std::cout<<"avt "<<nn<<" "<<(*__internals).edgesRef[e0]<<" "<<(*__internals).edgesRef[e1]<<" "<<(*__internals).edgesRef[e2]<<"\n";
1738                 inout(1,nn,tt[nn],(*__internals).edgesRef[e0],(*__internals).edgesRef[e1],(*__internals).edgesRef[e2],ii0,ii1,ii2);        //}
1739          (*__internals).__inout(1,nn,triangleedge[e0],triangleedge[e1],triangleedge[e2],
1740                                 (*__internals).edgesRef[e0],(*__internals).edgesRef[e1]
1741                                 ,(*__internals).edgesRef[e2],ii0,ii1,ii2);
1742          //if(tn0+tn1+tn2!=0) {
1743    //std::cout<<"apres "<<nn<<" "<<(*__internals).edgesRef[e0]<<" "<<(*__internals).edgesRef[e1]<<" //"<<(*__internals).edgesRef[e2]<<"\n";
1744          //}
1745          nn++;
1746    
1747  //std::cout<<"ii "<<ii0<<" "<<ii1<<" "<<ii2<<"\n";      
 //             std::cout<<"apres ref "<<(*__internals).edgesRef[e0]<<" "<<(*__internals).edgesRef[e1]<<" "<<(*__internals).edgesRef[e2]<<"\n";  
                nn++;  
1748      }      }
1749      nn=0;  ffout(0) << "fin premiere boucle\n";
         for(SurfaceMeshGenerator::Internals::EdgePairList::const_iterator currentTriangle  
      = (*(*__internals).triangleListes[0]).begin();  
      currentTriangle != (*(*__internals).triangleListes[0]).end();  
      ++currentTriangle)  
         {  
             const Edge& e0=(*currentTriangle)[0];  
             const Edge& e1=(*currentTriangle)[1];  
             const Edge& e2=(*currentTriangle)[2];  
   
   
             inout(2,nn,tt[nn],(*__internals).edgesRef[e0],(*__internals).edgesRef[e1],(*__internals).edgesRef[e2],ii0,ii1,ii2);  
   
   
             nn++;  
         }  
         nn=0;  
         for(SurfaceMeshGenerator::Internals::EdgePairList::const_iterator currentTriangle  
      = (*(*__internals).triangleListes[0]).begin();  
      currentTriangle != (*(*__internals).triangleListes[0]).end();  
      ++currentTriangle)  
         {  
             const Edge& e0=(*currentTriangle)[0];  
             const Edge& e1=(*currentTriangle)[1];  
             const Edge& e2=(*currentTriangle)[2];  
   
             inout(2,nn,tt[nn],(*__internals).edgesRef[e0],(*__internals).edgesRef[e1],(*__internals).edgesRef[e2],ii0,ii1,ii2);  
   
 //std::cout<<ii0<<" "<<ii1<<" "<<ii2<<"\n";  
               
             nn++;  
         }  
           
1750      nn=0;      nn=0;
1751      for(SurfaceMeshGenerator::Internals::EdgePairList::const_iterator currentTriangle      for(SurfaceMeshGenerator::Internals::EdgePairList::const_iterator currentTriangle
1752          = (*(*__internals).triangleListes[0]).begin();            = (*(*__internals).triangleListes[0]).begin();
1753          currentTriangle != (*(*__internals).triangleListes[0]).end();          currentTriangle != (*(*__internals).triangleListes[0]).end();
1754          ++currentTriangle)          ++currentTriangle) {
1755      {        const Edge& e0=(*currentTriangle)[0];
1756          //if the triangle is intersected        const Edge& e1=(*currentTriangle)[1];
1757  //std::cout<<"**************************************************"<<"\n";        const Edge& e2=(*currentTriangle)[2];
1758           const Edge& e0=(*currentTriangle)[0];        if(tt[nn][0]+tt[nn][1]+tt[nn][2]!=0) {
          const Edge& e1=(*currentTriangle)[1];  
        const Edge& e2=(*currentTriangle)[2];  
 //std::cout<<"  avt i123 "<<nn<<" "<<ii0<<" "<<ii1<<" "<<ii2<<" \n";  
        if(tt[nn][0]+tt[nn][1]+tt[nn][2]!=0)  
        {  
1759  //std::cout<<"-------------------------------------------------"<<"\n";  //std::cout<<"-------------------------------------------------"<<"\n";
1760  //std::cout<<"tt "<<nn<<" "<<tt[nn][0]<<" "<<tt[nn][1]<<" "<<tt[nn][2]<<"\n";  //std::cout<<"tt "<<nn<<" "<<tt[nn][0]<<" "<<tt[nn][1]<<" "<<tt[nn][2]<<"\n";
        }  
   
 //std::cout<<"ref "<<(*__internals).edgesRef[e0]<<" "<<(*__internals).edgesRef[e1]<<" "<<(*__internals).edgesRef[e2]<<"\n";  
 //std::cout<<"inter "<<nn<<" "<<(vv[0][0]).Vector3()[0]<<" "<<(vv[0][0]).Vector3()[1]<<" "<<(vv[0][0]).Vector3()[2]<<" \n";  
          // if(triangleedge[(*currentTriangle)[0]]+triangleedge[(*currentTriangle)[1]]+ triangleedge[(*currentTriangle)[2]]==2)  
 //std::cout<<"avt "<<(*__internals).edgesRef[e0]<<" "<<(*__internals).edgesRef[e1]<<" "<<(*__internals).edgesRef[e2]<<"\n";  
   
   //     if(nn==1254 or nn==1248 or nn==1252 or nn==1255 or nn==1249 )  
        {  
 //std::cout<<"-------------------------------------------------"<<"\n";  
   
1761  //std::cout<<"avt "<<nn<<" "<<(*__internals).edgesRef[e0]<<" "<<(*__internals).edgesRef[e1]<<" "<<(*__internals).edgesRef[e2]<<"\n";  //std::cout<<"avt "<<nn<<" "<<(*__internals).edgesRef[e0]<<" "<<(*__internals).edgesRef[e1]<<" "<<(*__internals).edgesRef[e2]<<"\n";
1762         }        }
1763         inout(3,nn,tt[nn],(*__internals).edgesRef[e0],(*__internals).edgesRef[e1],(*__internals).edgesRef[e2],ii0,ii1,ii2);        (*__internals).__inout(2,nn,triangleedge[e0],triangleedge[e1],triangleedge[e2],
1764                                 (*__internals).edgesRef[e0],(*__internals).edgesRef[e1]
1765      //   if(nn==1254 or nn==1248 or nn==1252 or nn==1255 or nn==1249 )                               ,(*__internals).edgesRef[e2],ii0,ii1,ii2);
1766         {        if(tt[nn][0]+tt[nn][1]+tt[nn][2]!=0) {
1767  //std::cout<<"apres "<<nn<<" "<<(*__internals).edgesRef[e0]<<" "<<(*__internals).edgesRef[e1]<<" "<<(*__internals).edgesRef[e2]<<"\n";  //std::cout<<"apres "<<nn<<" "<<(*__internals).edgesRef[e0]<<" "<<(*__internals).edgesRef[e1]<<" //"<<(*__internals).edgesRef[e2]<<"\n";
1768  //std::cout<<"ii "<<ii0<<" "<<ii1<<" "<<ii2<<"\n";        }
1769         }        
1770  //std::cout<<"apres "<<(*__internals).edgesRef[e0]<<" "<<(*__internals).edgesRef[e1]<<" "<<(*__internals).edgesRef[e2]<<"\n";        nn++;
1771  //std::cout<<" i123 "<<nn<<" "<<ii0<<" "<<ii1<<" "<<ii2<<" \n";      }
1772    ffout(0) << "fin deuxieme boucle\n";
1773    
1774  //std::cout<<e0.Ref()<<" "<<e1.Ref()<<" "<<e2.Ref()<<"\n";      nn=0;
1775          if(tt[nn][0]+tt[nn][1]+tt[nn][2]==2)      for(SurfaceMeshGenerator::Internals::EdgePairList::const_iterator currentTriangle
1776          {            = (*(*__internals).triangleListes[0]).begin();
1777            currentTriangle != (*(*__internals).triangleListes[0]).end();
1778            ++currentTriangle) {
1779          const Edge& e0=(*currentTriangle)[0];
1780          const Edge& e1=(*currentTriangle)[1];
1781          const Edge& e2=(*currentTriangle)[2];
1782          if(tt[nn][0]+tt[nn][1]+tt[nn][2]!=0) {
1783  //std::cout<<"-------------------------------------------------"<<"\n";  //std::cout<<"-------------------------------------------------"<<"\n";
1784              Vertex& V0 = *((*(*__internals).verticesListes[0])[(*currentTriangle)[0]]);  //std::cout<<"tt "<<nn<<" "<<tt[nn][0]<<" "<<tt[nn][1]<<" "<<tt[nn][2]<<"\n";
1785              Vertex& V1 = *((*(*__internals).verticesListes[0])[(*currentTriangle)[1]]);  //std::cout<<"avt "<<nn<<" "<<(*__internals).edgesRef[e0]<<" "<<(*__internals).edgesRef[e1]<<" "<<(*__internals).edgesRef[e2]<<"\n";
1786              Vertex& V2 = *((*(*__internals).verticesListes[0])[(*currentTriangle)[2]]);        }
1787  //std::cout<<"VO "<<V0.Vector3()[0]<<"  "<<V0.Vector3()[1]<<" "<<V0.Vector3()[2]<<"\n";        (*__internals).__inout(3,nn,triangleedge[e0],triangleedge[e1],triangleedge[e2],
1788  //std::cout<<"V1 "<<V1.Vector3()[0]<<"  "<<V1.Vector3()[1]<<" "<<V1.Vector3()[2]<<"\n";                               (*__internals).edgesRef[e0],(*__internals).edgesRef[e1]
1789  //std::cout<<"V2 "<<V2.Vector3()[0]<<"  "<<V2.Vector3()[1]<<" "<<V2.Vector3()[2]<<"\n";                               ,(*__internals).edgesRef[e2],ii0,ii1,ii2);
1790  //std::cout<<" "<<"\n";        if(tt[nn][0]+tt[nn][1]+tt[nn][2]!=0) {
1791  //     **************** DELETE************************************************  //std::cout<<"apres "<<nn<<" "<<(*__internals).edgesRef[e0]<<" "<<(*__internals).edgesRef[e1]<<" //"<<(*__internals).edgesRef[e2]<<"\n";
1792                 bool& i1=V0.Ref().insideDomain();        }
1793                 bool& i2=V1.Ref().insideDomain();        
1794                 bool& i3=V2.Ref().insideDomain();        nn++;
1795                 if(ii0==1)  }
1796                     i1=true;  ffout(0) << "fin troisieme boucle\n";
                else  
                    i1=false;  
                if(ii1==1)  
                    i2=true;  
                else  
                    i2=false;  
                if(ii2==1)  
                    i3=true;  
                else  
                    i3=false;  
   
 //std::cout<<"inside? "<<i1<<" "<<i2<<" "<<i3<<"\n";  
 //     ****************************************************************  
1797    
1798                //if(triangleedge[(*currentTriangle)[0]]+triangleedge[(*currentTriangle)[1]]==2)      nn=0;
1799                if(tt[nn][0]+tt[nn][1]==2)      for(SurfaceMeshGenerator::Internals::EdgePairList::const_iterator currentTriangle
1800                {            = (*(*__internals).triangleListes[0]).begin();
1801                    if(V0.Ref().insideDomain())          currentTriangle != (*(*__internals).triangleListes[0]).end();
1802                    {          ++currentTriangle) {
1803                        Vertex A=V0;        //if the triangle is intersected
1804                       // Vertex B=(listofVertex[(*currentTriangle)[0]]);        //std::cout<<"**************************************************"<<"\n";
1805                        //Vertex C=(listofVertex[(*currentTriangle)[1]]);        const Edge& e0=(*currentTriangle)[0];
1806                        Vertex B=(vv[nn][0]);        const Edge& e1=(*currentTriangle)[1];
1807                        Vertex C=(vv[nn][1]);        const Edge& e2=(*currentTriangle)[2];
1808                        vertexplot.push_back(A);        int tn0=triangleedge[e0];
1809                        vertexplot.push_back(B);        int tn1=triangleedge[e1];
1810                        vertexplot.push_back(C);        int tn2=triangleedge[e2];
1811  //std::cout<<"A "<<A.Vector3()[0]<<"  "<<A.Vector3()[1]<<" "<<A.Vector3()[2]<<"\n";        //std::cout<<"  avt i123 "<<nn<<" "<<ii0<<" "<<ii1<<" "<<ii2<<" \n";
1812  //std::cout<<"B "<<B.Vector3()[0]<<"  "<<B.Vector3()[1]<<" "<<B.Vector3()[2]<<"\n";        if(tn0+tn1+tn2!=0) {
1813  //std::cout<<"C "<<C.Vector3()[0]<<"  "<<C.Vector3()[1]<<" "<<C.Vector3()[2]<<"\n";          std::cout<<"-------------------------------------------------"<<"\n";
1814                    }          std::cout<<"tt "<<nn<<" "<<tn0<<" "<<tn1<<" "<<tn2<<"\n";
1815                    else if(V1.Ref().insideDomain() and V2.Ref().insideDomain())            std::cout<<"avt "<<nn<<" "<<(*__internals).edgesRef[e0]<<" "<<(*__internals).edgesRef[e1]<<" "<<(*__internals).edgesRef[e2]<<"\n";
1816                    {        }
1817                        Vertex A=V1;        
1818                        Vertex B=V2;  
1819                       // Vertex C=(listofVertex[(*currentTriangle)[0]]);        //std::cout<<"ref "<<(*__internals).edgesRef[e0]<<" "<<(*__internals).edgesRef[e1]<<" "<<(*__internals).edgesRef[e2]<<"\n";
1820                        //Vertex D=(listofVertex[(*currentTriangle)[1]]);        //std::cout<<"inter "<<nn<<" "<<(vv[0][0]).Vector3()[0]<<" "<<(vv[0][0]).Vector3()[1]<<" "<<(vv[0][0]).Vector3()[2]<<" \n";
1821                        Vertex C=(vv[nn][0]);        // if(triangleedge[(*currentTriangle)[0]]+triangleedge[(*currentTriangle)[1]]+ triangleedge[(*currentTriangle)[2]]==2)
1822                        Vertex D=(vv[nn][1]);        //std::cout<<"avt "<<(*__internals).edgesRef[e0]<<" "<<(*__internals).edgesRef[e1]<<" "<<(*__internals).edgesRef[e2]<<"\n";
1823                        vertexplot.push_back(A);  
1824                        vertexplot.push_back(C);        //std::cout<<"-------------------------------------------------"<<"\n";
1825                        vertexplot.push_back(D);  
1826                        vertexplot.push_back(A);        (*__internals).__inout(4,nn,tn0,tn1,tn2,
1827                        vertexplot.push_back(B);                               (*__internals).edgesRef[e0],(*__internals).edgesRef[e1]
1828                        vertexplot.push_back(D);                               ,(*__internals).edgesRef[e2],ii0,ii1,ii2);
1829                    }        if(tn0+tn1+tn2!=0) {
1830                    else  
1831                    {          std::cout<<"apres "<<nn<<" "<<(*__internals).edgesRef[e0]<<" "<<(*__internals).edgesRef[e1]<<" "<<(*__internals).edgesRef[e2]<<"\n";
1832  std::cout<<"erreur1\n";          std::cout<<"ii "<<ii0<<" "<<ii1<<" "<<ii2<<"\n";
1833                    }        }
1834                }        if(tn0+tn1+tn2==2) {
1835               // else if(triangleedge[(*currentTriangle)[0]]+triangleedge[(*currentTriangle)[2]]==2)          std::cout<<"-------------------------------------------------"<<"\n";
1836                else if(tt[nn][0]+tt[nn][2]==2)          Vertex& V0 = *((*(*__internals).verticesListes[0])[(*currentTriangle)[0]]);
1837                {          Vertex& V1 = *((*(*__internals).verticesListes[0])[(*currentTriangle)[1]]);
1838                    if(V1.Ref().insideDomain())          Vertex& V2 = *((*(*__internals).verticesListes[0])[(*currentTriangle)[2]]);
1839                    {          std::cout<<"VO "<<V0.Vector3()[0]<<"  "<<V0.Vector3()[1]<<" "<<V0.Vector3()[2]<<"\n";
1840                        Vertex A=V1;          std::cout<<"V1 "<<V1.Vector3()[0]<<"  "<<V1.Vector3()[1]<<" "<<V1.Vector3()[2]<<"\n";
1841                       // Vertex B=(listofVertex[(*currentTriangle)[0]]);          std::cout<<"V2 "<<V2.Vector3()[0]<<"  "<<V2.Vector3()[1]<<" "<<V2.Vector3()[2]<<"\n";
1842                       // Vertex C=(listofVertex[(*currentTriangle)[2]]);          std::cout<<" "<<"\n";
1843                        Vertex B=(vv[nn][0]);  
1844                        Vertex C=(vv[nn][2]);          // Pour Cécile :
1845                        vertexplot.push_back(A);          //  Cette partie m'a l'air assez dangereuse : tu es en train
1846                        vertexplot.push_back(B);          // de modifier un caractère essentiel de certains points du
1847                        vertexplot.push_back(C);          // maillage : s'ils sont ou non à l'intérieur du domaine !
1848            // (à cause de la référence bool&)
1849    
1850            bool i0=(ii0==1);
1851            bool i1 = (ii1==1);
1852            bool i2 = (ii2==1);
1853    
1854    std::cout<<"inside? "<<i0<<" "<<i1<<" "<<i2<<"\n";
1855    
1856            if(tn0+tn1==2) {
1857              if(i0) {
1858                Vertex A=V0;
1859  //std::cout<<"A "<<A.Vector3()[0]<<"  "<<A.Vector3()[1]<<" "<<A.Vector3()[2]<<"\n";  //std::cout<<"A "<<A.Vector3()[0]<<"  "<<A.Vector3()[1]<<" "<<A.Vector3()[2]<<"\n";
 //std::cout<<"B "<<B.Vector3()[0]<<"  "<<B.Vector3()[1]<<" "<<B.Vector3()[2]<<"\n";  
 //std::cout<<"C "<<C.Vector3()[0]<<"  "<<C.Vector3()[1]<<" "<<C.Vector3()[2]<<"\n";  
                   }  
                   else if(V0.Ref().insideDomain() and V2.Ref().insideDomain())  
                   {  
                       Vertex A=V0;  
                       Vertex B=V2;  
                    //   Vertex C=(listofVertex[(*currentTriangle)[0]]);  
                     //  Vertex D=(listofVertex[(*currentTriangle)[2]]);  
                       Vertex C=(vv[nn][0]);  
                       Vertex D=(vv[nn][2]);  
                       vertexplot.push_back(A);  
                       vertexplot.push_back(C);  
                       vertexplot.push_back(D);  
                       vertexplot.push_back(A);  
                       vertexplot.push_back(B);  
                       vertexplot.push_back(D);  
1860    
1861                    }                Vertex& B=(*(listofVertex[e0]));
1862                    else              Vertex& C=(*(listofVertex[e1]));
                   {  
 std::cout<<"erreur2\n";  
                   }  
               }  
               else  
               {  
                   if(V2.Ref().insideDomain())  
                   {  
                       Vertex A=V2;  
                    //   Vertex B=(listofVertex[(*currentTriangle)[1]]);  
                      // Vertex C=(listofVertex[(*currentTriangle)[2]]);  
                       Vertex B=(vv[nn][1]);  
                       Vertex C=(vv[nn][2]);  
                       vertexplot.push_back(A);  
                       vertexplot.push_back(B);  
                       vertexplot.push_back(C);  
 //std::cout<<"A "<<A.Vector3()[0]<<"  "<<A.Vector3()[1]<<" "<<A.Vector3()[2]<<"\n";  
1863  //std::cout<<"B "<<B.Vector3()[0]<<"  "<<B.Vector3()[1]<<" "<<B.Vector3()[2]<<"\n";  //std::cout<<"B "<<B.Vector3()[0]<<"  "<<B.Vector3()[1]<<" "<<B.Vector3()[2]<<"\n";
1864  //std::cout<<"C "<<C.Vector3()[0]<<"  "<<C.Vector3()[1]<<" "<<C.Vector3()[2]<<"\n";  //std::cout<<"C "<<C.Vector3()[0]<<"  "<<C.Vector3()[1]<<" "<<C.Vector3()[2]<<"\n";
1865                    }              ffout(0) << "----------\n";
1866                    else if(V1.Ref().insideDomain() and V0.Ref().insideDomain())              vertexplot.push_back(A);
1867                    {              vertexplot.push_back(B);
1868                        Vertex A=V0;              vertexplot.push_back(C);
1869                        Vertex B=V1;              
1870                    //    Vertex C=(listofVertex[(*currentTriangle)[1]]);            } else if(i1 and i2) {
1871                      //  Vertex D=(listofVertex[(*currentTriangle)[2]]);              Vertex A=V1;
1872                        Vertex C=(vv[nn][1]);              Vertex B=V2;
1873                        Vertex D=(vv[nn][2]);              Vertex& C=(*(listofVertex[e0]));
1874                        vertexplot.push_back(A);              Vertex& D=(*(listofVertex[e1]));
1875                        vertexplot.push_back(C);              vertexplot.push_back(A);
1876                        vertexplot.push_back(D);              vertexplot.push_back(C);
1877                        vertexplot.push_back(A);              vertexplot.push_back(D);
1878                        vertexplot.push_back(B);              vertexplot.push_back(A);
1879                        vertexplot.push_back(D);              vertexplot.push_back(B);
1880                    }              vertexplot.push_back(D);
1881                    else            } else {
1882                    {              std::cout<<"erreur1\n";
1883  std::cout<<"erreur3\n";            }
                   }  
               }  
1884          }          }
1885          if(tt[nn][0]+tt[nn][1]+tt[nn][2]==3)          // else if(triangleedge[(*currentTriangle)[0]]+triangleedge[(*currentTriangle)[2]]==2)
1886            else if(tn0+tn2==2) {
1887              if(i1) {
1888                Vertex A=V1;
1889                Vertex& B=(*(listofVertex[e0]));
1890                Vertex& C=(*(listofVertex[e2]));
1891                vertexplot.push_back(A);
1892                vertexplot.push_back(B);
1893                vertexplot.push_back(C);
1894                //std::cout<<"A "<<A.Vector3()[0]<<"  "<<A.Vector3()[1]<<" "<<A.Vector3()[2]<<"\n";
1895                //std::cout<<"B "<<B.Vector3()[0]<<"  "<<B.Vector3()[1]<<" "<<B.Vector3()[2]<<"\n";
1896                //std::cout<<"C "<<C.Vector3()[0]<<"  "<<C.Vector3()[1]<<" "<<C.Vector3()[2]<<"\n";
1897              } else if(i0 and i2) {
1898                Vertex A=V0;
1899                Vertex B=V2;
1900                Vertex& C=(*(listofVertex[e0]));
1901                Vertex& D=(*(listofVertex[e2]));
1902                vertexplot.push_back(A);
1903                vertexplot.push_back(C);
1904                vertexplot.push_back(D);
1905                vertexplot.push_back(A);
1906                vertexplot.push_back(B);
1907                vertexplot.push_back(D);
1908              } else {
1909                std::cout<<"erreur2\n";
1910              }
1911            } else {
1912              if(i2) {
1913                Vertex A=V2;
1914                Vertex& B=(*(listofVertex[e1]));
1915                Vertex& C=(*(listofVertex[e2]));
1916                vertexplot.push_back(A);
1917                vertexplot.push_back(B);
1918                vertexplot.push_back(C);
1919             //   std::cout<<"A "<<A.Vector3()[0]<<"  "<<A.Vector3()[1]<<" "<<A.Vector3()[2]<<"\n";
1920                //std::cout<<"B "<<B.Vector3()[0]<<"  "<<B.Vector3()[1]<<" "<<B.Vector3()[2]<<"\n";
1921                //std::cout<<"C "<<C.Vector3()[0]<<"  "<<C.Vector3()[1]<<" "<<C.Vector3()[2]<<"\n";
1922              } else if(i1 and i0) {
1923                Vertex A=V0;
1924                Vertex B=V1;
1925                std::cout<<"A "<<A.Vector3()[0]<<"  "<<A.Vector3()[1]<<" "<<A.Vector3()[2]<<"\n";
1926                std::cout<<"B "<<B.Vector3()[0]<<"  "<<B.Vector3()[1]<<" "<<B.Vector3()[2]<<"\n";
1927                Vertex& C=(*(listofVertex[e1]));
1928                Vertex& D=(*(listofVertex[e2]));
1929                std::cout<<"C "<<C.Vector3()[0]<<"  "<<C.Vector3()[1]<<" "<<C.Vector3()[2]<<"\n";
1930                vertexplot.push_back(A);
1931                vertexplot.push_back(C);
1932                vertexplot.push_back(D);
1933                vertexplot.push_back(A);
1934                vertexplot.push_back(B);
1935                vertexplot.push_back(D);
1936              } else {
1937                std::cout<<"erreur3\n";
1938              }
1939            }
1940          }
1941          if(tn0+tn1+tn2==3) {
1942            Vertex& V0 = *((*(*__internals).verticesListes[0])[(*currentTriangle)[0]]);
1943            Vertex& V1 = *((*(*__internals).verticesListes[0])[(*currentTriangle)[1]]);
1944            Vertex& V2 = *((*(*__internals).verticesListes[0])[(*currentTriangle)[2]]);
1945    
1946            // Pour Cécile : comme plus haut...
1947    
1948            bool i0=(ii0==1);
1949            bool i1 = (ii1==1);
1950            bool i2 = (ii2==1);
1951    
1952            std::cout<<"be careful triangle \n";
1953            //std::cout<<"inside? "<<V0.Ref().insideDomain()<<" "<<V1.Ref().insideDomain()<<" "<<V2.Ref().insideDomain()<<"\n";
1954            if(i0+i1+i2==1)
1955          {          {
1956              Vertex& V0 = *((*(*__internals).verticesListes[0])[(*currentTriangle)[0]]);          if(i0) {
1957              Vertex& V1 = *((*(*__internals).verticesListes[0])[(*currentTriangle)[1]]);            std::cout<<"be careful triangle 0\n";
1958              Vertex& V2 = *((*(*__internals).verticesListes[0])[(*currentTriangle)[2]]);            Vertex A=V0;
1959              bool& i1=V0.Ref().insideDomain();            Vertex& B=(*(listofVertex[e0]));
1960              bool& i2=V1.Ref().insideDomain();            Vertex& C=(*(listofVertex[e1]));
1961              bool& i3=V2.Ref().insideDomain();            Vertex& D=(*(listofVertex[e2]));
1962              if(ii0==1)            vertexplot.push_back(A);
1963                  i1=true;            vertexplot.push_back(B);
1964              else            vertexplot.push_back(D);
1965                  i1=false;            vertexplot.push_back(A);
1966              if(ii1==1)            vertexplot.push_back(C);
1967                  i2=true;            vertexplot.push_back(D);
1968              else            //std::cout<<"A "<<A.Vector3()[0]<<"  "<<A.Vector3()[1]<<" "<<A.Vector3()[2]<<"\n";
1969                  i2=false;            //std::cout<<"B "<<B.Vector3()[0]<<"  "<<B.Vector3()[1]<<" "<<B.Vector3()[2]<<"\n";
1970              if(ii2==1)            //std::cout<<"C "<<C.Vector3()[0]<<"  "<<C.Vector3()[1]<<" "<<C.Vector3()[2]<<"\n";
1971                  i3=true;            //std::cout<<"D "<<D.Vector3()[0]<<"  "<<D.Vector3()[1]<<" "<<D.Vector3()[2]<<"\n";
1972              else          } else if(i1) {
1973                  i3=false;            std::cout<<"be careful triangle 1 \n";
1974  //std::cout<<"be careful triangle \n";            Vertex A=V1;
1975  //std::cout<<"inside? "<<V0.Ref().insideDomain()<<" "<<V1.Ref().insideDomain()<<" "<<V2.Ref().insideDomain()<<"\n";            Vertex& B=(*(listofVertex[e0]));
1976              if(V0.Ref().insideDomain())            Vertex& C=(*(listofVertex[e1]));
1977              Vertex& D=(*(listofVertex[e2]));
1978              vertexplot.push_back(A);
1979              vertexplot.push_back(C);
1980              vertexplot.push_back(D);
1981              vertexplot.push_back(A);
1982              vertexplot.push_back(C);
1983              vertexplot.push_back(B);
1984            } else if(i2) {
1985              std::cout<<"be careful triangle 2 \n";
1986              Vertex A=V2;
1987              Vertex& B=(*(listofVertex[e0]));
1988              Vertex& C=(*(listofVertex[e1]));
1989              Vertex& D=(*(listofVertex[e2]));
1990              vertexplot.push_back(A);
1991              vertexplot.push_back(C);
1992              vertexplot.push_back(D);
1993              vertexplot.push_back(A);
1994              vertexplot.push_back(B);
1995              vertexplot.push_back(D);
1996              //std::cout<<"A "<<A.Vector3()[0]<<"  "<<A.Vector3()[1]<<" "<<A.Vector3()[2]<<"\n";
1997              //std::cout<<"B "<<B.Vector3()[0]<<"  "<<B.Vector3()[1]<<" "<<B.Vector3()[2]<<"\n";
1998              //std::cout<<"C "<<C.Vector3()[0]<<"  "<<C.Vector3()[1]<<" "<<C.Vector3()[2]<<"\n";
1999              //std::cout<<"D "<<D.Vector3()[0]<<"  "<<D.Vector3()[1]<<" "<<D.Vector3()[2]<<"\n";
2000            }
2001            }
2002            else
2003            {
2004                if(i0+i1==2)
2005              {              {
 std::cout<<"be careful triangle 0\n";  
2006                  Vertex A=V0;                  Vertex A=V0;
2007                  Vertex B=(vv[nn][0]);                  Vertex B=V1;
2008                  Vertex C=(vv[nn][1]);                  Vertex& C=(*(listofVertex[e1]));
2009                  Vertex D=(vv[nn][2]);                  Vertex& D=(*(listofVertex[e2]));
2010                  vertexplot.push_back(A);                  vertexplot.push_back(A);
2011                  vertexplot.push_back(B);                  vertexplot.push_back(C);
2012                  vertexplot.push_back(D);                  vertexplot.push_back(D);
2013                  vertexplot.push_back(A);                  vertexplot.push_back(A);
2014                  vertexplot.push_back(C);                  vertexplot.push_back(B);
2015                  vertexplot.push_back(D);                  vertexplot.push_back(D);
 //std::cout<<"A "<<A.Vector3()[0]<<"  "<<A.Vector3()[1]<<" "<<A.Vector3()[2]<<"\n";  
 //std::cout<<"B "<<B.Vector3()[0]<<"  "<<B.Vector3()[1]<<" "<<B.Vector3()[2]<<"\n";  
 //std::cout<<"C "<<C.Vector3()[0]<<"  "<<C.Vector3()[1]<<" "<<C.Vector3()[2]<<"\n";  
 //std::cout<<"D "<<D.Vector3()[0]<<"  "<<D.Vector3()[1]<<" "<<D.Vector3()[2]<<"\n";  
   
2016              }              }
2017              else if(V1.Ref().insideDomain())              else if(i0+i2==2)
2018              {              {
2019  std::cout<<"be careful triangle 1 \n";                  Vertex A=V0;
2020                  Vertex A=V1;                  Vertex B=V2;
2021                  Vertex B=(vv[nn][0]);                  Vertex& C=(*(listofVertex[e0]));
2022                  Vertex C=(vv[nn][1]);                  Vertex& D=(*(listofVertex[e2]));
                 Vertex D=(vv[nn][2]);  
2023                  vertexplot.push_back(A);                  vertexplot.push_back(A);
2024                  vertexplot.push_back(C);                  vertexplot.push_back(C);
2025                  vertexplot.push_back(D);                  vertexplot.push_back(D);
2026                  vertexplot.push_back(A);                  vertexplot.push_back(A);
                 vertexplot.push_back(C);  
2027                  vertexplot.push_back(B);                  vertexplot.push_back(B);
2028                                    vertexplot.push_back(D);
2029              }              }
2030              else if(V2.Ref().insideDomain())              else if(i1+i2==2)
2031              {              {
2032  std::cout<<"be careful triangle 2 \n";                  Vertex A=V1;
2033                  Vertex A=V2;                  Vertex B=V2;
2034                  Vertex B=(vv[nn][0]);                  Vertex& C=(*(listofVertex[e0]));
2035                  Vertex C=(vv[nn][1]);                  Vertex& D=(*(listofVertex[e1]));
                 Vertex D=(vv[nn][2]);  
2036                  vertexplot.push_back(A);                  vertexplot.push_back(A);
2037                  vertexplot.push_back(C);                  vertexplot.push_back(C);
2038                  vertexplot.push_back(D);                  vertexplot.push_back(D);
2039                  vertexplot.push_back(A);                  vertexplot.push_back(A);
2040                  vertexplot.push_back(B);                  vertexplot.push_back(B);
2041                  vertexplot.push_back(D);                  vertexplot.push_back(D);
 //std::cout<<"A "<<A.Vector3()[0]<<"  "<<A.Vector3()[1]<<" "<<A.Vector3()[2]<<"\n";  
 //std::cout<<"B "<<B.Vector3()[0]<<"  "<<B.Vector3()[1]<<" "<<B.Vector3()[2]<<"\n";  
 //std::cout<<"C "<<C.Vector3()[0]<<"  "<<C.Vector3()[1]<<" "<<C.Vector3()[2]<<"\n";  
 //std::cout<<"D "<<D.Vector3()[0]<<"  "<<D.Vector3()[1]<<" "<<D.Vector3()[2]<<"\n";  
                   
                   
2042              }              }
2043                else
2044    std::cout<<"3 points in\n";
2045          }          }
2046         nn++;        }
2047          nn++;
2048      }      }
2049          }*/
 }  
 */  
 size_t ntr=0;  
 for(size_t n = 0; n<scene.nbObjects(); ++n)  
 ntr+=(*(*__internals).verticesListes[n]).size();  
2050    
2051  s_mesh.setNumberOfVertices(ntr);    size_t ntr=0;
2052      for(size_t n = 0; n<scene.nbObjects(); ++n)
2053        ntr+=(*(*__internals).verticesListes[n]).size();
2054    
2055      s_mesh.setNumberOfVertices(ntr);
2056    
2057  // //  s_mesh.setNumberOfVertices((*__internals).verticesList.size());    // //  s_mesh.setNumberOfVertices((*__internals).verticesList.size());
2058    
2059  //std::cout<<"nb tot "<<s_mesh.numberOfVertices()<<"\n";    std::cout<<"nb tot "<<s_mesh.numberOfVertices()<<"\n";
2060    
2061    //! copies vertices into the s_mesh.    //! copies vertices into the s_mesh.
2062    int i = 0;    int i = 0;
2063  //    for(SurfaceMeshGenerator::Internals::VerticesList::iterator splitVertex    ////   for(SurfaceMeshGenerator::Internals::VerticesList::iterator splitVertex
2064  //      =  (*__internals).verticesList.begin();    ////  =  (*__internals).verticesList.begin();
2065  //      splitVertex != (*__internals).verticesList.end(); ++splitVertex) {    ////  splitVertex != (*__internals).verticesList.end(); ++splitVertex) {
2066  //  
2067  //  s_mesh.vertex(i) = (*(*splitVertex).second);    // //s_mesh.vertex(i) = (*(*splitVertex).second);
2068  //delete (*splitVertex).second;    ////delete (*splitVertex).second;
2069  // (*splitVertex).second = &s_mesh.vertex(i);    ////(*splitVertex).second = &s_mesh.vertex(i);
2070  // i++;    ////i++;
2071  //}    ////}
2072  for(size_t n = 0; n<scene.nbObjects(); ++n)    for(size_t n = 0; n<scene.nbObjects(); ++n)
2073  {      {
2074      for(SurfaceMeshGenerator::Internals::VerticesList::iterator splitVertex        for(SurfaceMeshGenerator::Internals::VerticesList::iterator splitVertex
2075          =  (*(*__internals).verticesListes[n]).begin();              =  (*(*__internals).verticesListes[n]).begin();
2076          splitVertex != (*(*__internals).verticesListes[n]).end(); ++splitVertex)            splitVertex != (*(*__internals).verticesListes[n]).end(); ++splitVertex)
2077      {          {
2078          s_mesh.vertex(i) = (*(*splitVertex).second);            s_mesh.vertex(i) = (*(*splitVertex).second);
2079          delete (*splitVertex).second;            delete (*splitVertex).second;
2080          (*splitVertex).second = &s_mesh.vertex(i);            (*splitVertex).second = &s_mesh.vertex(i);
2081          i++;            i++;
2082            }
2083      }      }
 }  
2084    
2085    
2086    
2087  //   std::set<const Cell*> toTreatHexahedra;    /*std::set<const Cell*> toTreatHexahedra;
2088  //   std::set_difference (cuttedHexahedra.begin(), cuttedHexahedra.end(),    std::set_difference (cuttedHexahedra.begin(), cuttedHexahedra.end(),
2089  //                     meshedHexahedra.begin(), meshedHexahedra.end(),                         meshedHexahedra.begin(), meshedHexahedra.end(),
2090  //                     inserter(toTreatHexahedra, toTreatHexahedra.begin()));                         inserter(toTreatHexahedra, toTreatHexahedra.begin()));*/
2091        
2092    i = 0;    i = 0;
2093    size_t ncell=0;    size_t ncell=0;
2094    
2095    for(size_t n = 0; n<scene.nbObjects(); ++n)    for(size_t n = 0; n<scene.nbObjects(); ++n)
2096            ncell+=(*(*__internals).triangleListes[n]).size();      ncell+=(*(*__internals).triangleListes[n]).size();
2097    
2098       s_mesh.setNumberOfCells(ncell);    s_mesh.setNumberOfCells(ncell);
2099  //std::cout<<"nb tot "<<s_mesh.numberOfCells()<<" \n";    std::cout<<"nb tot "<<s_mesh.numberOfCells()<<" \n";
2100    
2101    
2102  SurfaceMeshGenerator::Internals::MotherCellList::iterator motherCell    SurfaceMeshGenerator::Internals::MotherCellList::iterator motherCell
2103  = (*__internals).cellList.begin();      = (*__internals).cellList.begin();
2104    
2105    
2106  for(size_t n = 0; n<scene.nbObjects(); ++n)    for(size_t n = 0; n<scene.nbObjects(); ++n)
2107  {      {
2108    
2109      for(SurfaceMeshGenerator::Internals::EdgePairList::iterator currentTriangle        for(SurfaceMeshGenerator::Internals::EdgePairList::iterator currentTriangle
2110          = (*(*__internals).triangleListes[n]).begin();              = (*(*__internals).triangleListes[n]).begin();
2111          currentTriangle != (*(*__internals).triangleListes[n]).end();            currentTriangle != (*(*__internals).triangleListes[n]).end();
2112          ++currentTriangle)            ++currentTriangle)
2113      {          {
2114          //Vertex& V0 = *((*__internals).verticesList[(*currentTriangle)[0]]);            //Vertex& V0 = *((*__internals).verticesList[(*currentTriangle)[0]]);
2115          //Vertex& V1 = *((*__internals).verticesList[(*currentTriangle)[1]]);            //Vertex& V1 = *((*__internals).verticesList[(*currentTriangle)[1]]);
2116          //Vertex& V2 = *((*__internals).verticesList[(*currentTriangle)[2]]);            //Vertex& V2 = *((*__internals).verticesList[(*currentTriangle)[2]]);
2117          Vertex& V0 = *((*(*__internals).verticesListes[n])[(*currentTriangle)[0]]);            Vertex& V0 = *((*(*__internals).verticesListes[n])[(*currentTriangle)[0]]);
2118          Vertex& V1 = *((*(*__internals).verticesListes[n])[(*currentTriangle)[1]]);            Vertex& V1 = *((*(*__internals).verticesListes[n])[(*currentTriangle)[1]]);
2119          Vertex& V2 = *((*(*__internals).verticesListes[n])[(*currentTriangle)[2]]);            Vertex& V2 = *((*(*__internals).verticesListes[n])[(*currentTriangle)[2]]);
2120    
2121      s_mesh.cell(i) = Triangle(V0,V1,V2);            s_mesh.cell(i) = Triangle(V0,V1,V2);
2122  #warning SHOULD USE BC NUMBER (ie OBJECT NUMBER IN SCENE INSTEAD)  #warning SHOULD USE BC NUMBER (ie OBJECT NUMBER IN SCENE INSTEAD)
2123    
2124      s_mesh.cell(i).reference() = 1;            s_mesh.cell(i).reference() = 1;
2125      s_mesh.cell(i).setMother(*motherCell);            s_mesh.cell(i).setMother(*motherCell);
   
     i++;  
     motherCell++;  
   }  
 }  
2126    
2127    //   plot2(s_mesh, toTreatHexahedra, SMesh);            i++;
2128     //  plotmedit(s_mesh, toTreatHexahedra,nobj);            motherCell++;
2129            }
2130        }
2131    
2132  /*{    //  plot2(s_mesh, toTreatHexahedra, SMesh);
2133  std::ofstream sortie("plot.mesh");  /*ffout(0) << "plot medit \n";
2134      plotmedit(s_mesh, toTreatHexahedra,nobj);
2135    ffout(0) << "plot sortie \n";
2136      {
2137        std::ofstream sortie("plot.mesh");
2138    
2139      sortie<<"MeshVersionFormatted 1"<<" \n";      sortie<<"MeshVersionFormatted 1"<<" \n";
2140      sortie<<" \n";      sortie<<" \n";
# Line 1532  std::ofstream sortie("plot.mesh"); Line 2144  std::ofstream sortie("plot.mesh");
2144      sortie<<"Vertices"<<" \n";      sortie<<"Vertices"<<" \n";
2145      sortie<<vertexplot.size()<<"\n";      sortie<<vertexplot.size()<<"\n";
2146      for(size_t i=0;i<vertexplot.size();++i)      for(size_t i=0;i<vertexplot.size();++i)
2147      {        {
2148          sortie<<(*(&vertexplot[i])).Vector3()[0]<<" "<<(*(&vertexplot[i])).Vector3()[1]<<" "<<(*(&vertexplot[i])).Vector3()[2]<<" 0"<<" \n";          sortie<<(*(&vertexplot[i])).Vector3()[0]<<" "<<(*(&vertexplot[i])).Vector3()[1]<<" "<<(*(&vertexplot[i])).Vector3()[2]<<" 0"<<" \n";
2149      }        }
2150    
2151      sortie<<" \n";      sortie<<" \n";
2152      sortie<<"Triangles"<<"\n";      sortie<<"Triangles"<<"\n";
2153      sortie<<vertexplot.size()/3<<"\n";      sortie<<vertexplot.size()/3<<"\n";
2154    
2155      for (size_t i=0; i<vertexplot.size()/3; ++i)      for (size_t i=0; i<vertexplot.size()/3; ++i)
2156      {        {
2157          sortie<<3*i+1<<" "<<3*i+2<<" "<<3*i+3<<" "<<1<<"\n";          sortie<<3*i+1<<" "<<3*i+2<<" "<<3*i+3<<" "<<1<<"\n";
2158      }        }
2159      sortie<<" \n";      sortie<<" \n";
2160      sortie<<"End"<<"\n";      sortie<<"End"<<"\n";
2161  }*/    }*/
   
   
   
2162  }  }
2163    
2164    
# Line 1566  SurfaceMeshGenerator::~SurfaceMeshGenera Line 2175  SurfaceMeshGenerator::~SurfaceMeshGenera
2175  }  }
2176    
2177    
 void inout(size_t iter,const size_t nn,TinyVector<3,int> ttn,size_t& re0,size_t& re1,size_t& re2,size_t& ii0,size_t& ii1,size_t& ii2)  
 {  
   
        if(re0==1 and re1==1)  
        {  
            if(ttn[0]!=0 and iter>1)  
 std::cout<<"-- cas indetermine\n";  
        }  
        else if((re0+re1==1))  
        {  
            if(ttn[0]==0)  
            {  
 //std::cout<<"-- 1 0\n";  
                re0=0;  
                re1=0;  
            }  
            else  
            {  
 //std::cout<<"-- 1 in 0 out on sait traite\n";  
                ii0=re0;  
                ii1=re1;  
            }  
        }  
        else if((re0==2 and re1==0) or  
                (re0==0 and re1==2))  
        {  
            if(ttn[0]==0)  
            {  
                if(ttn[0]+ttn[1]+ttn[2]==0)  
                {  
                    re0=2;  
                    re1=2;  
                }  
                else  
                {  
 std::cout<<"-- 0 2 pas coupe\n";  
 std::cout<<ttn[0]<<" "<<ttn[1]<<" "<<ttn[2]<<"\n";  
 std::cout<<re0<<" "<<re1<<"\n";  
                }  
          /*      if(re0==2 and (ttn[1]==1 or ttn[2]==1))  
                {  
                    re0=-nn;  
                    ii0=-nn;  
                }  
                else  
                {  
                    re1=-nn;  
                    ii1=-nn;  
                }*/  
              //  re0=0;  
                //re1=0;  
            }  
            else  
            {  
 //std::cout<<"-- 0 2 on sait traite\n";  
                ii0=re0;  
                ii1=re1;  
            }  
        }  
        else if((re0+re1==3))  
        {  
            if(ttn[0]==0)  
            {  
 //std::cout<<"-- 1 2 pas coupe\n";  
                re0=2;  
                re1=2;  
            }  
            else  
            {  
 //std::cout<<"-- 1 2 on sait traite\n";  
                if(re0==1)  
                {  
                    ii0=0;  
                    ii1=1;  
                    re0=0;  
                }  
                else  
                {  
                    ii0=1;  
                    ii1=0;  
                    re1=0;  
                }  
            }  
        }  
   //     else  
 //std::cout<<"--0 out 0 out rien a faire\n";  
   
        if(re0==1 and re2==1)  
        {  
            if(ttn[1]!=0 and iter>1)  
 std::cout<<"--* cas indetermine\n";  
        }  
        else if((re0==1 and re2==0) or  
                (re0==0 and re2==1))  
        {  
            if(ttn[1]==0)  
            {  
 //std::cout<<"--* 1 O non inter\n";  
                re0=0;  
                re2=0;  
            }  
            else  
            {  
 //std::cout<<"--* 1 in 0 out on sait traite\n";  
                ii0=re0;  
                ii2=re2;  
            }  
        }  
        else if((re0==2 and re2==0) or  
                (re0==0 and re2==2))  
        {  
            if(ttn[1]==0)  
            {  
                if(ttn[0]+ttn[1]+ttn[2]==0)  
                {  
                    re0=2;  
                    re2=2;  
                }  
                else  
                {  
                     
 std::cout<<"--* 0 2 pas coupe\n";  
 std::cout<<ttn[0]<<" "<<ttn[1]<<" "<<ttn[2]<<"\n";  
 std::cout<<re0<<" "<<re2<<"\n";  
                }  
             /*   if(re0==2 and (ttn[0]==1 or ttn[2]==1))  
                {  
                    re0=-nn;  
                    ii0=-nn;  
                }  
                else  
                {  
                    re2=-nn;  
                    ii2=-nn;  
                }*/  
   
          //      re0=0;  
            //    re2=0;  
            }  
            else  
            {  
 //std::cout<<"--* 0 2 on sait traite\n";  
                ii0=re0;  
                ii2=re2;  
            }  
        }  
        else if((re0+re2==3))  
        {  
            if(ttn[1]==0)  
            {  
 //std::cout<<"--* 1 2 pas coupe\n";  
                re0=2;  
                re2=2;  
            }  
            else  
            {  
 //std::cout<<"--* 1 2 on sait traite\n";  
                if(re0==1)  
                {  
                    ii0=0;  
                    ii2=1;  
                    re0=0;  
                }  
                else  
                {  
                    ii0=1;  
                    ii2=0;  
                    re2=0;  
                }  
            }  
        }        
 //       else  
 //std::cout<<"--* 0 out 0 out rien a faire\n";  
   
        if(re2==1 and re1==1)  
        {  
            if(ttn[2]!=0 and iter>1)  
 std::cout<<"--+ cas indetermine\n";  
        }  
        else if((re2+re1==1))  
        {  
            if(ttn[2]==0)  
            {  
 //std::cout<<"--+ 1 0 non inter\n";  
                re2=0;  
                re1=0;  
            }  
            else  
            {  
 //std::cout<<"--+ 1 in 0 out on sait traite\n";  
                ii1=re1;  
                ii2=re2;  
            }  
        }  
        else if((re2==2 and re1==0) or  
                (re2==0 and re1==2))  
        {  
            if(ttn[2]==0)  
            {  
                if(ttn[0]+ttn[1]+ttn[2]==0)  
                {  
                    re1=2;  
                    re2=2;  
                }  
                else  
                {  
                     
 std::cout<<"--+ 0 2 pas coupe\n";  
 std::cout<<ttn[0]<<" "<<ttn[1]<<" "<<ttn[2]<<"\n";  
 std::cout<<re1<<" "<<re2<<"\n";  
                }  
             /*   if(re1==2 and (ttn[1]==1 or ttn[0]==1))  
                {  
                    re1=-nn;  
                    ii1=-nn;  
                }  
                else  
                {  
                    re2=-nn;  
                    ii2=-nn;  
                }*/  
                //re1=0;  
                //re2=0;  
            }  
            else  
            {  
 //std::cout<<"--+ 0 2 on sait traite\n";  
                ii1=re1;  
                ii2=re2;  
            }  
        }  
        else if((re2+re1==3))  
        {  
            if(ttn[2]==0)  
            {  
 //std::cout<<"--+ 1 2 pas coupe\n";  
                re2=2;  
                re1=2;  
            }  
            else  
            {  
 //std::cout<<"--+ 1 2 on sait traite\n";  
                if(re1==1)  
                {  
                    ii1=0;  
                    ii2=1;  
                    re1=0;  
                }  
                else  
                {  
                    ii1=1;  
                    ii2=0;  
                    re2=0;  
                }  
            }  
        }  
 //       else  
 //std::cout<<"--+ 0 out 0 out rien a faire\n";  
   
     if(iter>1 and (((re2==1 and re1==1 and re0==1) or re0+re1+re2==6)  
                    and ttn[0]+ttn[1]+ttn[2]==0))  
     {  
         re0=2;  
         re1=2;  
         re2=2;  
         ii0=1;  
         ii1=1;  
         ii2=1;  
     }  
     if(iter>2 and re0+re1+re2==0 and ttn[0]+ttn[1]+ttn[2]!=0)  
     {  
 std::cout<<"pbs de precision?? tr out cut\n";  
         ii0=0;  
         ii1=0;  
         ii2=0;  
     }  
     if(iter==1)  
     {  
         if(ttn[0]+ttn[1]+ttn[2]==0 and ((re0==0 and re1==1 and re2==1) or (re0==1 and re1==0 and re2==1)  
                                  or (re0==1 and re1==1 and re2==0) or (re0==0 and re1==0 and re2==1)  
                                  or (re0==0 and re1==1 and re2==0) or (re0==1 and re1==0 and re2==0)))  
         {  
             re0=0;  
             re1=0;  
             re2=0;  
         }  
     }  
   
     if(ii0==2)  
         ii0=1;  
     if(ii1==2)  
         ii1=1;  
     if(ii2==2)  
         ii2=1;  
2178    
 }  
2179    
2180  void Bary(const Vertex& P,TinyVector<3>& P0,TinyVector<3>& V12,TinyVector<3>& V02,size_t& compt)  void plotmedit(const SurfaceMeshOfTriangles& s_mesh,std::set<const Cell*>& toTreatHexahedra,size_t nobj)
2181  {  {
     //std::cout<<" tr vert "<<V02[0]<<" "<<V02[1]<<" "<<V02[2]<<"\n";  
     //std::cout<<" tr vert "<<V12[0]<<" "<<V12[1]<<" "<<V12[2]<<"\n";  
     //std::cout<<" P0 "<<P0[0]<<" "<<P0[1]<<" "<<P0[2]<<"\n";  
     //std::cout<<" P "<<P[0]<<" "<<P[1]<<" "<<P[2]<<"\n";  
   
2182    
2183      // **********************************************************************
2184        
2185    
2186      compt=0;    std::ofstream sortie("intersection.mesh");
2187      real_t det,eps=-1e-10;        
2188      det=V02[0]*V12[1]-V02[1]*V12[0];    sortie<<"MeshVersionFormatted 1"<<" \n";
2189      if(det<eps or det>-eps)    sortie<<" \n";
2190      sortie<<"Dimension"<<" \n";
2191      sortie<<"3"<<" \n";
2192      sortie<<" \n";
2193      sortie<<"Vertices"<<" \n";
2194      sortie<<8*toTreatHexahedra.size()+s_mesh.numberOfVertices()<<" \n";
2195      for (std::set<const Cell*>::iterator i = toTreatHexahedra.begin();
2196           i != toTreatHexahedra.end(); ++i)
2197      {      {
2198          //std::cout<<"det "<<det<<"\n";        const CartesianHexahedron& H = static_cast<const CartesianHexahedron&>(*(*i));
2199          TinyMatrix<2,2> mat;        for (size_t j=0; j<8; ++j)
         mat(0,0)=V02[0];  
         mat(1,0)=V02[1];  
         mat(0,1)=V12[0];  
         mat(1,1)=V12[1];  
         TinyVector<2> bary,secondmb;  
         secondmb[0]=P0[0];  
         secondmb[1]=P0[1];  
         bary=mat.invert()*secondmb;  
         //std::cout<<mat.invert()(0,0)<<" "<<mat.invert()(1,0)<<" "<<mat.invert()(0,1)<<" "<<mat.invert()(1,1)<<"\n";  
  //std::cout<<secondmb[0]<<" "<<secondmb[1]<<"\n";  
  //std::cout<<bary[0]<<" "<<bary[1]<<" "<<1.-bary[0]-bary[1]<<"\n";  
  //std::cout<<bary[0]>=eps<<" "<<bary[1]>=eps<<" "<<(1.-bary[0]-bary[1])>=eps<<"\n";  
         real_t temp=(1.-bary[0]-bary[1]);  
         if(bary[0]>=eps && bary[1]>=eps && temp>=eps)  
2200          {          {
2201              compt+=1;            sortie<<H(j)[0]<<" "<<H(j)[1]<<" "<<H(j)[2]<<" 0"<<" \n";
           //  std::cout<<bary[0]<<" "<<bary[1]<<" "<<1.-bary[0]-bary[1]<<" "<<compt<<"\n";  
      //std::cout<<verif[0]<<" "<<verif[1]<<" \n";  
      //std::cout<<P0[0]<<" "<<P0[2]<<" \n";  
2202          }          }
2203      }      }
2204      else    for(size_t i=0;i<s_mesh.numberOfVertices();++i)
         compt+=1;  
2205      {      {
2206          det=V02[0]*V12[2]-V02[2]*V12[0];        sortie<<s_mesh.vertex(i)[0]<<" "<<s_mesh.vertex(i)[1]<<" "<<s_mesh.vertex(i)[2]<<" 0"<<" \n";
         if(det<eps or det>-eps)  
         {  
             TinyMatrix<2,2> mat;  
             mat(0,0)=V02[0];  
             mat(1,0)=V02[2];  
             mat(1,1)=V12[2];  
             mat(0,1)=V12[0];  
             TinyVector<2> bary,secondmb;  
             secondmb[0]=P0[0];  
             secondmb[1]=P0[2];  
             bary=mat.invert()*secondmb;  
             //std::cout<<"------------------\n";  
      //std::cout<<mat.invert()(0,0)<<" "<<mat.invert()(1,0)<<" "<<mat.invert()(0,1)<<" "<<mat.invert()(1,1)<<"\n";  
      //std::cout<<"2 "<<bary[0]<<" "<<bary[1]<<" "<<1.-bary[0]-bary[1]<<"\n";  
   
             TinyVector<2> verif=mat*bary;  
             if(bary[0]>=eps && bary[1]>=eps && (1-bary[0]-bary[1])>=eps)  
             {  
                 compt+=1;  
                 //std::cout<<bary[0]<<" "<<bary[1]<<" "<<1.-bary[0]-bary[1]<<"\n";  
   //std::cout<<verif[0]<<" "<<verif[1]<<" \n";  
   //std::cout<<P0[0]<<" "<<P0[2]<<" \n";  
   
   
             }  
         }  
         else  
             compt+=1;  
2207      }      }
2208      {    //       for (size_t i=0; i<s_mesh.numberOfCells(); ++i)
2209          det=V02[1]*V12[2]-V02[2]*V12[1];    //       {
2210          if(det<eps or det>-eps)    //       const Triangle& T = s_mesh.cell(i);
2211          {    //       sortie<<T(0)[0]<<" "<<T(0)[1]<<" "<<T(0)[2]<<" 0"<<" \n";
2212              TinyMatrix<2,2> mat;    //       sortie<<T(1)[0]<<" "<<T(1)[1]<<" "<<T(1)[2]<<" 0"<<" \n";
2213              mat(0,0)=V02[1];    //       sortie<<T(2)[0]<<" "<<T(2)[1]<<" "<<T(2)[2]<<" 0"<<" \n";
2214              mat(1,0)=V02[2];    //
2215              mat(1,1)=V12[2];    //       }
2216              mat(0,1)=V12[1];    sortie<<" \n";
2217              TinyVector<2> bary,secondmb;    sortie<<"Triangles"<<"\n";
2218              secondmb[0]=P0[1];    sortie<<s_mesh.numberOfCells()<<"\n";
             secondmb[1]=P0[2];  
             bary=mat.invert()*secondmb;  
             //std::cout<<"------------------\n";  
   //std::cou    t<<mat.invert()(0,0)<<" "<<mat.invert()(1,0)<<" "<<mat.invert()(0,1)<<" "<<mat.invert()(1,1)<<"\n";  
   //std::cout<<"3 "<<bary[0]<<" "<<bary[1]<<" "<<1.-bary[0]-bary[1]<<"\n";  
2219    
2220              TinyVector<2> verif=mat*bary;    for (size_t i=0; i<s_mesh.numberOfCells(); ++i)
2221        {
2222          const Triangle& T = s_mesh.cell(i);
2223          sortie<<s_mesh.vertexNumber(T(0))+8*toTreatHexahedra.size()+1<<" "<<s_mesh.vertexNumber(T(1))+8*toTreatHexahedra.size()+1<<" "<<s_mesh.vertexNumber(T(2))+8*toTreatHexahedra.size()+1;
2224          //          sortie<<3*i+1+8*toTreatHexahedra.size()<<" "<<3*i+2+8*toTreatHexahedra.size()<<" //"<<3*i+3+8*toTreatHexahedra.size()<<" "<<1<<"\n";
2225          if(i<(nobj))
2226            sortie<<" 1"<<"\n";
2227          else
2228            sortie<<" 2"<<"\n";
2229        }
2230      sortie<<" \n";
2231      sortie<<"Hexahedra"<<"\n";
2232      sortie<<toTreatHexahedra.size()<<" \n";
2233    
2234              if(bary[0]>=eps && bary[1]>=eps && (1-bary[0]-bary[1])>=eps)    /*  for (std::set<const Cell*>::iterator i = toTreatHexahedra.begin();
2235              {        i != toTreatHexahedra.end(); ++i)
2236                  compt+=1;        {
2237                  //std::cout<<bary[0]<<" "<<bary[1]<<" "<<1.-bary[0]-bary[1]<<" "<<compt<<"\n";        const CartesianHexahedron& H = static_cast<const CartesianHexahedron&>(*(*i));
2238        //std::cout<<verif[0]<<" "<<verif[1]<<" \n";        for(size_t j=0; j<8; ++j)
2239        //std::cout<<P0[0]<<" "<<P0[1]<<" \n";        sortie<<SMesh.vertexNumber(H(j))<<" \n";
2240          }*/
2241      for(size_t i=0;i<toTreatHexahedra.size();i++)
2242        sortie<<8*i+1<<" "<<8*i+2<<" "<<8*i+3<<" "<<8*i+4<<" "<<8*i+5<<" "<<8*i+6<<" "<<8*i+7<<" "<<8*i+8<<" 3"<<"\n";
2243      sortie<<"End"<<" \n";
2244    
             }  
         }  
         else  
             compt+=1;  
     }  
       
     //    if(compt==3)  
     //  {  
     //std::cout<<P.Vector3()[0]<<" "<<P.Vector3()[1]<<" "<<P.Vector3()[2]<<" 3 point intersection \n";  
     //}  
2245  }  }
2246    
2247    // #include <vtkPolyDataMapper.h>
2248    // #include <vtkRenderWindow.h>
2249  void plotmedit(const SurfaceMeshOfTriangles& s_mesh,std::set<const Cell*>& toTreatHexahedra,size_t nobj)  // #include <vtkRenderWindowInteractor.h>
2250  {  // #include <vtkTriangle.h>
2251    // #include <vtkHexahedron.h>
2252      // **********************************************************************  // #include <vtkTetra.h>
2253        // #include <vtkUnstructuredGrid.h>
2254    // #include <vtkDataSetMapper.h>
2255        std::ofstream sortie("intersection.mesh");  // #include <vtkDataSetToPolyDataFilter.h>
2256          // #include <vtkDoubleArray.h>
2257        sortie<<"MeshVersionFormatted 1"<<" \n";  
2258        sortie<<" \n";  // void plot2(const SurfaceMeshOfTriangles& s,
2259        sortie<<"Dimension"<<" \n";  //         std::set<const Cell*>& hlist,
2260        sortie<<"3"<<" \n";  //         Structured3DMesh& m)
2261        sortie<<" \n";  // {
2262        sortie<<"Vertices"<<" \n";  //   vtkRenderer *ren1= vtkRenderer::New();
2263        sortie<<8*toTreatHexahedra.size()+s_mesh.numberOfVertices()<<" \n";  
2264      for (std::set<const Cell*>::iterator i = toTreatHexahedra.begin();  //   vtkActor *anActor1 = vtkActor::New();
2265      i != toTreatHexahedra.end(); ++i)  //   vtkActor *anActor2 = vtkActor::New();
2266      {  
2267          const CartesianHexahedron& H = static_cast<const CartesianHexahedron&>(*(*i));  //   {
2268          for (size_t j=0; j<8; ++j)  //     vtkPoints* vertices = vtkPoints::New();
2269          {  
2270              sortie<<H(j)[0]<<" "<<H(j)[1]<<" "<<H(j)[2]<<" 0"<<" \n";  //     vertices->SetNumberOfPoints(m.numberOfVertices());
2271          }  
2272      }  //     for (size_t i=0; i<m.numberOfVertices(); ++i) {
2273      for(size_t i=0;i<s_mesh.numberOfVertices();++i)  //       const Vertex& X = m.vertex(i);
2274      {  //       vertices->InsertPoint(i,X[0],X[1],X[2]);
2275          sortie<<s_mesh.vertex(i)[0]<<" "<<s_mesh.vertex(i)[1]<<" "<<s_mesh.vertex(i)[2]<<" 0"<<" \n";  //     }
2276      }  
2277  //       for (size_t i=0; i<s_mesh.numberOfCells(); ++i)  //     vtkHexahedron* hexahedron = vtkHexahedron::New();
2278  //       {  
2279  //         const Triangle& T = s_mesh.cell(i);  //     vtkUnstructuredGrid* anHexahedraGrid = vtkUnstructuredGrid::New();
2280  //         sortie<<T(0)[0]<<" "<<T(0)[1]<<" "<<T(0)[2]<<" 0"<<" \n";  //     anHexahedraGrid->Allocate(m.numberOfCells(),1);
2281  //         sortie<<T(1)[0]<<" "<<T(1)[1]<<" "<<T(1)[2]<<" 0"<<" \n";  
2282  //         sortie<<T(2)[0]<<" "<<T(2)[1]<<" "<<T(2)[2]<<" 0"<<" \n";  //     for (std::set<const Cell*>::iterator i = hlist.begin();
2283  //  //       i != hlist.end(); ++i) {
2284  //       }  //       const CartesianHexahedron& H = static_cast<const CartesianHexahedron&>(*(*i));
2285        sortie<<" \n";  //       for (size_t j=0; j<8; ++j)
2286        sortie<<"Triangles"<<"\n";  //      hexahedron->GetPointIds()->SetId(j,m.vertexNumber(H(j)));
2287        sortie<<s_mesh.numberOfCells()<<"\n";  
2288    //       anHexahedraGrid->InsertNextCell(hexahedron->GetCellType(),
2289        for (size_t i=0; i<s_mesh.numberOfCells(); ++i)  //                                    hexahedron->GetPointIds());
2290        {  //     }
2291            const Triangle& T = s_mesh.cell(i);  
2292            sortie<<s_mesh.vertexNumber(T(0))+8*toTreatHexahedra.size()+1<<" "<<s_mesh.vertexNumber(T(1))+8*toTreatHexahedra.size()+1<<" "<<s_mesh.vertexNumber(T(2))+8*toTreatHexahedra.size()+1;  //     anHexahedraGrid->SetPoints(vertices);
2293  //        sortie<<3*i+1+8*toTreatHexahedra.size()<<" "<<3*i+2+8*toTreatHexahedra.size()<<" //"<<3*i+3+8*toTreatHexahedra.size()<<" "<<1<<"\n";  
2294            if(i<(nobj))  //     vtkDoubleArray* values = vtkDoubleArray::New();
2295                sortie<<" 1"<<"\n";  //     values->Resize(m.numberOfVertices());
2296            else  //     for (size_t i=0; i<m.numberOfVertices(); ++i) {
2297                sortie<<" 2"<<"\n";  //       values->InsertNextValue(1);
2298        }  //     }
2299        sortie<<" \n";  
2300        sortie<<"Hexahedra"<<"\n";  //     anHexahedraGrid->GetPointData()->SetScalars(values);
2301      sortie<<toTreatHexahedra.size()<<" \n";  
2302    //     // We create an instance of vtkPolyDataMapper to map the polygonal data
2303      /*  for (std::set<const Cell*>::iterator i = toTreatHexahedra.begin();  //     // into graphics primitives. We connect the output of the cone souece
2304      i != toTreatHexahedra.end(); ++i)  //     // to the input of this mapper
2305  {  
2306          const CartesianHexahedron& H = static_cast<const CartesianHexahedron&>(*(*i));  //     vtkDataSetMapper* mapper = vtkDataSetMapper::New();
2307          for(size_t j=0; j<8; ++j)  //     mapper->SetInput(anHexahedraGrid);
2308              sortie<<SMesh.vertexNumber(H(j))<<" \n";  
2309  }*/  //     // create an actor to represent the cone. The actor coordinates rendering of
2310  for(size_t i=0;i<toTreatHexahedra.size();i++)  //     // the graphics primitives for a mapper. We set this actor's mapper to be
2311  sortie<<8*i+1<<" "<<8*i+2<<" "<<8*i+3<<" "<<8*i+4<<" "<<8*i+5<<" "<<8*i+6<<" "<<8*i+7<<" "<<8*i+8<<" 3"<<"\n";  //     // coneMapper which we created above.
2312      sortie<<"End"<<" \n";  //     //
2313    
2314    }  //     vtkProperty* prop = vtkProperty::New();
2315    //     prop->SetRepresentationToWireframe();
2316  /* #include <vtkPolyDataMapper.h>  
2317   #include <vtkRenderWindow.h>  //     anActor1->SetMapper( mapper );
2318   #include <vtkRenderWindowInteractor.h>  //     anActor1->SetProperty(prop);
2319   #include <vtkTriangle.h>  
2320   #include <vtkHexahedron.h>  //   }
2321   #include <vtkTetra.h>  
2322   #include <vtkUnstructuredGrid.h>  
2323   #include <vtkDataSetMapper.h>  //   {
2324   #include <vtkDataSetToPolyDataFilter.h>  //     vtkPoints* vertices = vtkPoints::New();
2325   #include <vtkDoubleArray.h>  
2326    //     vertices->SetNumberOfPoints(s.numberOfVertices());
2327    
2328    //     for (size_t i=0; i<s.numberOfVertices(); ++i) {
2329   void plot2(const SurfaceMeshOfTriangles& s,  //       const Vertex& X = s.vertex(i);
2330             std::set<const Cell*>& hlist,  //       vertices->InsertPoint(i,X[0],X[1],X[2]);
2331             Structured3DMesh& m)  //     }
2332   {  
2333     vtkRenderer *ren1= vtkRenderer::New();  //     vtkTriangle* triangle = vtkTriangle::New();
2334    
2335     vtkActor *anActor1 = vtkActor::New();  //     vtkUnstructuredGrid* aTriangleGrid = vtkUnstructuredGrid::New();
2336     vtkActor *anActor2 = vtkActor::New();  //     aTriangleGrid->Allocate(s.numberOfCells(),1);
2337    
2338     {  //     for (size_t i=0; i<s.numberOfCells(); ++i) {
2339     vtkPoints* vertices = vtkPoints::New();  //       const Triangle& T = s.cell(i);
2340    //       triangle->GetPointIds()->SetId(0,s.vertexNumber(T(0)));
2341     vertices->SetNumberOfPoints(m.numberOfVertices());  //       triangle->GetPointIds()->SetId(1,s.vertexNumber(T(1)));
2342    //       triangle->GetPointIds()->SetId(2,s.vertexNumber(T(2)));
2343     for (size_t i=0; i<m.numberOfVertices(); ++i) {  
2344       const Vertex& X = m.vertex(i);  //       aTriangleGrid->InsertNextCell(triangle->GetCellType(),
2345       vertices->InsertPoint(i,X[0],X[1],X[2]);  //                                  triangle->GetPointIds());
2346     }  //     }
2347    
2348     vtkHexahedron* hexahedron = vtkHexahedron::New();  //     aTriangleGrid->SetPoints(vertices);
2349    
2350     vtkUnstructuredGrid* anHexahedraGrid = vtkUnstructuredGrid::New();  //     vtkDoubleArray* values = vtkDoubleArray::New();
2351     anHexahedraGrid->Allocate(m.numberOfCells(),1);  //     values->Resize(s.numberOfVertices());
2352    //     for (size_t i=0; i<s.numberOfVertices(); ++i) {
2353     for (std::set<const Cell*>::iterator i = hlist.begin();  //       values->InsertNextValue(0);
2354          i != hlist.end(); ++i) {  //     }
2355       const CartesianHexahedron& H = static_cast<const CartesianHexahedron&>(*(*i));  
2356       for (size_t j=0; j<8; ++j)  //     aTriangleGrid->GetPointData()->SetScalars(values);
2357         hexahedron->GetPointIds()->SetId(j,m.vertexNumber(H(j)));  
2358    //     // We create an instance of vtkPolyDataMapper to map the polygonal data
2359       anHexahedraGrid->InsertNextCell(hexahedron->GetCellType(),  //     // into graphics primitives. We connect the output of the cone souece
2360                                      hexahedron->GetPointIds());  //     // to the input of this mapper
2361     }  
2362    //     vtkDataSetMapper* mapper = vtkDataSetMapper::New();
2363     anHexahedraGrid->SetPoints(vertices);  //     mapper->SetInput(aTriangleGrid);
2364    
2365     vtkDoubleArray* values = vtkDoubleArray::New();  
2366     values->Resize(m.numberOfVertices());  //     // create an actor to represent the cone. The actor coordinates rendering of
2367     for (size_t i=0; i<m.numberOfVertices(); ++i) {  //     // the graphics primitives for a mapper. We set this actor's mapper to be
2368       values->InsertNextValue(1);  //     // coneMapper which we created above.
2369     }  //     //
2370    //     anActor2->SetMapper( mapper );
2371     anHexahedraGrid->GetPointData()->SetScalars(values);  
2372    //   }
2373     // We create an instance of vtkPolyDataMapper to map the polygonal data  
2374     // into graphics primitives. We connect the output of the cone souece  //   // Create the Renderer and assign actors to it. A renderer is like a
2375     // to the input of this mapper  //   // viewport. It is part or all of a window on the screen and it is
2376    //   // responsible for drawing the actors it has.  We also set the background
2377     vtkDataSetMapper* mapper = vtkDataSetMapper::New();  //   // color here
2378     mapper->SetInput(anHexahedraGrid);  //   //
2379    //   //  ren1->AddActor( anActor1 );
2380     // create an actor to represent the cone. The actor coordinates rendering of  //   ren1->AddActor( anActor2 );
2381     // the graphics primitives for a mapper. We set this actor's mapper to be  //   ren1->SetBackground( 0.1, 0.2, 0.4 );
2382     // coneMapper which we created above.  
2383     //  //   //
2384    //   // Finally we create the render window which will show up on the screen
2385     vtkProperty* prop = vtkProperty::New();  //   // We put our renderer into the render window using AddRenderer. We also
2386     prop->SetRepresentationToWireframe();  //   // set the size to be 300 pixels by 300
2387    //   //
2388     anActor1->SetMapper( mapper );  //   vtkRenderWindow *renWin = vtkRenderWindow::New();
2389     anActor1->SetProperty(prop);  //   renWin->AddRenderer( ren1 );
2390    //   renWin->SetSize( 300, 300 );
2391     }  
2392    //   vtkRenderWindowInteractor* iren = vtkRenderWindowInteractor::New();
2393    //   iren->SetRenderWindow(renWin);
2394     {  
2395    vtkPoints* vertices = vtkPoints::New();  //   // now we loop over 360 degreeees and render the cone each time
2396    //   //
2397     vertices->SetNumberOfPoints(s.numberOfVertices());  //   // render the image
2398    //   renWin->Render();
2399     for (size_t i=0; i<s.numberOfVertices(); ++i) {  //   iren->Start();
2400       const Vertex& X = s.vertex(i);  //   //
2401       vertices->InsertPoint(i,X[0],X[1],X[2]);  //   // Free up any objects we created
2402     }  //   //
2403    //   //   cone->Delete();
2404     vtkTriangle* triangle = vtkTriangle::New();  //   //   coneMapper->Delete();
2405    //   //   coneActor->Delete();
2406     vtkUnstructuredGrid* aTriangleGrid = vtkUnstructuredGrid::New();  //   //   ren1->Delete();
2407     aTriangleGrid->Allocate(s.numberOfCells(),1);  // }
   
    for (size_t i=0; i<s.numberOfCells(); ++i) {  
      const Triangle& T = s.cell(i);  
      triangle->GetPointIds()->SetId(0,s.vertexNumber(T(0)));  
      triangle->GetPointIds()->SetId(1,s.vertexNumber(T(1)));  
      triangle->GetPointIds()->SetId(2,s.vertexNumber(T(2)));  
   
      aTriangleGrid->InsertNextCell(triangle->GetCellType(),  
                                   triangle->GetPointIds());  
    }  
   
    aTriangleGrid->SetPoints(vertices);  
   
    vtkDoubleArray* values = vtkDoubleArray::New();  
    values->Resize(s.numberOfVertices());  
    for (size_t i=0; i<s.numberOfVertices(); ++i) {  
      values->InsertNextValue(0);  
    }  
   
    aTriangleGrid->GetPointData()->SetScalars(values);  
   
    // We create an instance of vtkPolyDataMapper to map the polygonal data  
    // into graphics primitives. We connect the output of the cone souece  
    // to the input of this mapper  
   
    vtkDataSetMapper* mapper = vtkDataSetMapper::New();  
    mapper->SetInput(aTriangleGrid);  
   
   
    // create an actor to represent the cone. The actor coordinates rendering of  
    // the graphics primitives for a mapper. We set this actor's mapper to be  
    // coneMapper which we created above.  
    //  
    anActor2->SetMapper( mapper );  
   
    }  
   
    // Create the Renderer and assign actors to it. A renderer is like a  
    // viewport. It is part or all of a window on the screen and it is  
    // responsible for drawing the actors it has.  We also set the background  
    // color here  
    //  
    //  ren1->AddActor( anActor1 );  
    ren1->AddActor( anActor2 );  
    ren1->SetBackground( 0.1, 0.2, 0.4 );  
   
    //  
    // Finally we create the render window which will show up on the screen  
    // We put our renderer into the render window using AddRenderer. We also  
    // set the size to be 300 pixels by 300  
    //  
    vtkRenderWindow *renWin = vtkRenderWindow::New();  
    renWin->AddRenderer( ren1 );  
    renWin->SetSize( 300, 300 );  
   
    vtkRenderWindowInteractor* iren = vtkRenderWindowInteractor::New();  
    iren->SetRenderWindow(renWin);  
   
    // now we loop over 360 degreeees and render the cone each time  
    //  
    // render the image  
    renWin->Render();  
    iren->Start();  
    //  
    // Free up any objects we created  
    //  
  //   cone->Delete();  
  //   coneMapper->Delete();  
  //   coneActor->Delete();  
  //   ren1->Delete();  
  }*/  

Legend:
Removed from v.1.8  
changed lines
  Added in v.1.9

savannah-hackers-public@gnu.org
ViewVC Help
Powered by ViewVC 1.1.26