Commit b315066d authored by Michael Turner's avatar Michael Turner
Browse files

clean up of cad system

parent e4282905
......@@ -39,10 +39,13 @@
using namespace std;
namespace Nektar {
namespace NekMeshUtils {
namespace Nektar
{
namespace NekMeshUtils
{
CADSurf::CADSurf(int i, TopoDS_Shape in, vector<EdgeLoop> ein) : m_ID(i), m_edges(ein)
CADSurf::CADSurf(int i, TopoDS_Shape in, vector<EdgeLoop> ein)
: m_ID(i), m_edges(ein)
{
// this bit of code changes the units of the cad from mm opencascade
// defualt to m
......@@ -100,28 +103,35 @@ Array<OneD, NekDouble> CADSurf::locuv(Array<OneD, NekDouble> p)
if(projection.Distance(1) > 1.0)
{
stringstream ss;
cerr << "large locuv distance " << projection.Distance(1)/1000.0 << endl;
cerr << "large locuv distance "
<< projection.Distance(1)/1000.0 << endl;
}
}
//if the uv returned is slightly off the surface
//(which ShapeAnalysis_Surface can do sometimes)
if(uvr[0] < m_occSurface.FirstUParameter() ||
uvr[0] > m_occSurface.LastUParameter() ||
uvr[1] < m_occSurface.FirstVParameter() ||
uvr[1] > m_occSurface.LastVParameter())
{
if(uvr[0] < m_occSurface.FirstUParameter() && fabs(m_occSurface.FirstUParameter() - uvr[0]) < 1E-6)
if(uvr[0] < m_occSurface.FirstUParameter() &&
fabs(m_occSurface.FirstUParameter() - uvr[0]) < 1E-6)
{
uvr[0] = m_occSurface.FirstUParameter();
}
else if(uvr[0] > m_occSurface.LastUParameter() && fabs(m_occSurface.LastUParameter() - uvr[0]) < 1E-6)
else if(uvr[0] > m_occSurface.LastUParameter() &&
fabs(m_occSurface.LastUParameter() - uvr[0]) < 1E-6)
{
uvr[0] = m_occSurface.LastUParameter();
}
else if(uvr[1] < m_occSurface.FirstVParameter() && fabs(m_occSurface.FirstVParameter() - uvr[1]) < 1E-6)
else if(uvr[1] < m_occSurface.FirstVParameter() &&
fabs(m_occSurface.FirstVParameter() - uvr[1]) < 1E-6)
{
uvr[1] = m_occSurface.FirstVParameter();
}
else if(uvr[1] > m_occSurface.LastVParameter() && fabs(m_occSurface.LastVParameter() - uvr[1]) < 1E-6)
else if(uvr[1] > m_occSurface.LastVParameter() &&
fabs(m_occSurface.LastVParameter() - uvr[1]) < 1E-6)
{
uvr[1] = m_occSurface.LastVParameter();
}
......@@ -316,7 +326,8 @@ void CADSurf::Test(Array<OneD, NekDouble> uv)
{
if(fabs(uv[0]-m_occSurface.FirstUParameter()) > 1E-8)
{
error << "U(" << uv[0] << ") is less than Umin(" << m_occSurface.FirstUParameter() << ")";
error << "U(" << uv[0] << ") is less than Umin("
<< m_occSurface.FirstUParameter() << ")";
passed = false;
}
}
......@@ -324,7 +335,8 @@ void CADSurf::Test(Array<OneD, NekDouble> uv)
{
if(fabs(uv[0]-m_occSurface.LastUParameter()) > 1E-8)
{
error << "U(" << uv[0] << ") is greater than Umax(" << m_occSurface.LastUParameter() << ")";
error << "U(" << uv[0] << ") is greater than Umax("
<< m_occSurface.LastUParameter() << ")";
passed = false;
}
}
......@@ -332,7 +344,8 @@ void CADSurf::Test(Array<OneD, NekDouble> uv)
{
if(fabs(uv[1]-m_occSurface.FirstVParameter()) > 1E-8)
{
error << "V(" << uv[1] << ") is less than Vmin(" << m_occSurface.FirstVParameter() << ")";
error << "V(" << uv[1] << ") is less than Vmin("
<< m_occSurface.FirstVParameter() << ")";
passed = false;
}
}
......@@ -340,7 +353,8 @@ void CADSurf::Test(Array<OneD, NekDouble> uv)
{
if(fabs(uv[1]-m_occSurface.LastVParameter()) > 1E-8)
{
error << "V(" << uv[1] << ") is greater than Vmax(" << m_occSurface.LastVParameter() << ")";
error << "V(" << uv[1] << ") is greater than Vmax("
<< m_occSurface.LastVParameter() << ")";
passed = false;
}
}
......
......@@ -157,23 +157,38 @@ public:
m_correctNormal = false;
}
/**
* @brief surface needs to know if it is bounded by only two curves
*/
void SetTwoC()
{
m_hasTwoCurves = true;
}
/**
* @brief query two curves
*/
bool GetTwoC(){return m_hasTwoCurves;}
/**
* @brief return id
*/
int GetId(){return m_ID;}
/**
* @brief does unconstrained locuv to project point from anywhere
*/
NekDouble DistanceTo(Array<OneD, NekDouble> p);
/**
* @brief returns curvature at point uv
*/
NekDouble Curvature(Array<OneD, NekDouble> uv);
private:
/// Function which tests the the value of uv used is within the surface
void Test(Array<OneD, NekDouble> uv);
/// ID of surface.
int m_ID;
/// normal
......
......@@ -67,7 +67,9 @@ void CADSystem::SmallFeatureAnalysis(NekDouble min)
map<int, CADCurveSharedPtr>::iterator it;
for(it = m_curves.begin(); it != m_curves.end(); it++)
{
ASSERTL0(it->second->GetTotLength() > 1e-5,"curve too small for meshing");
ASSERTL0(it->second->GetTotLength() > 1e-5,
"curve too small for meshing");
if(it->second->GetTotLength() < 5.0 * min)
{
lens.push_back(pair<int,NekDouble>(
......@@ -81,7 +83,8 @@ void CADSystem::SmallFeatureAnalysis(NekDouble min)
cout << "\tWARNING" << endl;
for(int i = 0; i < lens.size(); i++)
{
cout << "\tCurve: " << lens[i].first << " has length: " << lens[i].second << endl;
cout << "\tCurve: " << lens[i].first << " has length: "
<< lens[i].second << endl;
}
if(lens.size() == 1)
cout << "\tThis curve's length is less than or close to minDelta" << endl;
......@@ -222,8 +225,9 @@ bool CADSystem::LoadCAD()
BRep_Tool::UVPoints(TopoDS::Edge(edge),TopoDS::Face(face),p1,p2);
m_verts[mapOfVerts.FindIndex(TopExp::FirstVertex(TopoDS::Edge(edge), Standard_True))]->
SetDegen(i, (p1.X() + p2.X())/2.0, (p1.Y() + p2.Y())/2.0);
m_verts[mapOfVerts.FindIndex(TopExp::FirstVertex(
TopoDS::Edge(edge), Standard_True))]->SetDegen(
i, (p1.X() + p2.X())/2.0, (p1.Y() + p2.Y())/2.0);
}
......@@ -274,11 +278,13 @@ bool CADSystem::LoadCAD()
if(wire != ow)
{
BRepBuilderAPI_MakeFace build(BRep_Tool::Surface(TopoDS::Face(face)),1e-7);
BRepBuilderAPI_MakeFace build(BRep_Tool::Surface(
TopoDS::Face(face)),1e-7);
build.Add(TopoDS::Wire(wire));
TopoDS_Shape newface = build.Shape();
wirefacecuts.push_back(newface);
BRepAdaptor_Surface b = BRepAdaptor_Surface(TopoDS::Face(newface));
BRepAdaptor_Surface b = BRepAdaptor_Surface(
TopoDS::Face(newface));
NekDouble u,v;
u = (b.LastUParameter()-b.FirstUParameter())/2.0;
v = (b.LastVParameter()-b.FirstVParameter())/2.0;
......@@ -291,6 +297,7 @@ bool CADSystem::LoadCAD()
{
if(j ==k) continue;
///TODO fix this test
BRepClass_FaceClassifier fc(TopoDS::Face(wirefacecuts[j]), centersofcutfaces[k], 1e-7);
//ASSERTL0(fc.State() == 1, "Internal face loops make this cad impossible to mesh");
}
......@@ -409,22 +416,6 @@ void CADSystem::AddSurf(int i, TopoDS_Shape in, vector<EdgeLoop> ein)
{
m_surfs[i]->SetTwoC();
}
/*//check the normals are interior
Array<OneD, NekDouble> bnds = newSurf->GetBounds();
Array<OneD, NekDouble> uv(2);
uv[0] = (bnds[1]+bnds[0])/2.0;
uv[1] = (bnds[3]+bnds[2])/2.0;
Array<OneD, NekDouble> N = newSurf->N(uv);
Array<OneD, NekDouble> P = newSurf->P(uv);
Array<OneD, NekDouble> tl(3);
tl[0] = P[0] + 0.01*N[0];
tl[1] = P[1] + 0.01*N[1];
tl[2] = P[2] + 0.01*N[2];
ASSERTL0(InsideShape(tl), "normal incorrectly orientated");*/
}
bool CADSystem::InsideShape(Array<OneD, NekDouble> loc)
......
......@@ -138,11 +138,17 @@ public:
return search->second;
}
/**
* @brief Gets map of all verts
*/
std::map<int, CADVertSharedPtr> GetVerts()
{
return m_verts;
}
/**
* @brief Gets number of vertices
*/
int GetNumVerts(){return m_verts.size();}
/**
......@@ -153,9 +159,13 @@ public:
*/
bool InsideShape(Array<OneD, NekDouble> loc);
/**
* @brief tests for features which are smaller that the specified mesh limit
*/
void SmallFeatureAnalysis(NekDouble min);
private:
/// Private function to add curve to CADSystem::m_verts.
void AddVert(int i, TopoDS_Shape in);
/// Private function to add curve to CADSystem::m_curves.
......@@ -172,8 +182,6 @@ private:
std::map<int, CADVertSharedPtr> m_verts;
/// occ master object
TopoDS_Shape shape;
bool smallfeatreAdjust;
};
typedef boost::shared_ptr<CADSystem> CADSystemSharedPtr;
......
......@@ -80,6 +80,9 @@ public:
degen = false;
}
/**
* @brief Get x,y,z location of the vertex
*/
Array<OneD, NekDouble> GetLoc()
{
Array<OneD, NekDouble> out(3);
......@@ -87,10 +90,19 @@ public:
return out;
}
/**
* @brief Return ID of the vertex
*/
int GetId(){return m_ID;}
/**
* @brief returns a node object of the cad vertex
*/
NodeSharedPtr GetNode(){return m_node;}
/**
* @brief if the vertex is degenerate manually set uv for that surface
*/
void SetDegen(int s, NekDouble u, NekDouble v)
{
degen = true;
......@@ -101,6 +113,9 @@ public:
m_node->SetCADSurf(s,uv);
}
/**
* @brief query is degenerate
*/
int IsDegen()
{
if(degen)
......
......@@ -212,6 +212,7 @@ void Octree::Build()
}
}
cout << Octants.size() << endl;
cout << minlimitedoct.size() << " octants contain min limited points" << endl;
vector<OctantSharedPtr> neighRevaluate;
......@@ -244,7 +245,7 @@ void Octree::Build()
}
}
if(maxdiff/mindiff > 1.5)
if(maxdiff/mindiff > 2.0)
{
vector<OctantSharedPtr> np = minlimitedoct[i]->GetNeighbourList();
neighRevaluate.insert(neighRevaluate.end(), np.begin(), np.end());
......@@ -252,17 +253,24 @@ void Octree::Build()
}
else
{
if(av < m_minDelta/5.0) av = m_minDelta/3.0;
minlimitedoct[i]->SetDelta(av);
}
}
cout << "done, neighbours" << endl;
cout << Octants.size() << endl;
cout << "done, neighbours " << neighRevaluate.size() << endl;
set<int> done;
for(int i = 0; i < neighRevaluate.size(); i++)
{
if(neighRevaluate[i]->IsLeaf())
{
AssignNeigbours(neighRevaluate[i]);
}
set<int>::iterator s = done.find(neighRevaluate[i]->GetId());
if(s == done.end())
if(neighRevaluate[i]->IsLeaf())
{
AssignNeigbours(neighRevaluate[i]);
done.insert(neighRevaluate[i]->GetId());
}
}
cout << "smoothing" << endl;
......@@ -372,6 +380,7 @@ void Octree::SubDivideMinLimited(OctantSharedPtr parent, vector<OctantSharedPtr>
}
if(children[i]->IsDeltaKnown())
{
if(av < m_minDelta/3.0) av = m_minDelta/3.0;
children[i]->SetDelta(av);
}
else
......@@ -379,12 +388,12 @@ void Octree::SubDivideMinLimited(OctantSharedPtr parent, vector<OctantSharedPtr>
children[i]->SetDelta(parent->GetDelta());
}
if(children[i]->GetDelta() < m_minDelta/10.0)
children[i]->SetDelta(m_minDelta/10.0);
if(children[i]->GetDelta() < m_minDelta/3.0)
children[i]->SetDelta(3.0);
if(maxdiff/mindiff > 1.5)
if(maxdiff/mindiff > 2.0)
{
if(children[i]->DX() / 2.0 > m_minDelta/10.0)
if(children[i]->DX() / 2.0 > m_minDelta/3.0)
{
SubDivideMinLimited(children[i], np);
}
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment