Changeset 2797


Ignore:
Timestamp:
01/12/10 11:17:20 (15 years ago)
Author:
Mathieu Morlighem
Message:

removed useless code

Location:
issm/trunk/src/c/Bamgx
Files:
4 edited

Legend:

Unmodified
Added
Removed
  • issm/trunk/src/c/Bamgx/Mesh2.h

    r2796 r2797  
    870870  GeometricalEdge * ProjectOnCurve( Edge & AB, Vertex &  A, Vertex & B,Real8 theta,
    871871                      Vertex & R,VertexOnEdge & BR,VertexOnGeom & GR);
    872    
    873  
    874   void WriteElements(std::ostream& f,Int4 * reft ,Int4 nbInT) const;
    875 
    876872 
    877873  Int4 Number(const Triangle & t) const  { return &t - triangles;}
     
    907903  int UnCrack();
    908904 
    909 #ifdef DEBUG
    910   void inline Check();
    911 #endif
    912 #ifdef DRAWING
    913   void  Draw() const ;
    914   void  InitDraw() const ;
    915   void   inquire()  ;
    916 #endif
    917905 friend std::ostream& operator <<(std::ostream& f,  const  Triangles & Th);
    918   void  Write(const char * filename);
    919906  void ConsGeometry(Real8 =-1.0,int *equiedges=0); // construct a geometry if no geo
    920907  void FillHoleInMesh() ;
     
    924911  void GeomToTriangles0(Int4 nbvx);// the  real constructor mesh generator
    925912  void PreInit(Int4,char * =0 );
    926   //
    927   void Write_nop5(OFortranUnFormattedFile * f,
    928                              Int4 &lnop5,Int4 &nef,Int4 &lgpdn,Int4 ndsr) const ;
    929 
    930913 
    931914}; // End Class Triangles
     
    970953  Real8 MinimalHmin() {return 2.0/coefIcoor;}
    971954  Real8 MaximalHmax() {return Max(pmax.x-pmin.x,pmax.y-pmin.y);}
    972   void ReadGeometry(const char * ) ;
    973955  void ReadGeometry(BamgGeom* bamggeom, BamgOpts* bamgopts);
    974   void ReadGeometry(MeshIstream & ,const char *)  ;
    975 
    976956  void EmptyGeometry();
    977957  Geometry() {EmptyGeometry();}// empty Geometry
    978958  void AfterRead();
    979959  Geometry(BamgGeom* bamggeom, BamgOpts* bamgopts) {EmptyGeometry();OnDisk=1;ReadGeometry(bamggeom,bamgopts);AfterRead();}
    980   Geometry(const char * filename) {EmptyGeometry();OnDisk=1;ReadGeometry(filename);AfterRead();}
    981 
    982   void ReadMetric(const char *,Real8 hmin,Real8 hmax,Real8 coef);
     960
    983961  const GeometricalVertex & operator[]  (Int4 i) const { return vertices[i];};
    984962  GeometricalVertex & operator[](Int4 i) { return vertices[i];};
     
    997975  GeometricalEdge *  Contening(const R2 P,  GeometricalEdge * start) const;
    998976  friend std::ostream& operator <<(std::ostream& f, const   Geometry & Gh);
    999   void Write(const char * filename);
    1000977  void WriteGeometry(BamgGeom* bamggeom, BamgOpts* bamgopts);
    1001  
    1002978}; // End Class Geometry
    1003979
  • issm/trunk/src/c/Bamgx/MeshRead.cpp

    r2795 r2797  
    529529        }
    530530}
    531 
    532 void Geometry::ReadGeometry(MeshIstream & f_in,const char * filename) 
    533 {
    534         long int verbosity=2;
    535   if(verbosity>1)
    536     cout << "  -- ReadGeometry " << filename << endl;
    537   assert(empty());
    538   nbiv=nbv=nbvx=0;
    539   nbe=nbt=nbtx=0;
    540   NbOfCurves=0;
    541  // BeginOfCurve=0;
    542   name=new char [strlen(filename)+1];
    543   strcpy(name,filename);
    544   Real8 Hmin = HUGE_VAL;// the infinie value
    545 //  Real8 MaximalAngleOfCorner = 10*Pi/180; ;
    546   Int4 hvertices =0;
    547   Int4 i;
    548   Int4 Version,dim=0;
    549   int field=0;
    550   int showfield=0;
    551   int NbErr=0;
    552 
    553   while (f_in.cm())
    554     {
    555       field=0;
    556       // warning ruse for on allocate fiedname at each time
    557       char fieldname[256];
    558       f_in.cm() >> fieldname ;
    559       f_in.err();
    560       if(f_in.eof()) break;
    561 //      cout <<  fieldname <<  " line " << LineNumber  <<endl;
    562       if (!strcmp(fieldname,"MeshVersionFormatted") )
    563         f_in  >> Version ;
    564       else if (!strcmp(fieldname,"End"))
    565         break;
    566       else if (!strcmp(fieldname,"end"))
    567         break;
    568       else if (!strcmp(fieldname,"Dimension"))
    569         {
    570           f_in   >>  dim ;
    571           if(verbosity>5)
    572             cout << "     Geom Record Dimension dim = " << dim << endl;       
    573           assert(dim ==2);
    574          }
    575        else if (!strcmp(fieldname,"hVertices"))
    576          {
    577            if (nbv <=0) {
    578              cerr<<"Error: the field Vertex is not found before hVertices " << filename<<endl;
    579              NbErr++;}       
    580            if(verbosity>5)
    581             cout << "     Geom Record hVertices nbv=" << nbv <<  endl;
    582            hvertices =1;
    583            for (i=0;i< nbv;i++)
    584              {
    585                Real4 h;
    586                f_in  >>  h ;
    587                vertices[i].m = Metric(h);
    588              }
    589          }
    590        else if (!strcmp(fieldname,"MetricVertices"))
    591          { hvertices =1;
    592            if (nbv <=0) {
    593              cerr<<"Error: the field Vertex is not found before MetricVertices " << filename<<endl;
    594              NbErr++;}       
    595            if(verbosity>5)
    596              cout << "     Geom Record MetricVertices nbv =" << nbv <<  endl;
    597            for (i=0;i< nbv;i++)
    598              {
    599                Real4 a11,a21,a22;
    600                f_in  >>  a11 >> a21 >> a22  ;
    601                vertices[i].m = Metric(a11,a21,a22);
    602              }
    603          }
    604        else if (!strcmp(fieldname,"h1h2VpVertices"))
    605          { hvertices =1;
    606            if (nbv <=0) {
    607              cerr<<"Error: the field Vertex is not found before h1h2VpVertices " << filename<<endl;
    608              NbErr++;}       
    609            if(verbosity>5)
    610              cout << "     Geom Record h1h2VpVertices nbv=" << nbv << endl;
    611 
    612            for (i=0;i< nbv;i++)
    613              {
    614                Real4 h1,h2,v1,v2;
    615                f_in  >> h1 >> h2 >>v1 >>v2 ;
    616                vertices[i].m = Metric(MatVVP2x2(1/(h1*h1),1/(h2*h2),D2(v1,v2)));
    617              }
    618          }
    619       else if (!strcmp(fieldname,"Vertices"))
    620         {
    621           assert(dim ==2);
    622           f_in   >>  nbv ;
    623           //          if(LineError) break;
    624           nbvx = nbv;
    625          
    626           vertices = new GeometricalVertex[nbvx];
    627           if(verbosity>5)
    628             cout << "     Geom Record Vertices nbv = " << nbv << "vertices = " << vertices<<endl;
    629           assert(nbvx >= nbv);
    630           nbiv = nbv;
    631           for (i=0;i<nbv;i++) {
    632             f_in  >> vertices[i].r.x  ;
    633             // if(LineError) break;
    634             f_in  >> vertices[i].r.y ;
    635             // if(LineError) break;
    636             f_in >>   vertices[i].ReferenceNumber   ;
    637             vertices[i].DirOfSearch=NoDirOfSearch;
    638             //            vertices[i].m.h = 0;
    639             vertices[i].color =0;
    640             vertices[i].Set();}
    641           // if(LineError) break;
    642             pmin =  vertices[0].r;
    643             pmax =  vertices[0].r;
    644             // recherche des extrema des vertices pmin,pmax
    645             for (i=0;i<nbv;i++) {
    646               pmin.x = Min(pmin.x,vertices[i].r.x);
    647               pmin.y = Min(pmin.y,vertices[i].r.y);
    648               pmax.x = Max(pmax.x,vertices[i].r.x);
    649               pmax.y = Max(pmax.y,vertices[i].r.y);
    650             }
    651            
    652               R2 DD05 = (pmax-pmin)*0.05;
    653               pmin -=  DD05;
    654               pmax +=  DD05;
    655            
    656             coefIcoor= (MaxICoor)/(Max(pmax.x-pmin.x,pmax.y-pmin.y));
    657             assert(coefIcoor >0);
    658             if (verbosity>5) {
    659               cout << "     Geom: min="<< pmin << "max ="<< pmax << " hmin = " << MinimalHmin() <<  endl;}
    660         }
    661       else if(!strcmp(fieldname,"MaximalAngleOfCorner")||!strcmp(fieldname,"AngleOfCornerBound"))
    662         {         
    663           f_in >> MaximalAngleOfCorner;
    664              
    665           if(verbosity>5)
    666             cout << "     Geom Record MaximalAngleOfCorner " << MaximalAngleOfCorner <<endl;
    667           MaximalAngleOfCorner *= Pi/180;
    668         }
    669       else if (!strcmp(fieldname,"Edges"))
    670         {
    671           if (nbv <=0) {
    672             cerr<<"Error: the field edges is not found before MetricVertices " << filename<<endl;
    673             NbErr++;}   
    674           else
    675             {
    676               int i1,i2;
    677               R2 zero2(0,0);
    678               f_in   >>  nbe ;
    679              
    680               edges = new GeometricalEdge[nbe];
    681               if(verbosity>5)
    682                 cout << "     Record Edges: Nb of Edge " << nbe <<endl;
    683               assert(edges);
    684               assert (nbv >0);
    685               Real4 *len =0;
    686               if (!hvertices)
    687                 {
    688                   len = new Real4[nbv];
    689                   for(i=0;i<nbv;i++)
    690                     len[i]=0;
    691                 }
    692              
    693               for (i=0;i<nbe;i++)
    694                 {
    695                   f_in  >> i1   >> i2 >>  edges[i].ref  ;
    696                  
    697                   i1--;i2--; // for C index
    698                   edges[i].v[0]=  vertices + i1;
    699                   edges[i].v[1]=  vertices + i2;
    700                   R2 x12 = vertices[i2].r-vertices[i1].r;
    701                   Real8 l12=sqrt((x12,x12));
    702                   edges[i].tg[0]=zero2;
    703                   edges[i].tg[1]=zero2;
    704                   edges[i].SensAdj[0] = edges[i].SensAdj[1] = -1;
    705                   edges[i].Adj[0] = edges[i].Adj[1] = 0;
    706                   edges[i].flag = 0;
    707                   if (!hvertices)
    708                     {
    709                       vertices[i1].color++;
    710                       vertices[i2].color++;
    711                       len[i1] += l12;
    712                       len[i2] += l12;
    713                     }
    714                  
    715                   Hmin = Min(Hmin,l12);
    716                 }
    717               // definition  the default of the given mesh size
    718               if (!hvertices)
    719                 {
    720                   for (i=0;i<nbv;i++)
    721                     if (vertices[i].color > 0)
    722                       vertices[i].m=  Metric(len[i] /(Real4) vertices[i].color);
    723                     else
    724                       vertices[i].m=  Metric(Hmin);
    725                   delete [] len;
    726                  
    727                   if(verbosity>3)
    728                     cout << "     Geom Hmin " << Hmin << endl;
    729                 }
    730              
    731             }
    732         }
    733       else if (!strcmp(fieldname,"EdgesTangence") ||!strcmp(fieldname,"TangentAtEdges")  )
    734         {
    735           int n,i,j,k;
    736           R2 tg;
    737           f_in  >> n ;
    738          
    739           if(verbosity>5)
    740             cout << "     Record TangentAtEdges: Nb of Edge " << n <<endl;
    741          
    742           for (k=0;k<n;k++)
    743             {
    744               f_in  >>  i  >> j ;
    745               f_in >> tg.x  >> tg.y ;
    746               assert( i <= nbe );
    747               assert( i > 0 );
    748               assert ( j == 1 || j==2 );
    749               i--;j--;// for C index
    750               edges[i].tg[j] = tg;
    751             }
    752         }
    753       else if (!strcmp(fieldname,"Corners"))
    754         {
    755           int i,j,n;
    756           f_in  >> n ;
    757           if(verbosity>5)
    758             cout << "     Record Corner: Nb of Corner " << n <<endl;
    759          
    760           for (i=0;i<n;i++) {     
    761             f_in  >>  j ;
    762             assert( j <= nbv );
    763             assert( j > 0 );
    764             j--;
    765             vertices[j].SetCorner();
    766             vertices[j].SetRequired();  }
    767         }
    768       else if (!strcmp(fieldname,"RequiredVertices"))
    769         {
    770           int i,j,n;
    771           f_in  >> n ;
    772 
    773           for (i=0;i<n;i++) {     
    774             f_in  >>  j ;
    775             assert( j <= nbv );
    776             assert( j > 0 );
    777             j--;
    778             vertices[j].SetRequired();  }
    779       }
    780       else if (!strcmp(fieldname,"RequiredEdges"))
    781         {
    782           int i,j,n;
    783           f_in  >> n ;
    784 
    785           for (i=0;i<n;i++) {     
    786             f_in  >>  j ;
    787             assert( j <= nbe );
    788             assert( j > 0 );
    789             j--;
    790             edges[j].SetRequired();  }
    791       }
    792     else if (!strcmp(fieldname,"SubDomain") || !strcmp(fieldname,"SubDomainFromGeom"))
    793       {
    794         f_in   >>  NbSubDomains ;
    795         if (NbSubDomains>0)
    796           {
    797             subdomains = new GeometricalSubDomain[  NbSubDomains];
    798             Int4 i0,i1,i2,i3;
    799             for (i=0;i<NbSubDomains;i++)
    800               {
    801                 f_in  >> i0  >>i1
    802                       >> i2  >>i3 ;
    803                
    804                 assert(i0 == 2);
    805                 assert(i1<=nbe && i1>0);
    806                 subdomains[i].edge=edges + (i1-1);
    807                 subdomains[i].sens = (int) i2;
    808                 subdomains[i].ref = i3;
    809               }
    810           }
    811       }
    812       else
    813         { // unkown field
    814           field = ++showfield;
    815           if(showfield==1) // just to show one time
    816             if (verbosity>3)
    817               cout << "    Warning we skip the field " << fieldname << " at line " << f_in.LineNumber << endl;
    818         }
    819       showfield=field; // just to show one time
    820     } // while !eof()
    821   // generation  de la geometrie
    822   // 1 construction des aretes
    823   // construire des aretes en chaque sommets
    824  
    825   if (nbv <=0) {
    826     cerr<<"Error: the field Vertex is not found in " << filename<<endl;
    827     NbErr++;}
    828   if(nbe <=0) {
    829     cerr <<"Error: the field Edges is not found in "<< filename<<endl
    830       ;NbErr++;}
    831   if(NbErr) MeshError(1);
    832 
    833  
    834 }
    835 
    836 
    837531}  // end of namespace bamg
  • issm/trunk/src/c/Bamgx/MeshWrite.cpp

    r2796 r2797  
    215215}
    216216
    217 void Geometry::Write(const char * filename)
    218 {
    219         long int verbosity=0;
    220 
    221   ofstream f(filename);
    222   if (f)
    223     {
    224       if(verbosity>1)
    225         cout << "  -- write geometry in file " << filename << endl;
    226        if (name) delete name;
    227        name = new char[strlen(filename)+1];
    228        strcpy(name,filename);
    229        OnDisk =1;
    230        f << *this;
    231     }
    232 }
    233 
    234217void Geometry::WriteGeometry(BamgGeom* bamggeom, BamgOpts* bamgopts){
    235218
     
    391374}
    392375
    393 ostream& operator <<(ostream& f, const   Geometry & Gh)
    394 {
    395    Int4  NbCorner=0;
    396    {
    397      f << "MeshVersionFormatted 0" <<endl;
    398      f << "\nDimension\n"  << 2 << endl;
    399 //     f << "\nIdentifier\n" ;
    400 //     WriteStr(f,Gh.identity);
    401 //     f <<endl;
    402    }
    403    int nbreqv=0;
    404    {
    405      
    406      f.precision(12);
    407      f << "\nVertices\n" << Gh.nbv <<endl;
    408      for (Int4 i=0;i<Gh.nbv;i++)
    409        {
    410          GeometricalVertex & v =  Gh.vertices[i];
    411          if (v.Required()) nbreqv++;
    412          f << v.r.x << " " << v.r.y << " " << v.ref() << endl;
    413          if (v.Corner()) NbCorner++;
    414        }
    415    }
    416    
    417    int nbcracked=0;
    418 
    419    {
    420      int nbreq=0;
    421      f << "\nEdges\n"<< Gh.nbe << endl;
    422      for(Int4 ie=0;ie<Gh.nbe;ie++)
    423        {
    424          
    425          GeometricalEdge & e = Gh.edges[ie];
    426          if (e.Required()) nbreq++;
    427          if (e.Cracked()) {
    428            Int4 ie1 = Gh.Number(e.link);
    429            if (ie <= ie1)  ++nbcracked;}
    430          f << Gh.Number(e[0])+1 << " " << Gh.Number(e[1])+1;
    431          f << " " << e.ref <<endl;
    432        }
    433      
    434      if (nbcracked)
    435        {
    436          f << "\nCrackedEdges\n"<< nbcracked<< endl;
    437          for(Int4 ie=0;ie<Gh.nbe;ie++)
    438            {
    439              GeometricalEdge & e = Gh.edges[ie];
    440              if (e.Cracked()) {
    441                Int4  ie1 = Gh.Number(e.link);
    442                if (ie <= ie1)  f << ie+1 << " " << ie1+1<< endl;
    443              }
    444            }
    445        }
    446      if(nbreq)
    447        {
    448          f << "\nRequiredEdges\n"<< nbreq<< endl;
    449          for(Int4 ie=0;ie<Gh.nbe;ie++)
    450            {
    451              GeometricalEdge & e = Gh.edges[ie];
    452              if (e.Required())
    453                f << ie+1 << endl;
    454            }
    455        }
    456      
    457      
    458      
    459    }
    460 
    461     f << "\nAngleOfCornerBound\n"
    462      << Gh.MaximalAngleOfCorner*180/Pi << endl;
    463     if (NbCorner)
    464       {
    465         f << "\nCorners\n" << NbCorner << endl;
    466         for (Int4 i=0,j=0;i<Gh.nbv;i++)
    467           {
    468             GeometricalVertex & v =  Gh.vertices[i];
    469             if (v.Corner())
    470               j++,f << Gh.Number(v)+1 << (j % 5 ? ' ' : '\n');
    471           }
    472        
    473      
    474       }
    475 
    476     if(nbreqv)
    477       {
    478         f << "\nRequiredVertices\n"<< nbreqv<< endl;
    479         for (Int4 j=0,i=0;i<Gh.nbv;i++)
    480           {
    481             GeometricalVertex & v =  Gh.vertices[i];
    482            
    483             if (v.Required())
    484               j++,f << i+1 << (j % 5 ? ' ' : '\n');
    485           }
    486         f << endl;
    487       }
    488    
    489     {
    490        Int4 i;
    491        f << "\nSubDomainFromGeom\n" ;
    492        f << Gh.NbSubDomains<< endl;
    493        for (i=0;i<Gh.NbSubDomains;i++)
    494          f << "2 " << Gh.Number(Gh.subdomains[i].edge)+1 << " " << Gh.subdomains[i].sens
    495            << " " << Gh.subdomains[i].ref << endl;       
    496      }
    497      {
    498        Int4 n=0,i;
    499 
    500        for(i=0;i< Gh.nbe;i++)
    501          {
    502            if(Gh.edges[i].TgA() && Gh.edges[i][0].Corner() )
    503              n++;
    504            if(Gh.edges[i].TgB() && Gh.edges[i][1].Corner() )
    505              n++;
    506          }
    507        if (n) {
    508          f << "TangentAtEdges " << n << endl;
    509          for(i=0;i< Gh.nbe;i++)
    510            {
    511              if (Gh.edges[i].TgA() && Gh.edges[i][0].Corner() )
    512                f << i+1 << " 1 " << Gh.edges[i].tg[0].x
    513                  << " " << Gh.edges[i].tg[0].y << endl;
    514              if (Gh.edges[i].TgB() && Gh.edges[i][1].Corner() )
    515                f << i+1 << " 2 " << Gh.edges[i].tg[1].x
    516                  << " " << Gh.edges[i].tg[1].y << endl;
    517            }
    518          
    519        }}
    520      //  f << " Not Yet Implemented" << endl;
    521      
    522      return f;
    523 }
    524 
    525 
    526376} // end of namespace bamg
  • issm/trunk/src/c/Bamgx/Metric.cpp

    r2790 r2797  
    10231023}
    10241024
    1025 void Geometry::ReadMetric(const char * fmetrix,Real8 hmin=1.0e-30,Real8 hmax=1.0e30,Real8 coef=1)
    1026 {
    1027         long int verbosity=0;
    1028 
    1029   hmin = Max(hmin,MinimalHmin());
    1030   MeshIstream f_metrix(fmetrix);
    1031   Int4 k,j;
    1032   f_metrix >>  k >> j ;
    1033   if(verbosity>1)
    1034     cout << "  -- ReadMetric  " << fmetrix
    1035          << ",  coef = " << coef
    1036          << ", hmin = " << hmin
    1037          << ", hmax = " << hmax
    1038          << (  (j == 1)? " Iso " : " AnIso " ) << endl;
    1039  
    1040   if (k != nbv ||  !(j == 1 || j == 3)) {
    1041     cerr << " Error Pb metrix " << k << " <> "
    1042          <<  nbv << " or  1 or 3  <> " << j << endl;
    1043     MeshError(1003);}
    1044          
    1045  
    1046   //  Int4 nberr = 0;
    1047   for (Int4 iv=0;iv<nbv;iv++)
    1048     {
    1049     Real8 h;
    1050     if (j == 1)
    1051       {
    1052       f_metrix >>  h ;
    1053       vertices[iv].m=Metric(Max(hmin,Min(hmax, h*coef)));
    1054       }
    1055     else if (j==3)
    1056       {
    1057         Real8 a,b,c;         
    1058         f_metrix >>  a >> b >> c  ;
    1059         MetricAnIso M(a,b,c);
    1060         MatVVP2x2 Vp(M/coef);
    1061         Vp.Maxh(hmin);
    1062       Vp.Minh(hmax);
    1063       vertices[iv].m = Vp;
    1064       }
    1065     }
    1066  
    1067 }
    1068 
    1069 
    10701025Real8 LengthInterpole(const MetricAnIso Ma,const  MetricAnIso Mb, R2 AB)
    10711026{
Note: See TracChangeset for help on using the changeset viewer.