mirror of
https://github.com/NGSolve/netgen.git
synced 2024-12-24 21:10:33 +05:00
Merge branch 'occ_closesurfaces' into 'master'
OCC - Identify shapes (for periodic/closesurface meshes) See merge request jschoeberl/netgen!462
This commit is contained in:
commit
55196514d2
@ -56,18 +56,18 @@ namespace netgen
|
||||
return false;
|
||||
auto & e = *other_ptr;
|
||||
|
||||
if(tol < Dist(GetCenter(), e.GetCenter()))
|
||||
if(tol < Dist(trafo(GetCenter()), e.GetCenter()))
|
||||
return false;
|
||||
|
||||
auto v0 = GetStartVertex().GetPoint();
|
||||
auto v1 = GetEndVertex().GetPoint();
|
||||
auto v0 = trafo(GetStartVertex().GetPoint());
|
||||
auto v1 = trafo(GetEndVertex().GetPoint());
|
||||
auto w0 = e.GetStartVertex().GetPoint();
|
||||
auto w1 = e.GetEndVertex().GetPoint();
|
||||
|
||||
// have two closed edges, use midpoints to compare
|
||||
if(Dist(v0,v1) < tol && Dist(w0,w1) < tol)
|
||||
{
|
||||
v1 = GetPoint(0.5);
|
||||
v1 = trafo(GetPoint(0.5));
|
||||
w1 = other_ptr->GetPoint(0.5);
|
||||
}
|
||||
|
||||
@ -75,6 +75,38 @@ namespace netgen
|
||||
(Dist(v0, w1) < tol && Dist(v1, w0) < tol) );
|
||||
}
|
||||
|
||||
bool GeometryFace :: IsMappedShape( const GeometryShape & other_, const Transformation<3> & trafo, double tol ) const
|
||||
{
|
||||
const auto other_ptr = dynamic_cast<const GeometryFace*>(&other_);
|
||||
if(!other_ptr)
|
||||
return false;
|
||||
auto & f = *other_ptr;
|
||||
|
||||
if(tol < Dist(GetCenter(), f.GetCenter()))
|
||||
return false;
|
||||
|
||||
// simple check: check if there is a bijective mapping of mapped edges
|
||||
auto & other_edges = f.edges;
|
||||
if(edges.Size() != other_edges.Size())
|
||||
return false;
|
||||
|
||||
auto nedges = edges.Size();
|
||||
Array<bool> is_mapped(nedges);
|
||||
is_mapped = false;
|
||||
|
||||
for(auto e : edges)
|
||||
{
|
||||
int found_mapping = 0;
|
||||
for(auto other_e : other_edges)
|
||||
if(e->IsMappedShape(*other_e, trafo, tol))
|
||||
found_mapping++;
|
||||
if(found_mapping != 1)
|
||||
return false;
|
||||
}
|
||||
|
||||
return true;
|
||||
}
|
||||
|
||||
void GeometryFace :: RestrictHTrig(Mesh& mesh,
|
||||
const PointGeomInfo& gi0,
|
||||
const PointGeomInfo& gi1,
|
||||
@ -169,6 +201,15 @@ namespace netgen
|
||||
|
||||
void NetgenGeometry :: ProcessIdentifications()
|
||||
{
|
||||
for(auto i : Range(vertices))
|
||||
vertices[i]->nr = i;
|
||||
for(auto i : Range(edges))
|
||||
edges[i]->nr = i;
|
||||
for(auto i : Range(faces))
|
||||
faces[i]->nr = i;
|
||||
for(auto i : Range(solids))
|
||||
solids[i]->nr = i;
|
||||
|
||||
auto mirror_identifications = [&] ( auto & shapes )
|
||||
{
|
||||
for(auto i : Range(shapes))
|
||||
@ -181,11 +222,39 @@ namespace netgen
|
||||
}
|
||||
};
|
||||
|
||||
auto tol = 1e-8 * bounding_box.Diam();
|
||||
for(auto & f : faces)
|
||||
for(auto & ident: f->identifications)
|
||||
for(auto e : static_cast<GeometryFace*>(ident.from)->edges)
|
||||
for(auto e_other : static_cast<GeometryFace*>(ident.to)->edges)
|
||||
if(e->IsMappedShape(*e_other, ident.trafo, tol))
|
||||
e->identifications.Append( {e, e_other, ident.trafo, ident.type, ident.name} );
|
||||
|
||||
for(auto & e : edges)
|
||||
for(auto & ident: e->identifications)
|
||||
{
|
||||
auto & from = static_cast<GeometryEdge&>(*ident.from);
|
||||
auto & to = static_cast<GeometryEdge&>(*ident.to);
|
||||
|
||||
GeometryVertex * pfrom[] = { &from.GetStartVertex(), &from.GetEndVertex() };
|
||||
GeometryVertex * pto[] = { &to.GetStartVertex(), &to.GetEndVertex() };
|
||||
|
||||
// swap points of other edge if necessary
|
||||
Point<3> p_from0 = ident.trafo(from.GetStartVertex().GetPoint());
|
||||
Point<3> p_from1 = ident.trafo(from.GetEndVertex().GetPoint());
|
||||
Point<3> p_to0 = ident.trafo(to.GetStartVertex().GetPoint());
|
||||
|
||||
if(Dist(p_from1, p_to0) < Dist(p_from1, p_to0))
|
||||
swap(pto[0], pto[1]);
|
||||
|
||||
for(auto i : Range(2))
|
||||
pfrom[i]->identifications.Append( {pfrom[i], pto[i], ident.trafo, ident.type, ident.name} );
|
||||
}
|
||||
|
||||
mirror_identifications(vertices);
|
||||
mirror_identifications(edges);
|
||||
mirror_identifications(faces);
|
||||
|
||||
// todo: propagate identifications faces -> edges -> vertices
|
||||
|
||||
auto find_primary = [&] (auto & shapes)
|
||||
{
|
||||
@ -729,6 +798,9 @@ namespace netgen
|
||||
sel[1] = s[1];
|
||||
sel[2] = tree.Find(trafo(mesh[s[1]]));
|
||||
sel[3] = tree.Find(trafo(mesh[s[0]]));
|
||||
for(auto i : Range(4))
|
||||
sel.GeomInfo()[i] = face.Project(mesh[sel[i]]);
|
||||
|
||||
sel.SetIndex(face.nr+1);
|
||||
mesh.AddSurfaceElement(sel);
|
||||
}
|
||||
|
@ -70,11 +70,19 @@ namespace netgen
|
||||
|
||||
class DLL_HEADER GeometryEdge : public GeometryShape
|
||||
{
|
||||
protected:
|
||||
GeometryVertex *start, *end;
|
||||
public:
|
||||
int domin=-1, domout=-1;
|
||||
|
||||
virtual const GeometryVertex& GetStartVertex() const = 0;
|
||||
virtual const GeometryVertex& GetEndVertex() const = 0;
|
||||
GeometryEdge( GeometryVertex &start_, GeometryVertex &end_ )
|
||||
: start(&start_), end(&end_)
|
||||
{}
|
||||
|
||||
virtual const GeometryVertex& GetStartVertex() const { return *start; }
|
||||
virtual const GeometryVertex& GetEndVertex() const { return *end; }
|
||||
virtual GeometryVertex& GetStartVertex() { return *start; }
|
||||
virtual GeometryVertex& GetEndVertex() { return *end; }
|
||||
virtual double GetLength() const = 0;
|
||||
virtual Point<3> GetCenter() const = 0;
|
||||
virtual Point<3> GetPoint(double t) const = 0;
|
||||
@ -100,11 +108,13 @@ namespace netgen
|
||||
class DLL_HEADER GeometryFace : public GeometryShape
|
||||
{
|
||||
public:
|
||||
Array<GeometryEdge*> edges;
|
||||
int domin=-1, domout=-1;
|
||||
|
||||
virtual Point<3> GetCenter() const = 0;
|
||||
virtual size_t GetNBoundaries() const = 0;
|
||||
virtual Array<Segment> GetBoundary(const Mesh& mesh) const = 0;
|
||||
|
||||
virtual PointGeomInfo Project(Point<3>& p) const = 0;
|
||||
// Project point using geo info. Fast if point is close to
|
||||
// parametrization in geo info.
|
||||
@ -143,6 +153,8 @@ namespace netgen
|
||||
newgi = Project(newp);
|
||||
}
|
||||
|
||||
virtual bool IsMappedShape( const GeometryShape & other, const Transformation<3> & trafo, double tolerance ) const override;
|
||||
|
||||
protected:
|
||||
void RestrictHTrig(Mesh& mesh,
|
||||
const PointGeomInfo& gi0,
|
||||
|
@ -156,8 +156,6 @@ namespace netgen
|
||||
if(!pi0.IsValid() || !pi1.IsValid())
|
||||
continue;
|
||||
|
||||
if(pi1<pi0)
|
||||
Swap(pi0,pi1);
|
||||
m_ident.Add(pi0, pi1, n);
|
||||
}
|
||||
m_ident.SetType( n, identifications.GetType(n) );
|
||||
@ -166,6 +164,59 @@ namespace netgen
|
||||
return ret;
|
||||
}
|
||||
|
||||
// Add between identified surface elements (only consider closesurface identifications)
|
||||
void FillCloseSurface( MeshingData & md)
|
||||
{
|
||||
static Timer timer("FillCloseSurface"); RegionTimer rtimer(timer);
|
||||
|
||||
auto & mesh = md.mesh;
|
||||
auto & identifications = mesh->GetIdentifications();
|
||||
auto nmax = identifications.GetMaxNr();
|
||||
|
||||
bool have_closesurfaces = false;
|
||||
for(auto i : Range(1,nmax+1))
|
||||
if(identifications.GetType(i) == Identifications::CLOSESURFACES)
|
||||
have_closesurfaces = true;
|
||||
if(!have_closesurfaces)
|
||||
return;
|
||||
|
||||
NgArray<int, PointIndex::BASE> map;
|
||||
for(auto identnr : Range(1,nmax+1))
|
||||
{
|
||||
if(identifications.GetType(identnr) != Identifications::CLOSESURFACES)
|
||||
continue;
|
||||
|
||||
identifications.GetMap(identnr, map);
|
||||
|
||||
for(auto & sel : mesh->SurfaceElements())
|
||||
{
|
||||
bool is_mapped = true;
|
||||
for(auto pi : sel.PNums())
|
||||
if(!PointIndex(map[pi]).IsValid())
|
||||
is_mapped = false;
|
||||
|
||||
if(!is_mapped)
|
||||
continue;
|
||||
|
||||
// in case we have symmetric mapping (used in csg), only map in one direction
|
||||
if(map[map[sel[0]]] == sel[0] && map[sel[0]] < sel[0])
|
||||
continue;
|
||||
|
||||
// insert prism
|
||||
auto np = sel.GetNP();
|
||||
Element el(2*np);
|
||||
for(auto i : Range(np))
|
||||
{
|
||||
el[i] = sel[i];
|
||||
el[i+np] = map[sel[i]];
|
||||
}
|
||||
|
||||
el.SetIndex(md.domain);
|
||||
mesh->AddVolumeElement(el);
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
void CloseOpenQuads( MeshingData & md)
|
||||
{
|
||||
auto & mesh = *md.mesh;
|
||||
@ -488,6 +539,9 @@ namespace netgen
|
||||
if (md[i].mesh->CheckOverlappingBoundary())
|
||||
throw NgException ("Stop meshing since boundary mesh is overlapping");
|
||||
|
||||
// TODO: FillCloseSurface is not working with CSG closesurfaces
|
||||
if(md[i].mesh->GetGeometry()->GetGeomType() == Mesh::GEOM_OCC)
|
||||
FillCloseSurface( md[i] );
|
||||
CloseOpenQuads( md[i] );
|
||||
MeshDomain(md[i]);
|
||||
});
|
||||
|
@ -7,8 +7,9 @@
|
||||
|
||||
namespace netgen
|
||||
{
|
||||
OCCEdge::OCCEdge(TopoDS_Shape edge_)
|
||||
: tedge(edge_.TShape()),
|
||||
OCCEdge::OCCEdge(TopoDS_Shape edge_, GeometryVertex & start_, GeometryVertex & end_)
|
||||
: GeometryEdge(start_, end_),
|
||||
tedge(edge_.TShape()),
|
||||
edge(TopoDS::Edge(edge_))
|
||||
{
|
||||
curve = BRep_Tool::Curve(edge, s0, s1);
|
||||
@ -18,26 +19,13 @@ namespace netgen
|
||||
if(verts.size() != 2)
|
||||
throw Exception("OCC edge does not have 2 vertices");
|
||||
|
||||
start = OCCVertex(verts[0]);
|
||||
end = OCCVertex(verts[1]);
|
||||
|
||||
// swap start/end if necessary
|
||||
double d00 = Dist(GetPoint(0), start.GetPoint());
|
||||
double d01 = Dist(GetPoint(0), end.GetPoint());
|
||||
double d00 = Dist(GetPoint(0), start->GetPoint());
|
||||
double d01 = Dist(GetPoint(0), end->GetPoint());
|
||||
if(d01 < d00)
|
||||
swap(start, end);
|
||||
}
|
||||
|
||||
const GeometryVertex& OCCEdge::GetStartVertex() const
|
||||
{
|
||||
return start;
|
||||
}
|
||||
|
||||
const GeometryVertex& OCCEdge::GetEndVertex() const
|
||||
{
|
||||
return end;
|
||||
}
|
||||
|
||||
double OCCEdge::GetLength() const
|
||||
{
|
||||
return props.Mass();
|
||||
|
@ -14,6 +14,7 @@ namespace netgen
|
||||
{
|
||||
class OCCEdge : public GeometryEdge
|
||||
{
|
||||
public:
|
||||
T_Shape tedge;
|
||||
TopoDS_Edge edge;
|
||||
Handle(Geom_Curve) curve;
|
||||
@ -21,16 +22,11 @@ namespace netgen
|
||||
GProp_GProps props;
|
||||
|
||||
public:
|
||||
OCCVertex start;
|
||||
OCCVertex end;
|
||||
|
||||
OCCEdge(TopoDS_Shape edge_);
|
||||
OCCEdge(TopoDS_Shape edge_, GeometryVertex & start_, GeometryVertex & end_);
|
||||
|
||||
auto Shape() const { return edge; }
|
||||
T_Shape TShape() const { return tedge; }
|
||||
|
||||
const GeometryVertex& GetStartVertex() const override;
|
||||
const GeometryVertex& GetEndVertex() const override;
|
||||
double GetLength() const override;
|
||||
Point<3> GetCenter() const override;
|
||||
Point<3> GetPoint(double t) const override;
|
||||
|
@ -98,7 +98,6 @@ namespace netgen
|
||||
|
||||
// auto cof = curve_on_face[ORIENTATION][edgenr];
|
||||
auto edge = edge_on_face[ORIENTATION][edgenr];
|
||||
OCCEdge gedge(edge);
|
||||
double s0, s1;
|
||||
auto cof = BRep_Tool::CurveOnSurface (edge, face, s0, s1);
|
||||
|
||||
|
@ -11,6 +11,22 @@ namespace netgen
|
||||
return occ2ng( Handle(BRep_TVertex)::DownCast(shape)->Pnt() );
|
||||
}
|
||||
|
||||
Transformation<3> occ2ng (const gp_Trsf & occ_trafo)
|
||||
{
|
||||
Transformation<3> trafo;
|
||||
auto v = occ_trafo.TranslationPart();
|
||||
auto m = occ_trafo.VectorialPart();
|
||||
auto & tv = trafo.GetVector();
|
||||
auto & tm = trafo.GetMatrix();
|
||||
for(auto i : Range(3))
|
||||
{
|
||||
tv[i] = v.Coord(i+1);
|
||||
for(auto k : Range(3))
|
||||
tm(i,k) = m(i+1,k+1);
|
||||
}
|
||||
return trafo;
|
||||
}
|
||||
|
||||
Box<3> GetBoundingBox( const TopoDS_Shape & shape )
|
||||
{
|
||||
Bnd_Box bb;
|
||||
|
@ -9,6 +9,7 @@
|
||||
#include <TopTools_IndexedMapOfShape.hxx>
|
||||
#include <TopoDS.hxx>
|
||||
#include <TopoDS_Vertex.hxx>
|
||||
#include <gp_Trsf.hxx>
|
||||
|
||||
#include "meshing.hpp"
|
||||
|
||||
@ -47,6 +48,8 @@ namespace netgen
|
||||
return occ2ng (BRep_Tool::Pnt (v));
|
||||
}
|
||||
|
||||
DLL_HEADER Transformation<3> occ2ng (const gp_Trsf & t);
|
||||
|
||||
inline gp_Pnt ng2occ (const Point<3> & p)
|
||||
{
|
||||
return gp_Pnt(p(0), p(1), p(2));
|
||||
|
@ -1108,20 +1108,24 @@ namespace netgen
|
||||
fsingular = esingular = vsingular = false;
|
||||
|
||||
// Add shapes
|
||||
for(auto v : GetVertices(shape))
|
||||
for(auto i1 : Range(1, vmap.Extent()+1))
|
||||
{
|
||||
auto v = vmap(i1);
|
||||
auto tshape = v.TShape();
|
||||
if(vertex_map.count(tshape)!=0)
|
||||
continue;
|
||||
vertex_map[tshape] = vertices.Size();
|
||||
auto occ_vertex = make_unique<OCCVertex>(TopoDS::Vertex(v));
|
||||
occ_vertex->nr = vertices.Size();
|
||||
vertex_map[tshape] = occ_vertex->nr;
|
||||
|
||||
if(global_shape_properties.count(tshape)>0)
|
||||
occ_vertex->properties = global_shape_properties[tshape];
|
||||
vertices.Append(std::move(occ_vertex));
|
||||
}
|
||||
|
||||
for(auto e : GetEdges(shape))
|
||||
for(auto i1 : Range(1, emap.Extent()+1))
|
||||
{
|
||||
auto e = emap(i1);
|
||||
auto tshape = e.TShape();
|
||||
auto edge = TopoDS::Edge(e);
|
||||
if(edge_map.count(tshape)!=0)
|
||||
@ -1129,15 +1133,15 @@ namespace netgen
|
||||
if(BRep_Tool::Degenerated(edge))
|
||||
continue;
|
||||
edge_map[tshape] = edges.Size();
|
||||
auto occ_edge = make_unique<OCCEdge>(edge);
|
||||
auto verts = GetVertices(e);
|
||||
auto occ_edge = make_unique<OCCEdge>(edge, *vertices[vertex_map[verts[0].TShape()]], *vertices[vertex_map[verts[1].TShape()]] );
|
||||
occ_edge->properties = global_shape_properties[tshape];
|
||||
occ_edge->start.nr = vertex_map[occ_edge->start.TShape()];
|
||||
occ_edge->end.nr = vertex_map[occ_edge->end.TShape()];
|
||||
edges.Append(std::move(occ_edge));
|
||||
}
|
||||
|
||||
for(auto f : GetFaces(shape))
|
||||
for(auto i1 : Range(1, fmap.Extent()+1))
|
||||
{
|
||||
auto f = fmap(i1);
|
||||
auto tshape = f.TShape();
|
||||
if(face_map.count(tshape)==0)
|
||||
{
|
||||
@ -1145,6 +1149,10 @@ namespace netgen
|
||||
auto k = faces.Size();
|
||||
face_map[tshape] = k;
|
||||
auto occ_face = make_unique<OCCFace>(f);
|
||||
|
||||
for(auto e : GetEdges(f))
|
||||
occ_face->edges.Append( edges[edge_map[e.TShape()]].get() );
|
||||
|
||||
if(global_shape_properties.count(tshape)>0)
|
||||
occ_face->properties = global_shape_properties[tshape];
|
||||
faces.Append(std::move(occ_face));
|
||||
@ -1162,8 +1170,9 @@ namespace netgen
|
||||
}
|
||||
|
||||
|
||||
for(auto s : GetSolids(shape))
|
||||
for(auto i1 : Range(1, somap.Extent()+1))
|
||||
{
|
||||
auto s = somap(i1);
|
||||
auto tshape = s.TShape();
|
||||
int k;
|
||||
if(solid_map.count(tshape)==0)
|
||||
@ -1194,6 +1203,9 @@ namespace netgen
|
||||
if(identifications.count(tshape))
|
||||
for(auto & ident : identifications[tshape])
|
||||
{
|
||||
if(shape_map.count(ident.from)==0 || shape_map.count(ident.to)==0)
|
||||
continue;
|
||||
|
||||
ShapeIdentification si{
|
||||
shapes[shape_map[ident.from]].get(),
|
||||
shapes[shape_map[ident.to]].get(),
|
||||
@ -1208,9 +1220,8 @@ namespace netgen
|
||||
add_identifications( edges, edge_map );
|
||||
add_identifications( faces, face_map );
|
||||
|
||||
ProcessIdentifications();
|
||||
|
||||
bounding_box = ::netgen::GetBoundingBox( shape );
|
||||
ProcessIdentifications();
|
||||
}
|
||||
|
||||
|
||||
@ -1940,25 +1951,6 @@ namespace netgen
|
||||
return occ2ng( props.CentreOfMass() );
|
||||
}
|
||||
|
||||
void OCCGeometry :: IdentifyEdges(const TopoDS_Shape & me, const TopoDS_Shape & you, string name, Identifications::ID_TYPE type)
|
||||
{
|
||||
auto cme = GetCenter(me);
|
||||
auto cyou = GetCenter(you);
|
||||
Transformation<3> trafo{cyou-cme};
|
||||
identifications[me.TShape()].push_back( {me.TShape(), you.TShape(), Transformation<3>(cyou - cme), name, type} );
|
||||
|
||||
auto vme = GetVertices(me);
|
||||
auto vyou = GetVertices(you);
|
||||
Point<3> pme0 = trafo(occ2ng(vme[0]));
|
||||
Point<3> pme1 = trafo(occ2ng(vme[1]));
|
||||
Point<3> pyou = occ2ng(vyou[0]);
|
||||
|
||||
bool do_swap = Dist(pme1, pyou) < Dist(pme0, pyou);
|
||||
|
||||
for(auto i : Range(2))
|
||||
identifications[vme[i].TShape()].push_back( {vme[i].TShape(), vyou[do_swap ? 1-i : i].TShape(), trafo, name, type} );
|
||||
}
|
||||
|
||||
bool IsMappedShape(const Transformation<3> & trafo, const TopoDS_Shape & me, const TopoDS_Shape & you)
|
||||
{
|
||||
if(me.ShapeType() != you.ShapeType()) return false;
|
||||
@ -1976,6 +1968,11 @@ namespace netgen
|
||||
std::map<T_Shape, T_Shape> vmap;
|
||||
|
||||
auto verts_me = GetVertices(me);
|
||||
auto verts_you = GetVertices(you);
|
||||
|
||||
if(verts_me.size() != verts_you.size())
|
||||
return false;
|
||||
|
||||
for (auto i : Range(verts_me.size()))
|
||||
{
|
||||
auto s = verts_me[i].TShape();
|
||||
@ -1986,8 +1983,7 @@ namespace netgen
|
||||
vmap[s] = nullptr;
|
||||
}
|
||||
|
||||
bool all_verts_mapped = true;
|
||||
for (auto vert : GetVertices(you))
|
||||
for (auto vert : verts_you)
|
||||
{
|
||||
auto s = vert.TShape();
|
||||
auto p = occ2ng(s);
|
||||
@ -1998,30 +1994,62 @@ namespace netgen
|
||||
return true;
|
||||
});
|
||||
if(!vert_mapped)
|
||||
{
|
||||
all_verts_mapped = false;
|
||||
break;
|
||||
}
|
||||
return false;
|
||||
}
|
||||
return all_verts_mapped;
|
||||
return true;
|
||||
}
|
||||
|
||||
void OCCGeometry :: IdentifyFaces(const TopoDS_Shape & solid, const TopoDS_Shape & me, const TopoDS_Shape & you, string name, Identifications::ID_TYPE type)
|
||||
void Identify(const TopoDS_Shape & me, const TopoDS_Shape & you, string name, Identifications::ID_TYPE type, std::optional<gp_Trsf> opt_trafo)
|
||||
{
|
||||
auto cme = GetCenter(me);
|
||||
auto cyou = GetCenter(you);
|
||||
Transformation<3> trafo(cyou-cme);
|
||||
gp_Trsf trafo;
|
||||
if(opt_trafo)
|
||||
{
|
||||
trafo = *opt_trafo;
|
||||
}
|
||||
else
|
||||
{
|
||||
auto v = GetCenter(you) - GetCenter(me);
|
||||
trafo.SetTranslation(gp_Vec(v[0], v[1], v[2]));
|
||||
}
|
||||
|
||||
identifications[me.TShape()].push_back
|
||||
(OCCIdentification { me.TShape(), you.TShape(), trafo, name, type });
|
||||
ListOfShapes list_me, list_you;
|
||||
list_me.push_back(me);
|
||||
list_you.push_back(you);
|
||||
Identify(list_me, list_you, name, type, trafo);
|
||||
}
|
||||
|
||||
auto edges_me = GetEdges(me);
|
||||
auto edges_you = GetEdges(you);
|
||||
void Identify(const ListOfShapes & me, const ListOfShapes & you, string name, Identifications::ID_TYPE type, gp_Trsf occ_trafo)
|
||||
{
|
||||
Transformation<3> trafo = occ2ng(occ_trafo);
|
||||
|
||||
for (auto e_me : edges_me)
|
||||
for (auto e_you : edges_you)
|
||||
if(IsMappedShape(trafo, e_me, e_you))
|
||||
IdentifyEdges(e_me, e_you, name, type);
|
||||
ListOfShapes id_me;
|
||||
ListOfShapes id_you;
|
||||
|
||||
if(auto faces_me = me.Faces(); faces_me.size()>0)
|
||||
{
|
||||
id_me = faces_me;
|
||||
id_you = you.Faces();
|
||||
}
|
||||
else if(auto edges_me = me.Edges(); edges_me.size()>0)
|
||||
{
|
||||
id_me = edges_me;
|
||||
id_you = you.Edges();
|
||||
}
|
||||
else
|
||||
{
|
||||
id_me = me.Vertices();
|
||||
id_you = you.Vertices();
|
||||
}
|
||||
|
||||
for(auto shape_me : id_me)
|
||||
for(auto shape_you : id_you)
|
||||
{
|
||||
if(!IsMappedShape(trafo, shape_me, shape_you))
|
||||
continue;
|
||||
|
||||
OCCGeometry::identifications[shape_me.TShape()].push_back
|
||||
(OCCIdentification { shape_me.TShape(), shape_you.TShape(), trafo, name, type });
|
||||
}
|
||||
}
|
||||
|
||||
void OCCParameters :: Print(ostream & ost) const
|
||||
|
@ -44,6 +44,7 @@ namespace netgen
|
||||
#define OCCGEOMETRYVISUALIZATIONFULLCHANGE 1 // Compute transformation matrices and redraw
|
||||
#define OCCGEOMETRYVISUALIZATIONHALFCHANGE 2 // Redraw
|
||||
|
||||
bool IsMappedShape(const Transformation<3> & trafo, const TopoDS_Shape & me, const TopoDS_Shape & you);
|
||||
|
||||
class EntityVisualizationCode
|
||||
{
|
||||
@ -341,12 +342,13 @@ namespace netgen
|
||||
bool ErrorInSurfaceMeshing ();
|
||||
|
||||
// void WriteOCC_STL(char * filename);
|
||||
static void IdentifyEdges(const TopoDS_Shape & me, const TopoDS_Shape & you, string name, Identifications::ID_TYPE type);
|
||||
static void IdentifyFaces(const TopoDS_Shape & solid,const TopoDS_Shape & me, const TopoDS_Shape & you, string name, Identifications::ID_TYPE type);
|
||||
|
||||
private:
|
||||
//bool FastProject (int surfi, Point<3> & ap, double& u, double& v) const;
|
||||
};
|
||||
|
||||
void Identify(const ListOfShapes & me, const ListOfShapes & you, string name, Identifications::ID_TYPE type, gp_Trsf occ_trafo);
|
||||
void Identify(const TopoDS_Shape & me, const TopoDS_Shape & you, string name, Identifications::ID_TYPE type, std::optional<gp_Trsf> opt_trafo);
|
||||
|
||||
|
||||
void PrintContents (OCCGeometry * geom);
|
||||
|
@ -199,14 +199,14 @@ py::object CastShape(const TopoDS_Shape & s)
|
||||
template <class TBuilder>
|
||||
void PropagateIdentifications (TBuilder & builder, TopoDS_Shape shape)
|
||||
{
|
||||
std::map<T_Shape, T_Shape> mod_map;
|
||||
std::map<T_Shape, std::set<T_Shape>> mod_map;
|
||||
std::map<T_Shape, bool> tshape_handled;
|
||||
|
||||
for (auto typ : { TopAbs_SOLID, TopAbs_FACE, TopAbs_EDGE, TopAbs_VERTEX })
|
||||
for (TopExp_Explorer e(shape, typ); e.More(); e.Next())
|
||||
{
|
||||
auto tshape = e.Current().TShape();
|
||||
mod_map[tshape] = tshape;
|
||||
mod_map[tshape].insert(tshape);
|
||||
tshape_handled[tshape] = false;
|
||||
}
|
||||
|
||||
@ -215,11 +215,8 @@ void PropagateIdentifications (TBuilder & builder, TopoDS_Shape shape)
|
||||
{
|
||||
auto tshape = e.Current().TShape();
|
||||
|
||||
auto & modified = builder.Modified(e.Current());
|
||||
if(modified.Size()!=1)
|
||||
continue;
|
||||
|
||||
mod_map[tshape] = modified.First().TShape();
|
||||
for (auto mods : builder.Modified(e.Current()))
|
||||
mod_map[tshape].insert(mods.TShape());
|
||||
}
|
||||
|
||||
for (auto typ : { TopAbs_SOLID, TopAbs_FACE, TopAbs_EDGE, TopAbs_VERTEX })
|
||||
@ -236,21 +233,33 @@ void PropagateIdentifications (TBuilder & builder, TopoDS_Shape shape)
|
||||
|
||||
auto tshape_mapped = mod_map[tshape];
|
||||
|
||||
for(auto & ident : OCCGeometry::identifications[tshape])
|
||||
for(auto ident : OCCGeometry::identifications[tshape])
|
||||
{
|
||||
// update existing identification
|
||||
if(tshape == tshape_mapped)
|
||||
{
|
||||
ident.to = mod_map[ident.to];
|
||||
ident.from = mod_map[ident.from];
|
||||
}
|
||||
else
|
||||
{
|
||||
OCCIdentification id_new = ident;
|
||||
id_new.to = mod_map[id_new.to];
|
||||
id_new.from = mod_map[id_new.from];
|
||||
OCCGeometry::identifications[mod_map[tshape]].push_back(id_new);
|
||||
}
|
||||
// nothing happened
|
||||
if(mod_map[ident.to].size()==1 && mod_map[ident.from].size() ==1)
|
||||
continue;
|
||||
|
||||
auto from = ident.from;
|
||||
auto to = ident.to;
|
||||
|
||||
for(auto from_mapped : mod_map[from])
|
||||
for(auto to_mapped : mod_map[to])
|
||||
{
|
||||
if(from==from_mapped && to==to_mapped)
|
||||
continue;
|
||||
|
||||
TopoDS_Shape s_from; s_from.TShape(from_mapped);
|
||||
TopoDS_Shape s_to; s_to.TShape(to_mapped);
|
||||
|
||||
if(!IsMappedShape(ident.trafo, s_from, s_to))
|
||||
continue;
|
||||
|
||||
OCCIdentification id_new = ident;
|
||||
id_new.to = to_mapped;
|
||||
id_new.from = from_mapped;
|
||||
auto id_owner = from == tshape ? from_mapped : to_mapped;
|
||||
OCCGeometry::identifications[id_owner].push_back(id_new);
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
@ -1071,27 +1080,11 @@ DLL_HEADER void ExportNgOCCShapes(py::module &m)
|
||||
BRepMesh_IncrementalMesh (shape, deflection, true);
|
||||
})
|
||||
|
||||
.def("Identify", [](const TopoDS_Shape & me, const TopoDS_Shape & you, string name, Identifications::ID_TYPE idtype) {
|
||||
// only edges supported, by now
|
||||
auto type = me.ShapeType();
|
||||
auto tyou = you.ShapeType();
|
||||
if(type != tyou)
|
||||
throw NgException ("Identify: cannot identify different shape types");
|
||||
|
||||
switch(type)
|
||||
{
|
||||
case TopAbs_VERTEX:
|
||||
case TopAbs_EDGE:
|
||||
OCCGeometry::IdentifyEdges(me, you, name, idtype);
|
||||
break;
|
||||
default:
|
||||
throw NgException ("Identify: unsupported shape type");
|
||||
break;
|
||||
}
|
||||
}, py::arg("other"), py::arg("name"), py::arg("type")=Identifications::PERIODIC, "Identify shapes for periodic meshing")
|
||||
|
||||
.def("Identify", OCCGeometry::IdentifyFaces, "Identify faces",
|
||||
py::arg("from"), py::arg("to"), py::arg("name"), py::arg("type")=Identifications::PERIODIC)
|
||||
.def("Identify", py::overload_cast<const TopoDS_Shape &, const TopoDS_Shape &, string, Identifications::ID_TYPE, std::optional<gp_Trsf>>(&Identify),
|
||||
py::arg("other"), py::arg("name"),
|
||||
py::arg("type")=Identifications::PERIODIC, py::arg("trafo")=nullopt,
|
||||
"Identify shapes for periodic meshing")
|
||||
|
||||
.def("Distance", [](const TopoDS_Shape& self,
|
||||
const TopoDS_Shape& other)
|
||||
@ -1634,6 +1627,11 @@ DLL_HEADER void ExportNgOCCShapes(py::module &m)
|
||||
}
|
||||
}, "set hpref for all elements of list")
|
||||
|
||||
.def("Identify", py::overload_cast<const ListOfShapes&, const ListOfShapes&, string, Identifications::ID_TYPE, gp_Trsf>(&Identify),
|
||||
py::arg("other"), py::arg("name"),
|
||||
py::arg("type")=Identifications::PERIODIC, py::arg("trafo"),
|
||||
"Identify shapes for periodic meshing")
|
||||
|
||||
;
|
||||
|
||||
|
||||
|
76
tests/pytest/test_occ_identifications.py
Normal file
76
tests/pytest/test_occ_identifications.py
Normal file
@ -0,0 +1,76 @@
|
||||
import pytest
|
||||
|
||||
from netgen.meshing import IdentificationType
|
||||
idtype = IdentificationType.CLOSESURFACES
|
||||
|
||||
def test_two_boxes():
|
||||
occ = pytest.importorskip("netgen.occ")
|
||||
inner = occ.Box((0,0,0), (1,1,1))
|
||||
trafo = occ.gp_Trsf().Scale(inner.center, 1.1)
|
||||
outer = trafo(inner)
|
||||
|
||||
inner.Identify(outer, "", idtype, trafo)
|
||||
shape = occ.Glue([outer-inner, inner])
|
||||
|
||||
geo = occ.OCCGeometry(shape)
|
||||
mesh = geo.GenerateMesh(maxh=0.3)
|
||||
have_prisms = False
|
||||
|
||||
for el in mesh.Elements3D():
|
||||
if len(el.vertices)==6:
|
||||
have_prisms = True
|
||||
break
|
||||
|
||||
assert have_prisms
|
||||
|
||||
def test_two_circles():
|
||||
occ = pytest.importorskip("netgen.occ")
|
||||
circ1 = occ.WorkPlane().Circle(1).Face()
|
||||
trafo = occ.gp_Trsf().Scale(circ1.center, 1.1)
|
||||
|
||||
circ2 = trafo(circ1)
|
||||
circ1.edges[0].Identify(circ2.edges[0], "", idtype, trafo)
|
||||
circ2 -= circ1
|
||||
shape = occ.Glue([circ1, circ2])
|
||||
|
||||
geo = occ.OCCGeometry(shape, 2)
|
||||
mesh = geo.GenerateMesh(maxh=0.2)
|
||||
have_quads = False
|
||||
|
||||
for el in mesh.Elements2D():
|
||||
if len(el.vertices)==4:
|
||||
have_quads = True
|
||||
break
|
||||
|
||||
assert have_quads
|
||||
|
||||
def test_cut_identified_face():
|
||||
occ = pytest.importorskip("netgen.occ")
|
||||
from netgen.occ import Z, Box, Cylinder, Glue, OCCGeometry
|
||||
box = Box((-1,-1,0), (1,1,1))
|
||||
cyl = Cylinder( (0,0,0), Z, 0.5, 1 )
|
||||
|
||||
box.faces.Min(Z).Identify(box.faces.Max(Z), "", idtype)
|
||||
shape = Glue([cyl, box])
|
||||
geo = OCCGeometry(shape)
|
||||
mesh = geo.GenerateMesh(maxh=0.5)
|
||||
|
||||
for el in mesh.Elements3D():
|
||||
assert len(el.vertices)==6
|
||||
|
||||
def test_identify_multiple_faces():
|
||||
occ = pytest.importorskip("netgen.occ")
|
||||
from netgen.occ import Z, Box, Cylinder, Glue, OCCGeometry, gp_Trsf
|
||||
box = Box((-1,-1,0), (1,1,1))
|
||||
cyl = Cylinder( (0,0,0), Z, 0.5, 1 )
|
||||
|
||||
shape = Glue([box, cyl])
|
||||
bot_faces = shape.faces[Z < 0.1]
|
||||
top_faces = shape.faces[Z > 0.1]
|
||||
bot_faces.Identify(top_faces, "", idtype, gp_Trsf.Translation((0,0,1)))
|
||||
|
||||
geo = OCCGeometry(shape)
|
||||
mesh = geo.GenerateMesh(maxh=0.3)
|
||||
|
||||
for el in mesh.Elements3D():
|
||||
assert len(el.vertices)==6
|
Loading…
Reference in New Issue
Block a user