revert pyramid-specific code in optimizations

This commit is contained in:
Matthias Hochsteger 2020-07-27 18:06:43 +02:00
parent c0b8b1c0cc
commit ae268637cf
3 changed files with 30 additions and 63 deletions

View File

@ -16,28 +16,16 @@ namespace netgen
// Calc badness of new element where pi1 and pi2 are replaced by pnew // Calc badness of new element where pi1 and pi2 are replaced by pnew
static double CalcBadReplacePoints (const Mesh::T_POINTS & points, const MeshingParameters & mp, const Element & elem, double h, PointIndex &pi1, PointIndex &pi2, MeshPoint &pnew) double CalcBadReplacePoints (const Mesh::T_POINTS & points, const MeshingParameters & mp, const Element & elem, double h, PointIndex &pi1, PointIndex &pi2, MeshPoint &pnew)
{ {
if (elem.GetType() != TET && elem.GetType() != PYRAMID) if (elem.GetType() != TET) return 0;
return 0;
MeshPoint p[5]; MeshPoint* p[] = {&points[elem[0]], &points[elem[1]], &points[elem[2]], &points[elem[3]]};
for (auto i : Range(elem.GetNP())) for (auto i : Range(4))
{ if(elem[i]==pi1 || elem[i]==pi2) p[i] = &pnew;
auto pi = elem[i];
if(pi==pi1 || pi==pi2)
p[i] = pnew;
else
p[i] = points[pi];
}
if (elem.GetType() == TET) return CalcTetBadness (*p[0], *p[1], *p[2], *p[3], h, mp);
return CalcTetBadness (p[0], p[1], p[2], p[3], h, mp);
if (elem.GetType() == PYRAMID)
return CalcTetBadness (p[0], p[1], p[2], p[4], h, mp)
+ CalcTetBadness (p[2], p[3], p[0], p[4], h, mp);
return 0;
} }
static ArrayMem<Element, 3> SplitElement (Element old, PointIndex pi0, PointIndex pi1, PointIndex pinew) static ArrayMem<Element, 3> SplitElement (Element old, PointIndex pi0, PointIndex pi1, PointIndex pinew)
@ -4107,11 +4095,15 @@ double MeshOptimize3d :: SplitImprove2Element (Mesh & mesh,
for (auto ei0 : has_both_points0) for (auto ei0 : has_both_points0)
{ {
if(mesh[ei0].GetType()!=TET)
return false;
badness_before += el_badness[ei0]; badness_before += el_badness[ei0];
badness_after += SplitElementBadness (mesh.Points(), mp, mesh[ei0], pi0, pi1, pnew); badness_after += SplitElementBadness (mesh.Points(), mp, mesh[ei0], pi0, pi1, pnew);
} }
for (auto ei1 : has_both_points1) for (auto ei1 : has_both_points1)
{ {
if(mesh[ei1].GetType()!=TET)
return false;
badness_before += el_badness[ei1]; badness_before += el_badness[ei1];
badness_after += SplitElementBadness (mesh.Points(), mp, mesh[ei1], pi2, pi3, pnew); badness_after += SplitElementBadness (mesh.Points(), mp, mesh[ei1], pi2, pi3, pnew);
} }

View File

@ -1,19 +1,6 @@
#ifndef FILE_IMPROVE3 #ifndef FILE_IMPROVE3
#define FILE_IMPROVE3 #define FILE_IMPROVE3
inline double
CalcBad (const Mesh::T_POINTS & points, const Element & elem, double h, const MeshingParameters & mp)
{
if (elem.GetType() == TET)
return CalcTetBadness (points[elem[0]], points[elem[1]],
points[elem[2]], points[elem[3]], h, mp);
if (elem.GetType() == PYRAMID)
return CalcTetBadness (points[elem[0]], points[elem[1]],
points[elem[2]], points[elem[4]], h, mp)
+ CalcTetBadness (points[elem[2]], points[elem[3]],
points[elem[0]], points[elem[4]], h, mp);
return 0;
}
extern double CalcTotalBad (const Mesh::T_POINTS & points, extern double CalcTotalBad (const Mesh::T_POINTS & points,
const Array<Element> & elements, const Array<Element> & elements,
@ -59,9 +46,13 @@ public:
double double
CalcBad (const Mesh::T_POINTS & points, const Element & elem, double h) CalcBad (const Mesh::T_POINTS & points, const Element & elem, double h)
{ {
return ::netgen::CalcBad(points, elem, h, mp); if (elem.GetType() == TET)
return CalcTetBadness (points[elem[0]], points[elem[1]],
points[elem[2]], points[elem[3]], h, mp);
return 0;
} }
double CalcTotalBad (const Mesh::T_POINTS & points, double CalcTotalBad (const Mesh::T_POINTS & points,
const Array<Element> & elements) const Array<Element> & elements)
{ {
@ -71,6 +62,16 @@ public:
}; };
inline double
CalcBad (const Mesh::T_POINTS & points, const Element & elem, double h, const MeshingParameters & mp)
{
if (elem.GetType() == TET)
return CalcTetBadness (points[elem[0]], points[elem[1]],
points[elem[2]], points[elem[3]], h, mp);
return 0;
}
extern int WrongOrientation (const Mesh::T_POINTS & points, const Element & el); extern int WrongOrientation (const Mesh::T_POINTS & points, const Element & el);

View File

@ -339,8 +339,8 @@ namespace netgen
{ {
static Timer tim("PointFunction - build elementsonpoint table"); RegionTimer reg(tim); static Timer tim("PointFunction - build elementsonpoint table"); RegionTimer reg(tim);
for (int i = 0; i < elements.Size(); i++) for (int i = 0; i < elements.Size(); i++)
if (!elements[i].IsDeleted()) if (elements[i].NP() == 4)
for (int j = 0; j < elements[i].GetNP(); j++) for (int j = 0; j < elements[i].NP(); j++)
elementsonpoint.Add (elements[i][j], i); elementsonpoint.Add (elements[i][j], i);
} }
@ -362,7 +362,8 @@ namespace netgen
for (int j = 0; j < elementsonpoint[actpind].Size(); j++) for (int j = 0; j < elementsonpoint[actpind].Size(); j++)
{ {
const Element & el = elements[elementsonpoint[actpind][j]]; const Element & el = elements[elementsonpoint[actpind][j]];
badness += CalcBad(points, el, -1, mp); badness += CalcTetBadness (points[el[0]], points[el[1]],
points[el[2]], points[el[3]], -1, mp);
} }
points[actpind] = Point<3> (hp); points[actpind] = Point<3> (hp);
@ -390,19 +391,6 @@ namespace netgen
vgrad += vgradi; vgrad += vgradi;
} }
if(el.GetType()==PYRAMID)
{
f += CalcTetBadnessGrad (points[el[0]],
points[el[1]],
points[el[2]],
points[el[4]], -1, 4, vgradi, mp);
vgrad += vgradi;
f += CalcTetBadnessGrad (points[el[2]],
points[el[3]],
points[el[0]],
points[el[4]], -1, 4, vgradi, mp);
vgrad += vgradi;
}
} }
points[actpind] = Point<3> (hp); points[actpind] = Point<3> (hp);
@ -435,20 +423,6 @@ namespace netgen
vgrad += vgradi; vgrad += vgradi;
} }
if(el.GetType()==PYRAMID)
{
f += CalcTetBadnessGrad (points[el[0]],
points[el[1]],
points[el[2]],
points[el[4]], -1, 4, vgradi, mp);
vgrad += vgradi;
f += CalcTetBadnessGrad (points[el[2]],
points[el[3]],
points[el[0]],
points[el[4]], -1, 4, vgradi, mp);
vgrad += vgradi;
}
} }
points[actpind] = Point<3> (hp); points[actpind] = Point<3> (hp);