rework of node and cfi to conform to new CAD system

library/NekMeshUtils/2DGenerator/2DGenerator.cpp

@@ -269,7 +269,7 @@ void Generator2D::MakeBL(int faceid)
             new Node(m_mesh->m_numNodes++, n[0], n[1], 0.0));
         CADSurfSharedPtr s = m_mesh->m_cad->GetSurf(faceid);
         Array<OneD, NekDouble> uv = s->locuv(n);
-        nn->SetCADSurf(faceid, s, uv);
+        nn->SetCADSurf(s, uv);
         nodeNormals[it->first] = nn;
     }
 

library/NekMeshUtils/CADSystem/CADSurf.cpp

@@ -107,7 +107,7 @@ void CADSurf::OrientateEdges(CADSurfSharedPtr surf,
 
         ein[i]->area = area;
 
-        point_xy cen;
+        point_xy cen(0.0,0.0);
         bg::centroid(polygon, cen);
 
         ein[i]->center    = Array<OneD, NekDouble>(2);

library/NekMeshUtils/CADSystem/CADVert.cpp

+////////////////////////////////////////////////////////////////////////////////
+//
+//  File: CADSystem.cpp
+//
+//  For more information, please see: http://www.nektar.info/
+//
+//  The MIT License
+//
+//  Copyright (c) 2006 Division of Applied Mathematics, Brown University (USA),
+//  Department of Aeronautics, Imperial College London (UK), and Scientific
+//  Computing and Imaging Institute, University of Utah (USA).
+//
+//  License for the specific language governing rights and limitations under
+//  Permission is hereby granted, free of charge, to any person obtaining a
+//  copy of this software and associated documentation files (the "Software"),
+//  to deal in the Software without restriction, including without limitation
+//  the rights to use, copy, modify, merge, publish, distribute, sublicense,
+//  and/or sell copies of the Software, and to permit persons to whom the
+//  Software is furnished to do so, subject to the following conditions:
+//
+//  The above copyright notice and this permission notice shall be included
+//  in all copies or substantial portions of the Software.
+//
+//  THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
+//  OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+//  FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
+//  THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+//  LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
+//  FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
+//  DEALINGS IN THE SOFTWARE.
+//
+//  Description: cad object methods.
+//
+////////////////////////////////////////////////////////////////////////////////
+#include "CADVert.h"
+#include "CADSurf.h"
+#include <NekMeshUtils/MeshElements/Node.h>
+
+using namespace std;
+
+namespace Nektar
+{
+namespace NekMeshUtils
+{
+
+void CADVert::SetDegen(int s, CADSurfSharedPtr su, NekDouble u, NekDouble v)
+{
+    degen     = true;
+    degensurf = s;
+    Array<OneD, NekDouble> uv(2);
+    uv[0] = u;
+    uv[1] = v;
+    m_node->SetCADSurf(su, uv);
+}
+
+Array<OneD, NekDouble> CADVert::GetLoc()
+{
+    Array<OneD, NekDouble> out(3);
+    out[0] = m_node->m_x;
+    out[1] = m_node->m_y;
+    out[2] = m_node->m_z;
+    return out;
+}
+}
+}

library/NekMeshUtils/CADSystem/CADVert.h

@@ -40,15 +40,21 @@
 
 #include <LibUtilities/BasicUtils/SharedArray.hpp>
 #include <LibUtilities/Memory/NekMemoryManager.hpp>
+#include <LibUtilities/BasicUtils/NekFactory.hpp>
 
 #include <NekMeshUtils/CADSystem/CADObject.h>
-#include <NekMeshUtils/MeshElements/Node.h>
 
 namespace Nektar
 {
 namespace NekMeshUtils
 {
 
+//forward decleration
+class Node;
+typedef boost::shared_ptr<Node> NodeSharedPtr;
+class CADSurf;
+typedef boost::shared_ptr<CADSurf> CADSurfSharedPtr;
+
 /**
  * @brief base class for CAD verticies.
  *
@@ -71,14 +77,7 @@ public:
     /**
      * @brief Get x,y,z location of the vertex
      */
-    Array<OneD, NekDouble> GetLoc()
-    {
-        Array<OneD, NekDouble> out(3);
-        out[0] = m_node->m_x;
-        out[1] = m_node->m_y;
-        out[2] = m_node->m_z;
-        return out;
-    }
+    Array<OneD, NekDouble> GetLoc();
 
     /**
      * @brief returns a node object of the cad vertex
@@ -91,15 +90,7 @@ public:
     /**
      * @brief if the vertex is degenerate manually set uv for that surface
      */
-    void SetDegen(int s, CADSurfSharedPtr su, NekDouble u, NekDouble v)
-    {
-        degen     = true;
-        degensurf = s;
-        Array<OneD, NekDouble> uv(2);
-        uv[0] = u;
-        uv[1] = v;
-        m_node->SetCADSurf(s, su, uv);
-    }
+    void SetDegen(int s, CADSurfSharedPtr su, NekDouble u, NekDouble v);
 
     /**
      * @brief query is degenerate

library/NekMeshUtils/CADSystem/CFI/CADCurveCFI.cpp

@@ -57,7 +57,7 @@ void CADCurveCFI::Initialise(int i, cfi::Line *in, NekDouble s)
 NekDouble CADCurveCFI::tAtArcLength(NekDouble s)
 {
     s /= m_scal;
-    Array<OneD, NekDouble> bds = Bounds();
+    Array<OneD, NekDouble> bds = GetBounds();
     NekDouble dt = (bds[1] - bds[0]) / 1000;
 
     NekDouble t   = bds[0];
@@ -96,7 +96,7 @@ NekDouble CADCurveCFI::loct(Array<OneD, NekDouble> xyz)
 
 NekDouble CADCurveCFI::Length(NekDouble ti, NekDouble tf)
 {
-    Array<OneD, NekDouble> bds = Bounds();
+    Array<OneD, NekDouble> bds = GetBounds();
     NekDouble dt = (bds[1] - bds[0]) / 1000;
 
     NekDouble t   = ti;
@@ -157,7 +157,7 @@ Array<OneD, NekDouble> CADCurveCFI::D2(NekDouble t)
     return out;
 }
 
-Array<OneD, NekDouble> CADCurveCFI::Bounds()
+Array<OneD, NekDouble> CADCurveCFI::GetBounds()
 {
     Array<OneD, NekDouble> t(2);
     t[0] = 0.0;
@@ -168,7 +168,7 @@ Array<OneD, NekDouble> CADCurveCFI::Bounds()
 
 Array<OneD, NekDouble> CADCurveCFI::GetMinMax()
 {
-    Array<OneD, NekDouble> bds = Bounds();
+    Array<OneD, NekDouble> bds = GetBounds();
 
     cfi::Position x1 = m_cfiEdge->calcXYZAtT(bds[0]);
     cfi::Position x2 = m_cfiEdge->calcXYZAtT(bds[1]);

library/NekMeshUtils/CADSystem/CFI/CADCurveCFI.h

@@ -53,29 +53,36 @@ public:
 
     static std::string key;
 
-    CADCurveCFI(){};
-
-    ~CADCurveCFI(){};
-
-    Array<OneD, NekDouble> Bounds();
-
-    NekDouble Length(NekDouble ti, NekDouble tf);
-
-    Array<OneD, NekDouble> P(NekDouble t);
+    CADCurveCFI()
+    {
+    }
 
-    Array<OneD, NekDouble> D2(NekDouble t);
+    ~CADCurveCFI()
+    {
+    }
 
-    NekDouble Curvature(NekDouble t)
+    virtual Array<OneD, NekDouble> GetBounds();
+    virtual NekDouble Length(NekDouble ti, NekDouble tf);
+    virtual Array<OneD, NekDouble> P(NekDouble t);
+    virtual Array<OneD, NekDouble> D2(NekDouble t);
+    virtual NekDouble Curvature(NekDouble t)
     {
         ASSERTL0(false, "not implemented in CFI engine");
         return 0;
     }
-
-    NekDouble tAtArcLength(NekDouble s);
-
-    Array<OneD, NekDouble> GetMinMax();
-
-    NekDouble loct(Array<OneD, NekDouble> xyz);
+    virtual Array<OneD, NekDouble> NormalWRT(NekDouble t, int surf)
+    {
+        ASSERTL0(false, "not implemented in CFI engine");
+        return Array<OneD, NekDouble>();
+    }
+    virtual Array<OneD, NekDouble> N(NekDouble t)
+    {
+        ASSERTL0(false, "not implemented in CFI engine");
+        return Array<OneD, NekDouble>();
+    }
+    virtual NekDouble tAtArcLength(NekDouble s);
+    virtual Array<OneD, NekDouble> GetMinMax();
+    virtual NekDouble loct(Array<OneD, NekDouble> xyz);
 
     void Initialise(int i, cfi::Line *in, NekDouble s);
 

library/NekMeshUtils/CADSystem/CFI/CADSurfCFI.cpp

@@ -45,12 +45,11 @@ namespace NekMeshUtils
 std::string CADSurfCFI::key = GetCADSurfFactory().RegisterCreatorFunction(
     "cfi", CADSurfCFI::create, "CADSurfCFI");
 
-void CADSurfCFI::Initialise(int i, cfi::Face *in, vector<EdgeLoopSharedPtr> ein)
+void CADSurfCFI::Initialise(int i, cfi::Face *in, NekDouble s)
 {
-    m_edges         = ein;
-    m_cfiSurface    = in;
-    m_correctNormal = true;
-    m_id            = i;
+    m_cfiSurface = in;
+    m_id         = i;
+    m_scal       = s;
 }
 
 Array<OneD, NekDouble> CADSurfCFI::GetBounds()

library/NekMeshUtils/CADSystem/CFI/CADSurfCFI.h

@@ -58,7 +58,7 @@ public:
 
     ~CADSurfCFI(){};
 
-    void Initialise(int i, cfi::Face *in, std::vector<EdgeLoopSharedPtr> ein);
+    void Initialise(int i, cfi::Face *in, NekDouble s);
     void SetScaling(NekDouble i)
     {
         m_scal = i;

library/NekMeshUtils/CADSystem/CFI/CADSystemCFI.cpp

@@ -82,8 +82,18 @@ bool CADSystemCFI::LoadCAD()
     model = cfiHandel.openModelFile(m_name.c_str());
 
     // make an assumption there are not multiple bodies in the solid
-    ASSERTL0(model->getEntityTotal(cfi::TYPE_BODY, cfi::SUBTYPE_ALL) == 1,
-             "cannot deal with multibodies");
+    if(model->getEntityTotal(cfi::TYPE_BODY, cfi::SUBTYPE_ALL) > 1)
+    {
+        cout << "NekMesh cannot deal with multiple CAD bodies"
+             << "it is going to assume the body you want is called W1" << endl;
+
+        body = static_cast<cfi::Body *>(model->getEntity("W1"));
+    }
+
+    if(model->getEntityTotal(cfi::TYPE_BODY, cfi::SUBTYPE_ALL) == 1)
+    {
+        body = model->getBodyEntity(1);
+    }
 
     // cfi doesnt mind stupid units so this scales everything back to meters
     // which is what nekmesh assumes its in
@@ -96,8 +106,6 @@ bool CADSystemCFI::LoadCAD()
         m_scal = 0.0254;
     }
 
-    body = model->getBodyEntity(1);
-
     // CFI does everything by string identifers
     // currently nekmesh cad system uses integer ids.
     // it really should use strings but doesnt currently
@@ -138,7 +146,6 @@ bool CADSystemCFI::LoadCAD()
     }
 
     // make the vertices and build a map of name to id
-    map<string, int> nameToVertId;
     map<string, cfi::Point *>::iterator vit;
     int i = 1; // from one to be consistent with oce
     for (vit = mapOfVerts.begin(); vit != mapOfVerts.end(); vit++, i++)
@@ -174,93 +181,25 @@ bool CADSystemCFI::LoadCAD()
         cfi::Oriented<cfi::TopoEntity *> orientatedFace = *it;
         cfi::Face *face = static_cast<cfi::Face *>(orientatedFace.entity);
         nameToFaceId[face->getName()] = i;
-        vector<cfi::Oriented<cfi::TopoEntity *> > *edgeList =
-            face->getChildList();
-
-        vector<EdgeLoopSharedPtr> edgeloops;
-        vector<vector<cfi::Oriented<cfi::TopoEntity *> > > cfiloops;
-        int done = 0;
-        while (done != edgeList->size())
-        {
-            EdgeLoopSharedPtr edgeloop = EdgeLoopSharedPtr(new EdgeLoop);
-            vector<cfi::Oriented<cfi::TopoEntity *> > cfiloop;
-            string firstVert;
-            vector<cfi::Oriented<cfi::TopoEntity *> > *vertList =
-                edgeList->at(done).entity->getChildList();
-            if (edgeList->at(done).orientation == 1)
-            {
-                firstVert = vertList->at(0).entity->getName();
-            }
-            else
-            {
-                firstVert = vertList->at(1).entity->getName();
-            }
 
-            edgeloop->edges.push_back(
-                m_curves[nameToCurveId[edgeList->at(done).entity->getName()]]);
-            cfiloop.push_back(edgeList->at(done));
-            adjsurfmap[nameToCurveId[edgeList->at(done).entity->getName()]]
-                .push_back(i);
-            edgeList->at(done).orientation == 1 ? edgeloop->edgeo.push_back(0)
-                                                : edgeloop->edgeo.push_back(1);
-
-            for (done++; done < edgeList->size(); done++)
-            {
-                bool end = false;
-                vertList = edgeList->at(done).entity->getChildList();
-                if (edgeList->at(done).orientation == 1)
-                {
-                    if (vertList->at(1).entity->getName() == firstVert)
-                    {
-                        end = true;
-                    }
-                }
-                else
-                {
-                    if (vertList->at(0).entity->getName() == firstVert)
-                    {
-                        end = true;
-                    }
-                }
-
-                edgeloop->edges.push_back(
-                    m_curves
-                        [nameToCurveId[edgeList->at(done).entity->getName()]]);
-                cfiloop.push_back(edgeList->at(done));
-                adjsurfmap[nameToCurveId[edgeList->at(done).entity->getName()]]
-                    .push_back(i);
-                edgeList->at(done).orientation == 1
-                    ? edgeloop->edgeo.push_back(0)
-                    : edgeloop->edgeo.push_back(1);
-
-                if (end)
-                {
-                    done++;
-                    break;
-                }
-            }
-            cfiloops.push_back(cfiloop);
-            edgeloops.push_back(edgeloop);
-        }
-
-        AddSurf(i, face, edgeloops);
+        AddSurf(i, face);
     }
 
+    exit(-1);
+
     // TODO identify Degenerated faces and setdegen on vertices accordinaly
 
     ASSERTL0(adjsurfmap.size() == m_curves.size(), "incorrect curve info");
 
     // This checks that all edges are bound by two surfaces, sanity check.
-    for (map<int, vector<int> >::iterator it = adjsurfmap.begin();
-         it != adjsurfmap.end(); it++)
+    if (!m_2d)
     {
-        ASSERTL0(it->second.size() == 2, "no three curve surfaces");
-        vector<CADSurfSharedPtr> sfs;
-        for (int i = 0; i < it->second.size(); i++)
+        map<int, CADCurveSharedPtr>::iterator it;
+        for (it = m_curves.begin(); it != m_curves.end(); it++)
         {
-            sfs.push_back(m_surfs[it->second[i]]);
+            ASSERTL0(it->second->GetAdjSurf().size() == 2,
+                     "curve is not joined to 2 surfaces");
         }
-        m_curves[it->first]->SetAdjSurf(sfs);
     }
 
     return true;
@@ -297,29 +236,107 @@ void CADSystemCFI::AddCurve(int i, cfi::Line *in)
     ASSERTL0(t[0] < t[1], "weirdness");
 
     vector<CADVertSharedPtr> vs;
-    vs.push_back(m_verts[nameToVertId[in->at(0)->getName()]]);
-    vs.push_back(m_verts[nameToVertId[in->at(1)->getName()]]);
+    vs.push_back(m_verts[nameToVertId[vertList->at(0).entity->getName()]]);
+    vs.push_back(m_verts[nameToVertId[vertList->at(1).entity->getName()]]);
     m_curves[i] = newCurve;
     m_curves[i]->SetVert(vs);
 }
 
-void CADSystemCFI::AddSurf(int i, cfi::Face *in,
-                           std::vector<EdgeLoopSharedPtr> ein)
+void CADSystemCFI::AddSurf(int i, cfi::Face *in)
 {
     CADSurfSharedPtr newSurf = GetCADSurfFactory().CreateInstance(key);
-    static_pointer_cast<CADSurfCFI>(newSurf)->Initialise(i, in, ein);
-    static_pointer_cast<CADSurfCFI>(newSurf)->SetScaling(m_scal);
-    m_surfs[i] = newSurf;
+    static_pointer_cast<CADSurfCFI>(newSurf)->Initialise(i, in, m_scal);
+
+    vector<cfi::Oriented<cfi::TopoEntity *> > *edgeList = in->getChildList();
+
+    vector<cfi::Oriented<cfi::TopoEntity *> >::iterator it;
+    for(it = edgeList->begin(); it != edgeList->end(); it++)
+    {
+        if( (*it).orientation == cfi::ORIENT_POSITIVE)
+        {
+            cout << "positive" << endl;
+        }
+        else if( (*it).orientation == cfi::ORIENT_NEGATIVE)
+        {
+            cout << "negative" << endl;
+        }
+        else
+        {
+            cout << "unsure" << endl;
+        }
+    }
+
+    vector<EdgeLoopSharedPtr> edgeloops;
+    /*int done = 0;
+    while (done != edgeList->size())
+    {
+        EdgeLoopSharedPtr edgeloop = EdgeLoopSharedPtr(new EdgeLoop);
+        string firstVert;
+        vector<cfi::Oriented<cfi::TopoEntity *> > *vertList =
+            edgeList->at(done).entity->getChildList();
+        if (edgeList->at(done).orientation == 1)
+        {
+            firstVert = vertList->at(0).entity->getName();
+        }
+        else
+        {
+            firstVert = vertList->at(1).entity->getName();
+        }
+
+        edgeloop->edges.push_back(
+            m_curves[nameToCurveId[edgeList->at(done).entity->getName()]]);
+        adjsurfmap[nameToCurveId[edgeList->at(done).entity->getName()]]
+            .push_back(i);
+        edgeList->at(done).orientation == 1 ? edgeloop->edgeo.push_back(0)
+                                            : edgeloop->edgeo.push_back(1);
+
+        for (done++; done < edgeList->size(); done++)
+        {
+            bool end = false;
+            vertList = edgeList->at(done).entity->getChildList();
+            if (edgeList->at(done).orientation == 1)
+            {
+                if (vertList->at(1).entity->getName() == firstVert)
+                {
+                    end = true;
+                }
+            }
+            else
+            {
+                if (vertList->at(0).entity->getName() == firstVert)
+                {
+                    end = true;
+                }
+            }
+
+            edgeloop->edges.push_back(
+                m_curves
+                    [nameToCurveId[edgeList->at(done).entity->getName()]]);
+            edgeList->at(done).orientation == 1
+                ? edgeloop->edgeo.push_back(0)
+                : edgeloop->edgeo.push_back(1);
+
+            if (end)
+            {
+                done++;
+                break;
+            }
+        }
+        edgeloops.push_back(edgeloop);
+    }*/
 
     // TODO find if surface has reversed normal or not
 
-    int tote = 0;
+    /*int tote = 0;
     for (int i = 0; i < ein.size(); i++)
     {
         tote += ein[i]->edges.size();
     }
 
     ASSERTL0(tote != 1, "cannot handle periodic curves");
+
+
+    m_surfs[i] = newSurf;*/
 }
 
 Array<OneD, NekDouble> CADSystemCFI::GetBoundingBox()

library/NekMeshUtils/CADSystem/CFI/CADSystemCFI.h

@@ -95,6 +95,7 @@ private:
     cfi::Cfi cfiHandel;
     cfi::Model *model;
     cfi::Body *body;
+    std::map<std::string, int> nameToVertId;
     std::map<std::string, int> nameToCurveId;
     std::map<std::string, int> nameToFaceId;
     std::map<std::string, std::vector<std::string> > mapVertToListEdge;

library/NekMeshUtils/CADSystem/CFI/CADVertCFI.cpp

@@ -36,6 +36,8 @@
 #include "CADSystemCFI.h"
 #include "CADVertCFI.h"
 
+#include <NekMeshUtils/MeshElements/Node.h>
+
 using namespace std;
 
 namespace Nektar

library/NekMeshUtils/CADSystem/OCE/CADVertOCE.cpp

@@ -36,6 +36,8 @@
 #include <NekMeshUtils/CADSystem/OCE/CADSystemOCE.h>
 #include <NekMeshUtils/CADSystem/OCE/CADVertOCE.h>
 
+#include <NekMeshUtils/MeshElements/Node.h>
+
 using namespace std;
 
 namespace Nektar

library/NekMeshUtils/CFIMesh/CFIMesh.cpp

@@ -67,7 +67,7 @@ void CFIMesh::Process()
     m_mesh->m_expDim   = 3;
     m_mesh->m_spaceDim = 3;
 
-    CADSystemCFISharedPtr cad =
+    /*CADSystemCFISharedPtr cad =
         boost::dynamic_pointer_cast<CADSystemCFI>(m_mesh->m_cad);
     map<string, int> nameToCurveId                  = cad->GetCFICurveId();
     map<string, int> nameToFaceId                   = cad->GetCFIFaceId();
@@ -104,11 +104,12 @@ void CFIMesh::Process()
             CADCurveSharedPtr c = m_mesh->m_cad->GetCurve(cid);
             NekDouble t         = c->loct(xyz);
             n->SetCADCurve(cid, c, t);
-            vector<CADSurfSharedPtr> ss = c->GetAdjSurf();
+            vector<pair<CADSurfSharedPtr, CADSystem::Orientation> > ss =
+                c->GetAdjSurf();
             for (int j = 0; j < ss.size(); j++)
             {
                 Array<OneD, NekDouble> uv = ss[j]->locuv(xyz);
-                n->SetCADSurf(ss[j]->GetId(), ss[j], uv);
+                n->SetCADSurf(ss[j].first->GetId(), ss[j].first, uv);
             }
         }
         else if (p->type == cfi::TYPE_FACE)
@@ -127,11 +128,12 @@ void CFIMesh::Process()
                 CADCurveSharedPtr c = m_mesh->m_cad->GetCurve(cid);
                 NekDouble t         = c->loct(xyz);
                 n->SetCADCurve(cid, c, t);
-                vector<CADSurfSharedPtr> ss = c->GetAdjSurf();
+                vector<pair<CADSurfSharedPtr, CADSystem::Orientation> > ss =
+                    c->GetAdjSurf();
                 for (int j = 0; j < ss.size(); j++)
                 {
                     Array<OneD, NekDouble> uv = ss[j]->locuv(xyz);
-                    n->SetCADSurf(ss[j]->GetId(), ss[j], uv);
+                    n->SetCADSurf(ss[j].first->GetId(), ss[j].first, uv);
                 }
             }
         }
@@ -169,7 +171,7 @@ void CFIMesh::Process()
         m_mesh->m_element[3].push_back(E);
     }
 
-    //there must be a surface mesh as well or the high-order doesnt work
+    // there must be a surface mesh as well or the high-order doesnt work
 
     vector<cfi::ElementDefinition> *tris =
         model->getElements(cfi::SUBTYPE_TR3, 3);
@@ -258,7 +260,7 @@ void CFIMesh::Process()
     ProcessEdges();
     ProcessFaces();
     ProcessElements();