Changeset 2847

Timestamp:

01/14/10 17:38:36 (15 years ago)

Author:

Mathieu Morlighem

Message:

added Bamgx/objects/Vertex.cpp

Location:

issm/trunk/src/c

Files:

• Bamgx/Mesh2.cpp (modified) (1 diff)
• Bamgx/objects/Vertex.cpp (added)
• Makefile.am (modified) (2 diffs)

issm/trunk/src/c/Bamgx/Mesh2.cpp

@@ -220 +220 @@
                 }
 
-
-
-
-                Real8  Vertex::Smoothing(Triangles & Th,const Triangles & BTh,Triangle  * & tstart ,Real8 omega)
-                {
-                        register Vertex * s  = this;
-                        Vertex &vP = *s,vPsave=vP;
-                        //  if (vP.on) return 0;// Don't move boundary vertex  
-
-                        register Triangle * tbegin= t , *tria = t , *ttc;
-
-                        register int k=0,kk=0,j = EdgesVertexTriangle[vint][0],jc;
-                        R2 P(s->r),PNew(0,0);
-                        //  cout << BTh.quadtree << " " <<  BTh.quadtree->root << endl;
-                        // assert(BTh.quadtree && BTh.quadtree->root);
-                        do {
-                                k++; 
-
-                                if (!tria->Hidden(j))
-                                {
-                                        Vertex &vQ = (*tria)[VerticesOfTriangularEdge[j][0]]; 
-
-                                        R2 Q = vQ,QP(P-Q);
-                                        Real8 lQP = LengthInterpole(vP,vQ,QP);
-                                        PNew += Q+QP/Max(lQP,1e-20);
-                                        kk ++;
-                                }
-                                ttc =  tria->TriangleAdj(j);
-                                jc = NextEdge[tria->NuEdgeTriangleAdj(j)];
-                                tria = ttc;
-                                j = NextEdge[jc];
-                                assert(k<2000);
-                        } while ( tbegin != tria); 
-                        if (kk<4) return 0;
-                        PNew = PNew/(Real8)kk;
-                        R2 Xmove((PNew-P)*omega);
-                        PNew = P+Xmove;
-                        Real8 delta=Norme2_2(Xmove); 
-
-
-                        // 
-                        Icoor2 deta[3];
-                        I2 IBTh  = BTh.toI2(PNew);
-
-                        tstart=BTh.FindTriangleContening(IBTh,deta,tstart);  
-
-                        if (tstart->det <0) 
-                        { // outside
-                                double ba,bb;
-                                TriangleAdjacent edge= CloseBoundaryEdge(IBTh,tstart,ba,bb) ;
-                                tstart = edge;
-                                vP.m= Metric(ba,*edge.EdgeVertex(0),bb,*edge.EdgeVertex(1));
-                        }
-                        else 
-                        { // inside
-                                Real8   aa[3];
-                                Real8 s = deta[0]+deta[1]+deta[2];
-                                aa[0]=deta[0]/s;
-                                aa[1]=deta[1]/s;
-                                aa[2]=deta[2]/s;
-                                vP.m = Metric(aa,(*tstart)[0],(*tstart)[1],(*tstart)[2]);
-                        }
-
-                        // recompute the det of the triangle
-                        vP.r = PNew;
-
-                        vP.i = Th.toI2(PNew);
-
-                        Vertex vPnew = vP;
-
-                        int ok=1;
-                        int loop=1;
-                        k=0;
-                        while (ok) 
-                        {
-                                ok =0;
-                                do {
-                                        k++; 
-                                        double detold = tria->det;
-                                        tria->det =  bamg::det( (*tria)[0],(*tria)[1]  ,(*tria)[2]);
-                                        if (loop) 
-                                        {
-                                                Vertex *v0,*v1,*v2,*v3;
-                                                if (tria->det<0) ok =1;                        
-                                                else if (tria->Quadrangle(v0,v1,v2,v3))
-                                                {
-                                                        vP = vPsave;
-                                                        Real8 qold =QuadQuality(*v0,*v1,*v2,*v3);
-                                                        vP = vPnew;
-                                                        Real8 qnew = QuadQuality(*v0,*v1,*v2,*v3);
-                                                        if (qnew<qold) ok = 1;
-                                                }
-                                                else if ( (double)tria->det < detold/2 ) ok=1;
-
-                                        }
-                                        tria->SetUnMarkUnSwap(0);
-                                        tria->SetUnMarkUnSwap(1);
-                                        tria->SetUnMarkUnSwap(2);
-                                        ttc =  tria->TriangleAdj(j);
-                                        jc = NextEdge[tria->NuEdgeTriangleAdj(j)];
-                                        tria = ttc;
-                                        j = NextEdge[jc];
-                                        assert(k<2000);
-                                } while ( tbegin != tria); 
-                                if (ok && loop) vP=vPsave; // no move 
-                                loop=0;
-                        }
-                        return delta;
-                }
-
-
-double QuadQuality(const Vertex & a,const Vertex &b,const Vertex &c,const Vertex &d) {
-        // calcul de 4 angles --
-        R2 A((R2)a),B((R2)b),C((R2)c),D((R2)d);
-        R2 AB(B-A),BC(C-B),CD(D-C),DA(A-D);
-        //  Move(A),Line(B),Line(C),Line(D),Line(A);
-        const Metric & Ma  = a;
-        const Metric & Mb  = b;
-        const Metric & Mc  = c;
-        const Metric & Md  = d;
-
-        double lAB=Norme2(AB);
-        double lBC=Norme2(BC);
-        double lCD=Norme2(CD);
-        double lDA=Norme2(DA);
-        AB /= lAB;  BC /= lBC;  CD /= lCD;  DA /= lDA;
-        // version aniso 
-        double cosDAB= Ma(DA,AB)/(Ma(DA)*Ma(AB)),sinDAB= Det(DA,AB);
-        double cosABC= Mb(AB,BC)/(Mb(AB)*Mb(BC)),sinABC= Det(AB,BC);
-        double cosBCD= Mc(BC,CD)/(Mc(BC)*Mc(CD)),sinBCD= Det(BC,CD);
-        double cosCDA= Md(CD,DA)/(Md(CD)*Md(DA)),sinCDA= Det(CD,DA);
-        double sinmin=Min(Min(sinDAB,sinABC),Min(sinBCD,sinCDA));
-        if (sinmin<=0) return sinmin;
-        return 1.0-Max(Max(Abs(cosDAB),Abs(cosABC)),Max(Abs(cosBCD),Abs(cosCDA)));
 }
 
-}
-

issm/trunk/src/c/Makefile.am

@@ -329 +329 @@
                                         ./Bamgx/objects/Triangles.cpp   \
                                         ./Bamgx/objects/Triangle.cpp    \
+                                        ./Bamgx/objects/Vertex.cpp      \
                                         ./Bamgx/objects/Geometry.cpp    \
                                         ./Bamgx/objects/QuadTree.cpp \
@@ -670 +671 @@
                                         ./Bamgx/objects/Triangles.cpp   \
                                         ./Bamgx/objects/Triangle.cpp    \
+                                        ./Bamgx/objects/Vertex.cpp      \
                                         ./Bamgx/objects/Geometry.cpp    \
                                         ./Bamgx/objects/QuadTree.cpp \