little occ code polish

This commit is contained in:
Joachim Schoeberl 2021-07-30 08:42:35 +02:00
parent 751f193d81
commit 62463b904e
3 changed files with 35 additions and 41 deletions

View File

@ -18,7 +18,7 @@ namespace netgen
#define VSMALL 1e-10 #define VSMALL 1e-10
DLL_HEADER bool merge_solids = 1; DLL_HEADER bool merge_solids = false;
// can you please explain what you intend to compute here (JS) !!! // can you please explain what you intend to compute here (JS) !!!
@ -317,9 +317,9 @@ namespace netgen
for (int i = 1; i <= nvertices; i++) for (int i = 1; i <= nvertices; i++)
{ {
gp_Pnt pnt = BRep_Tool::Pnt (TopoDS::Vertex(geom.vmap(i))); gp_Pnt pnt = BRep_Tool::Pnt (TopoDS::Vertex(geom.vmap(i)));
MeshPoint mp( Point<3>(pnt.X(), pnt.Y(), pnt.Z()) ); MeshPoint mp(occ2ng(pnt));
bool exists = 0; bool exists = false;
if (merge_solids) if (merge_solids)
for (PointIndex pi : mesh.Points().Range()) for (PointIndex pi : mesh.Points().Range())
if (Dist2 (mesh[pi], Point<3>(mp)) < eps*eps) if (Dist2 (mesh[pi], Point<3>(mp)) < eps*eps)
@ -496,7 +496,7 @@ namespace netgen
for (size_t i = 1; i <= mp.Size(); i++) for (size_t i = 1; i <= mp.Size(); i++)
{ {
bool exists = 0; bool exists = false;
tsearch.Start(); tsearch.Start();
// for (PointIndex j = first_ep; j < mesh.Points().End(); j++) // for (PointIndex j = first_ep; j < mesh.Points().End(); j++)

View File

@ -117,6 +117,11 @@ namespace netgen
return Point<3> (p.X(), p.Y(), p.Z()); return Point<3> (p.X(), p.Y(), p.Z());
} }
inline Vec<3> occ2ng (const gp_Vec & v)
{
return Vec<3> (v.X(), v.Y(), v.Z());
}
inline gp_Pnt ng2occ (const Point<3> & p) inline gp_Pnt ng2occ (const Point<3> & p)
{ {
return gp_Pnt(p(0), p(1), p(2)); return gp_Pnt(p(0), p(1), p(2));

View File

@ -19,34 +19,18 @@ namespace netgen
const PointGeomInfo & geominfo, const PointGeomInfo & geominfo,
Vec<3> & n) const Vec<3> & n) const
{ {
gp_Pnt pnt; GeomLProp_SLProps lprop(occface,geominfo.u,geominfo.v,1,1e-8);
gp_Vec du, dv;
/* if (lprop.IsNormalDefined())
double gu = geominfo.u;
double gv = geominfo.v;
if (fabs (gu) < 1e-3) gu = 0;
if (fabs (gv) < 1e-3) gv = 0;
occface->D1(gu,gv,pnt,du,dv);
*/
/*
occface->D1(geominfo.u,geominfo.v,pnt,du,dv);
n = Cross (Vec<3>(du.X(), du.Y(), du.Z()),
Vec<3>(dv.X(), dv.Y(), dv.Z()));
n.Normalize();
*/
GeomLProp_SLProps lprop(occface,geominfo.u,geominfo.v,1,1e-5);
double setu=geominfo.u,setv=geominfo.v;
if(lprop.D1U().Magnitude() < 1e-5 || lprop.D1V().Magnitude() < 1e-5)
{ {
n = occ2ng(lprop.Normal());
}
else
{
gp_Pnt pnt;
gp_Vec du, dv;
double setu=geominfo.u,setv=geominfo.v;
double ustep = 0.01*(umax-umin); double ustep = 0.01*(umax-umin);
double vstep = 0.01*(vmax-vmin); double vstep = 0.01*(vmax-vmin);
@ -57,9 +41,12 @@ namespace netgen
if(setu < umax) if(setu < umax)
{ {
lprop.SetParameters(setu,setv); lprop.SetParameters(setu,setv);
/*
n(0)+=lprop.Normal().X(); n(0)+=lprop.Normal().X();
n(1)+=lprop.Normal().Y(); n(1)+=lprop.Normal().Y();
n(2)+=lprop.Normal().Z(); n(2)+=lprop.Normal().Z();
*/
n += occ2ng(lprop.Normal());
} }
setu = geominfo.u; setu = geominfo.u;
while(setu > umin && (lprop.D1U().Magnitude() < 1e-5 || lprop.D1V().Magnitude() < 1e-5)) while(setu > umin && (lprop.D1U().Magnitude() < 1e-5 || lprop.D1V().Magnitude() < 1e-5))
@ -67,9 +54,12 @@ namespace netgen
if(setu > umin) if(setu > umin)
{ {
lprop.SetParameters(setu,setv); lprop.SetParameters(setu,setv);
/*
n(0)+=lprop.Normal().X(); n(0)+=lprop.Normal().X();
n(1)+=lprop.Normal().Y(); n(1)+=lprop.Normal().Y();
n(2)+=lprop.Normal().Z(); n(2)+=lprop.Normal().Z();
*/
n += occ2ng(lprop.Normal());
} }
setu = geominfo.u; setu = geominfo.u;
@ -78,9 +68,12 @@ namespace netgen
if(setv < vmax) if(setv < vmax)
{ {
lprop.SetParameters(setu,setv); lprop.SetParameters(setu,setv);
/*
n(0)+=lprop.Normal().X(); n(0)+=lprop.Normal().X();
n(1)+=lprop.Normal().Y(); n(1)+=lprop.Normal().Y();
n(2)+=lprop.Normal().Z(); n(2)+=lprop.Normal().Z();
*/
n += occ2ng(lprop.Normal());
} }
setv = geominfo.v; setv = geominfo.v;
while(setv > vmin && (lprop.D1U().Magnitude() < 1e-5 || lprop.D1V().Magnitude() < 1e-5)) while(setv > vmin && (lprop.D1U().Magnitude() < 1e-5 || lprop.D1V().Magnitude() < 1e-5))
@ -88,20 +81,17 @@ namespace netgen
if(setv > vmin) if(setv > vmin)
{ {
lprop.SetParameters(setu,setv); lprop.SetParameters(setu,setv);
/*
n(0)+=lprop.Normal().X(); n(0)+=lprop.Normal().X();
n(1)+=lprop.Normal().Y(); n(1)+=lprop.Normal().Y();
n(2)+=lprop.Normal().Z(); n(2)+=lprop.Normal().Z();
*/
n += occ2ng(lprop.Normal());
} }
setv = geominfo.v; setv = geominfo.v;
n.Normalize(); n.Normalize();
} }
else
{
n(0)=lprop.Normal().X();
n(1)=lprop.Normal().Y();
n(2)=lprop.Normal().Z();
}
if(glob_testout) if(glob_testout)
{ {
@ -112,7 +102,7 @@ namespace netgen
if (orient == TopAbs_REVERSED) n = -1*n; if (orient == TopAbs_REVERSED) n = -n;
// (*testout) << "GetNormalVector" << endl; // (*testout) << "GetNormalVector" << endl;
} }
@ -361,8 +351,7 @@ namespace netgen
gi.u = pspnew(0); gi.u = pspnew(0);
gi.v = pspnew(1); gi.v = pspnew(1);
gi.trignum = 1; gi.trignum = 1;
gp_Pnt val = occface->Value (gi.u, gi.v); p3d = occ2ng(occface->Value (gi.u, gi.v));
p3d = Point<3> (val.X(), val.Y(), val.Z());
}; };
} }
@ -376,7 +365,7 @@ namespace netgen
// try Newton's method ... // try Newton's method ...
gp_Pnt p(ap(0), ap(1), ap(2)); gp_Pnt p = ng2occ(ap);
double u = gi.u; double u = gi.u;
double v = gi.v; double v = gi.v;
@ -488,7 +477,7 @@ namespace netgen
gi.v = v; gi.v = v;
gi.trignum = 1; gi.trignum = 1;
ap = Point<3> (pnt.X(), pnt.Y(), pnt.Z()); ap = occ2ng(pnt); // Point<3> (pnt.X(), pnt.Y(), pnt.Z());
} }