EdgeSwapping() - some cleanup and parallelization of table building

This commit is contained in:
Matthias Hochsteger 2019-10-08 18:35:20 +02:00
parent 288bd2c3d8
commit a651a2d97e
2 changed files with 72 additions and 78 deletions

View File

@ -8,6 +8,7 @@ namespace netgen
{ {
using ngcore::ParallelForRange;
@ -28,7 +29,7 @@ namespace netgen
Array<bool> &swapped, Array<bool> &swapped,
const SurfaceElementIndex t1, const int o1, const SurfaceElementIndex t1, const int o1,
const int t, const int t,
NgArray<int,PointIndex::BASE> &pdef, Array<int,PointIndex> &pdef,
const bool check_only) const bool check_only)
{ {
bool should; bool should;
@ -175,6 +176,7 @@ namespace netgen
void MeshOptimize2d :: EdgeSwapping (Mesh & mesh, int usemetric) void MeshOptimize2d :: EdgeSwapping (Mesh & mesh, int usemetric)
{ {
static Timer timer("EdgeSwapping (2D)"); RegionTimer reg(timer); static Timer timer("EdgeSwapping (2D)"); RegionTimer reg(timer);
static Timer timer_nb("EdgeSwapping-Find neighbors");
if (usemetric) if (usemetric)
PrintMessage (3, "Edgeswapping, metric"); PrintMessage (3, "Edgeswapping, metric");
else else
@ -215,59 +217,68 @@ namespace netgen
Array<bool> swapped(mesh.GetNSE()); Array<bool> swapped(mesh.GetNSE());
NgArray<int,PointIndex::BASE> pdef(mesh.GetNP()); Array<int,PointIndex> pdef(mesh.GetNP());
NgArray<double,PointIndex::BASE> pangle(mesh.GetNP()); Array<double,PointIndex> pangle(mesh.GetNP());
// int e;
// double d;
// Vec3d nv1, nv2;
// double loch(-1);
static const double minangle[] = { 0, 1.481, 2.565, 3.627, 4.683, 5.736, 7, 9 }; static const double minangle[] = { 0, 1.481, 2.565, 3.627, 4.683, 5.736, 7, 9 };
for (int i = 0; i < seia.Size(); i++) if(faceindex == 0)
{ {
const Element2d & sel = mesh[seia[i]]; ParallelForRange( Range(pangle), [&] (auto myrange)
for (int j = 0; j < 3; j++) {
pangle[sel[j]] = 0.0; for (auto i : myrange)
pangle[i] = 0.0;
});
} }
// pangle = 0; else
for (int i = 0; i < seia.Size(); i++)
{ {
const Element2d & sel = mesh[seia[i]]; ParallelForRange( Range(seia), [&] (auto myrange)
for (int j = 0; j < 3; j++) {
{ for (auto i : myrange)
POINTTYPE typ = mesh[sel[j]].Type(); {
if (typ == FIXEDPOINT || typ == EDGEPOINT) const Element2d & sel = mesh[seia[i]];
{ for (int j = 0; j < 3; j++)
pangle[sel[j]] += pangle[sel[j]] = 0.0;
Angle (mesh[sel[(j+1)%3]] - mesh[sel[j]], }
mesh[sel[(j+2)%3]] - mesh[sel[j]]); });
}
}
} }
// for (PointIndex pi = PointIndex::BASE; ParallelForRange( Range(seia), [&] (auto myrange)
// pi < mesh.GetNP()+PointIndex::BASE; pi++) {
for (auto i : myrange)
// pdef = 0; {
for (int i = 0; i < seia.Size(); i++) const Element2d & sel = mesh[seia[i]];
{ for (int j = 0; j < 3; j++)
const Element2d & sel = mesh[seia[i]]; {
for (int j = 0; j < 3; j++) POINTTYPE typ = mesh[sel[j]].Type();
{ if (typ == FIXEDPOINT || typ == EDGEPOINT)
PointIndex pi = sel[j]; {
if (mesh[pi].Type() == INNERPOINT || mesh[pi].Type() == SURFACEPOINT) pangle[sel[j]] +=
pdef[pi] = -6; Angle (mesh[sel[(j+1)%3]] - mesh[sel[j]],
else mesh[sel[(j+2)%3]] - mesh[sel[j]]);
for (int j = 0; j < 8; j++) }
if (pangle[pi] >= minangle[j]) }
pdef[pi] = -1-j;
} }
} });
ParallelForRange( Range(seia), [&] (auto myrange)
{
for (auto i : myrange)
{
const Element2d & sel = mesh[seia[i]];
for (int j = 0; j < 3; j++)
{
PointIndex pi = sel[j];
if (mesh[pi].Type() == INNERPOINT || mesh[pi].Type() == SURFACEPOINT)
pdef[pi] = -6;
else
for (int j = 0; j < 8; j++)
if (pangle[pi] >= minangle[j])
pdef[pi] = -1-j;
}
}
});
/* /*
for (int i = 0; i < seia.Size(); i++) for (int i = 0; i < seia.Size(); i++)
@ -277,41 +288,23 @@ namespace netgen
pdef[sel[j]]++; pdef[sel[j]]++;
} }
*/ */
for (SurfaceElementIndex sei : seia) ParallelForRange( Range(seia), [&] (auto myrange)
for (PointIndex pi : mesh[sei].PNums<3>())
pdef[pi]++;
// for (int i = 0; i < seia.Size(); i++)
for (SurfaceElementIndex sei : seia)
for (int j = 0; j < 3; j++)
{ {
neighbors[sei].SetNr (j, -1); for (auto i : myrange)
neighbors[sei].SetOrientation (j, 0); {
} auto sei = seia[i];
for (PointIndex pi : mesh[sei].template PNums<3>())
/* AsAtomic(pdef[pi])++;
NgArray<Vec3d> normals(mesh.GetNP()); for (int j = 0; j < 3; j++)
for (i = 1; i <= mesh.GetNSE(); i++) {
{ neighbors[sei].SetNr (j, -1);
Element2d & hel = mesh.SurfaceElement(i); neighbors[sei].SetOrientation (j, 0);
if (hel.GetIndex() == faceindex) }
for (k = 1; k <= 3; k++) }
{ });
int pi = hel.PNum(k);
SelectSurfaceOfPoint (mesh.Point(pi), hel.GeomInfoPi(k));
int surfi = mesh.GetFaceDescriptor(faceindex).SurfNr();
GetNormalVector (surfi, mesh.Point(pi), normals.Elem(pi));
normals.Elem(pi) /= normals.Elem(pi).Length();
}
}
*/
/*
for (int i = 0; i < seia.Size(); i++)
{
const Element2d & sel = mesh[seia[i]];
*/
timer_nb.Start();
for (SurfaceElementIndex sei : seia) for (SurfaceElementIndex sei : seia)
{ {
const Element2d & sel = mesh[sei]; const Element2d & sel = mesh[sei];
@ -358,6 +351,7 @@ namespace netgen
} }
} }
} }
timer_nb.Stop();
for (SurfaceElementIndex sei : seia) for (SurfaceElementIndex sei : seia)
swapped[sei] = false; swapped[sei] = false;

View File

@ -37,7 +37,7 @@ public:
const NgArray<Point<3>* > & from, NgArray<Point<3>* > & dest); const NgArray<Point<3>* > & from, NgArray<Point<3>* > & dest);
bool EdgeSwapping (Mesh & mesh, const int usemetric, Array<Neighbour> &neighbors, Array<bool> &swapped, bool EdgeSwapping (Mesh & mesh, const int usemetric, Array<Neighbour> &neighbors, Array<bool> &swapped,
const SurfaceElementIndex t1, const int edge, const int t, NgArray<int,PointIndex::BASE> &pdef, const bool check_only=false); const SurfaceElementIndex t1, const int edge, const int t, Array<int,PointIndex> &pdef, const bool check_only=false);
void EdgeSwapping (Mesh & mesh, int usemetric); void EdgeSwapping (Mesh & mesh, int usemetric);
void CombineImprove (Mesh & mesh); void CombineImprove (Mesh & mesh);
void SplitImprove (Mesh & mesh); void SplitImprove (Mesh & mesh);