Commit 154e7e02 authored by Michael Turner's avatar Michael Turner
Browse files

fix 3d boundary layers

parent 271579a6
......@@ -165,7 +165,9 @@ Array<OneD, NekDouble> CADSurfOCE::locuv(Array<OneD, NekDouble> p)
NekDouble CADSurfOCE::Curvature(Array<OneD, NekDouble> uv)
{
#if defined(NEKTAR_DEBUG)
Test(uv);
#endif
Array<OneD, NekDouble> n = N(uv);
......@@ -245,7 +247,9 @@ void CADSurfOCE::ProjectTo(Array<OneD, NekDouble> &tp, Array<OneD, NekDouble> &u
Array<OneD, NekDouble> CADSurfOCE::P(Array<OneD, NekDouble> uv)
{
#if defined(NEKTAR_DEBUG)
Test(uv);
#endif
Array<OneD, NekDouble> location(3);
gp_Pnt loc;
......@@ -258,7 +262,9 @@ Array<OneD, NekDouble> CADSurfOCE::P(Array<OneD, NekDouble> uv)
Array<OneD, NekDouble> CADSurfOCE::N(Array<OneD, NekDouble> uv)
{
#if defined(NEKTAR_DEBUG)
Test(uv);
#endif
Array<OneD, NekDouble> normal(3);
gp_Pnt Loc;
......@@ -291,7 +297,9 @@ Array<OneD, NekDouble> CADSurfOCE::N(Array<OneD, NekDouble> uv)
Array<OneD, NekDouble> CADSurfOCE::D1(Array<OneD, NekDouble> uv)
{
#if defined(NEKTAR_DEBUG)
Test(uv);
#endif
Array<OneD, NekDouble> r(9);
gp_Pnt Loc;
......@@ -313,7 +321,9 @@ Array<OneD, NekDouble> CADSurfOCE::D1(Array<OneD, NekDouble> uv)
Array<OneD, NekDouble> CADSurfOCE::D2(Array<OneD, NekDouble> uv)
{
#if defined(NEKTAR_DEBUG)
Test(uv);
#endif
Array<OneD, NekDouble> r(18);
gp_Pnt Loc;
......@@ -388,10 +398,7 @@ void CADSurfOCE::Test(Array<OneD, NekDouble> uv)
}
error << " On Surface: " << GetId();
if (!passed)
{
cout << "Warning: " << error.str() << endl;
}
ASSERTL1(passed, "Warning: " + error.str());
}
}
}
......@@ -100,8 +100,6 @@ public:
unsigned int m_spaceDim;
/// a order tag to aid output, a bit of a hack
unsigned int m_nummode;
///
unsigned int m_numcomp;
/// List of mesh nodes.
std::vector<NodeSharedPtr> m_node;
/// Set of element vertices.
......
......@@ -430,8 +430,25 @@ void FaceMesh::Smoothing()
u0[0] -= (dF[0] * F[0] + dF[2] * F[1]);
u0[1] -= (dF[1] * F[0] + dF[3] * F[1]);
if (!(u0[0] < bounds[0] || u0[0] > bounds[1] || u0[1] < bounds[2] ||
u0[1] > bounds[3]))
bool inbounds = true;
if (u0[0] < bounds[0])
{
inbounds = false;
}
else if (u0[0] > bounds[1])
{
inbounds = false;
}
else if (u0[1] < bounds[2])
{
inbounds = false;
}
else if (u0[1] > bounds[3])
{
inbounds = false;
}
if (!inbounds)
{
continue;
}
......@@ -447,8 +464,7 @@ void FaceMesh::Smoothing()
NekDouble d = (di + dj) / 2.0;
NekDouble wij = sqrt((rj[0] - pu0[0]) * (rj[0] - pu0[0]) +
(rj[1] - pu0[1]) * (rj[1] - pu0[1]) +
(rj[2] - pu0[2]) * (rj[2] - pu0[2])) -
d;
(rj[2] - pu0[2]) * (rj[2] - pu0[2])) - d;
NekDouble umag = sqrt((uj[0] - u0[0]) * (uj[0] - u0[0]) +
(uj[1] - u0[1]) * (uj[1] - u0[1]));
......
......@@ -236,9 +236,9 @@ void BLMesh::GrowLayers()
bgi::rtree<boxI, bgi::quadratic<16> > TopTree;
map<int, bgi::rtree<boxI, bgi::quadratic<16> > > SubTrees;
ofstream file;
file.open("pts.3D");
file << "x y z value" << endl;
//ofstream file;
//file.open("pts.3D");
//file << "x y z value" << endl;
for(int l = 1; l < m_layer; l++)
{
......@@ -340,7 +340,7 @@ void BLMesh::GrowLayers()
bit->second->bl = l;
}
}
file.close();
//file.close();
}
inline bool sign(NekDouble a, NekDouble b)
......@@ -775,7 +775,7 @@ void BLMesh::BuildElements()
}
vector<int> tags;
tags.push_back(m_mesh->m_cad->GetNumSurf() < 100 ? 101 : 1001);
tags.push_back(m_id);
ElementSharedPtr E = GetElementFactory().
CreateInstance(LibUtilities::ePrism, pconf, pn, tags);
E->SetId(i);
......@@ -1008,7 +1008,7 @@ void BLMesh::Setup()
for(int i = 0; i < m_mesh->m_element[2].size(); i++)
{
//orientate the triangle
if(m_mesh->m_cad->GetSurf(m_mesh->m_element[2][i]->m_parentCAD->GetId())
if(!m_mesh->m_cad->GetSurf(m_mesh->m_element[2][i]->m_parentCAD->GetId())
->IsReversedNormal())
{
m_mesh->m_element[2][i]->Flip();
......
......@@ -54,12 +54,11 @@ public:
/**
*@brief default constructor
*/
BLMesh(MeshSharedPtr m, std::vector<unsigned int> bls,
NekDouble b,
int l,
NekDouble p) :
m_mesh(m), m_blsurfs(bls), m_bl(b), m_prog(p), m_layer(l)
BLMesh(MeshSharedPtr m, std::vector<unsigned int> bls, NekDouble b, int l,
NekDouble p, int id)
: m_mesh(m), m_blsurfs(bls), m_bl(b), m_prog(p), m_layer(l), m_id(id)
{
};
/**
......@@ -67,8 +66,14 @@ public:
*/
void Mesh();
std::vector<unsigned int> GetSymSurfs(){ return m_symSurfs;}
std::vector<unsigned int> GetBLSurfs(){ return m_blsurfs;}
std::vector<unsigned int> GetSymSurfs()
{
return m_symSurfs;
}
std::vector<unsigned int> GetBLSurfs()
{
return m_blsurfs;
}
std::map<NodeSharedPtr, NodeSharedPtr> GetSymNodes();
......@@ -100,7 +105,6 @@ public:
typedef boost::shared_ptr<blInfo> blInfoSharedPtr;
private:
void Setup();
void GrowLayers();
void Shrink();
......@@ -109,7 +113,8 @@ private:
bool IsPrismValid(ElementSharedPtr el);
NekDouble Proximity(NodeSharedPtr n, ElementSharedPtr el);
NekDouble Visability(std::vector<ElementSharedPtr> tris, Array<OneD, NekDouble> N);
NekDouble Visability(std::vector<ElementSharedPtr> tris,
Array<OneD, NekDouble> N);
Array<OneD, NekDouble> GetNormal(std::vector<ElementSharedPtr> tris);
/// mesh object containing surface mesh
......@@ -120,17 +125,17 @@ private:
NekDouble m_bl;
NekDouble m_prog;
int m_layer;
int m_id;
Array<OneD, NekDouble> m_layerT;
/// list of surfaces to be remeshed due to the boundary layer
std::vector<unsigned int> m_symSurfs;
/// data structure used to store and develop bl information
std::map<NodeSharedPtr, blInfoSharedPtr> m_blData;
std::map<NodeSharedPtr, std::vector<blInfoSharedPtr> > m_nToNInfo; //node to neighbouring information
std::map<ElementSharedPtr,ElementSharedPtr> m_priToTri;
std::map<NodeSharedPtr, std::vector<blInfoSharedPtr> >
m_nToNInfo; // node to neighbouring information
std::map<ElementSharedPtr, ElementSharedPtr> m_priToTri;
std::vector<ElementSharedPtr> m_psuedoSurface;
NekMatrix<NekDouble> m_deriv[3];
};
typedef boost::shared_ptr<BLMesh> BLMeshSharedPtr;
......
......@@ -145,7 +145,7 @@ void TetMesh::Mesh()
n.push_back(IdToNode[m_tetconnect[i][3]]);
ElmtConfig conf(LibUtilities::eTetrahedron, 1, false, false);
vector<int> tags;
tags.push_back(m_mesh->m_cad->GetNumSurf() < 100 ? 100 : 1000);
tags.push_back(m_id);
ElementSharedPtr E = GetElementFactory().CreateInstance(
LibUtilities::eTetrahedron, conf, n, tags);
......
......@@ -61,9 +61,9 @@ public:
/**
* @brief default constructor
*/
TetMesh(MeshSharedPtr m,
TetMesh(MeshSharedPtr m, int id,
std::vector<ElementSharedPtr> e = std::vector<ElementSharedPtr>())
: m_mesh(m), m_surface(e)
: m_mesh(m), m_surface(e), m_id(id)
{
};
......@@ -78,6 +78,7 @@ private:
std::vector<ElementSharedPtr> m_surface;
/// number of tetrahedra
int m_numtet;
int m_id;
/// conncetivity of the tets from the interface
std::vector<Array<OneD, int> > m_tetconnect;
......
......@@ -76,15 +76,13 @@ void VolumeMesh::Process()
if (m_config["blsurfs"].beenSet)
{
makeBL = true;
m_mesh->m_numcomp = 2;
makeBL = true;
ParseUtils::GenerateSeqVector(m_config["blsurfs"].as<string>().c_str(),
blSurfs);
}
else
{
makeBL = false;
m_mesh->m_numcomp = 1;
makeBL = false;
}
NekDouble prefix = 100;
......@@ -98,16 +96,10 @@ void VolumeMesh::Process()
{
BLMeshSharedPtr blmesh = MemoryManager<BLMesh>::AllocateSharedPtr(
m_mesh, blSurfs, m_config["blthick"].as<NekDouble>(),
m_config["bllayers"].as<int>(), m_config["blprog"].as<NekDouble>());
m_config["bllayers"].as<int>(), m_config["blprog"].as<NekDouble>(),
prefix + 1);
blmesh->Mesh();
/*ClearElementLinks();
ProcessVertices();
ProcessEdges();
ProcessFaces();
ProcessElements();
ProcessComposites();
return;*/
// remesh the correct surfaces
vector<unsigned int> symsurfs = blmesh->GetSymSurfs();
......@@ -280,7 +272,8 @@ void VolumeMesh::Process()
}
FaceMeshSharedPtr f = MemoryManager<FaceMesh>::AllocateSharedPtr(
symsurfs[i], m_mesh, cm, prefix + symsurfs[i]);
symsurfs[i], m_mesh, cm, symsurfs[i]);
f->OrientateCurves();
f->Mesh();
}
......@@ -306,11 +299,11 @@ void VolumeMesh::Process()
}
}
tet = MemoryManager<TetMesh>::AllocateSharedPtr(m_mesh, tetsurface);
tet = MemoryManager<TetMesh>::AllocateSharedPtr(m_mesh, prefix, tetsurface);
}
else
{
tet = MemoryManager<TetMesh>::AllocateSharedPtr(m_mesh);
tet = MemoryManager<TetMesh>::AllocateSharedPtr(m_mesh, prefix);
}
tet->Mesh();
......
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