ForkMaster/netgen
1
0
Fork 0
mirror of https://github.com/NGSolve/netgen.git synced 2025-02-11 21:03:07 +05:00
Code Issues Packages Projects Releases Wiki Activity

more PointIndex

Browse Source
This commit is contained in:
Joachim Schoeberl 2024-12-27 10:51:58 +01:00
parent 386edbf75e
commit c466fe8d07
5 changed files with 82 additions and 78 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

92
libsrc/meshing/classifyhpel.hpp
View File

@ -1,6 +1,6 @@
HPREF_ELEMENT_TYPE ClassifyTet(HPRefElement & el, INDEX_2_HASHTABLE<int> & edges, INDEX_2_HASHTABLE<int> & edgepoint_dom, HPREF_ELEMENT_TYPE ClassifyTet(HPRefElement & el, INDEX_2_HASHTABLE<int> & edges, INDEX_2_HASHTABLE<int> & edgepoint_dom,
TBitArray<PointIndex> & cornerpoint, TBitArray<PointIndex> & edgepoint, INDEX_3_HASHTABLE<int> & faces, INDEX_2_HASHTABLE<int> & face_edges, TBitArray<PointIndex> & cornerpoint, TBitArray<PointIndex> & edgepoint, INDEX_3_HASHTABLE<int> & faces, INDEX_2_HASHTABLE<int> & face_edges,
INDEX_2_HASHTABLE<int> & surf_edges, NgArray<int, PointIndex::BASE> & facepoint) INDEX_2_HASHTABLE<int> & surf_edges, Array<int, PointIndex> & facepoint)
{ {
int ep1(0), ep2(0), ep3(0), ep4(0), cp1(0), cp2(0), cp3(0), cp4(0), fp1, fp2, fp3, fp4; int ep1(0), ep2(0), ep3(0), ep4(0), cp1(0), cp2(0), cp3(0), cp4(0), fp1, fp2, fp3, fp4;
int isedge1(0), isedge2(0), isedge3(0), isedge4(0), isedge5(0), isedge6(0); int isedge1(0), isedge2(0), isedge3(0), isedge4(0), isedge5(0), isedge6(0);
@ -56,12 +56,12 @@ HPREF_ELEMENT_TYPE ClassifyTet(HPRefElement & el, INDEX_2_HASHTABLE<int> & edges
cp3 = cornerpoint.Test (el.pnums[pi3]); cp3 = cornerpoint.Test (el.pnums[pi3]);
cp4 = cornerpoint.Test (el.pnums[pi4]); cp4 = cornerpoint.Test (el.pnums[pi4]);
isedge1 = edges.Used (INDEX_2::Sort (el.pnums[j], el.pnums[k])); isedge1 = edges.Used (PointIndices<2>::Sort (el.pnums[j], el.pnums[k]));
isedge2 = edges.Used (INDEX_2::Sort (el.pnums[j], el.pnums[pi3])); isedge2 = edges.Used (PointIndices<2>::Sort (el.pnums[j], el.pnums[pi3]));
isedge3 = edges.Used (INDEX_2::Sort (el.pnums[j], el.pnums[pi4])); isedge3 = edges.Used (PointIndices<2>::Sort (el.pnums[j], el.pnums[pi4]));
isedge4 = edges.Used (INDEX_2::Sort (el.pnums[k], el.pnums[pi3])); isedge4 = edges.Used (PointIndices<2>::Sort (el.pnums[k], el.pnums[pi3]));
isedge5 = edges.Used (INDEX_2::Sort (el.pnums[k], el.pnums[pi4])); isedge5 = edges.Used (PointIndices<2>::Sort (el.pnums[k], el.pnums[pi4]));
isedge6 = edges.Used (INDEX_2::Sort (el.pnums[pi3], el.pnums[pi4])); isedge6 = edges.Used (PointIndices<2>::Sort (el.pnums[pi3], el.pnums[pi4]));
if (debug) if (debug)
{ {
@ -76,13 +76,13 @@ HPREF_ELEMENT_TYPE ClassifyTet(HPRefElement & el, INDEX_2_HASHTABLE<int> & edges
for (int j = 0; j < 4; j++) isface[j] = false; for (int j = 0; j < 4; j++) isface[j] = false;
for (int l = 0; l < 4; l++) for (int l = 0; l < 4; l++)
{ {
INDEX_3 i3(0,0,0); PointIndices<3> i3(PointIndex::INVALID, PointIndex::INVALID, PointIndex::INVALID);
switch (l) switch (l)
{ {
case 0: i3.I1() = el.pnums[k]; i3.I2() = el.pnums[pi3]; i3.I3() = el.pnums[pi4]; break; case 0: i3[0] = el.pnums[k]; i3[1] = el.pnums[pi3]; i3[2] = el.pnums[pi4]; break;
case 1: i3.I1() = el.pnums[j]; i3.I2() = el.pnums[pi3]; i3.I3() = el.pnums[pi4]; break; case 1: i3[0] = el.pnums[j]; i3[1] = el.pnums[pi3]; i3[2] = el.pnums[pi4]; break;
case 2: i3.I1() = el.pnums[j]; i3.I2() = el.pnums[k]; i3.I3() = el.pnums[pi4]; break; case 2: i3[0] = el.pnums[j]; i3[1] = el.pnums[k]; i3[2] = el.pnums[pi4]; break;
case 3: i3.I1() = el.pnums[j]; i3.I2() = el.pnums[k]; i3.I3() = el.pnums[pi3]; break; case 3: i3[0] = el.pnums[j]; i3[1] = el.pnums[k]; i3[2] = el.pnums[pi3]; break;
} }
i3.Sort(); i3.Sort();
if (faces.Used (i3)) if (faces.Used (i3))
@ -102,15 +102,15 @@ HPREF_ELEMENT_TYPE ClassifyTet(HPRefElement & el, INDEX_2_HASHTABLE<int> & edges
isfedge1 = isfedge2 = isfedge3 = isfedge4 = isfedge5 = isfedge6 = 0; isfedge1 = isfedge2 = isfedge3 = isfedge4 = isfedge5 = isfedge6 = 0;
for (int l = 0; l < 6; l++) for (int l = 0; l < 6; l++)
{ {
INDEX_2 i2(0,0); PointIndices<2> i2(PointIndex::INVALID, PointIndex::INVALID);
switch (l) switch (l)
{ {
case 0: i2.I1() = el.pnums[j]; i2.I2() = el[k]; break; case 0: i2[0] = el.pnums[j]; i2[1] = el[k]; break;
case 1: i2.I1() = el.pnums[j]; i2.I2() = el.pnums[pi3]; break; case 1: i2[0] = el.pnums[j]; i2[1] = el.pnums[pi3]; break;
case 2: i2.I1() = el.pnums[j]; i2.I2() = el.pnums[pi4]; break; case 2: i2[0] = el.pnums[j]; i2[1] = el.pnums[pi4]; break;
case 3: i2.I1() = el.pnums[k]; i2.I2() = el.pnums[pi3]; break; case 3: i2[0] = el.pnums[k]; i2[1] = el.pnums[pi3]; break;
case 4: i2.I1() = el.pnums[k]; i2.I2() = el.pnums[pi4]; break; case 4: i2[0] = el.pnums[k]; i2[1] = el.pnums[pi4]; break;
case 5: i2.I1() = el.pnums[pi3]; i2.I2() = el.pnums[pi4]; break; case 5: i2[0] = el.pnums[pi3]; i2[1] = el.pnums[pi4]; break;
} }
i2.Sort(); i2.Sort();
if (face_edges.Used (i2)) if (face_edges.Used (i2))
@ -142,7 +142,7 @@ HPREF_ELEMENT_TYPE ClassifyTet(HPRefElement & el, INDEX_2_HASHTABLE<int> & edges
fp1 = fp2 = fp3 = fp4 = 0; fp1 = fp2 = fp3 = fp4 = 0;
for (int l = 0; l < 4; l++) for (int l = 0; l < 4; l++)
{ {
int pti(0); PointIndex pti = PointIndex::INVALID;
switch (l) switch (l)
{ {
case 0: pti = el.pnums[j]; break; case 0: pti = el.pnums[j]; break;
@ -548,7 +548,7 @@ HPREF_ELEMENT_TYPE ClassifyTet(HPRefElement & el, INDEX_2_HASHTABLE<int> & edges
if (type != HP_NONE) if (type != HP_NONE)
{ {
int pnums[4]; PointIndex pnums[4];
pnums[0] = el.pnums[j]; pnums[0] = el.pnums[j];
pnums[1] = el.pnums[k]; pnums[1] = el.pnums[k];
pnums[2] = el.pnums[pi3]; pnums[2] = el.pnums[pi3];
@ -588,7 +588,7 @@ HPREF_ELEMENT_TYPE ClassifyTet(HPRefElement & el, INDEX_2_HASHTABLE<int> & edges
HPREF_ELEMENT_TYPE ClassifyPrism(HPRefElement & el, INDEX_2_HASHTABLE<int> & edges, INDEX_2_HASHTABLE<int> & edgepoint_dom, HPREF_ELEMENT_TYPE ClassifyPrism(HPRefElement & el, INDEX_2_HASHTABLE<int> & edges, INDEX_2_HASHTABLE<int> & edgepoint_dom,
TBitArray<PointIndex> & cornerpoint, TBitArray<PointIndex> & edgepoint, INDEX_3_HASHTABLE<int> & faces, INDEX_2_HASHTABLE<int> & face_edges, TBitArray<PointIndex> & cornerpoint, TBitArray<PointIndex> & edgepoint, INDEX_3_HASHTABLE<int> & faces, INDEX_2_HASHTABLE<int> & face_edges,
INDEX_2_HASHTABLE<int> & surf_edges, NgArray<int, PointIndex::BASE> & facepoint) INDEX_2_HASHTABLE<int> & surf_edges, Array<int, PointIndex> & facepoint)
{ {
HPREF_ELEMENT_TYPE type = HP_NONE; HPREF_ELEMENT_TYPE type = HP_NONE;
@ -622,7 +622,7 @@ HPREF_ELEMENT_TYPE ClassifyPrism(HPRefElement & el, INDEX_2_HASHTABLE<int> & edg
const ELEMENT_EDGE * eledges = MeshTopology::GetEdges1 (PRISM); const ELEMENT_EDGE * eledges = MeshTopology::GetEdges1 (PRISM);
for(int k=0;k<9;k++) for(int k=0;k<9;k++)
{ {
INDEX_2 i2 = INDEX_2 :: Sort(el.PNum(p[eledges[k][0]-1]),el.PNum(p[eledges[k][1]-1])); PointIndices<2> i2 = PointIndices<2> :: Sort(el.PNum(p[eledges[k][0]-1]),el.PNum(p[eledges[k][1]-1]));
if (edges.Used(i2)) edge_sing[k] = 2; if (edges.Used(i2)) edge_sing[k] = 2;
else edge_sing[k] = face_edges.Used(i2); else edge_sing[k] = face_edges.Used(i2);
} }
@ -630,16 +630,17 @@ HPREF_ELEMENT_TYPE ClassifyPrism(HPRefElement & el, INDEX_2_HASHTABLE<int> & edg
const ELEMENT_FACE * elfaces = MeshTopology::GetFaces1 (PRISM); const ELEMENT_FACE * elfaces = MeshTopology::GetFaces1 (PRISM);
for (int k=0;k<5;k++) for (int k=0;k<5;k++)
{ {
INDEX_3 i3; PointIndices<3> i3;
if(k<2) if(k<2)
i3 = INDEX_3::Sort(el.pnums[p[elfaces[k][0]-1]-1], el.pnums[p[elfaces[k][1]-1]-1], i3 = PointIndices<3>::Sort(el.pnums[p[elfaces[k][0]-1]-1], el.pnums[p[elfaces[k][1]-1]-1],
el.pnums[p[elfaces[k][2]-1]-1]); el.pnums[p[elfaces[k][2]-1]-1]);
else else
{ {
INDEX_4 i4 = INDEX_4(el.pnums[p[elfaces[k][0]-1]-1], el.pnums[p[elfaces[k][1]-1]-1], el.pnums[p[elfaces[k][2]-1]-1],el.pnums[p[elfaces[k][3]-1]-1]); PointIndices<4> i4 (el.pnums[p[elfaces[k][0]-1]-1], el.pnums[p[elfaces[k][1]-1]-1],
el.pnums[p[elfaces[k][2]-1]-1],el.pnums[p[elfaces[k][3]-1]-1]);
i4.Sort(); i4.Sort();
i3 = INDEX_3(i4.I1(), i4.I2(), i4.I3()); i3 = PointIndices<3>(i4[0], i4[1], i4[2]);
} }
if (faces.Used (i3)) if (faces.Used (i3))
@ -809,7 +810,7 @@ HPREF_ELEMENT_TYPE ClassifyPrism(HPRefElement & el, INDEX_2_HASHTABLE<int> & edg
if(type != HP_NONE) if(type != HP_NONE)
{ {
int pnums[6]; PointIndex pnums[6];
for(int j=0;j<6;j++) pnums[j] = el.PNum (p[j]); for(int j=0;j<6;j++) pnums[j] = el.PNum (p[j]);
for(int k=0;k<6;k++) el.pnums[k] = pnums[k]; for(int k=0;k<6;k++) el.pnums[k] = pnums[k];
} }
@ -825,15 +826,15 @@ HPREF_ELEMENT_TYPE ClassifyPrism(HPRefElement & el, INDEX_2_HASHTABLE<int> & edg
// #ifdef SABINE // #ifdef SABINE
HPREF_ELEMENT_TYPE ClassifyTrig(HPRefElement & el, INDEX_2_HASHTABLE<int> & edges, INDEX_2_HASHTABLE<int> & edgepoint_dom, HPREF_ELEMENT_TYPE ClassifyTrig(HPRefElement & el, INDEX_2_HASHTABLE<int> & edges, INDEX_2_HASHTABLE<int> & edgepoint_dom,
TBitArray<PointIndex> & cornerpoint, TBitArray<PointIndex> & edgepoint, INDEX_3_HASHTABLE<int> & faces, INDEX_2_HASHTABLE<int> & face_edges, TBitArray<PointIndex> & cornerpoint, TBitArray<PointIndex> & edgepoint, INDEX_3_HASHTABLE<int> & faces, INDEX_2_HASHTABLE<int> & face_edges,
INDEX_2_HASHTABLE<int> & surf_edges, NgArray<int, PointIndex::BASE> & facepoint, int dim, const FaceDescriptor & fd) INDEX_2_HASHTABLE<int> & surf_edges, Array<int, PointIndex> & facepoint, int dim, const FaceDescriptor & fd)
{ {
HPREF_ELEMENT_TYPE type = HP_NONE; HPREF_ELEMENT_TYPE type = HP_NONE;
int pnums[3]; PointIndex pnums[3];
int p[3]; int p[3];
INDEX_3 i3 (el.pnums[0], el.pnums[1], el.pnums[2]); PointIndices<3> i3 (el.pnums[0], el.pnums[1], el.pnums[2]);
i3.Sort(); i3.Sort();
bool sing_face = faces.Used (i3); bool sing_face = faces.Used (i3);
@ -882,7 +883,7 @@ HPREF_ELEMENT_TYPE ClassifyTrig(HPRefElement & el, INDEX_2_HASHTABLE<int> & edge
{ {
int ep1=p[eledges[k][0]-1]; int ep1=p[eledges[k][0]-1];
int ep2=p[eledges[k][1]-1]; int ep2=p[eledges[k][1]-1];
INDEX_2 i2(el.PNum(ep1),el.PNum(ep2)); PointIndices<2> i2(el.PNum(ep1),el.PNum(ep2));
if(edges.Used(i2)) if(edges.Used(i2))
{ {
@ -925,7 +926,7 @@ HPREF_ELEMENT_TYPE ClassifyTrig(HPRefElement & el, INDEX_2_HASHTABLE<int> & edge
int ep1=p[eledges[k][0]-1]; int ep1=p[eledges[k][0]-1];
int ep2=p[eledges[k][1]-1]; int ep2=p[eledges[k][1]-1];
INDEX_2 i2 = INDEX_2::Sort(el.PNum(ep1),el.PNum(ep2)); PointIndices<2> i2 = PointIndices<2>::Sort(el.PNum(ep1),el.PNum(ep2));
if(edges.Used(i2)) if(edges.Used(i2))
{ {
@ -950,8 +951,8 @@ HPREF_ELEMENT_TYPE ClassifyTrig(HPRefElement & el, INDEX_2_HASHTABLE<int> & edge
(dim==3 || edgepoint_dom.Used(INDEX_2(fd.SurfNr(),pnums[k])) || edgepoint_dom.Used(INDEX_2(-1,pnums[k])))) (dim==3 || edgepoint_dom.Used(INDEX_2(fd.SurfNr(),pnums[k])) || edgepoint_dom.Used(INDEX_2(-1,pnums[k]))))
//edgepoint, but not member of sing_edge on trig -> cp //edgepoint, but not member of sing_edge on trig -> cp
{ {
INDEX_2 i2a=INDEX_2::Sort(el.PNum(p[k]), el.PNum(p[(k+1)%3])); PointIndices<2> i2a = PointIndices<2>::Sort(el.PNum(p[k]), el.PNum(p[(k+1)%3]));
INDEX_2 i2b=INDEX_2::Sort(el.PNum(p[k]), el.PNum(p[(k+2)%3])); PointIndices<2> i2b = PointIndices<2>::Sort(el.PNum(p[k]), el.PNum(p[(k+2)%3]));
if(!edges.Used(i2a) && !edges.Used(i2b)) if(!edges.Used(i2a) && !edges.Used(i2b))
point_sing[p[k]-1] = 3; point_sing[p[k]-1] = 3;
@ -1302,7 +1303,7 @@ HPREF_ELEMENT_TYPE ClassifyTrig(HPRefElement & el, INDEX_2_HASHTABLE<int> & edge
#endif #endif
HPREF_ELEMENT_TYPE ClassifyQuad(HPRefElement & el, INDEX_2_HASHTABLE<int> & edges, INDEX_2_HASHTABLE<int> & edgepoint_dom, HPREF_ELEMENT_TYPE ClassifyQuad(HPRefElement & el, INDEX_2_HASHTABLE<int> & edges, INDEX_2_HASHTABLE<int> & edgepoint_dom,
TBitArray<PointIndex> & cornerpoint, TBitArray<PointIndex> & edgepoint, INDEX_3_HASHTABLE<int> & faces, INDEX_2_HASHTABLE<int> & face_edges, TBitArray<PointIndex> & cornerpoint, TBitArray<PointIndex> & edgepoint, INDEX_3_HASHTABLE<int> & faces, INDEX_2_HASHTABLE<int> & face_edges,
INDEX_2_HASHTABLE<int> & surf_edges, NgArray<int, PointIndex::BASE> & facepoint, int dim, const FaceDescriptor & fd) INDEX_2_HASHTABLE<int> & surf_edges, Array<int, PointIndex> & facepoint, int dim, const FaceDescriptor & fd)
{ {
HPREF_ELEMENT_TYPE type = HP_NONE; HPREF_ELEMENT_TYPE type = HP_NONE;
@ -1320,10 +1321,10 @@ HPREF_ELEMENT_TYPE ClassifyQuad(HPRefElement & el, INDEX_2_HASHTABLE<int> & edge
if (dim == 2) if (dim == 2)
{ {
ep1 = edgepoint_dom.Used (INDEX_2 (el.GetIndex(), el.PNumMod(j))); ep1 = edgepoint_dom.Used (PointIndices<2>(el.GetIndex(), el.PNumMod(j)));
ep2 = edgepoint_dom.Used (INDEX_2 (el.GetIndex(), el.PNumMod(j+1))); ep2 = edgepoint_dom.Used (PointIndices<2>(el.GetIndex(), el.PNumMod(j+1)));
ep3 = edgepoint_dom.Used (INDEX_2 (el.GetIndex(), el.PNumMod(j+2))); ep3 = edgepoint_dom.Used (PointIndices<2>(el.GetIndex(), el.PNumMod(j+2)));
ep4 = edgepoint_dom.Used (INDEX_2 (el.GetIndex(), el.PNumMod(j+3))); ep4 = edgepoint_dom.Used (PointIndices<2>(el.GetIndex(), el.PNumMod(j+3)));
} }
cp1 = cornerpoint.Test (el.PNumMod (j)); cp1 = cornerpoint.Test (el.PNumMod (j));
@ -1351,8 +1352,7 @@ HPREF_ELEMENT_TYPE ClassifyQuad(HPRefElement & el, INDEX_2_HASHTABLE<int> & edge
if(dim ==3) if(dim ==3)
{ {
INDEX_2 i2; PointIndices<2> i2 = PointIndices<2>(el.PNumMod (j), el.PNumMod (j+1));
i2 = INDEX_2(el.PNumMod (j), el.PNumMod (j+1));
// i2.Sort(); // i2.Sort();
isedge1 = edges.Used (i2); isedge1 = edges.Used (i2);
i2.Sort(); i2.Sort();
@ -1652,7 +1652,7 @@ HPREF_ELEMENT_TYPE ClassifyQuad(HPRefElement & el, INDEX_2_HASHTABLE<int> & edge
HPREF_ELEMENT_TYPE ClassifyHex(HPRefElement & el, INDEX_2_HASHTABLE<int> & edges, INDEX_2_HASHTABLE<int> & edgepoint_dom, HPREF_ELEMENT_TYPE ClassifyHex(HPRefElement & el, INDEX_2_HASHTABLE<int> & edges, INDEX_2_HASHTABLE<int> & edgepoint_dom,
TBitArray<PointIndex> & cornerpoint, TBitArray<PointIndex> & edgepoint, INDEX_3_HASHTABLE<int> & faces, INDEX_2_HASHTABLE<int> & face_edges, TBitArray<PointIndex> & cornerpoint, TBitArray<PointIndex> & edgepoint, INDEX_3_HASHTABLE<int> & faces, INDEX_2_HASHTABLE<int> & face_edges,
INDEX_2_HASHTABLE<int> & surf_edges, NgArray<int, PointIndex::BASE> & facepoint) INDEX_2_HASHTABLE<int> & surf_edges, Array<int, PointIndex> & facepoint)
{ {
HPREF_ELEMENT_TYPE type = HP_NONE; HPREF_ELEMENT_TYPE type = HP_NONE;
@ -1758,7 +1758,7 @@ HPREF_ELEMENT_TYPE ClassifyHex(HPRefElement & el, INDEX_2_HASHTABLE<int> & edges
HPREF_ELEMENT_TYPE ClassifyHex7 (HPRefElement & el, INDEX_2_HASHTABLE<int> & edges, INDEX_2_HASHTABLE<int> & edgepoint_dom, HPREF_ELEMENT_TYPE ClassifyHex7 (HPRefElement & el, INDEX_2_HASHTABLE<int> & edges, INDEX_2_HASHTABLE<int> & edgepoint_dom,
TBitArray<PointIndex> & cornerpoint, TBitArray<PointIndex> & edgepoint, TBitArray<PointIndex> & cornerpoint, TBitArray<PointIndex> & edgepoint,
INDEX_3_HASHTABLE<int> & faces, INDEX_2_HASHTABLE<int> & face_edges, INDEX_3_HASHTABLE<int> & faces, INDEX_2_HASHTABLE<int> & face_edges,
INDEX_2_HASHTABLE<int> & surf_edges, NgArray<int, PointIndex::BASE> & facepoint) INDEX_2_HASHTABLE<int> & surf_edges, Array<int, PointIndex> & facepoint)
{ {
// HPREF_ELEMENT_TYPE type = HP_NONE; // HPREF_ELEMENT_TYPE type = HP_NONE;
@ -1798,7 +1798,7 @@ HPREF_ELEMENT_TYPE ClassifyHex7 (HPRefElement & el, INDEX_2_HASHTABLE<int> & edg
HPREF_ELEMENT_TYPE ClassifySegm(HPRefElement & hpel, INDEX_2_HASHTABLE<int> & edges, INDEX_2_HASHTABLE<int> & edgepoint_dom, HPREF_ELEMENT_TYPE ClassifySegm(HPRefElement & hpel, INDEX_2_HASHTABLE<int> & edges, INDEX_2_HASHTABLE<int> & edgepoint_dom,
TBitArray<PointIndex> & cornerpoint, TBitArray<PointIndex> & edgepoint, TBitArray<PointIndex> & cornerpoint, TBitArray<PointIndex> & edgepoint,
INDEX_3_HASHTABLE<int> & faces, INDEX_2_HASHTABLE<int> & face_edges, INDEX_3_HASHTABLE<int> & faces, INDEX_2_HASHTABLE<int> & face_edges,
INDEX_2_HASHTABLE<int> & surf_edges, NgArray<int, PointIndex::BASE> & facepoint) INDEX_2_HASHTABLE<int> & surf_edges, Array<int, PointIndex> & facepoint)
{ {
int cp1 = cornerpoint.Test (hpel[0]); int cp1 = cornerpoint.Test (hpel[0]);
@ -1842,7 +1842,7 @@ HPREF_ELEMENT_TYPE ClassifySegm(HPRefElement & hpel, INDEX_2_HASHTABLE<int> & ed
HPREF_ELEMENT_TYPE ClassifyPyramid(HPRefElement & el, INDEX_2_HASHTABLE<int> & edges, INDEX_2_HASHTABLE<int> & edgepoint_dom, HPREF_ELEMENT_TYPE ClassifyPyramid(HPRefElement & el, INDEX_2_HASHTABLE<int> & edges, INDEX_2_HASHTABLE<int> & edgepoint_dom,
TBitArray<PointIndex> & cornerpoint, TBitArray<PointIndex> & edgepoint, TBitArray<PointIndex> & cornerpoint, TBitArray<PointIndex> & edgepoint,
INDEX_3_HASHTABLE<int> & faces, INDEX_2_HASHTABLE<int> & face_edges, INDEX_3_HASHTABLE<int> & faces, INDEX_2_HASHTABLE<int> & face_edges,
INDEX_2_HASHTABLE<int> & surf_edges, NgArray<int, PointIndex::BASE> & facepoint) INDEX_2_HASHTABLE<int> & surf_edges, Array<int, PointIndex> & facepoint)
{ {
// *testout << "classify pyramid, pnums = "; // *testout << "classify pyramid, pnums = ";
// for (int i = 0; i < 5; i++) *testout << el.pnums[i] << " "; // for (int i = 0; i < 5; i++) *testout << el.pnums[i] << " ";

47
libsrc/meshing/delaunay2d.cpp
View File

@ -38,7 +38,7 @@ namespace netgen
if(p1<p0) if(p1<p0)
Swap(p0,p1); Swap(p0,p1);
IVec<2> hash = {p0,p1}; PointIndices<2> hash = {p0,p1};
auto pos = edge_to_trig.Position(hash); auto pos = edge_to_trig.Position(hash);
if (pos == -1) return -1; if (pos == -1) return -1;
@ -53,7 +53,7 @@ namespace netgen
if(p1<p0) if(p1<p0)
Swap(p0,p1); Swap(p0,p1);
IVec<2> hash = {p0,p1}; PointIndices<2> hash = {p0,p1};
auto pos = edge_to_trig.Position(hash); auto pos = edge_to_trig.Position(hash);
if (pos == -1) if (pos == -1)
edge_to_trig[hash] = {eli, -1}; edge_to_trig[hash] = {eli, -1};
@ -80,7 +80,7 @@ namespace netgen
if(p1<p0) if(p1<p0)
Swap(p0,p1); Swap(p0,p1);
IVec<2> hash = {p0,p1}; PointIndices<2> hash = {p0,p1};
auto pos = edge_to_trig.Position(hash); auto pos = edge_to_trig.Position(hash);
auto i2 = edge_to_trig.GetData(pos); auto i2 = edge_to_trig.GetData(pos);
@ -93,7 +93,7 @@ namespace netgen
} }
void DelaunayMesh::AppendTrig( int pi0, int pi1, int pi2 ) void DelaunayMesh::AppendTrig( PointIndex pi0, PointIndex pi1, PointIndex pi2 )
{ {
DelaunayTrig el; DelaunayTrig el;
el[0] = pi0; el[0] = pi0;
@ -174,7 +174,8 @@ namespace netgen
{ {
const auto trig = trigs[i_trig]; const auto trig = trigs[i_trig];
if(trig[0]==-1) // if(trig[0]==-1)
if(!trig[0].IsValid())
continue; continue;
double rad2 = trig.Radius2(); double rad2 = trig.Radius2();
@ -254,9 +255,9 @@ namespace netgen
const DelaunayTrig & trig = trigs[j]; const DelaunayTrig & trig = trigs[j];
for (int k = 0; k < 3; k++) for (int k = 0; k < 3; k++)
{ {
int p1 = trig[k]; PointIndex p1 = trig[k];
int p2 = trig[(k+1)%3]; PointIndex p2 = trig[(k+1)%3];
IVec<2> edge{p1,p2}; PointIndices<2> edge{p1,p2};
edge.Sort(); edge.Sort();
bool found = false; bool found = false;
for (int l = 0; l < edges.Size(); l++) for (int l = 0; l < edges.Size(); l++)
@ -309,9 +310,9 @@ namespace netgen
for (int j : intersecting) for (int j : intersecting)
{ {
UnsetNeighbours(j); UnsetNeighbours(j);
trigs[j][0] = -1; trigs[j][0] = PointIndex::INVALID;
trigs[j][1] = -1; trigs[j][1] = PointIndex::INVALID;
trigs[j][2] = -1; trigs[j][2] = PointIndex::INVALID;
} }
for (auto edge : edges) for (auto edge : edges)
@ -331,7 +332,8 @@ namespace netgen
for (DelaunayTrig & trig : trigs) for (DelaunayTrig & trig : trigs)
{ {
if (trig[0] < 0) continue; // if (trig[0] < 0) continue;
if (!trig[0].IsValid()) continue;
Vec<3> n = Cross (P3(points[trig[1]])-P3(points[trig[0]]), Vec<3> n = Cross (P3(points[trig[1]])-P3(points[trig[0]]),
P3(points[trig[2]])-P3(points[trig[0]])); P3(points[trig[2]])-P3(points[trig[0]]));
@ -587,7 +589,7 @@ namespace netgen
t2.Start(); t2.Start();
Array<PointIndex> old_points; Array<PointIndex> old_points;
BitArray add_point(mesh.Points().Size()+1); TBitArray<PointIndex> add_point(mesh.Points().Size()+1);
Array<PointIndex> addpoints; Array<PointIndex> addpoints;
add_point.Clear(); add_point.Clear();
/* /*
@ -708,7 +710,8 @@ namespace netgen
for (auto & trig : dmesh.GetElements()) for (auto & trig : dmesh.GetElements())
{ {
if (trig[0] < 0) continue; // if (trig[0] < 0) continue;
if (!trig[0].IsValid()) continue;
Element2d el(trig[0], trig[1], trig[2]); Element2d el(trig[0], trig[1], trig[2]);
el.SetIndex (1); el.SetIndex (1);
@ -719,7 +722,7 @@ namespace netgen
while(!conforming) while(!conforming)
{ {
conforming = true; conforming = true;
BitArray marked_points(tempmesh.Points().Size()+1); TBitArray<PointIndex> marked_points(tempmesh.Points().Size()+1);
marked_points = false; marked_points = false;
// Check for trigs cutting a boundary edge (non-conforming mesh) // Check for trigs cutting a boundary edge (non-conforming mesh)
auto point_to_trigs = tempmesh.CreatePoint2SurfaceElementTable( 0 ); auto point_to_trigs = tempmesh.CreatePoint2SurfaceElementTable( 0 );
@ -776,7 +779,7 @@ namespace netgen
auto & el0 = tempmesh[cutting_trigs[0]]; auto & el0 = tempmesh[cutting_trigs[0]];
auto & el1 = tempmesh[cutting_trigs[1]]; auto & el1 = tempmesh[cutting_trigs[1]];
pi1 = el1[0]+el1[1]+el1[2] - pi2-pi3; pi1 = el1[0]-pi2+el1[1]-pi3+el1[2];
if(marked_points.Test(pi1)) continue; if(marked_points.Test(pi1)) continue;
@ -801,13 +804,13 @@ namespace netgen
// Mark edges and trigs as inside or outside, starting with boundary edges // Mark edges and trigs as inside or outside, starting with boundary edges
enum POSITION { UNKNOWN, BOUNDARY, INSIDE, OUTSIDE }; enum POSITION { UNKNOWN, BOUNDARY, INSIDE, OUTSIDE };
Array<POSITION, SurfaceElementIndex> trig_pos(tempmesh.SurfaceElements().Size()); Array<POSITION, SurfaceElementIndex> trig_pos(tempmesh.SurfaceElements().Size());
ngcore::ClosedHashTable<IVec<2>, POSITION> edge_pos(3*tempmesh.SurfaceElements().Size()); ngcore::ClosedHashTable<PointIndices<2>, POSITION> edge_pos(3*tempmesh.SurfaceElements().Size());
trig_pos = UNKNOWN; trig_pos = UNKNOWN;
for (auto & seg : tempmesh.LineSegments()) for (auto & seg : tempmesh.LineSegments())
{ {
ArrayMem<SurfaceElementIndex, 2> els; ArrayMem<SurfaceElementIndex, 2> els;
IVec<2> edge{seg[0], seg[1]}; PointIndices<2> edge{seg[0], seg[1]};
edge.Sort(); edge.Sort();
edge_pos[edge] = BOUNDARY; edge_pos[edge] = BOUNDARY;
@ -828,8 +831,8 @@ namespace netgen
else else
pos = OUTSIDE; pos = OUTSIDE;
IVec<2> e1{seg[0], pi2}; PointIndices<2> e1{seg[0], pi2};
IVec<2> e2{seg[1], pi2}; PointIndices<2> e2{seg[1], pi2};
e1.Sort(); e1.Sort();
e2.Sort(); e2.Sort();
if(!edge_pos.Used(e1)) if(!edge_pos.Used(e1))
@ -857,7 +860,7 @@ namespace netgen
// any edge of unknown trig already marked? // any edge of unknown trig already marked?
for(auto i : IntRange(3)) for(auto i : IntRange(3))
{ {
IVec<2> edge{el[(i+1)%3], el[(i+2)%3]}; PointIndices<2> edge{el[(i+1)%3], el[(i+2)%3]};
edge.Sort(); edge.Sort();
if(edge_pos.Used(edge) && edge_pos[edge]!=BOUNDARY) if(edge_pos.Used(edge) && edge_pos[edge]!=BOUNDARY)
{ {
@ -871,7 +874,7 @@ namespace netgen
if(trig_pos[sei] != UNKNOWN) if(trig_pos[sei] != UNKNOWN)
for(auto i : IntRange(3)) for(auto i : IntRange(3))
{ {
IVec<2> edge{el[(i+1)%3], el[(i+2)%3]}; PointIndices<2> edge{el[(i+1)%3], el[(i+2)%3]};
edge.Sort(); edge.Sort();
if(!edge_pos.Used(edge) || edge_pos[edge]==BOUNDARY) if(!edge_pos.Used(edge) || edge_pos[edge]==BOUNDARY)
edge_pos[edge] = trig_pos[sei]; edge_pos[edge] = trig_pos[sei];

10
libsrc/meshing/delaunay2d.hpp
View File

@ -22,7 +22,7 @@ namespace netgen
double r; double r;
double rad2; double rad2;
DelaunayTrig () = default; DelaunayTrig () = default;
DelaunayTrig (int p1, int p2, int p3) DelaunayTrig (PointIndex p1, PointIndex p2, PointIndex p3)
{ {
pnums[0] = p1; pnums[0] = p1;
pnums[1] = p2; pnums[1] = p2;
@ -38,19 +38,19 @@ namespace netgen
double Radius2() const { return rad2; } double Radius2() const { return rad2; }
Box<2> BoundingBox() const { return Box<2> (c-Vec<2>(r,r), c+Vec<2>(r,r)); } Box<2> BoundingBox() const { return Box<2> (c-Vec<2>(r,r), c+Vec<2>(r,r)); }
mutable PointIndex visited_pi = -1; mutable PointIndex visited_pi = PointIndex::INVALID; // -1;
}; };
class DelaunayMesh class DelaunayMesh
{ {
ngcore::ClosedHashTable<IVec<2>, IVec<2>> edge_to_trig; ngcore::ClosedHashTable<PointIndices<2>, IVec<2>> edge_to_trig;
Array<DelaunayTrig> trigs; Array<DelaunayTrig> trigs;
unique_ptr<DelaunayTree<2>> tree; unique_ptr<DelaunayTree<2>> tree;
Array<Point<2>, PointIndex> & points; Array<Point<2>, PointIndex> & points;
Array<int> closeels; Array<int> closeels;
Array<int> intersecting; Array<int> intersecting;
Array<IVec<2>> edges; Array<PointIndices<2>> edges;
int GetNeighbour( int eli, int edge ); int GetNeighbour( int eli, int edge );
@ -58,7 +58,7 @@ namespace netgen
void UnsetNeighbours( int eli ); void UnsetNeighbours( int eli );
void AppendTrig( int pi0, int pi1, int pi2 ); void AppendTrig( PointIndex pi0, PointIndex pi1, PointIndex pi2 );
public: public:
DelaunayMesh( Array<Point<2>, PointIndex> & points_, Box<2> box ); DelaunayMesh( Array<Point<2>, PointIndex> & points_, Box<2> box );

6
libsrc/meshing/hprefinement.cpp
View File

@ -605,7 +605,7 @@ namespace netgen
bool CheckSingularities(Mesh & mesh, INDEX_2_HASHTABLE<int> & edges, INDEX_2_HASHTABLE<int> & edgepoiclt_dom, bool CheckSingularities(Mesh & mesh, INDEX_2_HASHTABLE<int> & edges, INDEX_2_HASHTABLE<int> & edgepoiclt_dom,
NgBitArray & cornerpoint, NgBitArray & edgepoint, INDEX_3_HASHTABLE<int> & faces, INDEX_2_HASHTABLE<int> & face_edges, NgBitArray & cornerpoint, NgBitArray & edgepoint, INDEX_3_HASHTABLE<int> & faces, INDEX_2_HASHTABLE<int> & face_edges,
INDEX_2_HASHTABLE<int> & surf_edges, NgArray<int, PointIndex::BASE> & facepoint, int & levels, int & act_ref); INDEX_2_HASHTABLE<int> & surf_edges, Array<int, PointIndex> & facepoint, int & levels, int & act_ref);
bool ClassifyHPElements (Mesh & mesh, NgArray<HPRefElement> & elements, SplittingType split, int & act_ref, int & levels); bool ClassifyHPElements (Mesh & mesh, NgArray<HPRefElement> & elements, SplittingType split, int & act_ref, int & levels);
@ -1634,7 +1634,7 @@ namespace netgen
bool CheckSingularities(Mesh & mesh, INDEX_2_HASHTABLE<int> & edges, INDEX_2_HASHTABLE<int> & edgepoint_dom, bool CheckSingularities(Mesh & mesh, INDEX_2_HASHTABLE<int> & edges, INDEX_2_HASHTABLE<int> & edgepoint_dom,
TBitArray<PointIndex> & cornerpoint, TBitArray<PointIndex> & edgepoint, INDEX_3_HASHTABLE<int> & faces, INDEX_2_HASHTABLE<int> & face_edges, TBitArray<PointIndex> & cornerpoint, TBitArray<PointIndex> & edgepoint, INDEX_3_HASHTABLE<int> & faces, INDEX_2_HASHTABLE<int> & face_edges,
INDEX_2_HASHTABLE<int> & surf_edges, NgArray<int, PointIndex::BASE> & facepoint, int & levels, int & act_ref) INDEX_2_HASHTABLE<int> & surf_edges, Array<int, PointIndex> & facepoint, int & levels, int & act_ref)
{ {
bool sing = 0; bool sing = 0;
if (mesh.GetDimension() == 3) if (mesh.GetDimension() == 3)
@ -1902,7 +1902,7 @@ namespace netgen
INDEX_3_HASHTABLE<int> faces(mesh.GetNSE()+1); INDEX_3_HASHTABLE<int> faces(mesh.GetNSE()+1);
INDEX_2_HASHTABLE<int> face_edges(mesh.GetNSE()+1); INDEX_2_HASHTABLE<int> face_edges(mesh.GetNSE()+1);
INDEX_2_HASHTABLE<int> surf_edges(mesh.GetNSE()+1); INDEX_2_HASHTABLE<int> surf_edges(mesh.GetNSE()+1);
NgArray<int, PointIndex::BASE> facepoint(mesh.GetNP()); Array<int, PointIndex> facepoint(mesh.GetNP());
bool sing = CheckSingularities(mesh, edges, edgepoint_dom, bool sing = CheckSingularities(mesh, edges, edgepoint_dom,
cornerpoint, edgepoint, faces, face_edges, cornerpoint, edgepoint, faces, face_edges,

1
libsrc/meshing/meshtype.hpp
View File

@ -257,6 +257,7 @@ namespace netgen
void DoArchive (Archive & ar) { ar & i; } void DoArchive (Archive & ar) { ar & i; }
}; };
/* /*
inline PointIndex operator+ (PointIndex pi, int i) { return PointIndex(pi.i+i); } inline PointIndex operator+ (PointIndex pi, int i) { return PointIndex(pi.i+i); }
inline PointIndex operator+ (PointIndex pi, size_t i) { return PointIndex(pi.i+i); } inline PointIndex operator+ (PointIndex pi, size_t i) { return PointIndex(pi.i+i); }