Index: /issm/trunk/src/c/Bamgx/Mesh2.h =================================================================== --- /issm/trunk/src/c/Bamgx/Mesh2.h (revision 2944) +++ /issm/trunk/src/c/Bamgx/Mesh2.h (revision 2945) @@ -123,5 +123,5 @@ Int4 color; Vertex * to;// use in geometry Vertex to now the Mesh Vertex associed - VertexOnGeom * on; // if vint 8; // set with Triangles::SetVertexFieldOn() + VertexOnGeom * onGeometry; // if vint 8; // set with Triangles::SetVertexFieldOn() Vertex * onbv; // if vint == 16 on Background vertex Triangles::SetVertexFieldOnBTh() VertexOnEdge * onbe; // if vint == 32 on Background edge @@ -187,5 +187,5 @@ Vertex * v[2]; Int4 ref; - GeometricalEdge * on; + GeometricalEdge* onGeometry; Vertex & operator[](int i){return *v[i];}; Vertex * operator()(int i){return v[i];}; @@ -410,5 +410,4 @@ aa[a] |=8 ; } - void SetUnMarkUnSwap(int a){ @@ -567,5 +566,5 @@ operator const Real8 & () const {return abscisse;} int IsRequiredVertex(){ return this? (( abscisse<0 ? (gv?gv->Required():0):0 )) : 0;} - void SetOn(){mv->on=this;mv->vint=IsVertexOnGeom;} + void SetOn(){mv->onGeometry=this;mv->vint=IsVertexOnGeom;} inline void Set(const Triangles &,Int4,Triangles &); @@ -655,8 +654,4 @@ public: - enum TypeFileMesh { - AutoMesh=0,BDMesh=1,NOPOMesh=2,amMesh=3,am_fmtMesh=4,amdbaMesh=5, - ftqMesh=6,mshMesh=7}; - int static counter; // to kown the number of mesh in memory Geometry & Gh; // Geometry @@ -689,11 +684,8 @@ VertexOnEdge *VertexOnBThEdge; - Int4 NbCrackedVertices; - Int4 NbCrackedEdges; CrackedEdge *CrackedEdges; - R2 pmin,pmax; // extrema @@ -711,5 +703,4 @@ Triangles(Int4 i);//:BTh(*this),Gh(*new Geometry()){PreInit(i);} - ~Triangles(); Triangles(const char * ,Real8=-1) ; @@ -730,90 +721,85 @@ Triangles(const Triangles &,const int *flag,const int *bb); // truncature - void SetIntCoor(const char * from =0); - - // void RandomInit(); - // void CubeInit(int ,int); - - Real8 MinimalHmin() {return 2.0/coefIcoor;} - Real8 MaximalHmax() {return Max(pmax.x-pmin.x,pmax.y-pmin.y);} - const Vertex & operator[] (Int4 i) const { return vertices[i];}; - Vertex & operator[](Int4 i) { return vertices[i];}; - const Triangle & operator() (Int4 i) const { return triangles[i];}; - Triangle & operator()(Int4 i) { return triangles[i];}; - I2 toI2(const R2 & P) const { - return I2( (Icoor1) (coefIcoor*(P.x-pmin.x)) - ,(Icoor1) (coefIcoor*(P.y-pmin.y)) );} - R2 toR2(const I2 & P) const { - return R2( (double) P.x/coefIcoor+pmin.x, (double) P.y/coefIcoor+pmin.y);} - void Add( Vertex & s,Triangle * t,Icoor2 * =0) ; - void Insert(); - // void InsertOld(); - void ForceBoundary(); - void Heap(); - void FindSubDomain(int OutSide=0); - Int4 ConsRefTriangle(Int4 *) const; - void ShowHistogram() const; - void ShowRegulaty() const; // Add FH avril 2007 - // void ConsLinkTriangle(); - - void ReMakeTriangleContainingTheVertex(); - void UnMarkUnSwapTriangle(); - void SmoothMetric(Real8 raisonmax) ; - void BoundAnisotropy(Real8 anisomax,double hminaniso= 1e-100) ; - void MaxSubDivision(Real8 maxsubdiv); - Edge** MakeGeometricalEdgeToEdge(); - void SetVertexFieldOn(); - void SetVertexFieldOnBTh(); - Int4 SplitInternalEdgeWithBorderVertices(); - void MakeQuadrangles(double costheta); - int SplitElement(int choice); - void MakeQuadTree(); - void NewPoints( Triangles &,int KeepVertices =1 ); - Int4 InsertNewPoints(Int4 nbvold,Int4 & NbTSwap) ; - void NewPoints(int KeepVertices=1){ NewPoints(*this,KeepVertices);} - void ReNumberingTheTriangleBySubDomain(bool justcompress=false); - void ReNumberingVertex(Int4 * renu); - void SmoothingVertex(int =3,Real8=0.3); - Metric MetricAt (const R2 &) const; - GeometricalEdge * ProjectOnCurve( Edge & AB, Vertex & A, Vertex & B,Real8 theta, - Vertex & R,VertexOnEdge & BR,VertexOnGeom & GR); - - Int4 Number(const Triangle & t) const { return &t - triangles;} - Int4 Number(const Triangle * t) const { return t - triangles;} - Int4 Number(const Vertex & t) const { return &t - vertices;} - Int4 Number(const Vertex * t) const { return t - vertices;} - Int4 Number(const Edge & t) const { return &t - edges;} - Int4 Number(const Edge * t) const { return t - edges;} - Int4 Number2(const Triangle * t) const { - // if(t>= triangles && t < triangles + nbt ) - return t - triangles; - // else return t - OutSidesTriangles; - } - - Vertex * NearestVertex(Icoor1 i,Icoor1 j) ; - Triangle * FindTriangleContening(const I2 & ,Icoor2 [3],Triangle *tstart=0) const; - - void ReadMesh(BamgMesh* bamgmesh, BamgOpts* bamgopts); - void WriteMesh(BamgMesh* bamgmesh,BamgOpts* bamgopts); - - void ReadMetric(BamgOpts* bamgopts,const Real8 hmin,const Real8 hmax,const Real8 coef); - void WriteMetric(BamgOpts* bamgopts); - void IntersectConsMetric(BamgOpts* bamgopts); - void BuildMetric0(BamgOpts* bamgopts); - void BuildMetric1(BamgOpts* bamgopts); - void IntersectGeomMetric(BamgOpts* bamgopts); - - - int isCracked() const {return NbCrackedVertices != 0;} - int Crack(); - int UnCrack(); - - void ConsGeometry(Real8 =-1.0,int *equiedges=0); // construct a geometry if no geo - void FillHoleInMesh() ; - int CrackMesh(); + void SetIntCoor(const char * from =0); + + Real8 MinimalHmin() {return 2.0/coefIcoor;} + Real8 MaximalHmax() {return Max(pmax.x-pmin.x,pmax.y-pmin.y);} + const Vertex & operator[] (Int4 i) const { return vertices[i];}; + Vertex & operator[](Int4 i) { return vertices[i];}; + const Triangle & operator() (Int4 i) const { return triangles[i];}; + Triangle & operator()(Int4 i) { return triangles[i];}; + I2 toI2(const R2 & P) const { + return I2( (Icoor1) (coefIcoor*(P.x-pmin.x)) + ,(Icoor1) (coefIcoor*(P.y-pmin.y)) );} + R2 toR2(const I2 & P) const { + return R2( (double) P.x/coefIcoor+pmin.x, (double) P.y/coefIcoor+pmin.y); + } + void Add( Vertex & s,Triangle * t,Icoor2 * =0) ; + void Insert(); + void ForceBoundary(); + void Heap(); + void FindSubDomain(int OutSide=0); + Int4 ConsRefTriangle(Int4 *) const; + void ShowHistogram() const; + void ShowRegulaty() const; // Add FH avril 2007 + + void ReMakeTriangleContainingTheVertex(); + void UnMarkUnSwapTriangle(); + void SmoothMetric(Real8 raisonmax) ; + void BoundAnisotropy(Real8 anisomax,double hminaniso= 1e-100) ; + void MaxSubDivision(Real8 maxsubdiv); + Edge** MakeGeometricalEdgeToEdge(); + void SetVertexFieldOn(); + void SetVertexFieldOnBTh(); + Int4 SplitInternalEdgeWithBorderVertices(); + void MakeQuadrangles(double costheta); + int SplitElement(int choice); + void MakeQuadTree(); + void NewPoints( Triangles &,int KeepVertices =1 ); + Int4 InsertNewPoints(Int4 nbvold,Int4 & NbTSwap) ; + void NewPoints(int KeepVertices=1){ NewPoints(*this,KeepVertices);} + void ReNumberingTheTriangleBySubDomain(bool justcompress=false); + void ReNumberingVertex(Int4 * renu); + void SmoothingVertex(int =3,Real8=0.3); + Metric MetricAt (const R2 &) const; + GeometricalEdge * ProjectOnCurve( Edge & AB, Vertex & A, Vertex & B,Real8 theta, + Vertex & R,VertexOnEdge & BR,VertexOnGeom & GR); + + Int4 Number(const Triangle & t) const { return &t - triangles;} + Int4 Number(const Triangle * t) const { return t - triangles;} + Int4 Number(const Vertex & t) const { return &t - vertices;} + Int4 Number(const Vertex * t) const { return t - vertices;} + Int4 Number(const Edge & t) const { return &t - edges;} + Int4 Number(const Edge * t) const { return t - edges;} + Int4 Number2(const Triangle * t) const { + // if(t>= triangles && t < triangles + nbt ) + return t - triangles; + // else return t - OutSidesTriangles; + } + + Vertex * NearestVertex(Icoor1 i,Icoor1 j) ; + Triangle * FindTriangleContening(const I2 & ,Icoor2 [3],Triangle *tstart=0) const; + + void ReadMesh(BamgMesh* bamgmesh, BamgOpts* bamgopts); + void WriteMesh(BamgMesh* bamgmesh,BamgOpts* bamgopts); + + void ReadMetric(BamgOpts* bamgopts,const Real8 hmin,const Real8 hmax,const Real8 coef); + void WriteMetric(BamgOpts* bamgopts); + void IntersectConsMetric(BamgOpts* bamgopts); + void BuildMetric0(BamgOpts* bamgopts); + void BuildMetric1(BamgOpts* bamgopts); + void IntersectGeomMetric(BamgOpts* bamgopts); + + int isCracked() const {return NbCrackedVertices != 0;} + int Crack(); + int UnCrack(); + + void ConsGeometry(Real8 =-1.0,int *equiedges=0); // construct a geometry if no geo + void FillHoleInMesh() ; + int CrackMesh(); private: - void GeomToTriangles1(Int4 nbvx,int KeepVertices=1);// the real constructor mesh adaption - void GeomToTriangles0(Int4 nbvx);// the real constructor mesh generator - void PreInit(Int4,char * =0 ); + void GeomToTriangles1(Int4 nbvx,int KeepVertices=1);// the real constructor mesh adaption + void GeomToTriangles0(Int4 nbvx); // the real constructor mesh generator + void PreInit(Int4,char* =NULL ); }; @@ -978,6 +964,6 @@ v[0] = ThNew.vertices + Th.Number(v[0]); v[1] = ThNew.vertices + Th.Number(v[1]); - if (on) - on = ThNew.Gh.edges+Th.Gh.Number(on); + if (onGeometry) + onGeometry = ThNew.Gh.edges+Th.Gh.Number(onGeometry); if (adj[0]) adj[0] = ThNew.edges + Th.Number(adj[0]); if (adj[1]) adj[1] = ThNew.edges + Th.Number(adj[1]); @@ -1080,5 +1066,5 @@ inline void Triangles::SetVertexFieldOn(){ for (Int4 i=0;iCrackedEdges){ + if(verbose>3) printf(" processing CrackedEdges\n"); + for (i=0;iNumCrackedEdges;i++){ + throw ErrorException(__FUNCT__,exprintf("Not implemented yet")); + } + } + else{ + if(verbose>3) printf(" no CrackedEdges found\n"); } @@ -912,5 +923,5 @@ } /*}}}1*/ - /*FUNCTION Geometry::Geometry::ProjectOnCurve {{{1*/ + /*FUNCTION Geometry::ProjectOnCurve {{{1*/ GeometricalEdge* Geometry::ProjectOnCurve(const Edge & e,Real8 s,Vertex &V,VertexOnGeom &GV ) const { Real8 save_s=s; @@ -918,9 +929,9 @@ retry: s=save_s; - GeometricalEdge* on = e.on; + GeometricalEdge* on = e.onGeometry; if (!on){ throw ErrorException(__FUNCT__,exprintf("!on")); } - if (!e[0].on || !e[1].on){ + if (!e[0].onGeometry || !e[1].onGeometry){ throw ErrorException(__FUNCT__,exprintf("!e[0].on || !e[1].on")); } @@ -932,9 +943,9 @@ Real8 lge[mxe+1]; int bge=mxe/2,tge=bge; - ge[bge] = e.on; + ge[bge] = e.onGeometry; sensge[bge]=1; R2 V0 = v0,V1=v1,V01=V1-V0; - VertexOnGeom vg0= *v0.on, vg1=*v1.on; + VertexOnGeom vg0= *v0.onGeometry, vg1=*v1.onGeometry; // GeometricalEdge * eg0 = e.on,* eg1 = e.on, *eg=NULL; Index: /issm/trunk/src/c/Bamgx/objects/Triangles.cpp =================================================================== --- /issm/trunk/src/c/Bamgx/objects/Triangles.cpp (revision 2944) +++ /issm/trunk/src/c/Bamgx/objects/Triangles.cpp (revision 2945) @@ -436,5 +436,5 @@ throw ErrorException(__FUNCT__,exprintf("EdgesOnGeometricEdge error: We must have : (i>0 && j >0 && i <= nbe && j <= Gh.nbe)")); } - edges[i-1].on = Gh.edges + j-1; + edges[i-1].onGeometry = Gh.edges + j-1; } } @@ -652,5 +652,5 @@ k=0; for (i=0;iNumEdgesOnGeometricEdge=k; @@ -660,7 +660,7 @@ int count=0; for (i=0;iEdgesOnGeometricEdge[count*2+0]=(double)i+1; //back to Matlab indexing - bamgmesh->EdgesOnGeometricEdge[count*2+1]=(double)Gh.Number(edges[i].on)+1; //back to Matlab indexing + bamgmesh->EdgesOnGeometricEdge[count*2+1]=(double)Gh.Number(edges[i].onGeometry)+1; //back to Matlab indexing count=count+1; } @@ -1498,9 +1498,9 @@ edges[add].v[0] = &triangles[it][VerticesOfTriangularEdge[j][0]]; edges[add].v[1] = &triangles[it][VerticesOfTriangularEdge[j][1]]; - edges[add].on=0; + edges[add].onGeometry=0; if (iCracked()) k++; + if(edges[i].onGeometry->Cracked()) k++; if( k==0) return 0; CurrentTh = this; @@ -1822,5 +1822,5 @@ const int nbe0 = nbe; for (k=i=0;iCracked()) + if(edges[i].onGeometry->Cracked()) { e[nbe] = e[i]; @@ -1828,5 +1828,5 @@ e[nbe].v[0] = e[i].v[1]; e[nbe].v[1] = e[i].v[0]; - e[nbe].on = e[i].on->link ; // fqux + e[nbe].onGeometry = e[i].onGeometry->link ; // faux CrackedEdges[k++]=CrackedEdge(edges,i,nbe); nbe++; @@ -1966,5 +1966,5 @@ throw ErrorException(__FUNCT__,exprintf("Missing Edge %i, v0=%i,v1=%i",i ,Number(edges[i][0]),Number(edges[i][1]))); } - if ( nbswp >=0 && edges[i].on->Cracked()) + if ( nbswp >=0 && edges[i].onGeometry->Cracked()) ta.SetCracked(); } @@ -2260,8 +2260,8 @@ for (i=0;iIsRequiredVertex()) { + if(!edges[i][j].onGeometry->IsRequiredVertex()) { throw ErrorException(__FUNCT__,exprintf("adj and vertex required esge(?)")); } @@ -2723,5 +2723,5 @@ edges[nbe].v[1]=b->to;; edges[nbe].ref = e->ref; - edges[nbe].on = e; + edges[nbe].onGeometry = e; edges[nbe].adj[0] = 0; edges[nbe].adj[1] = 0; @@ -2817,5 +2817,5 @@ edges[nbe].v[1]=vb; edges[nbe].ref = e->ref; - edges[nbe].on = e; + edges[nbe].onGeometry = e; edges[nbe].adj[0] = PreviousNewEdge; if(PreviousNewEdge) @@ -2850,5 +2850,5 @@ edges[nbe].v[1]=vb; edges[nbe].ref = e->ref; - edges[nbe].on = e; + edges[nbe].onGeometry = e; edges[nbe].adj[0] = PreviousNewEdge; edges[nbe].adj[1] = 0; @@ -2991,11 +2991,11 @@ Edge & ei = BTh.edges[iedge]; for(int je=0;je<2;je++) // for the 2 extremites - if (!ei.on->Mark() && ei[je].on->IsRequiredVertex() ) + if (!ei.onGeometry->Mark() && ei[je].onGeometry->IsRequiredVertex() ) { // a begin of curve - int nc = ei.on->CurveNumber; + int nc = ei.onGeometry->CurveNumber; if( - ei.on==Gh.curves[nc].be && - (GeometricalVertex *) *ei[je].on == &(*Gh.curves[nc].be)[Gh.curves[nc].kb] // same extremity + ei.onGeometry==Gh.curves[nc].be && + (GeometricalVertex *) *ei[je].onGeometry == &(*Gh.curves[nc].be)[Gh.curves[nc].kb] // same extremity ) { @@ -3052,10 +3052,10 @@ int k0equi=jedgeequi,k1equi; Edge * peequi= BTh.edges+iedgeequi; - GeometricalEdge *ongequi = peequi->on; + GeometricalEdge *ongequi = peequi->onGeometry; Real8 sNew=Lstep;// abcisse of the new points (phase==1) L=0;// length of the curve Int4 i=0;// index of new points on the curve - register GeometricalVertex * GA0 = *(*peequi)[k0equi].on; + register GeometricalVertex * GA0 = *(*peequi)[k0equi].onGeometry; Vertex *A0; A0 = GA0->to; // the vertex in new mesh @@ -3069,5 +3069,5 @@ if(ongequi->Required() ) { - GeometricalVertex *GA1 = *(*peequi)[1-k0equi].on; + GeometricalVertex *GA1 = *(*peequi)[1-k0equi].onGeometry; A1 = GA1->to; // } @@ -3079,8 +3079,8 @@ k1equi= 1 - k0equi; - if (!pe || !ee.on){ + if (!pe || !ee.onGeometry){ throw ErrorException(__FUNCT__,exprintf("!pe || !ee.on")); } - ee.on->SetMark(); + ee.onGeometry->SetMark(); Vertex & v0=ee[0], & v1=ee[1]; R2 AB= (R2) v1 - (R2) v0; @@ -3124,5 +3124,5 @@ A1->ReferenceNumber = eeequi.ref; A1->DirOfSearch =NoDirOfSearch; - e->on = ongequi; + e->onGeometry = ongequi; e->v[0]= A0; e->v[1]= A1; @@ -3139,12 +3139,12 @@ } - if (ee.on->CurveNumber!=ei.on->CurveNumber){ + if (ee.onGeometry->CurveNumber!=ei.onGeometry->CurveNumber){ throw ErrorException(__FUNCT__,exprintf("ee.on->CurveNumber!=ei.on->CurveNumber")); } - if ( ee[k1].on->IsRequiredVertex()) { - if (!eeequi[k1equi].on->IsRequiredVertex()){ + if ( ee[k1].onGeometry->IsRequiredVertex()) { + if (!eeequi[k1equi].onGeometry->IsRequiredVertex()){ throw ErrorException(__FUNCT__,exprintf("!eeequi[k1equi].on->IsRequiredVertex()")); } - register GeometricalVertex * GA1 = *eeequi[k1equi].on; + register GeometricalVertex * GA1 = *eeequi[k1equi].onGeometry; A1=GA1->to;// the vertex in new mesh if (A1-vertices<0 || A1-vertices>=nbv){ @@ -3165,5 +3165,5 @@ { Edge *e = edges + nbe++; - e->on = ongequi; + e->onGeometry = ongequi; e->v[0]= A0; e->v[1]= A1; @@ -3507,17 +3507,17 @@ { Edge * ei = edges+i; - GeometricalEdge *on = ei->on; - e[Gh.Number(on)] = ei; + GeometricalEdge *onGeometry = ei->onGeometry; + e[Gh.Number(onGeometry)] = ei; } for ( i=0;ion; + GeometricalEdge *onGeometry = ei->onGeometry; int j= ii; - while (!(*on)[j].Required()) { - Adj(on,j); // next geom edge + while (!(*onGeometry)[j].Required()) { + Adj(onGeometry,j); // next geom edge j=1-j; - if (e[Gh.Number(on)]) break; // optimisation - e[Gh.Number(on)] = ei; + if (e[Gh.Number(onGeometry)]) break; // optimisation + e[Gh.Number(onGeometry)] = ei; } } @@ -3681,5 +3681,5 @@ { Vertex & bv = Bh[i]; - if (!bv.on) { + if (!bv.onGeometry) { vertices[nbv].r = bv.r; vertices[nbv++].m = bv.m;} @@ -4474,5 +4474,5 @@ Triangle * lastT = triangles + nbt; for (i=0;i 10000000) { - throw ErrorException(__FUNCT__,exprintf("Loop in forcing Egde, nb de swap=%i, nb of try swap (%i)too big",NbSwap,l)); + throw ErrorException(__FUNCT__,exprintf("Loop in forcing Egde, nb de swap=%i, nb of try swap (%i) too big",NbSwap,l)); } Vertex *aa = tc.EdgeVertex(0), *bb = tc.EdgeVertex(1);