Changeset 2782
- Timestamp:
- 01/07/10 14:31:05 (15 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
issm/trunk/src/c/Bamgx/Mesh2.cpp
r2740 r2782 4128 4128 } 4129 4129 4130 void Triangles::FillHoleInMesh(){ 4131 4130 void Triangles::FillHoleInMesh() { 4132 4131 Triangles * OldCurrentTh =CurrentTh; 4133 long int verbosity=200;4134 4135 4132 CurrentTh=this; 4136 4133 // Int4 NbTold = nbt; 4137 4134 // generation of the integer coor 4138 4139 // double coef = coefIcoor; 4140 // recherche des extrema des vertices pmin,pmax 4141 Int4 i; 4142 if(verbosity>2) 4143 cout << " -- FillHoleInMesh: Nb of vertices =" << nbv 4135 { 4136 4137 // double coef = coefIcoor; 4138 // recherche des extrema des vertices pmin,pmax 4139 Int4 i; 4140 if(verbosity>2) 4141 cout << " -- FillHoleInMesh: Nb of vertices =" << nbv 4144 4142 << " Pmin = "<< pmin << " Pmax = "<< pmax << endl; 4145 4143 4146 assert(ordre);4147 for (i=0;i<nbv;i++)4148 ordre[i]= 0 ;4149 4150 4151 NbSubDomains =0;4152 4153 // generation of the adjacence of the triangles4154 SetOfEdges4 * edge4= new SetOfEdges4(nbt*3,nbv);4155 Int4 * st = new Int4[nbt*3];4156 for (i=0;i<nbt*3;i++)4157 st[i]=-1;4158 Int4 kk =0;4159 for (i=0;i<nbe;i++)4160 kk += (i == edge4->addtrie(Number(edges[i][0]),Number(edges[i][1])));4161 if (kk != nbe)4162 {4163 cerr << " Some Double edge in the mesh, the number is " << kk-nbe << endl;4164 MeshError(1002,this);4165 }4166 for (i=0;i<nbt;i++)4167 for (int j=0;j<3;j++)4168 {4169 // Int4 i0,i1;4170 Int4 k =edge4->addtrie(Number(triangles[i][VerticesOfTriangularEdge[j][0]]),4171 Number(triangles[i][VerticesOfTriangularEdge[j][1]]));4172 Int4 invisible = triangles[i].Hidden(j);4173 if(st[k]==-1)4174 4175 else if(st[k]>=0) {4176 assert( ! triangles[i].TriangleAdj(j) && !triangles[st[k] / 3].TriangleAdj((int) (st[k]%3)));4177 4178 triangles[i].SetAdj2(j,triangles + st[k] / 3,(int) (st[k]%3));4179 if (invisible) triangles[i].SetHidden(j);4180 if (k<nbe) {4181 triangles[i].SetLocked(j);4182 }4183 st[k]=-2-st[k]; }4184 else {4185 cerr << " The edge ("4144 assert(ordre); 4145 for (i=0;i<nbv;i++) 4146 ordre[i]= 0 ; 4147 4148 4149 NbSubDomains =0; 4150 4151 // generation of the adjacence of the triangles 4152 SetOfEdges4 * edge4= new SetOfEdges4(nbt*3,nbv); 4153 Int4 * st = new Int4[nbt*3]; 4154 for (i=0;i<nbt*3;i++) 4155 st[i]=-1; 4156 Int4 kk =0; 4157 for (i=0;i<nbe;i++) 4158 kk += (i == edge4->addtrie(Number(edges[i][0]),Number(edges[i][1]))); 4159 if (kk != nbe) 4160 { 4161 cerr << " Some Double edge in the mesh, the number is " << kk-nbe << endl; 4162 MeshError(1002,this); 4163 } 4164 for (i=0;i<nbt;i++) 4165 for (int j=0;j<3;j++) 4166 { 4167 // Int4 i0,i1; 4168 Int4 k =edge4->addtrie(Number(triangles[i][VerticesOfTriangularEdge[j][0]]), 4169 Number(triangles[i][VerticesOfTriangularEdge[j][1]])); 4170 Int4 invisible = triangles[i].Hidden(j); 4171 if(st[k]==-1) 4172 st[k]=3*i+j; 4173 else if(st[k]>=0) { 4174 assert( ! triangles[i].TriangleAdj(j) && !triangles[st[k] / 3].TriangleAdj((int) (st[k]%3))); 4175 4176 triangles[i].SetAdj2(j,triangles + st[k] / 3,(int) (st[k]%3)); 4177 if (invisible) triangles[i].SetHidden(j); 4178 if (k<nbe) { 4179 triangles[i].SetLocked(j); 4180 } 4181 st[k]=-2-st[k]; } 4182 else { 4183 cerr << " The edge (" 4186 4184 << Number(triangles[i][VerticesOfTriangularEdge[j][0]]) 4187 4185 << " , " 4188 4186 << Number(triangles[i][VerticesOfTriangularEdge[j][1]]) 4189 4187 << " ) is in more than 2 triangles " <<k <<endl; 4190 cerr << " Edge " << j << " Of Triangle " << i << endl;4191 cerr << " Edge " << (-st[k]+2)%3 << " Of Triangle " << (-st[k]+2)/3 << endl;4192 cerr << " Edge " << triangles[(-st[k]+2)/3].NuEdgeTriangleAdj((int)((-st[k]+2)%3))4188 cerr << " Edge " << j << " Of Triangle " << i << endl; 4189 cerr << " Edge " << (-st[k]+2)%3 << " Of Triangle " << (-st[k]+2)/3 << endl; 4190 cerr << " Edge " << triangles[(-st[k]+2)/3].NuEdgeTriangleAdj((int)((-st[k]+2)%3)) 4193 4191 << " Of Triangle " << Number(triangles[(-st[k]+2)/3].TriangleAdj((int)((-st[k]+2)%3))) << endl; 4194 MeshError(9999,this);} 4195 4196 4197 } 4198 4199 if(verbosity>5) { 4200 cout << " On Mesh " << name << endl; 4201 cout << " - The number of Vertices = " << nbv << endl; 4202 cout << " - The number of Triangles = " << nbt << endl; 4203 cout << " - The number of given edge = " << nbe << endl; 4204 cout << " - The number of all edges = " << edge4->nb() << endl; 4205 cout << " - The Euler number = 1-Nb Of Hole = " << nbt-edge4->nb()+nbv << endl; 4206 } 4207 4208 4209 // check the consistant of edge[].adj and the geometrical required vertex 4210 Int4 k=0; 4211 for (i=0;i<edge4->nb();i++) 4212 if (st[i] >=0) // edge alone 4213 if (i < nbe) 4214 { 4215 Int4 i0=edge4->i(i);ordre[i0] = vertices+i0; 4216 Int4 i1=edge4->j(i);ordre[i1] = vertices+i1; 4192 MeshError(9999,this);} 4193 4194 4217 4195 } 4218 else { 4219 k++; 4220 if (verbosity>20 && k <20) 4221 { 4222 Int4 i0=edge4->i(i); 4223 Int4 i1=edge4->j(i); 4224 cerr << " Lose boundary edges " << i << " : " << i0 << " " << i1 << endl; 4225 } 4196 if(verbosity>5) { 4197 cout << " On Mesh " << name << endl; 4198 cout << " - The number of Vertices = " << nbv << endl; 4199 cout << " - The number of Triangles = " << nbt << endl; 4200 cout << " - The number of given edge = " << nbe << endl; 4201 cout << " - The number of all edges = " << edge4->nb() << endl; 4202 cout << " - The Euler number = 1-Nb Of Hole = " << nbt-edge4->nb()+nbv << endl; } 4203 4204 4205 // check the consistant of edge[].adj and the geometrical required vertex 4206 Int4 k=0; 4207 for (i=0;i<edge4->nb();i++) 4208 if (st[i] >=0) // edge alone 4209 if (i < nbe) 4210 { 4211 Int4 i0=edge4->i(i);ordre[i0] = vertices+i0; 4212 Int4 i1=edge4->j(i);ordre[i1] = vertices+i1; 4213 } 4214 else { 4215 k++; 4216 if (verbosity>20 && k <20) 4217 { 4218 Int4 i0=edge4->i(i); 4219 Int4 i1=edge4->j(i); 4220 cerr << " Lose boundary edges " << i << " : " << i0 << " " << i1 << endl; 4221 } 4222 } 4223 4224 if(k != 0) { 4225 if (verbosity>20) 4226 { 4227 cout << " The given edge are " << endl; 4228 for (int i=0;i< nbe;i++) 4229 cout << " Edge " << i << " : " << Number(edges[i][0]) << " " << Number(edges[i][1]) 4230 << " " << edges[i].ref << endl; 4231 } 4232 cerr << k << " boundary edges are not defined as edges " << endl; 4233 MeshError(9998,this); 4226 4234 } 4227 4228 if(k != 0) { 4229 if (verbosity>20) 4230 { 4231 cout << " The given edge are " << endl; 4232 for (int i=0;i< nbe;i++) 4233 cout << " Edge " << i << " : " << Number(edges[i][0]) << " " << Number(edges[i][1]) 4234 << " " << edges[i].ref << endl; 4235 } 4236 cerr << k << " boundary edges are not defined as edges " << endl; 4237 MeshError(9998,this); 4238 } 4239 printf("salope1\n"); 4240 4241 // generation of the mesh with boundary points 4242 Int4 nbvb = 0; 4243 for (i=0;i<nbv;i++) 4244 { 4245 vertices[i].t=0; 4246 vertices[i].vint=0; 4247 if (ordre[i]) 4248 ordre[nbvb++] = ordre[i]; 4249 } 4250 4251 Triangle *savetriangles= triangles; 4252 Int4 savenbt=nbt; 4253 Int4 savenbtx=nbtx; 4254 SubDomain * savesubdomains = subdomains; 4255 subdomains = 0; 4256 4257 Int4 Nbtriafillhole = 2*nbvb; 4258 Triangle * triafillhole =new Triangle[Nbtriafillhole]; 4259 if (verbosity>9) 4260 cout << " Nbtriafillhole triafillhole*" << triafillhole << endl; 4261 triangles = triafillhole; 4262 4263 nbt=2; 4264 nbtx= Nbtriafillhole; 4265 4266 for (i=2 ; det( ordre[0]->i, ordre[1]->i, ordre[i]->i ) == 0;){ 4267 if ( ++i >= nbvb) { 4268 cerr << "FillHoleInMesh: All the vertices are aline " << nbvb << endl; 4269 MeshError(998,this); 4270 } 4271 } 4272 Exchange( ordre[2], ordre[i]); 4273 4274 Vertex * v0=ordre[0], *v1=ordre[1]; 4275 4276 4277 triangles[0](0) = 0; // sommet pour infini 4278 triangles[0](1) = v0; 4279 triangles[0](2) = v1; 4280 4281 triangles[1](0) = 0;// sommet pour infini 4282 triangles[1](2) = v0; 4283 triangles[1](1) = v1; 4284 const int e0 = OppositeEdge[0]; 4285 const int e1 = NextEdge[e0]; 4286 const int e2 = PreviousEdge[e0]; 4287 triangles[0].SetAdj2(e0, &triangles[1] ,e0); 4288 triangles[0].SetAdj2(e1, &triangles[1] ,e2); 4289 triangles[0].SetAdj2(e2, &triangles[1] ,e1); 4290 4291 triangles[0].det = -1; // faux triangles 4292 triangles[1].det = -1; // faux triangles 4293 4294 triangles[0].SetTriangleContainingTheVertex(); 4295 triangles[1].SetTriangleContainingTheVertex(); 4296 4297 triangles[0].link=&triangles[1]; 4298 triangles[1].link=&triangles[0]; 4299 4300 printf("salope2\n"); 4301 // nbtf = 2; 4302 if ( !quadtree ) delete quadtree; // ->ReInitialise(); 4303 4304 quadtree = new QuadTree(this,0); 4305 quadtree->Add(*v0); 4306 quadtree->Add(*v1); 4307 4308 // on ajoute les sommets un a un 4309 Int4 NbSwap=0; 4310 for (Int4 icount=2; icount<nbvb; icount++) { 4311 4312 Vertex *vi = ordre[icount]; 4313 // cout << " Add vertex " << Number(vi) << endl; 4314 Icoor2 dete[3]; 4315 Triangle *tcvi = FindTriangleContening(vi->i,dete); 4316 quadtree->Add(*vi); 4317 Add(*vi,tcvi,dete); 4318 NbSwap += vi->Optim(1,1); 4319 4320 }// end loop on icount 4321 printf("salope3\n"); 4322 4323 //Int4 nbtfillhole = nbt; 4324 // inforce the boundary 4325 TriangleAdjacent ta(0,0); 4326 Int4 nbloss = 0,knbe=0; 4327 for ( i = 0; i < nbe; i++) 4328 if (st[i] >=0) // edge alone => on border ... FH oct 2009 4329 { 4330 Vertex & a=edges[i][0], & b = edges[i][1]; 4331 if (a.t && b.t) // le bug est la si maillage avec des bod non raffine 1. 4332 { 4333 knbe++; 4334 if (ForceEdge(a,b,ta)<0) 4335 nbloss++; 4336 } 4337 } 4338 if(nbloss) 4339 { 4340 cerr << " we loss some " << nbloss << " " << " edges other " << knbe << endl; 4341 MeshError(1100,this); 4342 } 4343 FindSubDomain(1); 4344 // remove all the hole 4345 // remove all the good sub domain 4346 Int4 krm =0; 4347 for (i=0;i<nbt;i++) 4348 if (triangles[i].link) // remove triangles 4349 { 4350 krm++; 4351 for (int j=0;j<3;j++) 4352 { 4353 TriangleAdjacent ta = triangles[i].Adj(j); 4354 Triangle & tta = * (Triangle *) ta; 4355 if(! tta.link) // edge between remove and not remove 4356 { // change the link of ta; 4357 int ja = ta; 4358 Vertex *v0= ta.EdgeVertex(0); 4359 Vertex *v1= ta.EdgeVertex(1); 4360 Int4 k =edge4->addtrie(v0?Number(v0):nbv,v1? Number(v1):nbv); 4361 assert(st[k] >=0); 4362 tta.SetAdj2(ja,savetriangles + st[k] / 3,(int) (st[k]%3)); 4363 ta.SetLock(); 4364 st[k]=-2-st[k]; 4365 } 4366 } 4367 } 4368 Int4 NbTfillHoll =0; 4369 for (i=0;i<nbt;i++) 4370 if (triangles[i].link) { 4371 triangles[i]=Triangle((Vertex *) NULL,(Vertex *) NULL,(Vertex *) NULL); 4372 triangles[i].color=-1; 4373 } 4374 else 4375 { 4376 triangles[i].color= savenbt+ NbTfillHoll++; 4377 } 4378 // cout << savenbt+NbTfillHoll << " " << savenbtx << endl; 4379 assert(savenbt+NbTfillHoll <= savenbtx ); 4380 // copy of the outside triangles in saveTriangles 4381 for (i=0;i<nbt;i++) 4382 if(triangles[i].color>=0) 4383 { 4384 savetriangles[savenbt]=triangles[i]; 4385 savetriangles[savenbt].link=0; 4386 savenbt++; 4387 } 4388 // gestion of the adj 4389 k =0; 4390 Triangle * tmax = triangles + nbt; 4391 for (i=0;i<savenbt;i++) 4392 { 4393 Triangle & ti = savetriangles[i]; 4394 for (int j=0;j<3;j++) 4395 { 4396 Triangle * ta = ti.TriangleAdj(j); 4397 int aa = ti.NuEdgeTriangleAdj(j); 4398 int lck = ti.Locked(j); 4399 if (!ta) k++; // bug 4400 else if ( ta >= triangles && ta < tmax) 4401 { 4402 ta= savetriangles + ta->color; 4403 ti.SetAdj2(j,ta,aa); 4404 if(lck) ti.SetLocked(j); 4405 } 4406 } 4407 } 4408 // OutSidesTriangles = triangles; 4409 // Int4 NbOutSidesTriangles = nbt; 4410 4411 // restore triangles; 4412 nbt=savenbt; 4413 nbtx=savenbtx; 4414 delete [] triangles; 4415 delete [] subdomains; 4416 triangles = savetriangles; 4417 subdomains = savesubdomains; 4418 // cout << triangles << " <> " << OutSidesTriangles << endl; 4419 /* k=0; 4420 for (i=0;i<nbt;i++) 4421 for (int j=0;j<3;j++) 4422 if (!triangles[i].TriangleAdj(j)) 4423 k++; 4424 4425 */ 4426 4427 printf("salope4\n"); 4428 if (k) { 4429 cerr << "Error Nb of triangles edge alone = " << k << endl; 4430 MeshError(9997,this); 4431 } 4432 printf("salope4a\n"); 4433 FindSubDomain(); 4434 // cout << " NbTOld = " << NbTold << " == " << nbt - NbOutT << " " << nbt << endl; 4435 4436 // 4437 4438 printf("salope4b\n"); 4439 delete edge4; 4440 printf("salope4ba\n"); 4441 delete [] st; 4442 printf("salope4bb %i\n",nbv); 4443 for (i=0;i<nbv;i++){ 4444 printf("%i\n",i); 4445 quadtree->Add(vertices[i]); 4446 } 4447 printf("salope4c\n"); 4448 4449 SetVertexFieldOn(); 4450 printf("salope5a\n"); 4451 4452 for (i=0;i<nbe;i++){ 4453 if(edges[i].on){ 4454 for(int j=0;j<2;j++){ 4455 if (!edges[i].adj[j]){ 4456 if(!edges[i][j].on->IsRequiredVertex()) { 4457 cerr << " Erreur adj et sommet requis edges [" << i << "][ " << j << "]= " 4458 << Number(edges[i][j]) << " : " << " on = " << Gh.Number(edges[i].on) ; 4459 if (edges[i][j].on->OnGeomVertex()) 4235 // generation of the mesh with boundary points 4236 Int4 nbvb = 0; 4237 for (i=0;i<nbv;i++) 4238 { 4239 vertices[i].t=0; 4240 vertices[i].vint=0; 4241 if (ordre[i]) 4242 ordre[nbvb++] = ordre[i]; 4243 } 4244 4245 Triangle *savetriangles= triangles; 4246 Int4 savenbt=nbt; 4247 Int4 savenbtx=nbtx; 4248 SubDomain * savesubdomains = subdomains; 4249 subdomains = 0; 4250 4251 Int4 Nbtriafillhole = 2*nbvb; 4252 Triangle * triafillhole =new Triangle[Nbtriafillhole]; 4253 if (verbosity>9) 4254 cout << " Nbtriafillhole triafillhole*" << triafillhole << endl; 4255 triangles = triafillhole; 4256 4257 nbt=2; 4258 nbtx= Nbtriafillhole; 4259 4260 for (i=2 ; det( ordre[0]->i, ordre[1]->i, ordre[i]->i ) == 0;) 4261 if ( ++i >= nbvb) { 4262 cerr << "FillHoleInMesh: All the vertices are aline " << nbvb << endl; 4263 MeshError(998,this); } 4264 Exchange( ordre[2], ordre[i]); 4265 4266 Vertex * v0=ordre[0], *v1=ordre[1]; 4267 4268 4269 triangles[0](0) = 0; // sommet pour infini 4270 triangles[0](1) = v0; 4271 triangles[0](2) = v1; 4272 4273 triangles[1](0) = 0;// sommet pour infini 4274 triangles[1](2) = v0; 4275 triangles[1](1) = v1; 4276 const int e0 = OppositeEdge[0]; 4277 const int e1 = NextEdge[e0]; 4278 const int e2 = PreviousEdge[e0]; 4279 triangles[0].SetAdj2(e0, &triangles[1] ,e0); 4280 triangles[0].SetAdj2(e1, &triangles[1] ,e2); 4281 triangles[0].SetAdj2(e2, &triangles[1] ,e1); 4282 4283 triangles[0].det = -1; // faux triangles 4284 triangles[1].det = -1; // faux triangles 4285 4286 triangles[0].SetTriangleContainingTheVertex(); 4287 triangles[1].SetTriangleContainingTheVertex(); 4288 4289 triangles[0].link=&triangles[1]; 4290 triangles[1].link=&triangles[0]; 4291 4292 #ifdef DEBUG 4293 triangles[0].check(); 4294 triangles[1].check(); 4295 #endif 4296 // nbtf = 2; 4297 if ( !quadtree ) 4298 delete quadtree; // ->ReInitialise(); 4299 4300 quadtree = new QuadTree(this,0); 4301 quadtree->Add(*v0); 4302 quadtree->Add(*v1); 4303 4304 // on ajoute les sommets un a un 4305 Int4 NbSwap=0; 4306 for (Int4 icount=2; icount<nbvb; icount++) { 4307 4308 Vertex *vi = ordre[icount]; 4309 // cout << " Add vertex " << Number(vi) << endl; 4310 Icoor2 dete[3]; 4311 Triangle *tcvi = FindTriangleContening(vi->i,dete); 4312 quadtree->Add(*vi); 4313 Add(*vi,tcvi,dete); 4314 NbSwap += vi->Optim(1,1); 4315 4316 #ifdef DRAWING2 4317 cout << Number(vi) << " " << NbSwap << endl; 4318 reffecran(); 4319 Draw(); 4320 vi->Draw(); 4321 inquire(); 4322 #endif 4323 }// end loop on icount 4324 #ifdef DRAWING1 4325 inquire(); 4326 #endif 4327 4328 //Int4 nbtfillhole = nbt; 4329 // inforce the boundary 4330 TriangleAdjacent ta(0,0); 4331 Int4 nbloss = 0,knbe=0; 4332 for ( i = 0; i < nbe; i++) 4333 if (st[i] >=0) // edge alone => on border ... FH oct 2009 4334 { 4335 Vertex & a=edges[i][0], & b = edges[i][1]; 4336 if (a.t && b.t) // le bug est la si maillage avec des bod non raffine 1. 4337 { 4338 knbe++; 4339 if (ForceEdge(a,b,ta)<0) 4340 nbloss++; 4341 } 4342 } 4343 if(nbloss) 4344 { 4345 cerr << " we loss some " << nbloss << " " << " edges other " << knbe << endl; 4346 MeshError(1100,this); 4347 } 4348 FindSubDomain(1); 4349 // remove all the hole 4350 // remove all the good sub domain 4351 Int4 krm =0; 4352 for (i=0;i<nbt;i++) 4353 if (triangles[i].link) // remove triangles 4354 { 4355 krm++; 4356 for (int j=0;j<3;j++) 4357 { 4358 TriangleAdjacent ta = triangles[i].Adj(j); 4359 Triangle & tta = * (Triangle *) ta; 4360 if(! tta.link) // edge between remove and not remove 4361 { // change the link of ta; 4362 int ja = ta; 4363 Vertex *v0= ta.EdgeVertex(0); 4364 Vertex *v1= ta.EdgeVertex(1); 4365 Int4 k =edge4->addtrie(v0?Number(v0):nbv,v1? Number(v1):nbv); 4366 assert(st[k] >=0); 4367 tta.SetAdj2(ja,savetriangles + st[k] / 3,(int) (st[k]%3)); 4368 ta.SetLock(); 4369 st[k]=-2-st[k]; 4370 } 4371 } 4372 } 4373 Int4 NbTfillHoll =0; 4374 for (i=0;i<nbt;i++) 4375 if (triangles[i].link) { 4376 triangles[i]=Triangle((Vertex *) NULL,(Vertex *) NULL,(Vertex *) NULL); 4377 triangles[i].color=-1; 4378 } 4379 else 4380 { 4381 triangles[i].color= savenbt+ NbTfillHoll++; 4382 #ifdef DEBUG 4383 triangles[i].check(); 4384 #endif 4385 } 4386 // cout << savenbt+NbTfillHoll << " " << savenbtx << endl; 4387 assert(savenbt+NbTfillHoll <= savenbtx ); 4388 // copy of the outside triangles in saveTriangles 4389 for (i=0;i<nbt;i++) 4390 if(triangles[i].color>=0) 4391 { 4392 savetriangles[savenbt]=triangles[i]; 4393 savetriangles[savenbt].link=0; 4394 savenbt++; 4395 } 4396 // gestion of the adj 4397 k =0; 4398 Triangle * tmax = triangles + nbt; 4399 for (i=0;i<savenbt;i++) 4400 { 4401 Triangle & ti = savetriangles[i]; 4402 for (int j=0;j<3;j++) 4403 { 4404 Triangle * ta = ti.TriangleAdj(j); 4405 int aa = ti.NuEdgeTriangleAdj(j); 4406 int lck = ti.Locked(j); 4407 if (!ta) k++; // bug 4408 else if ( ta >= triangles && ta < tmax) 4409 { 4410 ta= savetriangles + ta->color; 4411 ti.SetAdj2(j,ta,aa); 4412 if(lck) ti.SetLocked(j); 4413 } 4414 } 4415 } 4416 // OutSidesTriangles = triangles; 4417 // Int4 NbOutSidesTriangles = nbt; 4418 4419 // restore triangles; 4420 nbt=savenbt; 4421 nbtx=savenbtx; 4422 delete [] triangles; 4423 delete [] subdomains; 4424 triangles = savetriangles; 4425 subdomains = savesubdomains; 4426 // cout << triangles << " <> " << OutSidesTriangles << endl; 4427 /* k=0; 4428 for (i=0;i<nbt;i++) 4429 for (int j=0;j<3;j++) 4430 if (!triangles[i].TriangleAdj(j)) 4431 k++; 4432 */ 4433 if (k) { 4434 cerr << "Error Nb of triangles edge alone = " << k << endl; 4435 MeshError(9997,this); 4436 } 4437 FindSubDomain(); 4438 // cout << " NbTOld = " << NbTold << " == " << nbt - NbOutT << " " << nbt << endl; 4439 4440 // 4441 4442 delete edge4; 4443 delete [] st; 4444 for (i=0;i<nbv;i++) 4445 quadtree->Add(vertices[i]); 4446 4447 SetVertexFieldOn(); 4448 4449 for (i=0;i<nbe;i++) 4450 if(edges[i].on) 4451 for(int j=0;j<2;j++) 4452 if (!edges[i].adj[j]) 4453 if(!edges[i][j].on->IsRequiredVertex()) { 4454 cerr << " Erreur adj et sommet requis edges [" << i << "][ " << j << "]= " 4455 << Number(edges[i][j]) << " : " << " on = " << Gh.Number(edges[i].on) ; 4456 if (edges[i][j].on->OnGeomVertex()) 4460 4457 cerr << " vertex " << Gh.Number(edges[i][j].on->gv); 4461 else if (edges[i][j].on->OnGeomEdge())4458 else if (edges[i][j].on->OnGeomEdge()) 4462 4459 cerr << "Edges " << Gh.Number(edges[i][j].on->ge); 4463 else4460 else 4464 4461 cerr << " = " << edges[i][j].on ; 4465 cerr << endl; 4466 } 4467 } 4468 } 4469 } 4470 } 4471 4472 printf("salope5\n"); 4462 cerr << endl; 4463 } 4464 } 4473 4465 CurrentTh=OldCurrentTh; 4474 4466 }
Note:
See TracChangeset
for help on using the changeset viewer.