more PointIndex

This commit is contained in:
Joachim Schoeberl 2024-12-26 16:32:50 +01:00
parent 9c9b4ea880
commit a675c42d89
8 changed files with 111 additions and 44 deletions

View File

@ -329,8 +329,8 @@ namespace ngcore
case 1:
{
size_t oldval = nd;
while (blocknr+1>nd) {
nd.compare_exchange_weak (oldval, blocknr+1);
while (blocknr-IndexBASE<IndexType>()+1>nd) {
nd.compare_exchange_weak (oldval, blocknr-IndexBASE<IndexType>()+1);
oldval = nd;
}
break;

View File

@ -115,14 +115,19 @@ class INDEX_2
public:
///
INDEX_2 () { }
INDEX_2 (const INDEX_2&) = default;
INDEX_2 (INDEX_2&&) = default;
INDEX_2 & operator= (const INDEX_2&) = default;
INDEX_2 & operator= (INDEX_2&&) = default;
///
constexpr INDEX_2 (INDEX ai1, INDEX ai2)
: i{ai1, ai2} { }
// { i[0] = ai1; i[1] = ai2; }
///
constexpr INDEX_2 (const INDEX_2 & in2)
: i{in2.i[0], in2.i[1]} { }
// constexpr INDEX_2 (const INDEX_2 & in2)
// : i{in2.i[0], in2.i[1]} { }
// { i[0] = in2.i[0]; i[1] = in2.i[1]; }
@ -206,12 +211,14 @@ public:
///
INDEX_3 () { }
///
INDEX_3 (INDEX ai1, INDEX ai2, INDEX ai3)
{ i[0] = ai1; i[1] = ai2; i[2] = ai3; }
constexpr INDEX_3 (INDEX ai1, INDEX ai2, INDEX ai3)
: i{ai1, ai2, ai3} { }
// { i[0] = ai1; i[1] = ai2; i[2] = ai3; }
///
INDEX_3 (const INDEX_3 & in2)
{ i[0] = in2.i[0]; i[1] = in2.i[1]; i[2] = in2.i[2]; }
constexpr INDEX_3 (const INDEX_3 & in2)
: i{in2.i[0], in2.i[1], in2.i[2]} { }
// { i[0] = in2.i[0]; i[1] = in2.i[1]; i[2] = in2.i[2]; }
static INDEX_3 Sort (INDEX_3 i3)
@ -479,5 +486,12 @@ namespace netgen
{
return HashValue2(IVec<2,netgen::INDEX>(ind[0], ind[1]), mask);
}
constexpr inline size_t HashValue2 (const netgen::INDEX_3 & ind, size_t mask)
{
return HashValue2(IVec<3,netgen::INDEX>(ind[0], ind[1], ind[2]), mask);
}
}
#endif

View File

@ -698,7 +698,7 @@ namespace netgen
// int i, j;
int cnt = 0;
int found;
bool found;
double len2, maxlen2;
INDEX_2 ep;
@ -707,7 +707,7 @@ namespace netgen
do
{
found = 0;
found = false;
maxlen2 = 1e30;
// for (int i = 1; i <= mesh.GetNE(); i++)
@ -775,7 +775,7 @@ namespace netgen
{
maxlen2 = len2;
ep = i2;
found = 1;
found = true;
}
}
}
@ -833,8 +833,8 @@ namespace netgen
e1.Sort();
e2.Sort();
int used1 = edgenumber.Used (e1);
int used2 = edgenumber.Used (e2);
bool used1 = edgenumber.Used (e1);
bool used2 = edgenumber.Used (e2);
if (used1 && !used2)
{
@ -2906,7 +2906,8 @@ namespace netgen
// INDEX_2_HASHTABLE<int> cutedges(10 + 5 * (mtets.Size()+mprisms.Size()+mtris.Size()+mquads.Size()));
// INDEX_2_CLOSED_HASHTABLE<PointIndex> cutedges(10 + 9 * (mtets.Size()+mprisms.Size()+mtris.Size()+mquads.Size()));
ClosedHashTable<INDEX_2, PointIndex> cutedges(10 + 9 * (mtets.Size()+mprisms.Size()+mtris.Size()+mquads.Size()));
// ClosedHashTable<INDEX_2, PointIndex> cutedges(10 + 9 * (mtets.Size()+mprisms.Size()+mtris.Size()+mquads.Size()));
ClosedHashTable<PointIndices<2>, PointIndex> cutedges(10 + 9 * (mtets.Size()+mprisms.Size()+mtris.Size()+mquads.Size()));
bool noprojection = false;
NgProfiler::StopTimer (timer1a);

View File

@ -851,16 +851,17 @@ namespace netgen
tets_with_3_bnd_points.SetSize(cnt);
static Timer t1("Build face table"); t1.Start();
ngcore::ClosedHashTable< ngcore::IVec<3>, int > face_table( 4*cnt + 3 );
for(auto ei : tets_with_3_bnd_points)
for(auto j : Range(4))
// ngcore::ClosedHashTable< ngcore::IVec<3>, int > face_table( 4*cnt + 3 );
ngcore::ClosedHashTable< PointIndices<3>, int > face_table( 4*cnt + 3 );
for (auto ei : tets_with_3_bnd_points)
for (auto j : Range(4))
{
PointIndices<3> i3_ = tempels[ei].GetFace (j);
ngcore::IVec<3> i3 = {i3_[0], i3_[1], i3_[2]};
if(bnd_points[i3[0]] && bnd_points[i3[1]] && bnd_points[i3[2]])
PointIndices<3> i3 = tempels[ei].GetFace (j);
// ngcore::IVec<3> i3 = {i3_[0], i3_[1], i3_[2]};
if(bnd_points[i3[0]] && bnd_points[i3[1]] && bnd_points[i3[2]])
{
i3.Sort();
face_table.Set( i3, true );
i3.Sort();
face_table.Set( i3, true );
}
}
t1.Stop();
@ -871,7 +872,8 @@ namespace netgen
for (int i = 1; i <= mesh.GetNOpenElements(); i++)
{
const Element2d & tri = mesh.OpenElement(i);
ngcore::IVec<3,PointIndex> i3(tri[0], tri[1], tri[2]);
// ngcore::IVec<3,PointIndex> i3(tri[0], tri[1], tri[2]);
PointIndices<3> i3(tri[0], tri[1], tri[2]);
i3.Sort();
if(!face_table.Used(i3))
openels.Append(i);
@ -1016,7 +1018,7 @@ namespace netgen
for (int j = 0; j < 4; j++)
{
pp[j] = &mesh.Point(el[j]);
tetpi[j] = el[j];
tetpi[j] = el[j]-IndexBASE<PointIndex>()+1;
}
Point3d tetpmin(*pp[0]);
@ -1045,7 +1047,7 @@ namespace netgen
for (int k = 1; k <= 3; k++)
{
tripp[k-1] = &mesh.Point (tri.PNum(k));
tripi[k-1] = tri.PNum(k);
tripi[k-1] = tri.PNum(k)-IndexBASE<PointIndex>()+1;
}
if (IntersectTetTriangle (&pp[0], &tripp[0], tetpi, tripi))
@ -1133,13 +1135,10 @@ namespace netgen
for (auto & el : tempels)
for (int j = 0; j < 4; j++)
el.NB(j) = 0;
TABLE<int,PointIndex::BASE> elsonpoint(mesh.GetNP());
/*
for (int i = 0; i < tempels.Size(); i++)
{
const DelaunayTet & el = tempels[i];
*/
TABLE<int,PointIndex::BASE> elsonpoint(mesh.GetNP());
for (const DelaunayTet & el : tempels)
{
PointIndices<4> i4(el[0], el[1], el[2], el[3]);
@ -1158,7 +1157,25 @@ namespace netgen
elsonpoint.Add (i4.I1(), i+1);
elsonpoint.Add (i4.I2(), i+1);
}
*/
TableCreator<int, PointIndex> creator(mesh.GetNP());
while (!creator.Done())
{
for (int i = 0; i < tempels.Size(); i++)
{
const DelaunayTet & el = tempels[i];
PointIndices<4> i4(el[0], el[1], el[2], el[3]);
i4.Sort();
creator.Add (i4[0], i+1);
creator.Add (i4[1], i+1);
}
creator++;
}
auto elsonpoint = creator.MoveTable();
// cout << "elsonpoint mem: ";
// elsonpoint.PrintMemInfo(cout);

View File

@ -713,7 +713,7 @@ namespace netgen
if (mesh.IsSegment (pi1, pi2)) continue;
INDEX_2 ii2 (pi1, pi2);
PointIndices<2> ii2 (pi1, pi2);
ii2.Sort();
if (els_on_edge.Used (ii2))
{
@ -739,7 +739,7 @@ namespace netgen
if (mesh.LegalTrig(sel)) continue;
// find longest edge
INDEX_2 edge;
PointIndices<2> edge;
double edge_len = 0;
PointIndex pi1, pi2, pi3, pi4;
PointGeomInfo gi1, gi2, gi3, gi4;

View File

@ -452,8 +452,10 @@ namespace netgen
lock.Lock();
timestamp = NextTimeStamp();
int maxn = max2 (s[0], s[1]);
maxn += 1-PointIndex::BASE;
// int maxn = max2 (s[0], s[1]);
// maxn += 1-PointIndex::BASE;
int maxn = max2 (s[0]-IndexBASE<PointIndex>()+1,
s[1]-IndexBASE<PointIndex>()+1);
/*
if (maxn > ptyps.Size())
@ -540,12 +542,19 @@ namespace netgen
void Mesh :: SetSurfaceElement (SurfaceElementIndex sei, const Element2d & el)
{
/*
int maxn = el[0];
for (int i = 1; i < el.GetNP(); i++)
if (el[i] > maxn) maxn = el[i];
maxn += 1-PointIndex::BASE;
*/
PointIndex maxpi = el[0];
for (int i = 1; i < el.GetNP(); i++)
if (el[i] > maxpi) maxpi = el[i];
int maxn = maxpi-IndexBASE<PointIndex>()+1;
if (maxn <= points.Size())
{
for (int i = 0; i < el.GetNP(); i++)
@ -843,9 +852,12 @@ namespace netgen
outfile.setf (ios::fixed, ios::floatfield);
outfile.setf (ios::showpoint);
/*
PointIndex pi;
for (pi = PointIndex::BASE;
pi < GetNP()+PointIndex::BASE; pi++)
*/
for (PointIndex pi : (*this).Points().Range())
{
outfile.width(22);
outfile << (*this)[pi](0)/scale << " ";
@ -1743,7 +1755,7 @@ namespace netgen
}
maxglob = comm.AllReduce (maxglob, NG_MPI_MAX);
int numglob = maxglob+1-PointIndex::BASE;
int numglob = maxglob+1-IndexBASE<PointIndex>();
if (comm.Rank() > 0)
{
comm.Send (globnum, 0, 200);

View File

@ -299,7 +299,7 @@ namespace netgen
Element2d :: Element2d (int pi1, int pi2, int pi3)
Element2d :: Element2d (PointIndex pi1, PointIndex pi2, PointIndex pi3)
{
pnum[0] = pi1;
pnum[1] = pi2;
@ -322,7 +322,7 @@ namespace netgen
is_curved = false;
}
Element2d :: Element2d (int pi1, int pi2, int pi3, int pi4)
Element2d :: Element2d (PointIndex pi1, PointIndex pi2, PointIndex pi3, PointIndex pi4)
{
pnum[0] = pi1;
pnum[1] = pi2;

View File

@ -297,8 +297,13 @@ namespace netgen
{
public:
PointIndices () = default;
PointIndices (INDEX_2 i2) : INDEX_2(i2) { ; }
PointIndices (PointIndex i1, PointIndex i2) : INDEX_2(i1,i2) { ; }
PointIndices (const PointIndices&) = default;
PointIndices (PointIndices&&) = default;
PointIndices & operator= (const PointIndices&) = default;
PointIndices & operator= (PointIndices&&) = default;
constexpr PointIndices (INDEX_2 i2) : INDEX_2(i2) { ; }
constexpr PointIndices (PointIndex i1, PointIndex i2) : INDEX_2(i1,i2) { ; }
PointIndex operator[] (int i) const { return PointIndex(INDEX_2::operator[](i)); }
PointIndex & operator[] (int i) { return reinterpret_cast<PointIndex&>(INDEX_2::operator[](i)); }
using INDEX_2::Sort;
@ -310,8 +315,12 @@ namespace netgen
{
public:
PointIndices () = default;
PointIndices (INDEX_3 i3) : INDEX_3(i3) { ; }
PointIndices (PointIndex i1, PointIndex i2, PointIndex i3) : INDEX_3(i1,i2,i3) { ; }
PointIndices (const PointIndices&) = default;
PointIndices (PointIndices&&) = default;
PointIndices & operator= (const PointIndices&) = default;
PointIndices & operator= (PointIndices&&) = default;
constexpr PointIndices (INDEX_3 i3) : INDEX_3(i3) { ; }
constexpr PointIndices (PointIndex i1, PointIndex i2, PointIndex i3) : INDEX_3(i1,i2,i3) { ; }
PointIndex operator[] (int i) const { return PointIndex(INDEX_3::operator[](i)); }
PointIndex & operator[] (int i) { return reinterpret_cast<PointIndex&>(INDEX_3::operator[](i)); }
using INDEX_3::Sort;
@ -336,6 +345,20 @@ namespace netgen
}
namespace ngcore
{
template <typename T> constexpr inline T InvalidHash();
template <>
constexpr inline netgen::PointIndices<2> InvalidHash<netgen::PointIndices<2>> ()
{ return netgen::PointIndices<2>{netgen::PointIndex::INVALID, netgen::PointIndex::INVALID}; }
template <>
constexpr inline netgen::PointIndices<3> InvalidHash<netgen::PointIndices<3>> ()
{ return netgen::PointIndices<3>{netgen::PointIndex::INVALID, netgen::PointIndex::INVALID, netgen::PointIndex::INVALID}; }
}
namespace std
{
// structured binding support
@ -584,9 +607,9 @@ namespace netgen
///
DLL_HEADER Element2d (ELEMENT_TYPE type);
///
DLL_HEADER Element2d (int pi1, int pi2, int pi3);
DLL_HEADER Element2d (PointIndex pi1, PointIndex pi2, PointIndex pi3);
///
DLL_HEADER Element2d (int pi1, int pi2, int pi3, int pi4);
DLL_HEADER Element2d (PointIndex pi1, PointIndex pi2, PointIndex pi3, PointIndex pi4);
///
ELEMENT_TYPE GetType () const { return typ; }
///