![(please configure the [header_logo] section in trac.ini)](/trac/issm/chrome/common/trac_banner.png)
Changeset 2801
- Timestamp:
- 01/12/10 14:02:07 (15 years ago)
- Location:
- issm/trunk/src/c/Bamgx
- Files:
-
- 2 edited
Legend:
- Unmodified
- Added
- Removed
-
issm/trunk/src/c/Bamgx/Mesh2.cpp
r2789 r2801 146 146 int ap1 = PreviousEdge[a1]; 147 147 int ap2 = PreviousEdge[a2]; 148 #ifdef DRAWING1149 couleur(0);150 t1->Draw();151 t2->Draw();152 #endif153 148 #ifdef DEBUG1 154 149 t1->check(); … … 185 180 t2->check(); 186 181 #endif 187 #ifdef DRAWING1188 couleur(1);189 t1->Draw();190 t2->Draw();191 #endif192 #ifdef DRAWING1193 if( CurrentTh)194 CurrentTh->inquire();195 #endif196 197 182 } // end swap 198 183 … … 574 559 NewItem(A,Metric(ba,v0,bb,v1)); 575 560 return; 576 /*577 cerr << nbegin << det(vi,vj,b) << " deti= " << deti << " detj=" <<detj << endl;578 cerr << "SplitEdge on boucle A" << A << " B = " << B << endl;579 580 #ifdef DRAWING581 reffecran();582 Bh.Draw();583 penthickness(5);584 Move(A);585 Line(B);586 penthickness(1);587 588 Bh.inquire();589 penthickness(5);590 Move(A);591 Line(B);592 penthickness(1);593 Bh.inquire();594 #endif595 MeshError(997);*/596 561 } 597 562 } // end while (t->det <0) … … 961 926 cout << " detvasasb =" << detvasasb << "detvbsasb = " << detvbsasb 962 927 << " " << pva << " " << pvb << " " <<CurrentTh <<endl; 963 #ifdef DRAWING1964 reffecran();965 CurrentTh->Draw();966 penthickness(10);967 pva->MoveTo();pvb->LineTo();968 penthickness(1);969 CurrentTh->inquire();970 #endif971 928 } 972 929 assert( ( (detvasasb <= 0) && (detvbsasb >= 0)) || ( (detvasasb >= 0) && (detvbsasb <= 0))); … … 993 950 << " if1 = " << ((ndet1 >0) && (ndet2 >0)) 994 951 << " if2 = " << ((dets1 <=0 && dets2 <=0) || (dets2 >=0 && detsb >=0)) << endl; 995 #ifdef DRAWING996 couleur(0);997 t1->Draw();998 t2->Draw();999 #endif1000 952 } 1001 953 #endif … … 1005 957 #ifdef DEBUG 1006 958 if (ForDebugging) { 1007 #ifdef DRAWING1008 couleur(4);1009 t1->Draw();1010 t2->Draw();1011 rattente(1);1012 #endif1013 959 } 1014 960 #endif … … 1081 1027 } 1082 1028 1083 #ifdef DRAWING11084 a.MoveTo();b.LineTo();1085 #endif1086 1087 1029 while (v2 != &b) { 1088 1030 TriangleAdjacent tc = Previous(Adj(tta)); … … 1128 1070 ForDebugging = 1; 1129 1071 #endif 1130 #ifdef DRAWING11131 if ( CurrentTh ) {1132 reffecran();1133 couleur(6);1134 CurrentTh->Draw();1135 couleur(1);1136 penthickness(10);1137 a.MoveTo();b.LineTo();1138 penthickness(1);1139 CurrentTh->inquire();1140 couleur(6);1141 l=0;1142 reffecran();1143 while (ks=SwapForForcingEdge( va, vb, tc, detss, det1,det2,NbSwap) && (l++ < 1000))1144 cerr << " " << CurrentTh->Number(tc.EdgeVertex(0))<<" " <<CurrentTh->Number(tc.EdgeVertex(1)) << " ";1145 }1146 #endif1147 1072 MeshError(990); 1148 1073 } … … 1577 1502 { 1578 1503 cout << " Pb swap the point s is on a edge =>swap " << iedge << " " << *tt[izerodet] << endl; 1579 #ifdef DRAWING1580 if( CurrentTh && withrgraphique)1581 {1582 reffecran();1583 1584 DrawMark(s.r);1585 CurrentTh->inquire();1586 DrawMark(s.r);1587 rattente(1);1588 }1589 #endif1590 1504 } 1591 1505 assert(rswap); … … 1597 1511 tt[2]->check(); 1598 1512 #endif 1599 #ifdef DRAWING11600 tt[0]->Draw();1601 tt[1]->Draw();1602 tt[2]->Draw();1603 #endif1604 1605 1513 } 1606 1514 … … 1665 1573 tcvi = FindTriangleContening(vi.i,dete); 1666 1574 cout << (*tcvi)[1] << (*tcvi)[2] << endl; 1667 #ifdef DRAWING11668 inquire();1669 penthickness(5);1670 DrawMark(vi.r);1671 penthickness(1);1672 inquire();1673 #endif1674 1675 1575 MeshError(1001,this); 1676 1576 } … … 1678 1578 1679 1579 quadtree->Add(vi); 1680 #ifdef DRAWING11681 DrawMark(vi.r);1682 #endif1683 1580 assert (tcvi && tcvi->det >= 0) ;// internal 1684 1581 Add(vi,tcvi,dete); … … 1758 1655 tcvj = FindTriangleContening(vj.i,dete); 1759 1656 cout << (*tcvj)[1] << (*tcvj)[2] << endl; 1760 #ifdef DRAWING11761 inquire();1762 penthickness(5);1763 DrawMark(vj.r);1764 penthickness(1);1765 1766 inquire();1767 #endif1768 1769 1657 MeshError(1001,this); 1770 1658 } … … 1772 1660 1773 1661 quadtree->Add(vj); 1774 #ifdef DRAWING11775 DrawMark(vj.r);1776 #endif1777 1662 assert (tcvj && tcvj->det >= 0) ;// internal 1778 1663 Add(vj,tcvj,dete); … … 1787 1672 nbv = iv; 1788 1673 } 1789 1790 #ifdef DRAWING11791 inquire();1792 #endif1793 1674 1794 1675 for (i=nbvold;i<nbv;i++) … … 1823 1704 } 1824 1705 if(NbErr) { 1825 #ifdef DRAWING1826 Int4 kkk=0;1827 // UnMarkUnSwapTriangle();1828 // for (i=0;i<nbv;i++)1829 // kkk += vertices[i].Optim(0);1830 if(verbosity>3)1831 cout << " Nb of swap louche " << kkk << endl;1832 if(kkk) {1833 for (i=0;i<nbt;i++)1834 if (triangles[i].link)1835 {1836 double dd =Det(triangles[i][1].r-triangles[i][0].r,triangles[i][2].r-triangles[i][0].r);1837 if(dd <=0)1838 {1839 NbErr++;1840 cerr << " xxxdet triangle i " << i << " = " << triangles[i].det ;1841 cerr << " xxxdet triangle " << dd ;1842 cerr << " xxxLes trois sommets " ;1843 cerr << Number(triangles[i][0]) << " " << Number(triangles[i][1]) << " "1844 << Number(triangles[i][2]) << endl;1845 cerr << triangles[i][0].r << triangles[i][1].r << triangles[i][2].r << endl;1846 cerr << triangles[i][0].i << triangles[i][1].i << triangles[i][2].i << endl;1847 }1848 } }1849 inquire();1850 #endif1851 // MeshError(11);1852 1706 } 1853 1707 … … 1906 1760 nbtold = nbt; 1907 1761 nbvold = nbv; 1908 #ifdef DRAWING11909 inquire();1910 #endif1911 1762 1912 1763 // default size of IntersectionTriangle … … 1949 1800 1950 1801 }// for triangle 1951 1952 #ifdef DRAWING11953 cout << " -------------------------------------------- " << endl;1954 inquire();1955 reffecran();1956 Draw();1957 penthickness(2);1958 #endif1959 1802 1960 1803 if (!InsertNewPoints(nbvold,NbTSwap)) … … 2043 1886 nbtold = nbt; 2044 1887 nbvold = nbv; 2045 #ifdef DRAWING12046 inquire();2047 #endif2048 2049 1888 // default size of IntersectionTriangle 2050 1889 … … 2129 1968 2130 1969 2131 #ifdef DRAWING12132 penthickness(2); Move(A);Line(B); penthickness(1);2133 #endif2134 1970 const Int4 NbvOld = nbv; 2135 1971 lIntTria.SplitEdge(Bh,A,B); … … 2213 2049 R2 VP = V10*s + V11*(1-s); 2214 2050 int sss = (c11-c10) >0 ? 1 : -1; 2215 #ifdef DRAWING12216 penthickness(2);2217 Move((R2) vi0);2218 Line(VP);2219 penthickness(1);2220 #endif2221 2051 // find the 2 point by dichotomie 2222 2052 //cout << " t =" << Number(t) << " c0 " << c0 ; … … 2260 2090 2261 2091 2262 #ifdef DRAWING12263 Move(vertices[i0].r);2264 Line(vertices[i1].r);2265 #endif2266 2092 #ifdef TRACETRIANGLE 2267 2093 if(trace) { … … 2304 2130 } 2305 2131 #endif 2306 #ifdef DRAWING12307 Move(vertices[i0].r);2308 Line(vertices[i1].r);2309 DrawMark(C);2310 #endif2311 2132 // update the new point points of the list 2312 2133 for (ip0=i0;i0 != (ipp0 = vertices[ip0].ReferenceNumber);ip0=ipp0) … … 2317 2138 } // for (i0= .... 2318 2139 }// for triangle 2319 2320 #ifdef DRAWING12321 cout << " -------------------------------------------- " << endl;2322 inquire();2323 reffecran();2324 Draw();2325 penthickness(2);2326 #endif2327 2140 2328 2141 // remove of all the double points … … 2338 2151 vertices[kkk] = vertices[i]; 2339 2152 vertices[kkk].i = toI2(vertices[kkk].r); 2340 #ifdef DRAWING12341 DrawMark(vertices[kkk]);2342 #endif2343 2153 vertices[kkk++].ReferenceNumber = 0; 2344 2154 2345 2155 } 2346 #ifdef DRAWING12347 penthickness(1);2348 #endif2349 2156 2350 2157 // insertion part --- … … 2377 2184 2378 2185 quadtree->Add(*vi); // 2379 #ifdef DRAWING12380 DrawMark(vi->r);2381 #endif2382 2186 assert (tcvi->det >= 0) ;// internal 2383 2187 Add(*vi,tcvi,dete),NbSwap += vi->Optim(1); … … 2385 2189 } 2386 2190 cout << " Nb swap = " << NbSwap << " after " ; 2387 #ifdef DRAWING12388 inquire();2389 #endif2390 2191 2391 2192 for (i=nbvold;i<nbv;i++) … … 2514 2315 Add(*vi,tcvi,dete); 2515 2316 NbSwap += vi->Optim(1,0); 2516 #ifdef DRAWING12517 inquire();2518 #endif2519 2317 2520 2318 }// fin de boucle en icount … … 2569 2367 << " Init Total Cpu Time = " << time3 - time0 << "s " << endl; 2570 2368 2571 #ifdef DRAWING12572 inquire();2573 #endif2574 2369 CurrentTh=OldCurrentTh; 2575 2370 } … … 2630 2425 { 2631 2426 long int verbosity=20; 2632 //#define DRAWING12633 2427 2634 2428 if (verbosity >2) … … 2654 2448 for (Int4 itt=0;itt<nbt;itt++) 2655 2449 triangles[itt].link=0; // par defaut pas de couleur 2656 #ifdef DRAWING12657 reffecran();2658 #endif2659 2450 2660 2451 Int4 NbSubDomTot =0; … … 2665 2456 Int4 i = 0; // niveau de la pile 2666 2457 t->link = t ; // sd forme d'un triangle cicular link 2667 #ifdef DRAWING12668 t->Draw(NbSubDomTot-1);2669 #endif2670 2671 2458 2672 2459 HeapTriangle[i] =t ; … … 2684 2471 { 2685 2472 i++; 2686 #ifdef DRAWING12687 ta->Draw(NbSubDomTot-1);2688 #endif2689 2473 ta->link = t->link ; // on chaine les triangles 2690 2474 t->link = ta ; // d'un meme sous domaine … … 2751 2535 t = t1->link; 2752 2536 mark[it]=k; 2753 #ifdef DRAWING12754 t1->Draw(k);2755 #endif2756 2537 subdomains[k].head = t1; 2757 2538 // cout << " New -- " << Number(t1) << " " << it << endl; 2758 2539 do {// cout << " k " << k << " " << Number(t) << endl; 2759 2540 mark[Number(t)]=k; 2760 #ifdef DRAWING12761 t->Draw(k);2762 #endif2763 2541 t=t->link; 2764 2542 } while (t!=t1); 2765 #ifdef DRAWING12766 t1->Draw(k);2767 #endif2768 2543 mark[it]=k++;} 2769 2544 // else if(mark[it] == -2 ) triangles[it].Draw(999); … … 2902 2677 kkk++; 2903 2678 assert(mark[Number(tt)]<0); 2904 #ifdef DRAWING12905 tt->Draw(i);2906 #endif2907 2679 mark[Number(tt)]=i; 2908 2680 tt=tt->link; … … 2937 2709 2938 2710 } 2939 2940 #ifdef DRAWING12941 inquire();2942 #endif2943 2711 NbOutT=0; 2944 2712 for (it=0;it<nbt;it++) … … 2948 2716 if (verbosity> 2) 2949 2717 cout << " Nb of Sub borned Domain = " << NbSubDomTot << " NbOutTriangles = " << NbOutT <<endl; 2950 #ifdef DRAWING12951 inquire();2952 #endif2953 2954 //#undef DRAWING12955 2718 2956 2719 … … 3267 3030 void Triangles::GeomToTriangles1(Int4 inbvx,int KeepBackVertices) 3268 3031 { 3269 //#define DRAWING13270 3032 Gh.NbRef++;// add a ref to Gh 3271 3033 … … 3296 3058 PreInit(inbvx); 3297 3059 BTh.SetVertexFieldOn(); 3298 #ifdef DRAWING3299 if (withrgraphique)3300 { BTh.InitDraw();3301 reffecran();3302 CurrentTh = this;}3303 #endif3304 3060 int * bcurve = new int[Gh.NbOfCurves]; // 3305 3061 … … 3542 3298 PreviousNewEdge->adj[1] = e; 3543 3299 PreviousNewEdge = e; 3544 #ifdef DRAWING13545 e->Draw();3546 // A0->Draw();3547 A1->m.Draw(*A1);3548 A1->Draw(Number(A1));3549 #endif3550 3300 A0=A1; 3551 3301 sNew += Lstep; … … 3592 3342 PreviousNewEdge->adj[1] = e; 3593 3343 PreviousNewEdge = e; 3594 // cout << "Last new edge " << nbe << " " << " on " << Gh.Number(pe->on) 3595 // << " of curve =" <<pe->on->CurveNumber <<endl; 3596 3597 3598 #ifdef DRAWING1 3599 e->Draw(); 3600 A1->Draw(); 3601 A0->Draw(); 3602 // inquire(); 3603 #endif 3344 3604 3345 assert(i==NbCreatePointOnCurve); 3605 3346 … … 3657 3398 3658 3399 3659 #ifdef DRAWING13660 reffecran();3661 InitDraw();3662 Draw();3663 inquire();3664 #endif3665 3666 3400 Insert(); 3667 3401 ForceBoundary(); 3668 3402 FindSubDomain(); 3669 3403 3670 #ifdef DRAWING13671 reffecran();3672 Draw();3673 inquire();3674 #endif3675 // NewPointsOld(*this) ;3676 // BTh.ReMakeTriangleContainingTheVertex(); // FH change => put in NewPoints3677 // for (Int4 iv=0;iv<BTh.nbv;iv++)3678 // BTh[iv].i = toI2(BTh[iv].r);3679 3404 NewPoints(BTh,KeepBackVertices) ; 3680 3405 CurrentTh = 0; 3681 //#undef DRAWING13682 3406 } 3683 3407 … … 3689 3413 Int4 i,NbOfCurves=0,NbNewPoints,NbEdgeCurve; 3690 3414 Real8 lcurve, lstep,s; 3691 #ifdef DRAWING3692 if (withrgraphique)3693 {3694 Gh.InitDraw() ;3695 CurrentTh = this; }3696 #endif3697 3415 3698 3416 R2 AB; … … 3731 3449 // generation of the curves 3732 3450 assert(! edges); 3733 #ifdef DRAWING13734 reffecran();3735 #endif3736 3451 // 2 step 3737 3452 // --step=0 to compute the number of edges + alloc at end … … 3774 3489 edges[nbe].adj[0] = 0; 3775 3490 edges[nbe].adj[1] = 0; 3776 #ifdef DRAWING13777 edges[nbe].Draw();3778 #endif3779 3491 nbe++;} 3780 3492 } … … 3879 3591 if(PreviousNewEdge) 3880 3592 PreviousNewEdge->adj[1] = &edges[nbe]; 3881 #ifdef DRAWING13882 vb->Draw();3883 edges[nbe].Draw();3884 #endif3885 3593 PreviousNewEdge = edges + nbe; 3886 3594 nbe++; … … 3936 3644 PreviousNewEdge->adj[1] = & edges[nbe]; 3937 3645 3938 3939 #ifdef DRAWING13940 edges[nbe].Draw();3941 #endif3942 3646 nbe++;} 3943 3647 else … … 3965 3669 } // for (step=0;step<2;step++) 3966 3670 3967 #ifdef DRAWING13968 reffecran();3969 InitDraw();3970 Draw();3971 inquire();3972 #endif3973 3974 3671 Insert(); 3975 3672 ForceBoundary(); 3976 3673 FindSubDomain(); 3977 3674 3978 #ifdef DRAWING13979 reffecran();3980 Draw();3981 inquire();3982 #endif3983 3675 // NewPointsOld(*this) ; 3984 3676 NewPoints(*this,0) ; … … 4497 4189 } 4498 4190 4499 /** -- old with a bug we loss some time last swap4500 4501 Int4 Triangle::Optim(Int2 i,int koption)4502 {4503 // turn in the positif sens around vertex s4504 register Int4 NbSwap =0;4505 register Vertex * s = ns[i];4506 register Triangle * tbegin=0 , *t = this , *ttc;4507 register int k=0,j = EdgesVertexTriangle[i][0],jc;4508 tbegin=t;4509 do {4510 k++;4511 #ifdef DEBUG4512 assert( s == & (*t)[VerticesOfTriangularEdge[j][1]] );4513 #endif4514 #ifdef DRAWING14515 t->Draw();4516 DrawMark( s->r);4517 #endif4518 ttc = t->at[j];4519 jc = NextEdge[t->aa[j]&3];4520 cout << *t << " " << VerticesOfTriangularEdge[j][1] << "\n\t try swap " << * ttc << " " << jc ;4521 while ( ttc->swap(jc,koption)) {4522 NbSwap++,assert(k++<20000);4523 ttc = t->at[j];4524 jc = NextEdge[t->aa[j]&3];4525 cout << "\n\t s " << *ttc << " " << jc << endl;4526 }4527 cout << endl;4528 t = ttc;4529 j = NextEdge[jc];4530 assert(k<20000);4531 } while ( (tbegin != t));4532 4533 return NbSwap;4534 }4535 */4536 4191 Int4 Triangle::Optim(Int2 i,int koption) 4537 4192 { … … 4626 4281 cout << " MakeQuadTree" << endl; 4627 4282 if ( !quadtree ) quadtree = new QuadTree(this); 4628 4629 4630 #ifdef DRAWING14631 quadtree->Draw();4632 rattente(1);4633 reffecran();4634 quadtree->Draw();4635 rattente(1);4636 #endif4637 4283 4638 4284 } … … 4731 4377 R2 PQ = vQ.r - vP.r; 4732 4378 Real8 l = log(LengthInterpole(vP,vQ,PQ)); 4733 #ifdef DRAWING24734 if (l>1.4) {4735 penthickness(3);4736 vP.MoveTo(),vQ.LineTo();4737 penthickness(1);4738 cout << " l = " << l << Number(vP) << " edge = " << Number(vQ) << endl;4739 }4740 #endif4741 4379 nbedges++; 4742 4380 k = (int) ((l - lmin)*delta); … … 5052 4690 } 5053 4691 assert(a>= vertices && a < vertices+nbv); 5054 #ifdef DRAWING15055 a->Draw();5056 #endif5057 4692 // int k=0; 5058 4693 t = a->t; … … 5083 4718 assert( kkkk++ < 2000 ); 5084 4719 j= OppositeVertex[jj]; 5085 5086 #ifdef DRAWING15087 t->Draw();5088 #endif5089 4720 dete[j] = detop; //det(*b,*s1,*s2); 5090 4721 jn = NextVertex[j]; -
issm/trunk/src/c/Bamgx/Mesh2.h
r2799 r2801 40 40 #include <unistd.h> 41 41 #endif 42 #ifdef DRAWING43 #include "rgraph.hpp"44 #endif45 42 46 43 extern int SHOW; … … 90 87 {return a<0 ? Pi + a :a - Pi ;} 91 88 92 #ifdef DRAWING93 extern Real4 xGrafCoef,yGrafCoef,xGrafOffSet,yGrafOffSet;94 extern R2 GrafPMin,GrafPMax;95 extern Real8 Grafh;96 #endif97 98 89 Icoor2 inline det(const I2 &a,const I2 & b,const I2 &c) 99 90 { … … 157 148 r = (x*i + y*j) >=0;} 158 149 return r;} 159 #ifdef DRAWING160 void Draw() {161 if (dir!= MaxICoor) {162 Icoor2 x(dir/2),y1(MaxICoor/2-Abs(x)),y(dir%2?-y1:y1);163 R2 D(x,y);164 double eps = Grafh/Norme2(D)/20;165 D = D*eps;166 rmoveto(D.x,D.y);167 }168 }169 #endif170 171 172 173 174 150 }; 175 151 ///////////////////////////////////////////////////////////////////////////////////// … … 203 179 {f << "(" << v.i << "," << v.r << MatVVP2x2(v.m) << ")" ; return f;} 204 180 inline void Set(const Vertex & rec,const Triangles &,Triangles &); 205 206 #ifdef DRAWING207 void Draw(Int4 =-1) const ;208 void MoveTo() const { rmoveto(r.x,r.y); }209 void LineTo() const { rlineto(r.x,r.y); }210 #endif211 181 }; 212 182 … … 282 252 inline void Set(const Triangles &,Int4,Triangles &); 283 253 284 #ifdef DRAWING285 void Draw(Int4 = -1) const ;286 #endif287 254 }; // end of Edge class 288 255 … … 352 319 inline void Set(const GeometricalEdge & rec,const Geometry & Th ,Geometry & ThNew); 353 320 354 #ifdef DRAWING355 void Draw(Int4 =-1);356 #endif357 358 321 }; 359 322 … … 493 456 void inline checka(Int1 a); 494 457 void inline check(); 495 #endif496 497 #ifdef DRAWING498 void Draw(Int4 i=-1) const;499 int swapDRAW(Int2 a1);500 501 458 #endif 502 459 … … 1428 1385 #endif 1429 1386 1430 1431 1432 1433 #ifdef DRAWING 1434 extern Real4 xGrafCoef,yGrafCoef,xGrafOffSet,yGrafOffSet; // R2 -> I2 transform 1435 extern R2 Gpmin,Gpmax; 1436 //extern Real8 Gh; 1437 // cf routine ILineTo IMoveto 1438 1439 extern void IMoveTo(long i,long j); 1440 extern void ILineTo(long i,long j); 1441 extern char Getxyc(long &i,long &j); 1442 extern void Draw(float ,float ); 1443 extern void Draw(long ,long ); 1444 extern void DrawMark(R2 r); 1445 //inline void DrawMark(D2 r) {DrawMark(R2(r.x,r.y));} 1446 inline void Move(I2 x) {IMoveTo(x.x,x.y);} 1447 inline void Move(R2 x) {rmoveto(x.x,x.y);} 1448 //inline void Move(D2 x) {rmoveto(x.x,x.y);} 1449 inline void Line(I2 x){ILineTo(x.x,x.y);} 1450 inline void Line(R2 x) {rlineto(x.x,x.y);} 1451 //inline void Line(D2 x) {rlineto(x.x,x.y);} 1387 } 1452 1388 #endif 1453 1454 }1455 1456 1457 #endif1458
Note:
See TracChangeset
for help on using the changeset viewer.
