diff --git a/libsrc/meshing/basegeom.cpp b/libsrc/meshing/basegeom.cpp index 6df0cf6e..36782786 100644 --- a/libsrc/meshing/basegeom.cpp +++ b/libsrc/meshing/basegeom.cpp @@ -27,29 +27,59 @@ namespace netgen int NetgenGeometry :: GenerateMesh (shared_ptr & mesh, MeshingParameters & mparam) { - if (!mesh) return 1; + multithread.percent = 0; - if (mparam.perfstepsstart <= MESHCONST_MESHVOLUME) + if(mparam.perfstepsstart <= MESHCONST_ANALYSE) { - multithread.task = "Volume meshing"; - - MESHING3_RESULT res = - MeshVolume (mparam, *mesh); - - if (res != MESHING3_OK) return 1; - - if (multithread.terminate) return 0; - - RemoveIllegalElements (*mesh); - if (multithread.terminate) return 0; - - MeshQuality3d (*mesh); + if(!mesh) + mesh = make_shared(); + mesh->geomtype = GetGeomType(); + Analyse(*mesh, mparam); } - - if (multithread.terminate || mparam.perfstepsend <= MESHCONST_MESHVOLUME) + if(multithread.terminate || mparam.perfstepsend <= MESHCONST_ANALYSE) return 0; + if(mparam.perfstepsstart <= MESHCONST_MESHEDGES) + FindEdges(*mesh, mparam); + + if(multithread.terminate || mparam.perfstepsend <= MESHCONST_MESHEDGES) + return 0; + + if (mparam.perfstepsstart <= MESHCONST_MESHSURFACE) + { + MeshSurface(*mesh, mparam); + mesh->CalcSurfacesOfNode(); + } + + if (multithread.terminate || mparam.perfstepsend <= MESHCONST_MESHSURFACE) + return 0; + + if (mparam.perfstepsstart <= MESHCONST_OPTSURFACE) + OptimizeSurface(*mesh, mparam); + + if (multithread.terminate || mparam.perfstepsend <= MESHCONST_OPTSURFACE) + return 0; + + + if(mparam.perfstepsstart <= MESHCONST_MESHVOLUME) + { + multithread.task = "Volume meshing"; + + MESHING3_RESULT res = MeshVolume (mparam, *mesh); + + if (res != MESHING3_OK) return 1; + if (multithread.terminate) return 0; + + RemoveIllegalElements (*mesh); + if (multithread.terminate) return 0; + + MeshQuality3d (*mesh); + } + + if (multithread.terminate || mparam.perfstepsend <= MESHCONST_MESHVOLUME) + return 0; + if (mparam.perfstepsstart <= MESHCONST_OPTVOLUME) { @@ -58,9 +88,9 @@ namespace netgen OptimizeVolume (mparam, *mesh); if (multithread.terminate) return 0; } - + FinalizeMesh(*mesh); return 0; - } + } const Refinement & NetgenGeometry :: GetRefinement () const diff --git a/libsrc/meshing/basegeom.hpp b/libsrc/meshing/basegeom.hpp index f6e63046..1c2223b3 100644 --- a/libsrc/meshing/basegeom.hpp +++ b/libsrc/meshing/basegeom.hpp @@ -22,9 +22,16 @@ namespace netgen virtual const Refinement & GetRefinement () const; - virtual void DoArchive(Archive&) + virtual void DoArchive(Archive&) { throw NgException("DoArchive not implemented for " + Demangle(typeid(*this).name())); } + virtual Mesh::GEOM_TYPE GetGeomType() const { return Mesh::NO_GEOM; } + virtual void Analyse(Mesh& mesh, + const MeshingParameters& mparam) {} + virtual void FindEdges(Mesh& mesh, const MeshingParameters& mparam) {} + virtual void MeshSurface(Mesh& mesh, const MeshingParameters& mparam) {} + virtual void OptimizeSurface(Mesh& mesh, const MeshingParameters& mparam) {} + virtual void FinalizeMesh(Mesh& mesh) const {} virtual void Save (string filename) const; virtual void SaveToMeshFile (ostream & /* ost */) const { ; } }; diff --git a/libsrc/meshing/improve2.cpp b/libsrc/meshing/improve2.cpp index 2fc9f5f4..9fcfc25a 100644 --- a/libsrc/meshing/improve2.cpp +++ b/libsrc/meshing/improve2.cpp @@ -862,70 +862,81 @@ namespace netgen // TODO: split also bad trigs, nut just illegal ones if (mesh.LegalTrig(sel)) continue; - for (int j = 0; j < 3; j++) + // find longest edge + INDEX_2 edge; + double edge_len = 0; + PointIndex pi1, pi2, pi3, pi4; + PointGeomInfo gi1, gi2, gi3, gi4; + for(auto j : Range(1,4)) { - PointIndex pi1 = sel.PNumMod(j+2); - PointIndex pi2 = sel.PNumMod(j+3); - PointIndex pi3 = sel.PNumMod(j+1); - PointIndex pi4; - PointGeomInfo gi1 = sel.GeomInfoPiMod(j+2); - PointGeomInfo gi2 = sel.GeomInfoPiMod(j+3); - PointGeomInfo gi3 = sel.GeomInfoPiMod(j+1); - PointGeomInfo gi4; - - if (mesh.IsSegment (pi1, pi2)) continue; - - // get neighbor element - INDEX_2 ii2 (pi1, pi2); - ii2.Sort(); - auto els = els_on_edge.Get(ii2); - SurfaceElementIndex other_i = get<0>(els); - if(other_i==sei) other_i = get<1>(els); - auto & other = mesh[other_i]; - - // find opposite point of neighbor element - for (int j = 0; j < 3; j++) - if(other[j]!=pi1 && other[j]!=pi2) - { - pi4 = other[j]; - gi4 = other.GeomInfoPi(j); - break; - } - - // split edge pi1,pi2 - Point<3> p5; - PointIndex pi5; - PointGeomInfo gi5; - - mesh.GetGeometry()->GetRefinement().PointBetween (mesh[pi1], mesh[pi2], 0.5, - faceindex, - gi1, gi2, p5, gi5); - - pi5 = mesh.AddPoint(p5); - - Element2d e1(3); - e1.SetIndex(faceindex); - e1={ {pi1,gi1}, {pi5,gi5}, {pi3,gi3} }; - mesh.AddSurfaceElement( e1 ); - - Element2d e2(3); - e2.SetIndex(faceindex); - e2 ={ {pi5,gi5}, {pi2,gi2}, {pi3,gi3} }; - mesh.AddSurfaceElement( e2 ); - - Element2d e3(3); - e3.SetIndex(faceindex); - e3 ={ {pi1,gi1}, {pi4,gi4}, {pi5,gi5} }; - mesh.AddSurfaceElement( e3 ); - - Element2d e4(3); - e4.SetIndex(faceindex); - e4 ={ {pi4,gi4}, {pi2,gi2}, {pi5,gi5} }; - mesh.AddSurfaceElement( e4 ); - - sel.Delete(); - other.Delete(); + auto test_pi1 = sel.PNumMod(j); + auto test_pi2 = sel.PNumMod(j+1); + if (mesh.IsSegment(test_pi1, test_pi2)) + continue; + auto len = (mesh[test_pi2]-mesh[test_pi1]).Length(); + if(len > edge_len) + { + edge = {test_pi1, test_pi2}; + edge.Sort(); + edge_len = len; + pi1 = test_pi1; + pi2 = test_pi2; + pi3 = sel.PNumMod(j+2); + gi1 = sel.GeomInfoPiMod(j); + gi2 = sel.GeomInfoPiMod(j+1); + gi3 = sel.GeomInfoPiMod(j+2); + } } + if(!edge_len) + throw Exception("Couldn't find edge to split, something is wrong"); + // get neighbor element + auto els = els_on_edge.Get(edge); + SurfaceElementIndex other_i = get<0>(els); + if(other_i==sei) other_i = get<1>(els); + auto & other = mesh[other_i]; + + // find opposite point of neighbor element + for (int j = 0; j < 3; j++) + if(other[j]!=pi1 && other[j]!=pi2) + { + pi4 = other[j]; + gi4 = other.GeomInfoPi(j); + break; + } + + // split edge pi1,pi2 + Point<3> p5; + PointIndex pi5; + PointGeomInfo gi5; + + mesh.GetGeometry()->GetRefinement().PointBetween (mesh[pi1], mesh[pi2], 0.5, + faceindex, + gi1, gi2, p5, gi5); + + pi5 = mesh.AddPoint(p5); + + Element2d e1(3); + e1.SetIndex(faceindex); + e1={ {pi1,gi1}, {pi5,gi5}, {pi3,gi3} }; + mesh.AddSurfaceElement( e1 ); + + Element2d e2(3); + e2.SetIndex(faceindex); + e2 ={ {pi5,gi5}, {pi2,gi2}, {pi3,gi3} }; + mesh.AddSurfaceElement( e2 ); + + Element2d e3(3); + e3.SetIndex(faceindex); + e3 ={ {pi1,gi1}, {pi4,gi4}, {pi5,gi5} }; + mesh.AddSurfaceElement( e3 ); + + Element2d e4(3); + e4.SetIndex(faceindex); + e4 ={ {pi4,gi4}, {pi2,gi2}, {pi5,gi5} }; + mesh.AddSurfaceElement( e4 ); + + sel.Delete(); + other.Delete(); } mesh.SetNextTimeStamp(); diff --git a/libsrc/meshing/meshclass.hpp b/libsrc/meshing/meshclass.hpp index 372430ee..db9c96c3 100644 --- a/libsrc/meshing/meshclass.hpp +++ b/libsrc/meshing/meshclass.hpp @@ -864,7 +864,7 @@ namespace netgen /// friend class Meshing3; - + // only for saving the geometry enum GEOM_TYPE { NO_GEOM = 0, GEOM_2D = 1, GEOM_CSG = 10, GEOM_STL = 11, GEOM_OCC = 12, GEOM_ACIS = 13 }; GEOM_TYPE geomtype; diff --git a/libsrc/occ/occgenmesh.cpp b/libsrc/occ/occgenmesh.cpp index 45721490..c8642400 100644 --- a/libsrc/occ/occgenmesh.cpp +++ b/libsrc/occ/occgenmesh.cpp @@ -602,7 +602,7 @@ namespace netgen void OCCMeshSurface (OCCGeometry & geom, Mesh & mesh, - MeshingParameters & mparam) + const MeshingParameters & mparam) { static Timer t("OCCMeshSurface"); RegionTimer r(t); @@ -796,7 +796,7 @@ namespace netgen // Philippose - 15/01/2009 double maxh = geom.face_maxh[k-1]; //double maxh = mparam.maxh; - mparam.checkoverlap = 0; + // mparam.checkoverlap = 0; // int noldpoints = mesh->GetNP(); int noldsurfel = mesh.GetNSE(); @@ -916,7 +916,7 @@ namespace netgen } void OCCOptimizeSurface(OCCGeometry & geom, Mesh & mesh, - MeshingParameters & mparam) + const MeshingParameters & mparam) { const char * savetask = multithread.task; multithread.task = "Optimizing surface"; @@ -991,7 +991,7 @@ namespace netgen - void OCCSetLocalMeshSize(OCCGeometry & geom, Mesh & mesh, + void OCCSetLocalMeshSize(const OCCGeometry & geom, Mesh & mesh, const MeshingParameters & mparam, const OCCParameters& occparam) { static Timer t1("OCCSetLocalMeshSize"); @@ -1279,197 +1279,6 @@ namespace netgen mesh.LoadLocalMeshSize (mparam.meshsizefilename); } - - - - int OCCGenerateMesh (OCCGeometry & geom, shared_ptr & mesh, MeshingParameters & mparam, - const OCCParameters& occparam) - { - multithread.percent = 0; - - if (mparam.perfstepsstart <= MESHCONST_ANALYSE) - { - if(mesh.get() == nullptr) - mesh = make_shared(); - mesh->geomtype = Mesh::GEOM_OCC; - - OCCSetLocalMeshSize(geom,*mesh, mparam, occparam); - } - - if (multithread.terminate || mparam.perfstepsend <= MESHCONST_ANALYSE) - return TCL_OK; - - if (mparam.perfstepsstart <= MESHCONST_MESHEDGES) - { - OCCFindEdges (geom, *mesh, mparam); - - /* - cout << "Removing redundant points" << endl; - - int i, j; - int np = mesh->GetNP(); - NgArray equalto; - - equalto.SetSize (np); - equalto = 0; - - for (i = 1; i <= np; i++) - { - for (j = i+1; j <= np; j++) - { - if (!equalto[j-1] && (Dist2 (mesh->Point(i), mesh->Point(j)) < 1e-12)) - equalto[j-1] = i; - } - } - - for (i = 1; i <= np; i++) - if (equalto[i-1]) - { - cout << "Point " << i << " is equal to Point " << equalto[i-1] << endl; - for (j = 1; j <= mesh->GetNSeg(); j++) - { - Segment & seg = mesh->LineSegment(j); - if (seg[0] == i) seg[0] = equalto[i-1]; - if (seg[1] == i) seg[1] = equalto[i-1]; - } - } - - cout << "Removing degenerated segments" << endl; - for (j = 1; j <= mesh->GetNSeg(); j++) - { - Segment & seg = mesh->LineSegment(j); - if (seg[0] == seg[1]) - { - mesh->DeleteSegment(j); - cout << "Deleting Segment " << j << endl; - } - } - - mesh->Compress(); - */ - - /* - for (int i = 1; i <= geom.fmap.Extent(); i++) - { - Handle(Geom_Surface) hf1 = - BRep_Tool::Surface(TopoDS::Face(geom.fmap(i))); - for (int j = i+1; j <= geom.fmap.Extent(); j++) - { - Handle(Geom_Surface) hf2 = - BRep_Tool::Surface(TopoDS::Face(geom.fmap(j))); - if (hf1 == hf2) cout << "face " << i << " and face " << j << " lie on same surface" << endl; - } - } - */ - -#ifdef LOG_STREAM - (*logout) << "Edges meshed" << endl - << "time = " << GetTime() << " sec" << endl - << "points: " << mesh->GetNP() << endl; -#endif - } - - if (multithread.terminate || mparam.perfstepsend <= MESHCONST_MESHEDGES) - return TCL_OK; - - if (mparam.perfstepsstart <= MESHCONST_MESHSURFACE) - { - OCCMeshSurface (geom, *mesh, mparam); - if (multithread.terminate) return TCL_OK; - -#ifdef LOG_STREAM - (*logout) << "Surfaces meshed" << endl - << "time = " << GetTime() << " sec" << endl - << "points: " << mesh->GetNP() << endl; -#endif - -#ifdef STAT_STREAM - (*statout) << mesh->GetNSeg() << " & " - << mesh->GetNSE() << " & - &" - << GetTime() << " & " << endl; -#endif - - // MeshQuality2d (*mesh); - mesh->CalcSurfacesOfNode(); - } - - if (multithread.terminate || mparam.perfstepsend <= MESHCONST_MESHSURFACE) - return TCL_OK; - - if (mparam.perfstepsstart <= MESHCONST_OPTSURFACE) - { - OCCOptimizeSurface(geom, *mesh, mparam); - } - - if (multithread.terminate || mparam.perfstepsend <= MESHCONST_OPTSURFACE) - return TCL_OK; - - if (mparam.perfstepsstart <= MESHCONST_MESHVOLUME) - { - multithread.task = "Volume meshing"; - - MESHING3_RESULT res = MeshVolume (mparam, *mesh); - - if (res != MESHING3_OK) return TCL_ERROR; - if (multithread.terminate) return TCL_OK; - - RemoveIllegalElements (*mesh); - if (multithread.terminate) return TCL_OK; - - MeshQuality3d (*mesh); - -#ifdef STAT_STREAM - (*statout) << GetTime() << " & "; -#endif - -#ifdef LOG_STREAM - (*logout) << "Volume meshed" << endl - << "time = " << GetTime() << " sec" << endl - << "points: " << mesh->GetNP() << endl; -#endif - } - - if (multithread.terminate || mparam.perfstepsend <= MESHCONST_MESHVOLUME) - return TCL_OK; - - if (mparam.perfstepsstart <= MESHCONST_OPTVOLUME) - { - multithread.task = "Volume optimization"; - - OptimizeVolume (mparam, *mesh); - if (multithread.terminate) return TCL_OK; - -#ifdef STAT_STREAM - (*statout) << GetTime() << " & " - << mesh->GetNE() << " & " - << mesh->GetNP() << " " << '\\' << '\\' << " \\" << "hline" << endl; -#endif - -#ifdef LOG_STREAM - (*logout) << "Volume optimized" << endl - << "time = " << GetTime() << " sec" << endl - << "points: " << mesh->GetNP() << endl; -#endif - - // cout << "Optimization complete" << endl; - - } - - /* - (*testout) << "NP: " << mesh->GetNP() << endl; - for (int i = 1; i <= mesh->GetNP(); i++) - (*testout) << mesh->Point(i) << endl; - - (*testout) << endl << "NSegments: " << mesh->GetNSeg() << endl; - for (int i = 1; i <= mesh->GetNSeg(); i++) - (*testout) << mesh->LineSegment(i) << endl; - */ - - for (int i = 0; i < mesh->GetNDomains(); i++) - if (geom.snames.Size()) - mesh->SetMaterial (i+1, geom.snames[i]); - return TCL_OK; - } } #endif diff --git a/libsrc/occ/occgeom.cpp b/libsrc/occ/occgeom.cpp index 0faf7f0a..49c9742e 100644 --- a/libsrc/occ/occgeom.cpp +++ b/libsrc/occ/occgeom.cpp @@ -74,6 +74,36 @@ void STEP_GetEntityName(const TopoDS_Shape & theShape, STEPCAFControl_Reader * a } } + void OCCGeometry :: Analyse(Mesh& mesh, + const MeshingParameters& mparam) + { + OCCSetLocalMeshSize(*this, mesh, mparam, occparam); + } + + void OCCGeometry :: FindEdges(Mesh& mesh, + const MeshingParameters& mparam) + { + OCCFindEdges(*this, mesh, mparam); + } + + void OCCGeometry :: MeshSurface(Mesh& mesh, + const MeshingParameters& mparam) + { + OCCMeshSurface(*this, mesh, mparam); + } + + void OCCGeometry :: OptimizeSurface(Mesh& mesh, + const MeshingParameters& mparam) + { + OCCOptimizeSurface(*this, mesh, mparam); + } + + void OCCGeometry :: FinalizeMesh(Mesh& mesh) const + { + for (int i = 0; i < mesh.GetNDomains(); i++) + if (snames.Size()) + mesh.SetMaterial (i+1, snames[i]); + } void OCCGeometry :: PrintNrShapes () { @@ -1703,10 +1733,10 @@ void STEP_GetEntityName(const TopoDS_Shape & theShape, STEPCAFControl_Reader * a DLL_HEADER extern OCCParameters occparam; OCCParameters occparam; - int OCCGeometry :: GenerateMesh (shared_ptr & mesh, MeshingParameters & mparam) - { - return OCCGenerateMesh (*this, mesh, mparam, occparam); - } + // int OCCGeometry :: GenerateMesh (shared_ptr & mesh, MeshingParameters & mparam) + // { + // return OCCGenerateMesh (*this, mesh, mparam, occparam); + // } } diff --git a/libsrc/occ/occgeom.hpp b/libsrc/occ/occgeom.hpp index 750c6e29..6b0415d0 100644 --- a/libsrc/occ/occgeom.hpp +++ b/libsrc/occ/occgeom.hpp @@ -183,12 +183,33 @@ namespace netgen return a00*a11*a22 + a01*a12*a20 + a10*a21*a02 - a20*a11*a02 - a10*a01*a22 - a21*a12*a00; } + class DLL_HEADER OCCParameters + { + public: + /// Factor for meshing close edges + double resthcloseedgefac = 2.; + /// Enable / Disable detection of close edges + int resthcloseedgeenable = true; + + /// Minimum edge length to be used for dividing edges to mesh points + double resthminedgelen = 0.001; + + /// Enable / Disable use of the minimum edge length (by default use 1e-4) + int resthminedgelenenable = true; + + /*! + Dump all the OpenCascade specific meshing parameters + to console + */ + void Print (ostream & ost) const; + }; class OCCGeometry : public NetgenGeometry { Point<3> center; + OCCParameters occparam; public: TopoDS_Shape shape; @@ -239,11 +260,26 @@ namespace netgen emap.Clear(); vmap.Clear(); } + + Mesh::GEOM_TYPE GetGeomType() const override + { return Mesh::GEOM_OCC; } + + void SetOCCParameters(const OCCParameters& par) + { cout << "set occ pars to = " << par.resthcloseedgefac << endl; occparam = par; } + + void Analyse(Mesh& mesh, + const MeshingParameters& mparam) override; + void FindEdges(Mesh& mesh, + const MeshingParameters& mparam) override; + void MeshSurface(Mesh& mesh, + const MeshingParameters& mparam) override; + void OptimizeSurface(Mesh& mesh, + const MeshingParameters& mparam) override; + void FinalizeMesh(Mesh& mesh) const override; + DLL_HEADER void Save (string filename) const override; - DLL_HEADER virtual void Save (string filename) const; - - void DoArchive(Archive& ar); + void DoArchive(Archive& ar) override; DLL_HEADER void BuildFMap(); @@ -391,37 +427,9 @@ namespace netgen // void WriteOCC_STL(char * filename); - DLL_HEADER virtual int GenerateMesh (shared_ptr & mesh, MeshingParameters & mparam); + // DLL_HEADER virtual int GenerateMesh (shared_ptr & mesh, MeshingParameters & mparam); - DLL_HEADER virtual const Refinement & GetRefinement () const; - }; - - - - class DLL_HEADER OCCParameters - { - public: - - /// Factor for meshing close edges - double resthcloseedgefac = 2.; - - - /// Enable / Disable detection of close edges - int resthcloseedgeenable = true; - - - /// Minimum edge length to be used for dividing edges to mesh points - double resthminedgelen = 0.001; - - - /// Enable / Disable use of the minimum edge length (by default use 1e-4) - int resthminedgelenenable = true; - - /*! - Dump all the OpenCascade specific meshing parameters - to console - */ - void Print (ostream & ost) const; + DLL_HEADER const Refinement & GetRefinement () const override; }; @@ -434,15 +442,12 @@ namespace netgen // Philippose - 31.09.2009 // External access to the mesh generation functions within the OCC // subsystem (Not sure if this is the best way to implement this....!!) - DLL_HEADER extern int OCCGenerateMesh (OCCGeometry & occgeometry, shared_ptr & mesh, - MeshingParameters & mparam, const OCCParameters& occparam); - - DLL_HEADER extern void OCCSetLocalMeshSize(OCCGeometry & geom, Mesh & mesh, const MeshingParameters & mparam, + DLL_HEADER extern void OCCSetLocalMeshSize(const OCCGeometry & geom, Mesh & mesh, const MeshingParameters & mparam, const OCCParameters& occparam); - DLL_HEADER extern void OCCMeshSurface (OCCGeometry & geom, Mesh & mesh, MeshingParameters & mparam); + DLL_HEADER extern void OCCMeshSurface (OCCGeometry & geom, Mesh & mesh, const MeshingParameters & mparam); - DLL_HEADER extern void OCCOptimizeSurface (OCCGeometry & geom, Mesh & mesh, MeshingParameters & mparam); + DLL_HEADER extern void OCCOptimizeSurface (OCCGeometry & geom, Mesh & mesh, const MeshingParameters & mparam); DLL_HEADER extern void OCCFindEdges (OCCGeometry & geom, Mesh & mesh, const MeshingParameters & mparam); } diff --git a/libsrc/occ/occpkg.cpp b/libsrc/occ/occpkg.cpp index de9367d7..7ebc13d3 100644 --- a/libsrc/occ/occpkg.cpp +++ b/libsrc/occ/occpkg.cpp @@ -53,6 +53,8 @@ namespace netgen atof (Tcl_GetVar (interp, "::stloptions.resthminedgelen", 0)); occparam.resthminedgelenenable = atoi (Tcl_GetVar (interp, "::stloptions.resthminedgelenenable", 0)); + if(auto geo = dynamic_pointer_cast(ng_geometry); geo) + geo->SetOCCParameters(occparam); } }; diff --git a/libsrc/occ/python_occ.cpp b/libsrc/occ/python_occ.cpp index 6a847e14..d3a38c15 100644 --- a/libsrc/occ/python_occ.cpp +++ b/libsrc/occ/python_occ.cpp @@ -183,11 +183,12 @@ DLL_HEADER void ExportNgOCC(py::module &m) CreateOCCParametersFromKwargs(occparam, kwargs); CreateMPfromKwargs(mp, kwargs); } + geo->SetOCCParameters(occparam); auto mesh = make_shared(); - SetGlobalMesh(mesh); + geo->GenerateMesh(mesh, mp); mesh->SetGeometry(geo); + SetGlobalMesh(mesh); ng_geometry = geo; - OCCGenerateMesh(*geo, mesh, mp, occparam); return mesh; }, py::arg("mp") = nullptr, py::call_guard(), diff --git a/libsrc/stlgeom/meshstlsurface.cpp b/libsrc/stlgeom/meshstlsurface.cpp index d7b25e8f..25932abe 100644 --- a/libsrc/stlgeom/meshstlsurface.cpp +++ b/libsrc/stlgeom/meshstlsurface.cpp @@ -867,6 +867,22 @@ void STLSurfaceOptimization (STLGeometry & geom, break; } } + while(mesh.CheckOverlappingBoundary()) + { + for(const auto & el : mesh.SurfaceElements()) + { + if(el.BadElement()) + { + cout << "Restrict localh at el nr " << el << endl; + for(const auto& p : el.PNums()) + { + const auto& pnt = mesh[p]; + mesh.RestrictLocalH(pnt, 0.5*mesh.GetH(pnt)); + } + } + } + optmesh.SplitImprove(mesh); + } //(*testout) << "optimize, after, step = " << meshparam.optimize2d[j-1] << mesh.Point (3679) << endl; } diff --git a/ng/ngpkg.cpp b/ng/ngpkg.cpp index 5c0b8638..85052353 100644 --- a/ng/ngpkg.cpp +++ b/ng/ngpkg.cpp @@ -537,6 +537,7 @@ namespace netgen // delete ng_geometry; // ng_geometry = hgeom; ng_geometry = shared_ptr (hgeom); + geometryregister[i]->SetParameters(interp); mesh.reset(); return TCL_OK;