Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
Public Member Functions | Private Member Functions | Private Attributes | Friends | List of all members
Nektar::NekMeshUtils::FaceMesh Class Reference

class for surface meshes on individual surfaces (paramter plane meshes) More...

#include <FaceMesh.h>

Collaboration diagram for Nektar::NekMeshUtils::FaceMesh:
Collaboration graph
[legend]

Public Member Functions

 FaceMesh (const int id, MeshSharedPtr m, const std::map< int, CurveMeshSharedPtr > &cmeshes, const int comp)
 Default constructor. More...
 
void Mesh ()
 mesh exectuation command More...
 
bool ValidateCurves ()
 validate the curve meshes More...
 

Private Member Functions

void OrientateCurves ()
 Get the boundries of the surface and extracts the nodes from the curve meshes in the correct order. More...
 
void Stretching ()
 Calculate the paramter plane streching factor. More...
 
void Smoothing ()
 performs node smoothing on face More...
 
void DiagonalSwap ()
 performs diagonal swapping of edges More...
 
void BuildLocalMesh ()
 build a local version of mesh elements More...
 
void OptimiseLocalMesh ()
 function which calls the optimisation routines More...
 
bool Validate ()
 Validate the surface mesh base on the octree and real dimensions of the edges. More...
 
void AddNewPoint (Array< OneD, NekDouble > uv)
 adds a new stiener point to the triangulation for meshing More...
 
void MakeBL ()
 adds a quad layer around any interior loops More...
 

Private Attributes

MeshSharedPtr m_mesh
 mesh pointer More...
 
CADSurfSharedPtr m_cadsurf
 CAD surface. More...
 
std::map< int, CurveMeshSharedPtrm_curvemeshes
 Map of the curve meshes which bound the surfaces. More...
 
std::vector
< CADSystem::EdgeLoopSharedPtr
m_edgeloops
 data structure containing the edges, their order and oreientation for the surface More...
 
int m_id
 id of the surface mesh More...
 
std::vector< std::vector
< NodeSharedPtr > > 
orderedLoops
 list of boundary nodes in their order loops More...
 
std::vector< NodeSharedPtrm_stienerpoints
 list of stiener points in the triangulation More...
 
NekDouble m_str
 pplane stretching More...
 
std::vector< std::vector
< NodeSharedPtr > > 
m_connec
 triangle connectiviities More...
 
NodeSet m_localNodes
 local set of nodes More...
 
EdgeSet m_localEdges
 local set of edges More...
 
std::vector< ElementSharedPtrm_localElements
 local list of elements More...
 
NodeSet m_inBoundary
 set of nodes which are in the boundary (easier to identify conflicts with) More...
 
int m_compId
 identity to put into element tags More...
 

Friends

class MemoryManager< FaceMesh >
 

Detailed Description

class for surface meshes on individual surfaces (paramter plane meshes)

Definition at line 52 of file FaceMesh.h.

Constructor & Destructor Documentation

Nektar::NekMeshUtils::FaceMesh::FaceMesh ( const int  id,
MeshSharedPtr  m,
const std::map< int, CurveMeshSharedPtr > &  cmeshes,
const int  comp 
)
inline

Default constructor.

Definition at line 60 of file FaceMesh.h.

References m_cadsurf, m_edgeloops, m_id, and m_mesh.

64  : m_mesh(m), m_curvemeshes(cmeshes), m_id(id), m_compId(comp)
65 
66  {
67  m_cadsurf = m_mesh->m_cad->GetSurf(m_id);
68  m_edgeloops = m_cadsurf->GetEdges();
69  };
std::map< int, CurveMeshSharedPtr > m_curvemeshes
Map of the curve meshes which bound the surfaces.
Definition: FaceMesh.h:135
int m_compId
identity to put into element tags
Definition: FaceMesh.h:158
std::vector< CADSystem::EdgeLoopSharedPtr > m_edgeloops
data structure containing the edges, their order and oreientation for the surface ...
Definition: FaceMesh.h:138
int m_id
id of the surface mesh
Definition: FaceMesh.h:140
CADSurfSharedPtr m_cadsurf
CAD surface.
Definition: FaceMesh.h:133
MeshSharedPtr m_mesh
mesh pointer
Definition: FaceMesh.h:131

Member Function Documentation

void Nektar::NekMeshUtils::FaceMesh::AddNewPoint ( Array< OneD, NekDouble uv)
private

adds a new stiener point to the triangulation for meshing

Definition at line 1109 of file FaceMesh.cpp.

References class_topology::Node.

1110 {
1111  // adds a new point but checks that there are no other points nearby first
1112  Array<OneD, NekDouble> np = m_cadsurf->P(uv);
1113  NekDouble npDelta = m_mesh->m_octree->Query(np);
1114 
1115  NodeSharedPtr n = boost::shared_ptr<Node>(
1116  new Node(m_mesh->m_numNodes++, np[0], np[1], np[2]));
1117 
1118  bool add = true;
1119 
1120  for (int i = 0; i < orderedLoops.size(); i++)
1121  {
1122  for (int j = 0; j < orderedLoops[i].size(); j++)
1123  {
1124  NekDouble r = orderedLoops[i][j]->Distance(n);
1125 
1126  if (r < npDelta / 2.0)
1127  {
1128  add = false;
1129  break;
1130  }
1131  }
1132  }
1133 
1134  if (add)
1135  {
1136  for (int i = 0; i < m_stienerpoints.size(); i++)
1137  {
1138  NekDouble r = m_stienerpoints[i]->Distance(n);
1139 
1140  if (r < npDelta / 2.0)
1141  {
1142  add = false;
1143  break;
1144  }
1145  }
1146  }
1147 
1148  if (add)
1149  {
1150  n->SetCADSurf(m_id, m_cadsurf, uv);
1151  m_stienerpoints.push_back(n);
1152  }
1153 }
std::vector< std::vector< NodeSharedPtr > > orderedLoops
list of boundary nodes in their order loops
Definition: FaceMesh.h:142
int m_id
id of the surface mesh
Definition: FaceMesh.h:140
boost::shared_ptr< Node > NodeSharedPtr
Definition: Node.h:50
double NekDouble
std::vector< NodeSharedPtr > m_stienerpoints
list of stiener points in the triangulation
Definition: FaceMesh.h:144
CADSurfSharedPtr m_cadsurf
CAD surface.
Definition: FaceMesh.h:133
MeshSharedPtr m_mesh
mesh pointer
Definition: FaceMesh.h:131
void Nektar::NekMeshUtils::FaceMesh::BuildLocalMesh ( )
private

build a local version of mesh elements

Definition at line 883 of file FaceMesh.cpp.

References Nektar::LibUtilities::NekFactory< tKey, tBase, >::CreateInstance(), Nektar::LibUtilities::eTriangle, Nektar::NekMeshUtils::GetElementFactory(), and Nektar::iterator.

884 {
885  /*************************
886  // build a local set of nodes edges and elemenets for optimstaion prior to
887  putting them into m_mesh
888  */
889 
890  for (int i = 0; i < m_connec.size(); i++)
891  {
892  ElmtConfig conf(LibUtilities::eTriangle, 1, false, false);
893 
894  vector<int> tags;
895  tags.push_back(m_compId);
897  LibUtilities::eTriangle, conf, m_connec[i], tags);
898  E->m_parentCAD = m_cadsurf;
899 
900  vector<NodeSharedPtr> nods = E->GetVertexList();
901  for (int j = 0; j < nods.size(); j++)
902  {
903  // nodes are already unique some will insert some wont
904  m_localNodes.insert(nods[j]);
905  }
906 
907  E->SetId(m_localElements.size());
908  m_localElements.push_back(E);
909  }
910 
911  for (int i = 0; i < m_localElements.size(); ++i)
912  {
913  for (int j = 0; j < m_localElements[i]->GetEdgeCount(); ++j)
914  {
915  pair<EdgeSet::iterator, bool> testIns;
916  EdgeSharedPtr ed = m_localElements[i]->GetEdge(j);
917  // look for edge in m_mesh edgeset from curves
918  EdgeSet::iterator s = m_mesh->m_edgeSet.find(ed);
919  if (!(s == m_mesh->m_edgeSet.end()))
920  {
921  ed = *s;
922  m_localElements[i]->SetEdge(j, *s);
923  }
924 
925  testIns = m_localEdges.insert(ed);
926 
927  if (testIns.second)
928  {
929  EdgeSharedPtr ed2 = *testIns.first;
930  ed2->m_elLink.push_back(
931  pair<ElementSharedPtr, int>(m_localElements[i], j));
932  }
933  else
934  {
935  EdgeSharedPtr e2 = *(testIns.first);
936  m_localElements[i]->SetEdge(j, e2);
937  e2->m_elLink.push_back(
938  pair<ElementSharedPtr, int>(m_localElements[i], j));
939  }
940  }
941  }
942 }
tBaseSharedPtr CreateInstance(tKey idKey BOOST_PP_COMMA_IF(MAX_PARAM) BOOST_PP_ENUM_BINARY_PARAMS(MAX_PARAM, tParam, x))
Create an instance of the class referred to by idKey.
Definition: NekFactory.hpp:162
std::vector< ElementSharedPtr > m_localElements
local list of elements
Definition: FaceMesh.h:154
ElementFactory & GetElementFactory()
Definition: Element.cpp:47
int m_compId
identity to put into element tags
Definition: FaceMesh.h:158
std::vector< std::vector< NodeSharedPtr > > m_connec
triangle connectiviities
Definition: FaceMesh.h:148
EdgeSet m_localEdges
local set of edges
Definition: FaceMesh.h:152
boost::shared_ptr< Edge > EdgeSharedPtr
Shared pointer to an edge.
Definition: Edge.h:135
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
boost::shared_ptr< Element > ElementSharedPtr
Definition: Edge.h:49
CADSurfSharedPtr m_cadsurf
CAD surface.
Definition: FaceMesh.h:133
MeshSharedPtr m_mesh
mesh pointer
Definition: FaceMesh.h:131
NodeSet m_localNodes
local set of nodes
Definition: FaceMesh.h:150
void Nektar::NekMeshUtils::FaceMesh::DiagonalSwap ( )
private

performs diagonal swapping of edges

Definition at line 488 of file FaceMesh.cpp.

References ASSERTL0, Nektar::LibUtilities::NekFactory< tKey, tBase, >::CreateInstance(), Nektar::LibUtilities::eTriangle, Nektar::NekMeshUtils::GetElementFactory(), and Nektar::iterator.

489 {
490  map<int, int> idealConnec;
491  map<int, int> actualConnec;
492  map<int, vector<EdgeSharedPtr> > nodetoedge;
493  // figure out ideal node count and actual node count
494  EdgeSet::iterator eit;
495  for (eit = m_localEdges.begin(); eit != m_localEdges.end(); eit++)
496  {
497  nodetoedge[(*eit)->m_n1->m_id].push_back(*eit);
498  nodetoedge[(*eit)->m_n2->m_id].push_back(*eit);
499  }
500  NodeSet::iterator nit;
501  for (nit = m_localNodes.begin(); nit != m_localNodes.end(); nit++)
502  {
503  // this routine is broken and needs looking at
504  NodeSet::iterator f = m_inBoundary.find((*nit));
505  if (f == m_inBoundary.end()) // node isnt on curve so skip
506  {
507  // node is interior
508  idealConnec[(*nit)->m_id] = 6;
509  }
510  else
511  {
512  // need to identify the two other nodes on the boundary to find
513  // interior angle
514  /*vector<NodeSharedPtr> ns;
515  vector<EdgeSharedPtr> e = nodetoedge[(*nit)->m_id];
516  for (int i = 0; i < e.size(); i++)
517  {
518  if (!e[i]->onCurve)
519  continue; // the linking nodes are not going to exist on
520  // interior edges
521 
522  if (e[i]->m_n1 == (*nit))
523  ns.push_back(e[i]->m_n2);
524  else
525  ns.push_back(e[i]->m_n1);
526  }
527  ASSERTL0(ns.size() == 2,
528  "failed to find 2 nodes in the angle system");
529 
530  idealConnec[(*nit)->m_id] =
531  ceil((*nit)->Angle(ns[0], ns[1]) / 3.142 * 3) + 1;*/
532  idealConnec[(*nit)->m_id] = 4;
533  }
534  }
535  for (nit = m_localNodes.begin(); nit != m_localNodes.end(); nit++)
536  {
537  actualConnec[(*nit)->m_id] = nodetoedge[(*nit)->m_id].size();
538  }
539 
540  // edgeswapping fun times
541  // perfrom edge swap based on node defect and then angle
542  for (int q = 0; q < 4; q++)
543  {
544  int edgesStart = m_localEdges.size();
545  EdgeSet edges = m_localEdges;
546  m_localEdges.clear();
547 
548  int swappedEdges = 0;
549 
551 
552  for (it = edges.begin(); it != edges.end(); it++)
553  {
554  EdgeSharedPtr e = *it;
555 
556  NodeSet::iterator f1 = m_inBoundary.find((*it)->m_n1);
557  NodeSet::iterator f2 = m_inBoundary.find((*it)->m_n2);
558  if (f1 != m_inBoundary.end() && f2 != m_inBoundary.end())
559  {
560  m_localEdges.insert(e);
561  continue;
562  }
563 
564  ElementSharedPtr tri1 = e->m_elLink[0].first;
565  ElementSharedPtr tri2 = e->m_elLink[1].first;
566 
567  NodeSharedPtr n1 = e->m_n1;
568  NodeSharedPtr n2 = e->m_n2;
569 
570  vector<NodeSharedPtr> nt = tri1->GetVertexList();
571 
572  // identify node a,b,c,d of the swapping
573  NodeSharedPtr A, B, C, D;
574  if (nt[0] != n1 && nt[0] != n2)
575  {
576  C = nt[0];
577  B = nt[1];
578  A = nt[2];
579  }
580  else if (nt[1] != n1 && nt[1] != n2)
581  {
582  C = nt[1];
583  B = nt[2];
584  A = nt[0];
585  }
586  else if (nt[2] != n1 && nt[2] != n2)
587  {
588  C = nt[2];
589  B = nt[0];
590  A = nt[1];
591  }
592  else
593  {
594  ASSERTL0(false, "failed to identify verticies in tri1");
595  }
596 
597  nt = tri2->GetVertexList();
598 
599  if (nt[0] != n1 && nt[0] != n2)
600  {
601  D = nt[0];
602  }
603  else if (nt[1] != n1 && nt[1] != n2)
604  {
605  D = nt[1];
606  }
607  else if (nt[2] != n1 && nt[2] != n2)
608  {
609  D = nt[2];
610  }
611  else
612  {
613  ASSERTL0(false, "failed to identify verticies in tri2");
614  }
615 
616  // determine signed area of alternate config
617  // //cout<< A->GetNumCADSurf() << " " << B->GetNumCADSurf() << " "
618  // // << C->GetNumCADSurf() << " " << D->GetNumCADSurf() << endl;
619  // ofstream file;
620  // file.open("pts.3D");
621  // file << "x y z value" << endl;
622  // file << A->m_x << " " << A->m_y << " " << A->m_z << endl;
623  // file << B->m_x << " " << B->m_y << " " << B->m_z << endl;
624  // file << C->m_x << " " << C->m_y << " " << C->m_z << endl;
625  // file << D->m_x << " " << D->m_y << " " << D->m_z << endl;
626  // file.close();
627  Array<OneD, NekDouble> ai, bi, ci, di;
628  ai = A->GetCADSurfInfo(m_id);
629  bi = B->GetCADSurfInfo(m_id);
630  ci = C->GetCADSurfInfo(m_id);
631  di = D->GetCADSurfInfo(m_id);
632 
633  NekDouble CDA, CBD;
634 
635  CDA = 0.5 * (-di[0] * ci[1] + ai[0] * ci[1] + ci[0] * di[1] -
636  ai[0] * di[1] - ci[0] * ai[1] + di[0] * ai[1]);
637 
638  CBD = 0.5 * (-bi[0] * ci[1] + di[0] * ci[1] + ci[0] * bi[1] -
639  di[0] * bi[1] - ci[0] * di[1] + bi[0] * di[1]);
640 
641  // if signed area of the swapping triangles is less than zero
642  // that configuration is invalid and swap cannot be performed
643  if (!(CDA > 0.001 && CBD > 0.001))
644  {
645  m_localEdges.insert(e);
646  continue;
647  }
648 
649  bool swap = false; // assume do not swap
650 
651  if (q < 2)
652  {
653  int nodedefectbefore = 0;
654  nodedefectbefore +=
655  abs(actualConnec[A->m_id] - idealConnec[A->m_id]);
656  nodedefectbefore +=
657  abs(actualConnec[B->m_id] - idealConnec[B->m_id]);
658  nodedefectbefore +=
659  abs(actualConnec[C->m_id] - idealConnec[C->m_id]);
660  nodedefectbefore +=
661  abs(actualConnec[D->m_id] - idealConnec[D->m_id]);
662 
663  int nodedefectafter = 0;
664  nodedefectafter +=
665  abs(actualConnec[A->m_id] - 1 - idealConnec[A->m_id]);
666  nodedefectafter +=
667  abs(actualConnec[B->m_id] - 1 - idealConnec[B->m_id]);
668  nodedefectafter +=
669  abs(actualConnec[C->m_id] + 1 - idealConnec[C->m_id]);
670  nodedefectafter +=
671  abs(actualConnec[D->m_id] + 1 - idealConnec[D->m_id]);
672 
673  if (nodedefectafter < nodedefectbefore)
674  {
675  swap = true;
676  }
677  }
678  else
679  {
680  NekDouble minanglebefore = C->Angle(A, B);
681  minanglebefore = min(minanglebefore, A->Angle(B, C));
682  minanglebefore = min(minanglebefore, B->Angle(A, C));
683  minanglebefore = min(minanglebefore, B->Angle(A, D));
684  minanglebefore = min(minanglebefore, A->Angle(B, D));
685  minanglebefore = min(minanglebefore, D->Angle(A, B));
686 
687  NekDouble minangleafter = C->Angle(B, D);
688  minangleafter = min(minangleafter, D->Angle(B, C));
689  minangleafter = min(minangleafter, B->Angle(C, D));
690  minangleafter = min(minangleafter, C->Angle(A, D));
691  minangleafter = min(minangleafter, A->Angle(C, D));
692  minangleafter = min(minangleafter, D->Angle(A, C));
693 
694  if (minangleafter > minanglebefore)
695  {
696  swap = true;
697  }
698  }
699 
700  if (swap)
701  {
702  actualConnec[A->m_id]--;
703  actualConnec[B->m_id]--;
704  actualConnec[C->m_id]++;
705  actualConnec[D->m_id]++;
706 
707  // make the 4 other edges
708  EdgeSharedPtr CA, AD, DB, BC, CAt, ADt, DBt, BCt;
709  CAt = boost::shared_ptr<Edge>(new Edge(C, A));
710  ADt = boost::shared_ptr<Edge>(new Edge(A, D));
711  DBt = boost::shared_ptr<Edge>(new Edge(D, B));
712  BCt = boost::shared_ptr<Edge>(new Edge(B, C));
713 
714  vector<EdgeSharedPtr> es = tri1->GetEdgeList();
715  for (int i = 0; i < 3; i++)
716  {
717  if (es[i] == CAt)
718  {
719  CA = es[i];
720  }
721  if (es[i] == BCt)
722  {
723  BC = es[i];
724  }
725  }
726  es = tri2->GetEdgeList();
727  for (int i = 0; i < 3; i++)
728  {
729  if (es[i] == DBt)
730  {
731  DB = es[i];
732  }
733  if (es[i] == ADt)
734  {
735  AD = es[i];
736  }
737  }
738 
739  // now sort out links for the 4 edges surrounding the patch
740  vector<pair<ElementSharedPtr, int> > links;
741 
742  links = CA->m_elLink;
743  CA->m_elLink.clear();
744  for (int i = 0; i < links.size(); i++)
745  {
746  if (links[i].first->GetId() == tri1->GetId())
747  {
748  continue;
749  }
750  CA->m_elLink.push_back(links[i]);
751  }
752 
753  links = BC->m_elLink;
754  BC->m_elLink.clear();
755  for (int i = 0; i < links.size(); i++)
756  {
757  if (links[i].first->GetId() == tri1->GetId())
758  {
759  continue;
760  }
761  BC->m_elLink.push_back(links[i]);
762  }
763 
764  links = AD->m_elLink;
765  AD->m_elLink.clear();
766  for (int i = 0; i < links.size(); i++)
767  {
768  if (links[i].first->GetId() == tri2->GetId())
769  {
770  continue;
771  }
772  AD->m_elLink.push_back(links[i]);
773  }
774 
775  links = DB->m_elLink;
776  DB->m_elLink.clear();
777  for (int i = 0; i < links.size(); i++)
778  {
779  if (links[i].first->GetId() == tri2->GetId())
780  {
781  continue;
782  }
783  DB->m_elLink.push_back(links[i]);
784  }
785 
786  EdgeSharedPtr newe = boost::shared_ptr<Edge>(new Edge(C, D));
787 
788  vector<NodeSharedPtr> t1, t2;
789  t1.push_back(B);
790  t1.push_back(D);
791  t1.push_back(C);
792  t2.push_back(A);
793  t2.push_back(C);
794  t2.push_back(D);
795 
796  ElmtConfig conf(LibUtilities::eTriangle, 1, false, false);
797  vector<int> tags = tri1->GetTagList();
798 
799  int id1 = tri1->GetId();
800  int id2 = tri2->GetId();
801 
803  LibUtilities::eTriangle, conf, t1, tags);
804  tags = tri2->GetTagList();
806  LibUtilities::eTriangle, conf, t2, tags);
807 
808  ntri1->SetId(id1);
809  ntri2->SetId(id2);
810  ntri1->m_parentCAD = m_cadsurf;
811  ntri2->m_parentCAD = m_cadsurf;
812 
813  vector<EdgeSharedPtr> t1es = ntri1->GetEdgeList();
814  for (int i = 0; i < 3; i++)
815  {
816  if (t1es[i] == DB)
817  {
818  ntri1->SetEdge(i, DB);
819  DB->m_elLink.push_back(
820  pair<ElementSharedPtr, int>(ntri1, i));
821  }
822  else if (t1es[i] == BC)
823  {
824  ntri1->SetEdge(i, BC);
825  BC->m_elLink.push_back(
826  pair<ElementSharedPtr, int>(ntri1, i));
827  }
828  else if (t1es[i] == newe)
829  {
830  ntri1->SetEdge(i, newe);
831  newe->m_elLink.push_back(
832  pair<ElementSharedPtr, int>(ntri1, i));
833  }
834  else
835  {
836  ASSERTL0(false, "weird edge in new tri 1");
837  }
838  }
839  vector<EdgeSharedPtr> t2es = ntri2->GetEdgeList();
840  for (int i = 0; i < 3; i++)
841  {
842  if (t2es[i] == CA)
843  {
844  ntri2->SetEdge(i, CA);
845  CA->m_elLink.push_back(
846  pair<ElementSharedPtr, int>(ntri2, i));
847  }
848  else if (t2es[i] == AD)
849  {
850  ntri2->SetEdge(i, AD);
851  AD->m_elLink.push_back(
852  pair<ElementSharedPtr, int>(ntri2, i));
853  }
854  else if (t2es[i] == newe)
855  {
856  ntri2->SetEdge(i, newe);
857  newe->m_elLink.push_back(
858  pair<ElementSharedPtr, int>(ntri2, i));
859  }
860  else
861  {
862  ASSERTL0(false, "weird edge in new tri 2");
863  }
864  }
865 
866  m_localEdges.insert(newe);
867 
868  m_localElements[id1] = ntri1;
869  m_localElements[id2] = ntri2;
870 
871  swappedEdges++;
872  }
873  else
874  {
875  m_localEdges.insert(e);
876  }
877  }
878 
879  ASSERTL0(m_localEdges.size() == edgesStart, "mismatch edge count");
880  }
881 }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
tBaseSharedPtr CreateInstance(tKey idKey BOOST_PP_COMMA_IF(MAX_PARAM) BOOST_PP_ENUM_BINARY_PARAMS(MAX_PARAM, tParam, x))
Create an instance of the class referred to by idKey.
Definition: NekFactory.hpp:162
std::vector< ElementSharedPtr > m_localElements
local list of elements
Definition: FaceMesh.h:154
ElementFactory & GetElementFactory()
Definition: Element.cpp:47
EdgeSet m_localEdges
local set of edges
Definition: FaceMesh.h:152
int m_id
id of the surface mesh
Definition: FaceMesh.h:140
boost::shared_ptr< Node > NodeSharedPtr
Definition: Node.h:50
double NekDouble
boost::shared_ptr< Edge > EdgeSharedPtr
Shared pointer to an edge.
Definition: Edge.h:135
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
boost::shared_ptr< Element > ElementSharedPtr
Definition: Edge.h:49
boost::unordered_set< EdgeSharedPtr, EdgeHash > EdgeSet
Definition: Edge.h:161
CADSurfSharedPtr m_cadsurf
CAD surface.
Definition: FaceMesh.h:133
NodeSet m_localNodes
local set of nodes
Definition: FaceMesh.h:150
NodeSet m_inBoundary
set of nodes which are in the boundary (easier to identify conflicts with)
Definition: FaceMesh.h:156
void Nektar::NekMeshUtils::FaceMesh::MakeBL ( )
private

adds a quad layer around any interior loops

void Nektar::NekMeshUtils::FaceMesh::Mesh ( )

mesh exectuation command

Definition at line 127 of file FaceMesh.cpp.

References ASSERTL0.

128 {
129  Stretching();
130  OrientateCurves();
131 
132  int numPoints = 0;
133  for (int i = 0; i < orderedLoops.size(); i++)
134  {
135  numPoints += orderedLoops[i].size();
136  for (int j = 0; j < orderedLoops[i].size(); j++)
137  {
138  m_inBoundary.insert(orderedLoops[i][j]);
139  }
140  }
141 
142  stringstream ss;
143  ss << "3 points required for triangulation, " << numPoints << " in loop"
144  << endl;
145  ss << "curves: ";
146  for (int i = 0; i < m_edgeloops.size(); i++)
147  {
148  for (int j = 0; j < m_edgeloops[i]->edges.size(); j++)
149  {
150  ss << m_edgeloops[i]->edges[j]->GetId() << " ";
151  }
152  }
153 
154  ASSERTL0(numPoints > 2, ss.str());
155 
156  // create interface to triangle thirdparty library
157  TriangleInterfaceSharedPtr pplanemesh =
159 
160  vector<Array<OneD, NekDouble> > centers;
161  for (int i = 0; i < m_edgeloops.size(); i++)
162  {
163  centers.push_back(m_edgeloops[i]->center);
164  }
165 
166  pplanemesh->Assign(orderedLoops, centers, m_id, m_str);
167 
168  pplanemesh->Mesh();
169 
170  pplanemesh->Extract(m_connec);
171 
172  bool repeat = true;
173  int meshcounter = 1;
174 
175  // continuously remesh until all triangles conform to the spacing in the
176  // octree
177  while (repeat)
178  {
179  repeat = Validate();
180  if (!repeat)
181  {
182  break;
183  }
184  m_connec.clear();
185  pplanemesh->AssignStiener(m_stienerpoints);
186  pplanemesh->Mesh();
187  pplanemesh->Extract(m_connec);
188  meshcounter++;
189  // break;
190  }
191 
192  // build a local version of the mesh (one set of triangles). this is done
193  // so edge connectivity infomration can be used for optimisation
194  BuildLocalMesh();
195 
197 
198  // make new elements and add to list from list of nodes and connectivity
199  // from triangle removing unnesercary infomration from the elements
200  for (int i = 0; i < m_localElements.size(); i++)
201  {
202  vector<EdgeSharedPtr> e = m_localElements[i]->GetEdgeList();
203  for (int j = 0; j < e.size(); j++)
204  {
205  e[j]->m_elLink.clear();
206  }
207  m_mesh->m_element[2].push_back(m_localElements[i]);
208  }
209 
210  if (m_mesh->m_verbose)
211  {
212  cout << "\r "
213  " "
214  " ";
215  cout << scientific << "\r\t\tFace " << m_id << endl
216  << "\t\t\tNodes: " << m_localNodes.size() << endl
217  << "\t\t\tEdges: " << m_localEdges.size() << endl
218  << "\t\t\tTriangles: " << m_localElements.size() << endl
219  << "\t\t\tLoops: " << m_edgeloops.size() << endl
220  << "\t\t\tSTR: " << m_str << endl
221  << endl;
222  }
223 }
std::vector< std::vector< NodeSharedPtr > > orderedLoops
list of boundary nodes in their order loops
Definition: FaceMesh.h:142
bool Validate()
Validate the surface mesh base on the octree and real dimensions of the edges.
Definition: FaceMesh.cpp:998
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
boost::shared_ptr< TriangleInterface > TriangleInterfaceSharedPtr
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
void Stretching()
Calculate the paramter plane streching factor.
Definition: FaceMesh.cpp:944
std::vector< ElementSharedPtr > m_localElements
local list of elements
Definition: FaceMesh.h:154
std::vector< std::vector< NodeSharedPtr > > m_connec
triangle connectiviities
Definition: FaceMesh.h:148
std::vector< CADSystem::EdgeLoopSharedPtr > m_edgeloops
data structure containing the edges, their order and oreientation for the surface ...
Definition: FaceMesh.h:138
EdgeSet m_localEdges
local set of edges
Definition: FaceMesh.h:152
int m_id
id of the surface mesh
Definition: FaceMesh.h:140
void OrientateCurves()
Get the boundries of the surface and extracts the nodes from the curve meshes in the correct order...
Definition: FaceMesh.cpp:1155
std::vector< NodeSharedPtr > m_stienerpoints
list of stiener points in the triangulation
Definition: FaceMesh.h:144
NekDouble m_str
pplane stretching
Definition: FaceMesh.h:146
void BuildLocalMesh()
build a local version of mesh elements
Definition: FaceMesh.cpp:883
void OptimiseLocalMesh()
function which calls the optimisation routines
Definition: FaceMesh.cpp:225
MeshSharedPtr m_mesh
mesh pointer
Definition: FaceMesh.h:131
NodeSet m_localNodes
local set of nodes
Definition: FaceMesh.h:150
NodeSet m_inBoundary
set of nodes which are in the boundary (easier to identify conflicts with)
Definition: FaceMesh.h:156
void Nektar::NekMeshUtils::FaceMesh::OptimiseLocalMesh ( )
private

function which calls the optimisation routines

Definition at line 225 of file FaceMesh.cpp.

226 {
227  // each optimisation algorithm is based on the work in chapter 19
228  DiagonalSwap();
229 
230  Smoothing();
231 
232  DiagonalSwap();
233 
234  Smoothing();
235 }
void DiagonalSwap()
performs diagonal swapping of edges
Definition: FaceMesh.cpp:488
void Smoothing()
performs node smoothing on face
Definition: FaceMesh.cpp:237
void Nektar::NekMeshUtils::FaceMesh::OrientateCurves ( )
private

Get the boundries of the surface and extracts the nodes from the curve meshes in the correct order.

Definition at line 1155 of file FaceMesh.cpp.

References Nektar::StdRegions::eForwards.

1156 {
1157  // this could be a second run on orentate so clear some info
1158  orderedLoops.clear();
1159 
1160  // create list of bounding loop nodes
1161  for (int i = 0; i < m_edgeloops.size(); i++)
1162  {
1163  vector<NodeSharedPtr> cE;
1164  for (int j = 0; j < m_edgeloops[i]->edges.size(); j++)
1165  {
1166  int cid = m_edgeloops[i]->edges[j]->GetId();
1167  vector<NodeSharedPtr> edgePoints =
1168  m_curvemeshes[cid]->GetMeshPoints();
1169 
1170  int numPoints = m_curvemeshes[cid]->GetNumPoints();
1171 
1172  if (m_edgeloops[i]->edgeo[j] == CADOrientation::eForwards)
1173  {
1174  for (int k = 0; k < numPoints - 1; k++)
1175  {
1176  cE.push_back(edgePoints[k]);
1177  }
1178  }
1179  else
1180  {
1181  for (int k = numPoints - 1; k > 0; k--)
1182  {
1183  cE.push_back(edgePoints[k]);
1184  }
1185  }
1186  }
1187  orderedLoops.push_back(cE);
1188  }
1189 }
std::vector< std::vector< NodeSharedPtr > > orderedLoops
list of boundary nodes in their order loops
Definition: FaceMesh.h:142
std::map< int, CurveMeshSharedPtr > m_curvemeshes
Map of the curve meshes which bound the surfaces.
Definition: FaceMesh.h:135
std::vector< CADSystem::EdgeLoopSharedPtr > m_edgeloops
data structure containing the edges, their order and oreientation for the surface ...
Definition: FaceMesh.h:138
void Nektar::NekMeshUtils::FaceMesh::Smoothing ( )
private

performs node smoothing on face

Definition at line 237 of file FaceMesh.cpp.

References ASSERTL0, Nektar::iterator, and class_topology::Node.

238 {
239  EdgeSet::iterator eit;
240  NodeSet::iterator nit;
241 
242  Array<OneD, NekDouble> bounds = m_cadsurf->GetBounds();
243 
244  map<int, vector<EdgeSharedPtr> > connectingedges;
245 
246  map<int, vector<ElementSharedPtr> > connectingelements;
247 
248  for (eit = m_localEdges.begin(); eit != m_localEdges.end(); eit++)
249  {
250  connectingedges[(*eit)->m_n1->m_id].push_back(*eit);
251  connectingedges[(*eit)->m_n2->m_id].push_back(*eit);
252  }
253 
254  for (int i = 0; i < m_localElements.size(); i++)
255  {
256  vector<NodeSharedPtr> v = m_localElements[i]->GetVertexList();
257  for (int j = 0; j < 3; j++)
258  {
259  connectingelements[v[j]->m_id].push_back(m_localElements[i]);
260  }
261  }
262 
263  // perform 4 runs of elastic relaxation based on the octree
264  for (int q = 0; q < 4; q++)
265  {
266  for (nit = m_localNodes.begin(); nit != m_localNodes.end(); nit++)
267  {
268  NodeSet::iterator f = m_inBoundary.find((*nit));
269 
270  // node is on curve so skip
271  if (f != m_inBoundary.end())
272  {
273  continue;
274  }
275 
276  // this can be real nodes or dummy nodes depending on the system
277  vector<NodeSharedPtr> connodes;
278 
279  vector<EdgeSharedPtr> edges = connectingedges[(*nit)->m_id];
280  vector<ElementSharedPtr> els = connectingelements[(*nit)->m_id];
281 
282  vector<NodeSharedPtr> nodesystem;
283  vector<NekDouble> lamp;
284 
285  for (int i = 0; i < edges.size(); i++)
286  {
287  vector<NekDouble> lambda;
288 
289  NodeSharedPtr J;
290  if (*nit == edges[i]->m_n1)
291  {
292  J = edges[i]->m_n2;
293  }
294  else if (*nit == edges[i]->m_n2)
295  {
296  J = edges[i]->m_n1;
297  }
298  else
299  {
300  ASSERTL0(false, "could not find node");
301  }
302 
303  Array<OneD, NekDouble> ui = (*nit)->GetCADSurfInfo(m_id);
304  Array<OneD, NekDouble> uj = J->GetCADSurfInfo(m_id);
305 
306  for (int j = 0; j < els.size(); j++)
307  {
308  vector<NodeSharedPtr> v = els[j]->GetVertexList();
309 
310  // elememt is adjacent to J therefore no intersection on IJ
311  if (v[0] == J || v[1] == J || v[2] == J)
312  {
313  continue;
314  }
315 
316  // need to find other edge
317  EdgeSharedPtr AtoB;
318  bool found = false;
319  vector<EdgeSharedPtr> es = els[j]->GetEdgeList();
320  for (int k = 0; k < es.size(); k++)
321  {
322  if (!(es[k]->m_n1 == *nit || es[k]->m_n2 == *nit))
323  {
324  found = true;
325  AtoB = es[k];
326  break;
327  }
328  }
329  ASSERTL0(found, "failed to find edge to test");
330 
331  Array<OneD, NekDouble> A = AtoB->m_n1->GetCADSurfInfo(m_id);
332  Array<OneD, NekDouble> B = AtoB->m_n2->GetCADSurfInfo(m_id);
333 
334  NekDouble lam = ((A[0] - uj[0]) * (B[1] - A[1]) -
335  (A[1] - uj[1]) * (B[0] - A[0])) /
336  ((ui[0] - uj[0]) * (B[1] - A[1]) -
337  (ui[1] - uj[1]) * (B[0] - A[0]));
338 
339  if (!(lam < 0) && !(lam > 1))
340  lambda.push_back(lam);
341  }
342 
343  if (lambda.size() > 0)
344  {
345  sort(lambda.begin(), lambda.end());
346  // make a new dummy node based on the system
347  Array<OneD, NekDouble> ud(2);
348  ud[0] = uj[0] + lambda[0] * (ui[0] - uj[0]);
349  ud[1] = uj[1] + lambda[0] * (ui[1] - uj[1]);
350  Array<OneD, NekDouble> locd = m_cadsurf->P(ud);
351  NodeSharedPtr dn = boost::shared_ptr<Node>(
352  new Node(0, locd[0], locd[1], locd[2]));
353  dn->SetCADSurf(m_id, m_cadsurf, ud);
354 
355  nodesystem.push_back(dn);
356  lamp.push_back(lambda[0]);
357  }
358  else
359  {
360  nodesystem.push_back(J);
361  lamp.push_back(1.0);
362  }
363  }
364 
365  Array<OneD, NekDouble> u0(2);
366  u0[0] = 0.0;
367  u0[1] = 0.0;
368 
369  for (int i = 0; i < nodesystem.size(); i++)
370  {
371  Array<OneD, NekDouble> uj = nodesystem[i]->GetCADSurfInfo(m_id);
372  u0[0] += uj[0] / nodesystem.size();
373  u0[1] += uj[1] / nodesystem.size();
374  }
375 
376  /*Array<OneD, NekDouble> pu0 = m_cadsurf->P(u0);
377  NekDouble di = m_mesh->m_octree->Query(pu0);
378  Array<OneD, NekDouble> F(2, 0.0), dF(4, 0.0);
379  for (int i = 0; i < nodesystem.size(); i++)
380  {
381  Array<OneD, NekDouble> rj = nodesystem[i]->GetLoc();
382  Array<OneD, NekDouble> uj = nodesystem[i]->GetCADSurfInfo(m_id);
383  NekDouble dj = m_mesh->m_octree->Query(rj);
384  NekDouble d = (di + dj) / 2.0;
385  NekDouble wij = sqrt((rj[0] - pu0[0]) * (rj[0] - pu0[0]) +
386  (rj[1] - pu0[1]) * (rj[1] - pu0[1]) +
387  (rj[2] - pu0[2]) * (rj[2] - pu0[2])) -
388  d;
389 
390  NekDouble umag = sqrt((uj[0] - u0[0]) * (uj[0] - u0[0]) +
391  (uj[1] - u0[1]) * (uj[1] - u0[1]));
392 
393  F[0] += wij * (uj[0] - u0[0]) / umag;
394  F[1] += wij * (uj[1] - u0[1]) / umag;
395 
396  Array<OneD, NekDouble> d1 = m_cadsurf->D1(u0);
397 
398  Array<OneD, NekDouble> dw(2, 0.0);
399  dw[0] = -2.0 *
400  ((rj[0] - pu0[0]) * d1[3] + (rj[1] - pu0[1]) * d1[4] +
401  (rj[2] - pu0[2]) * d1[5]);
402  dw[1] = -2.0 *
403  ((rj[0] - pu0[0]) * d1[6] + (rj[1] - pu0[1]) * d1[7] +
404  (rj[2] - pu0[2]) * d1[8]);
405 
406  Array<OneD, NekDouble> drhs(4, 0.0);
407  drhs[0] = 2.0 * ((uj[0] - u0[0]) * (uj[0] - u0[0]) - umag) /
408  umag / umag;
409  drhs[1] =
410  2.0 * ((uj[0] - u0[0]) * (uj[1] - u0[1])) / umag / umag;
411  drhs[2] = drhs[1];
412  drhs[3] = 2.0 * ((uj[1] - u0[1]) * (uj[1] - u0[1]) - umag) /
413  umag / umag;
414 
415  dF[0] += dw[0] * (uj[0] - u0[0]) / umag + wij * drhs[0];
416 
417  dF[1] += dw[0] * (uj[1] - u0[1]) / umag + wij * drhs[1];
418 
419  dF[2] += dw[1] * (uj[0] - u0[0]) / umag + wij * drhs[2];
420 
421  dF[3] += dw[1] * (uj[1] - u0[1]) / umag + wij * drhs[3];
422  }
423 
424  NekDouble det = dF[0] * dF[3] - dF[1] * dF[2];
425  NekDouble tmp = dF[0];
426  dF[0] = dF[3] / det;
427  dF[3] = tmp / det;
428  dF[1] *= -1.0 / det;
429  dF[2] *= -1.0 / det;
430 
431  u0[0] -= (dF[0] * F[0] + dF[2] * F[1]);
432  u0[1] -= (dF[1] * F[0] + dF[3] * F[1]);*/
433 
434  bool inbounds = true;
435  if (u0[0] < bounds[0])
436  {
437  inbounds = false;
438  }
439  else if (u0[0] > bounds[1])
440  {
441  inbounds = false;
442  }
443  else if (u0[1] < bounds[2])
444  {
445  inbounds = false;
446  }
447  else if (u0[1] > bounds[3])
448  {
449  inbounds = false;
450  }
451 
452  if (!inbounds)
453  {
454  continue;
455  }
456 
457  /*Array<OneD, NekDouble> FN(2, 0.0);
458  pu0 = m_cadsurf->P(u0);
459  di = m_mesh->m_octree->Query(pu0);
460  for (int i = 0; i < nodesystem.size(); i++)
461  {
462  Array<OneD, NekDouble> rj = nodesystem[i]->GetLoc();
463  Array<OneD, NekDouble> uj = nodesystem[i]->GetCADSurfInfo(m_id);
464  NekDouble dj = m_mesh->m_octree->Query(rj);
465  NekDouble d = (di + dj) / 2.0;
466  NekDouble wij = sqrt((rj[0] - pu0[0]) * (rj[0] - pu0[0]) +
467  (rj[1] - pu0[1]) * (rj[1] - pu0[1]) +
468  (rj[2] - pu0[2]) * (rj[2] - pu0[2])) - d;
469 
470  NekDouble umag = sqrt((uj[0] - u0[0]) * (uj[0] - u0[0]) +
471  (uj[1] - u0[1]) * (uj[1] - u0[1]));
472 
473  FN[0] += wij * (uj[0] - u0[0]) / umag;
474  FN[1] += wij * (uj[1] - u0[1]) / umag;
475  }
476 
477  if (F[0] * F[0] + F[1] * F[1] < FN[0] * FN[0] + FN[1] * FN[1])
478  {
479  continue;
480  }*/
481 
482  Array<OneD, NekDouble> l2 = m_cadsurf->P(u0);
483  (*nit)->Move(l2, m_id, u0);
484  }
485  }
486 }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
std::vector< ElementSharedPtr > m_localElements
local list of elements
Definition: FaceMesh.h:154
EdgeSet m_localEdges
local set of edges
Definition: FaceMesh.h:152
int m_id
id of the surface mesh
Definition: FaceMesh.h:140
boost::shared_ptr< Node > NodeSharedPtr
Definition: Node.h:50
double NekDouble
boost::shared_ptr< Edge > EdgeSharedPtr
Shared pointer to an edge.
Definition: Edge.h:135
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
CADSurfSharedPtr m_cadsurf
CAD surface.
Definition: FaceMesh.h:133
NodeSet m_localNodes
local set of nodes
Definition: FaceMesh.h:150
NodeSet m_inBoundary
set of nodes which are in the boundary (easier to identify conflicts with)
Definition: FaceMesh.h:156
void Nektar::NekMeshUtils::FaceMesh::Stretching ( )
private

Calculate the paramter plane streching factor.

Definition at line 944 of file FaceMesh.cpp.

945 {
946  // define a sampling and calculate the aspect ratio of the paramter plane
947  m_str = 0.0;
948  Array<OneD, NekDouble> bnds = m_cadsurf->GetBounds();
949 
950  NekDouble dxu = int(bnds[1] - bnds[0] < bnds[3] - bnds[2]
951  ? 40
952  : (bnds[1] - bnds[0]) / (bnds[3] - bnds[2]) * 40);
953  NekDouble dxv = int(bnds[3] - bnds[2] < bnds[1] - bnds[0]
954  ? 40
955  : (bnds[3] - bnds[2]) / (bnds[1] - bnds[0]) * 40);
956 
957  NekDouble du = (bnds[1] - bnds[0]) / dxu;
958  NekDouble dv = (bnds[3] - bnds[2]) / dxv;
959 
960  int ct = 0;
961 
962  for (int i = 0; i < dxu; i++)
963  {
964  for (int j = 0; j < dxv; j++)
965  {
966  Array<OneD, NekDouble> uv(2);
967  uv[0] = bnds[0] + i * du;
968  uv[1] = bnds[2] + j * dv;
969  if (i == dxu - 1)
970  {
971  uv[0] = bnds[1];
972  }
973  if (j == dxv - 1)
974  {
975  uv[1] = bnds[3];
976  }
977  Array<OneD, NekDouble> r = m_cadsurf->D1(uv);
978 
979  NekDouble ru = sqrt(r[3] * r[3] + r[4] * r[4] + r[5] * r[5]);
980  NekDouble rv = sqrt(r[6] * r[6] + r[7] * r[7] + r[8] * r[8]);
981 
982  ru *= du;
983  rv *= dv;
984 
985  if (rv < 1E-8)
986  {
987  continue;
988  }
989 
990  m_str += ru / rv;
991  ct++;
992  }
993  }
994 
995  m_str /= ct;
996 }
double NekDouble
NekDouble m_str
pplane stretching
Definition: FaceMesh.h:146
CADSurfSharedPtr m_cadsurf
CAD surface.
Definition: FaceMesh.h:133
bool Nektar::NekMeshUtils::FaceMesh::Validate ( )
private

Validate the surface mesh base on the octree and real dimensions of the edges.

Definition at line 998 of file FaceMesh.cpp.

999 {
1000  // check all edges in the current mesh for length against the octree
1001  // if the octree is not conformed to add a new point inside the triangle
1002  // if no new points are added meshing can stop
1003  int pointBefore = m_stienerpoints.size();
1004  for (int i = 0; i < m_connec.size(); i++)
1005  {
1006  Array<OneD, NekDouble> r(3), a(3);
1007 
1008  vector<Array<OneD, NekDouble> > info;
1009 
1010  for (int j = 0; j < 3; j++)
1011  {
1012  info.push_back(m_connec[i][j]->GetCADSurfInfo(m_id));
1013  }
1014 
1015  r[0] = m_connec[i][0]->Distance(m_connec[i][1]);
1016  r[1] = m_connec[i][1]->Distance(m_connec[i][2]);
1017  r[2] = m_connec[i][2]->Distance(m_connec[i][0]);
1018 
1019  a[0] = m_connec[i][0]->Angle(m_connec[i][1]->GetLoc(),
1020  m_connec[i][2]->GetLoc(),
1021  m_cadsurf->N(info[0]));
1022  a[1] = m_connec[i][1]->Angle(m_connec[i][2]->GetLoc(),
1023  m_connec[i][0]->GetLoc(),
1024  m_cadsurf->N(info[1]));
1025  a[2] = m_connec[i][2]->Angle(m_connec[i][0]->GetLoc(),
1026  m_connec[i][1]->GetLoc(),
1027  m_cadsurf->N(info[2]));
1028 
1029  NekDouble d1 = m_mesh->m_octree->Query(m_connec[i][0]->GetLoc());
1030  NekDouble d2 = m_mesh->m_octree->Query(m_connec[i][1]->GetLoc());
1031  NekDouble d3 = m_mesh->m_octree->Query(m_connec[i][2]->GetLoc());
1032 
1033  Array<OneD, NekDouble> uvc(2);
1034  uvc[0] = (info[0][0] + info[1][0] + info[2][0]) / 3.0;
1035  uvc[1] = (info[0][1] + info[1][1] + info[2][1]) / 3.0;
1036 
1037  Array<OneD, NekDouble> locc = m_cadsurf->P(uvc);
1038  NekDouble d4 = m_mesh->m_octree->Query(locc);
1039 
1040  NekDouble d = (d1 + d2 + d3 + d4) / 4.0;
1041 
1042  vector<bool> valid(3);
1043  valid[0] = r[0] < d * 1.5;
1044  valid[1] = r[1] < d * 1.5;
1045  valid[2] = r[2] < d * 1.5;
1046 
1047  vector<bool> angValid(3);
1048  angValid[0] = a[0] / M_PI * 180.0 > 20.0 && a[0] / M_PI * 180.0 < 120.0;
1049  angValid[1] = a[1] / M_PI * 180.0 > 20.0 && a[1] / M_PI * 180.0 < 120.0;
1050  angValid[2] = a[2] / M_PI * 180.0 > 20.0 && a[2] / M_PI * 180.0 < 120.0;
1051 
1052  int numValid = 0;
1053  int numAngValid = 0;
1054  for (int j = 0; j < 3; j++)
1055  {
1056  if (valid[j])
1057  {
1058  numValid++;
1059  }
1060  if (angValid[j])
1061  {
1062  numAngValid++;
1063  }
1064  }
1065 
1066  // if numvalid is zero no work to be done
1067  /*if (numValid != 3)
1068  {
1069  AddNewPoint(uvc);
1070  }*/
1071 
1072  if (numValid != 3 || numAngValid != 3)
1073  {
1074  // break the bad edge
1075  /*int a=0, b=0;
1076  if(!valid[0])
1077  {
1078  a = 0;
1079  b = 1;
1080  }
1081  else if(!valid[1])
1082  {
1083  a = 1;
1084  b = 2;
1085  }
1086  else if(!valid[2])
1087  {
1088  a = 2;
1089  b = 0;
1090  }
1091 
1092  Array<OneD, NekDouble> uvn(2);
1093  uvn[0] = (info[a][0] + info[b][0]) / 2.0;
1094  uvn[1] = (info[a][1] + info[b][1]) / 2.0;*/
1095  AddNewPoint(uvc);
1096  }
1097  }
1098 
1099  if (m_stienerpoints.size() == pointBefore)
1100  {
1101  return false;
1102  }
1103  else
1104  {
1105  return true;
1106  }
1107 }
std::vector< std::vector< NodeSharedPtr > > m_connec
triangle connectiviities
Definition: FaceMesh.h:148
int m_id
id of the surface mesh
Definition: FaceMesh.h:140
double NekDouble
std::vector< NodeSharedPtr > m_stienerpoints
list of stiener points in the triangulation
Definition: FaceMesh.h:144
void AddNewPoint(Array< OneD, NekDouble > uv)
adds a new stiener point to the triangulation for meshing
Definition: FaceMesh.cpp:1109
CADSurfSharedPtr m_cadsurf
CAD surface.
Definition: FaceMesh.h:133
MeshSharedPtr m_mesh
mesh pointer
Definition: FaceMesh.h:131
bool Nektar::NekMeshUtils::FaceMesh::ValidateCurves ( )

validate the curve meshes

Definition at line 47 of file FaceMesh.cpp.

48 {
49  vector<int> curvesInSurface;
50  for (int i = 0; i < m_edgeloops.size(); i++)
51  {
52  for (int j = 0; j < m_edgeloops[i]->edges.size(); j++)
53  {
54  curvesInSurface.push_back(m_edgeloops[i]->edges[j]->GetId());
55  }
56  }
57 
58  bool error = false;
59 
60  for (int i = 0; i < curvesInSurface.size(); i++)
61  {
62  vector<EdgeSharedPtr> es =
63  m_curvemeshes[curvesInSurface[i]]->GetMeshEdges();
64 
65  for (int j = i; j < curvesInSurface.size(); j++)
66  {
67  if (i == j)
68  {
69  continue;
70  }
71 
72  vector<EdgeSharedPtr> es2 =
73  m_curvemeshes[curvesInSurface[j]]->GetMeshEdges();
74 
75  for (int l = 0; l < es.size(); l++)
76  {
77  Array<OneD, NekDouble> P1 = es[l]->m_n1->GetCADSurfInfo(m_id);
78  Array<OneD, NekDouble> P2 = es[l]->m_n2->GetCADSurfInfo(m_id);
79  for (int k = 0; k < es2.size(); k++)
80  {
81  if (es[l]->m_n1 == es2[k]->m_n1 ||
82  es[l]->m_n1 == es2[k]->m_n2 ||
83  es[l]->m_n2 == es2[k]->m_n1 ||
84  es[l]->m_n2 == es2[k]->m_n2)
85  {
86  continue;
87  }
88 
89  Array<OneD, NekDouble> P3 =
90  es2[k]->m_n1->GetCADSurfInfo(m_id);
91  Array<OneD, NekDouble> P4 =
92  es2[k]->m_n2->GetCADSurfInfo(m_id);
93 
94  NekDouble den = (P4[0] - P3[0]) * (P2[1] - P1[1]) -
95  (P2[0] - P1[0]) * (P4[1] - P3[1]);
96 
97  if (fabs(den) < 1e-8)
98  {
99  continue;
100  }
101 
102  NekDouble t = ((P1[0] - P3[0]) * (P4[1] - P3[1]) -
103  (P4[0] - P3[0]) * (P1[1] - P3[1])) /
104  den;
105  NekDouble u =
106  (P1[0] - P3[0] + t * (P2[0] - P1[0])) / (P4[0] - P3[0]);
107 
108  if (t < 1.0 && t > 0.0 && u < 1.0 && u > 0.0)
109  {
110  Array<OneD, NekDouble> uv(2);
111  uv[0] = P1[0] + t * (P2[0] - P1[0]);
112  uv[1] = P1[1] + t * (P2[1] - P1[1]);
113  Array<OneD, NekDouble> loc = m_cadsurf->P(uv);
114  cout << endl
115  << "Curve mesh error at " << loc[0] << " "
116  << loc[1] << " " << loc[2] << " on face " << m_id
117  << endl;
118  error = true;
119  }
120  }
121  }
122  }
123  }
124  return error;
125 }
std::map< int, CurveMeshSharedPtr > m_curvemeshes
Map of the curve meshes which bound the surfaces.
Definition: FaceMesh.h:135
std::vector< CADSystem::EdgeLoopSharedPtr > m_edgeloops
data structure containing the edges, their order and oreientation for the surface ...
Definition: FaceMesh.h:138
int m_id
id of the surface mesh
Definition: FaceMesh.h:140
double NekDouble
CADSurfSharedPtr m_cadsurf
CAD surface.
Definition: FaceMesh.h:133

Friends And Related Function Documentation

friend class MemoryManager< FaceMesh >
friend

Definition at line 55 of file FaceMesh.h.

Member Data Documentation

CADSurfSharedPtr Nektar::NekMeshUtils::FaceMesh::m_cadsurf
private

CAD surface.

Definition at line 133 of file FaceMesh.h.

Referenced by FaceMesh().

int Nektar::NekMeshUtils::FaceMesh::m_compId
private

identity to put into element tags

Definition at line 158 of file FaceMesh.h.

std::vector<std::vector<NodeSharedPtr> > Nektar::NekMeshUtils::FaceMesh::m_connec
private

triangle connectiviities

Definition at line 148 of file FaceMesh.h.

std::map<int, CurveMeshSharedPtr> Nektar::NekMeshUtils::FaceMesh::m_curvemeshes
private

Map of the curve meshes which bound the surfaces.

Definition at line 135 of file FaceMesh.h.

std::vector<CADSystem::EdgeLoopSharedPtr> Nektar::NekMeshUtils::FaceMesh::m_edgeloops
private

data structure containing the edges, their order and oreientation for the surface

Definition at line 138 of file FaceMesh.h.

Referenced by FaceMesh().

int Nektar::NekMeshUtils::FaceMesh::m_id
private

id of the surface mesh

Definition at line 140 of file FaceMesh.h.

Referenced by FaceMesh().

NodeSet Nektar::NekMeshUtils::FaceMesh::m_inBoundary
private

set of nodes which are in the boundary (easier to identify conflicts with)

Definition at line 156 of file FaceMesh.h.

EdgeSet Nektar::NekMeshUtils::FaceMesh::m_localEdges
private

local set of edges

Definition at line 152 of file FaceMesh.h.

std::vector<ElementSharedPtr> Nektar::NekMeshUtils::FaceMesh::m_localElements
private

local list of elements

Definition at line 154 of file FaceMesh.h.

NodeSet Nektar::NekMeshUtils::FaceMesh::m_localNodes
private

local set of nodes

Definition at line 150 of file FaceMesh.h.

MeshSharedPtr Nektar::NekMeshUtils::FaceMesh::m_mesh
private

mesh pointer

Definition at line 131 of file FaceMesh.h.

Referenced by FaceMesh().

std::vector<NodeSharedPtr> Nektar::NekMeshUtils::FaceMesh::m_stienerpoints
private

list of stiener points in the triangulation

Definition at line 144 of file FaceMesh.h.

NekDouble Nektar::NekMeshUtils::FaceMesh::m_str
private

pplane stretching

Definition at line 146 of file FaceMesh.h.

std::vector<std::vector<NodeSharedPtr> > Nektar::NekMeshUtils::FaceMesh::orderedLoops
private

list of boundary nodes in their order loops

Definition at line 142 of file FaceMesh.h.