Commit 45fdad16 authored by Dave Moxey's avatar Dave Moxey

Merge branch 'feature/CADfix-export' into 'master'

CADfix developments in NekMesh

Closes #162 and #164

See merge request !938
parents 24feb9db ea52836a
Pipeline #1660 passed with stages
in 191 minutes and 26 seconds
......@@ -9,6 +9,10 @@ v5.1.0
**CardiacEPSolver**
- Added additional parameter sets to Fenton-Karma model (!1119)
**NekMesh**
- Improved boundary layer splitting and output to CADfix (!938)
v5.0.1
------
**Library**
......
......@@ -24,7 +24,7 @@ IF(NEKTAR_USE_MESHGEN)
IF(CFI_INCLUDE_DIR)
MESSAGE(STATUS "cfi libraries: ${CFI_LIBRARY_API}")
MESSAGE(STATUS "Found CFI Libraries: ${CFI_LIBRARY_API}")
INCLUDE_DIRECTORIES(NekMeshUtils ${CFI_INCLUDE_DIR_HXX})
INCLUDE_DIRECTORIES(NekMeshUtils ${CFI_INCLUDE_DIR})
......@@ -37,3 +37,5 @@ IF(NEKTAR_USE_MESHGEN)
ENDIF()
ENDIF()
ENDIF()
INCLUDE_DIRECTORIES(SYSTEM ${CFI_INCLUDE_DIR_HXX})
......@@ -77,12 +77,13 @@ public:
*
* @return Array of two entries, min and max parametric coordinate.
*/
virtual Array<OneD, NekDouble> GetBounds() = 0;
NEKMESHUTILS_EXPORT virtual Array<OneD, NekDouble> GetBounds() = 0;
/**
* @brief Returns the minimum and maximum parametric coords t of the curve.
*/
virtual void GetBounds(NekDouble &tmin, NekDouble &tmax) = 0;
NEKMESHUTILS_EXPORT virtual void GetBounds(
NekDouble &tmin, NekDouble &tmax) = 0;
/**
* @brief Calculates the arclength between the two paremetric points \p ti
......@@ -92,7 +93,8 @@ public:
* @param tf Second parametric coordinate.
* @return Arc length between \p ti and \p tf.
*/
virtual NekDouble Length(NekDouble ti, NekDouble tf) = 0;
NEKMESHUTILS_EXPORT virtual NekDouble Length(
NekDouble ti, NekDouble tf) = 0;
/**
* @brief Gets the location (x,y,z) in an array out of the curve at
......@@ -101,7 +103,7 @@ public:
* @param t Parametric coordinate
* @return Array of x,y,z
*/
virtual Array<OneD, NekDouble> P(NekDouble t) = 0;
NEKMESHUTILS_EXPORT virtual Array<OneD, NekDouble> P(NekDouble t) = 0;
/**
* @brief Gets the location (x,y,z) in an array out of the curve at
......@@ -109,17 +111,18 @@ public:
*
* @param t Parametric coordinate
*/
virtual void P(NekDouble t, NekDouble &x, NekDouble &y, NekDouble &z) = 0;
NEKMESHUTILS_EXPORT virtual void P(
NekDouble t, NekDouble &x, NekDouble &y, NekDouble &z) = 0;
/**
* @brief Gets the second derivatives at t
*/
virtual Array<OneD, NekDouble> D2(NekDouble t) = 0;
NEKMESHUTILS_EXPORT virtual Array<OneD, NekDouble> D2(NekDouble t) = 0;
/**
* @brief Calculates the radius of curvature of the curve at point t
*/
virtual NekDouble Curvature(NekDouble t) = 0;
NEKMESHUTILS_EXPORT virtual NekDouble Curvature(NekDouble t) = 0;
/**
* @brief Calculates the parametric coordinate and arclength location
......@@ -130,14 +133,14 @@ public:
*
* @todo This really needs improving for accuracy.
*/
virtual NekDouble tAtArcLength(NekDouble s) = 0;
NEKMESHUTILS_EXPORT virtual NekDouble tAtArcLength(NekDouble s) = 0;
/**
* @brief Gets the start and end of the curve.
*
* @return Array with 6 entries of endpoints x1,y1,z1,x2,y2,z2.
*/
virtual Array<OneD, NekDouble> GetMinMax() = 0;
NEKMESHUTILS_EXPORT virtual Array<OneD, NekDouble> GetMinMax() = 0;
/**
* @brief set the ids of the surfaces either side of the curve
......@@ -186,25 +189,26 @@ public:
* @brief locates a point in the parametric space. returns the
* distance to the point and passes t by reference and updates it
*/
virtual NekDouble loct(Array<OneD, NekDouble> xyz, NekDouble &t) = 0;
NEKMESHUTILS_EXPORT virtual NekDouble loct(
Array<OneD, NekDouble> xyz, NekDouble &t) = 0;
/**
* @brief Returns the orientation of the curve with respect to a given
* surface by id surf
*/
CADOrientation::Orientation GetOrienationWRT(int surf);
NEKMESHUTILS_EXPORT CADOrientation::Orientation GetOrienationWRT(int surf);
/**
* @brief Returns the normal to the curve which is orientate with respect
* to the surface surf
*/
Array<OneD, NekDouble> NormalWRT(NekDouble t, int surf);
NEKMESHUTILS_EXPORT Array<OneD, NekDouble> NormalWRT(NekDouble t, int surf);
/**
* @brief Returns the normal to a curve, it will always point in the concave
* direction
*/
virtual Array<OneD, NekDouble> N(NekDouble t) = 0;
NEKMESHUTILS_EXPORT virtual Array<OneD, NekDouble> N(NekDouble t) = 0;
protected:
/// Length of edge
......
......@@ -36,6 +36,7 @@
#define NEKMESHUTILS_CADSYSTEM_CADOBJ
#include <LibUtilities/Memory/NekMemoryManager.hpp>
#include <NekMeshUtils/NekMeshUtilsDeclspec.h>
namespace Nektar
{
......@@ -48,7 +49,8 @@ enum cadType
{
eVert,
eCurve,
eSurf
eSurf,
eOther
};
}
......
......@@ -49,16 +49,6 @@ namespace Nektar
namespace NekMeshUtils
{
Array<OneD, NekDouble> CADSurf::locuv(Array<OneD, NekDouble> p)
{
NekDouble dist;
Array<OneD, NekDouble> uv = locuv(p, dist);
WARNINGL1(dist < 1e-3, "large locuv distance");
return uv;
}
void CADSurf::OrientateEdges(CADSurfSharedPtr surf,
vector<EdgeLoopSharedPtr> &ein)
{
......
......@@ -113,13 +113,13 @@ public:
*
* @return Array of 4 entries with parametric umin,umax,vmin,vmax.
*/
virtual Array<OneD, NekDouble> GetBounds() = 0;
NEKMESHUTILS_EXPORT virtual Array<OneD, NekDouble> GetBounds() = 0;
/**
* @brief Get the limits of the parametric space for the surface.
*/
virtual void GetBounds(NekDouble &umin, NekDouble &umax, NekDouble &vmin,
NekDouble &vmax) = 0;
NEKMESHUTILS_EXPORT virtual void GetBounds(
NekDouble &umin, NekDouble &umax, NekDouble &vmin, NekDouble &vmax) = 0;
/**
* @brief Get the normal vector at parametric point u,v.
......@@ -127,7 +127,8 @@ public:
* @param uv Array of u and v parametric coords.
* @return Array of xyz components of normal vector.
*/
virtual Array<OneD, NekDouble> N(Array<OneD, NekDouble> uv) = 0;
NEKMESHUTILS_EXPORT virtual Array<OneD, NekDouble> N(
Array<OneD, NekDouble> uv) = 0;
/**
* @brief Get the set of first derivatives at parametric point u,v
......@@ -135,7 +136,8 @@ public:
* @param uv Array of u and v parametric coords.
* @return Array of xyz copmonents of first derivatives.
*/
virtual Array<OneD, NekDouble> D1(Array<OneD, NekDouble> uv) = 0;
NEKMESHUTILS_EXPORT virtual Array<OneD, NekDouble> D1(
Array<OneD, NekDouble> uv) = 0;
/**
* @brief Get the set of second derivatives at parametric point u,v
......@@ -143,7 +145,8 @@ public:
* @param uv array of u and v parametric coords
* @return array of xyz copmonents of second derivatives
*/
virtual Array<OneD, NekDouble> D2(Array<OneD, NekDouble> uv) = 0;
NEKMESHUTILS_EXPORT virtual Array<OneD, NekDouble> D2(
Array<OneD, NekDouble> uv) = 0;
/**
* @brief Get the x,y,z at parametric point u,v.
......@@ -151,15 +154,17 @@ public:
* @param uv Array of u and v parametric coords.
* @return Array of x,y,z location.
*/
virtual Array<OneD, NekDouble> P(Array<OneD, NekDouble> uv) = 0;
NEKMESHUTILS_EXPORT virtual Array<OneD, NekDouble> P(
Array<OneD, NekDouble> uv) = 0;
/**
* @brief Get the x,y,z at parametric point u,v.
*
* @param uv Array of u and v parametric coords.
*/
virtual void P(Array<OneD, NekDouble> uv, NekDouble &x, NekDouble &y,
NekDouble &z) = 0;
NEKMESHUTILS_EXPORT virtual void P(
Array<OneD, NekDouble> uv, NekDouble &x, NekDouble &y,
NekDouble &z) = 0;
/**
* @brief Performs a reverse look up to find u,v and x,y,z.
......@@ -168,8 +173,8 @@ public:
* @param p Array of xyz location
* @return The parametric location of xyz on this surface
*/
virtual Array<OneD, NekDouble> locuv(Array<OneD, NekDouble> p,
NekDouble &dist) = 0;
NEKMESHUTILS_EXPORT virtual Array<OneD, NekDouble> locuv(
Array<OneD, NekDouble> p, NekDouble &dist) = 0;
/**
* @brief overload function of locuv ommiting the dist parameter
......@@ -177,23 +182,30 @@ public:
* it will produce a warning. To do large distance projection use the other
* locuv method
*/
Array<OneD, NekDouble> locuv(Array<OneD, NekDouble> p);
Array<OneD, NekDouble> locuv(Array<OneD, NekDouble> p)
{
NekDouble dist;
Array<OneD, NekDouble> uv = locuv(p, dist);
WARNINGL1(dist < 1e-3, "large locuv distance");
return uv;
}
/**
* @brief Returns the bounding box of the surface
*/
virtual Array<OneD, NekDouble> BoundingBox() = 0;
NEKMESHUTILS_EXPORT virtual Array<OneD, NekDouble> BoundingBox() = 0;
/**
* @brief returns curvature at point uv
*/
virtual NekDouble Curvature(Array<OneD, NekDouble> uv) = 0;
NEKMESHUTILS_EXPORT virtual NekDouble Curvature(
Array<OneD, NekDouble> uv) = 0;
/**
* @brief Is the surface defined by a planar surface (i.e not nurbs and is
* flat)
*/
virtual bool IsPlanar() = 0;
NEKMESHUTILS_EXPORT virtual bool IsPlanar() = 0;
/**
* @brief query reversed normal
......
......@@ -71,11 +71,9 @@ public:
/**
* @brief Default constructor.
*/
CADSystem(std::string name) : m_name(name)
CADSystem(std::string name, std::string engine)
: m_name(name), m_engine(engine)
{
m_2d = false;
m_cfiMesh = false;
m_verbose = false;
}
virtual ~CADSystem()
......@@ -90,29 +88,52 @@ public:
return m_name;
}
/**
* @brief Set the CAD engine up to handle 2D geometries explicitly.
*/
void Set2D()
{
m_2d = true;
}
/**
* @brief Returns whether the CAD engine is set in 2D mode.
*/
bool Is2D()
{
return m_2d;
}
void SetNACA(std::string i)
/**
* @brief Gets a configuration option.
*
* @param key Configuration key.
*
* @return The configuration value.
*/
const std::string &GetConfig(const std::string &key) const
{
m_naca = i;
auto it = m_config.find(key);
return it->second;
}
void SetCFIMesh()
/**
* @brief Sets a configuration option @p key to @p value.
*
* @param key Configuration key.
* @param value The configuration value.
*/
void SetConfig(const std::string &key, const std::string &value)
{
m_cfiMesh = true;
m_config[key] = value;
}
void SetVerbose()
/**
* @brief Sets the verbose flag for additional output.
*/
void SetVerbose(bool verbose = true)
{
m_verbose = true;
m_verbose = verbose;
}
/**
......@@ -120,7 +141,7 @@ public:
*
* @return true if completed successfully
*/
virtual bool LoadCAD() = 0;
NEKMESHUTILS_EXPORT virtual bool LoadCAD() = 0;
/**
* @brief Returns bounding box of the domain.
......@@ -130,7 +151,7 @@ public:
*
* @return Array with 6 entries: xmin, xmax, ymin, ymax, zmin and zmax.
*/
virtual Array<OneD, NekDouble> GetBoundingBox() = 0;
NEKMESHUTILS_EXPORT virtual Array<OneD, NekDouble> GetBoundingBox() = 0;
/**
* @brief Get the number of surfaces.
......@@ -182,7 +203,7 @@ public:
}
/**
* @brief Gets map of all vertices
* @brief Gets map of all vertices.
*/
std::map<int, CADVertSharedPtr> GetVerts()
{
......@@ -192,7 +213,7 @@ public:
/**
* @brief Gets number of vertices
*/
int GetNumVerts()
int GetNumVerts() const
{
return m_verts.size();
}
......@@ -205,9 +226,19 @@ public:
NEKMESHUTILS_EXPORT Array<OneD, NekDouble> GetPeriodicTranslationVector(
int first, int second);
/**
* @brief Return the name of the CAD system engine (e.g. OCE).
*/
std::string GetEngine() const
{
return m_engine;
}
protected:
/// Name of cad file
/// Name of CAD file
std::string m_name;
/// Name of CAD engine.
std::string m_engine;
/// Map of curves
std::map<int, CADCurveSharedPtr> m_curves;
/// Map of surfaces
......@@ -215,13 +246,12 @@ protected:
/// Map of vertices
std::map<int, CADVertSharedPtr> m_verts;
/// Verbosity
bool m_verbose;
bool m_verbose = false;
/// 2D cad flag
bool m_2d;
/// Will the CAD be used with a CFI mesh flag
bool m_cfiMesh;
/// string of 4 digit NACA code to be created
std::string m_naca;
bool m_2d = false;
/// Configuration options which might be used per-engine, serialised in
/// strings.
std::map<std::string, std::string> m_config;
/**
* @brief Reports basic properties to screen.
......
......@@ -133,7 +133,7 @@ Array<OneD, NekDouble> CADCurveCFI::P(NekDouble t)
return out;
}
void P(NekDouble t, NekDouble &x, NekDouble &y, NekDouble &z)
void CADCurveCFI::P(NekDouble t, NekDouble &x, NekDouble &y, NekDouble &z)
{
cfi::Position p = m_cfiEdge->calcXYZAtT(t);
......
......@@ -68,11 +68,13 @@ public:
virtual Array<OneD, NekDouble> D2(NekDouble t);
virtual NekDouble Curvature(NekDouble t)
{
boost::ignore_unused(t);
ASSERTL0(false, "Function: Curvature not implemented in CFI engine");
return 0;
}
virtual Array<OneD, NekDouble> N(NekDouble t)
{
boost::ignore_unused(t);
ASSERTL0(false, "Function: N not implemented in CFI engine");
return Array<OneD, NekDouble>();
}
......@@ -82,6 +84,11 @@ public:
void Initialise(int i, cfi::Line *in, NekDouble s);
cfi::Line *GetCfiPointer()
{
return m_cfiEdge;
}
private:
/// cfi object
cfi::Line *m_cfiEdge;
......
////////////////////////////////////////////////////////////////////////////////
//
// File: CADElementCFI.h
//
// 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).
//
// 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: Storage for CFI element handle.
//
////////////////////////////////////////////////////////////////////////////////
#ifndef NEKMESHUTILS_CADSYSTEM_CFI_CADELEMENTCFI
#define NEKMESHUTILS_CADSYSTEM_CFI_CADELEMENTCFI
#include "CADSystemCFI.h"
#include "../CADObject.h"
namespace Nektar
{
namespace NekMeshUtils
{
class CADElementCFI : public CADObject
{
public:
CADElementCFI(cfi::MeshableEntity *cfiEntity) : m_cfiEntity(cfiEntity)
{
m_type = CADType::eOther;
}
~CADElementCFI() = default;
cfi::MeshableEntity *GetCfiPointer()
{
return m_cfiEntity;
}
private:
/// CFI object for surface.
cfi::MeshableEntity *m_cfiEntity;
};
typedef std::shared_ptr<CADElementCFI> CADElementCFISharedPtr;
}
}
#endif
......@@ -56,22 +56,22 @@ Array<OneD, NekDouble> CADSurfCFI::GetBounds()
Array<OneD, NekDouble> b(4);
cfi::UVBox bx = m_cfiSurface->calcUVBox();
b[0] = bx.uLower * m_scal;
b[1] = bx.uUpper * m_scal;
b[2] = bx.vLower * m_scal;
b[3] = bx.vUpper * m_scal;
b[0] = bx.uLower;
b[1] = bx.uUpper;
b[2] = bx.vLower;
b[3] = bx.vUpper;
return b;
}
virtual void CADSurfCFI::GetBounds(NekDouble &umin, NekDouble &umax,
NekDouble &vmin, NekDouble &vmax)
void CADSurfCFI::GetBounds(NekDouble &umin, NekDouble &umax,
NekDouble &vmin, NekDouble &vmax)
{
cfi::UVBox bx = m_cfiSurface->calcUVBox();
umin = b[0];
umax = b[1];
vmin = b[2];
vmax = b[3];
umin = bx.uLower;
umax = bx.uUpper;
vmin = bx.vLower;
vmax = bx.vUpper;
}
Array<OneD, NekDouble> CADSurfCFI::locuv(Array<OneD, NekDouble> p,
......@@ -118,7 +118,7 @@ Array<OneD, NekDouble> CADSurfCFI::P(Array<OneD, NekDouble> uv)
return out;
}
void P(Array<OneD, NekDouble> uv, NekDouble &x, NekDouble &y, NekDouble &z)
void CADSurfCFI::P(Array<OneD, NekDouble> uv, NekDouble &x, NekDouble &y, NekDouble &z)
{
cfi::UVPosition uvp(uv[0], uv[1]);
cfi::Position p = m_cfiSurface->calcXYZAtUV(uvp);
......
......@@ -65,7 +65,7 @@ public:
Array<OneD, NekDouble> GetBounds();
virtual void GetBounds(NekDouble &umin, NekDouble &umax, NekDouble &vmin,
void GetBounds(NekDouble &umin, NekDouble &umax, NekDouble &vmin,
NekDouble &vmax);
Array<OneD, NekDouble> N(Array<OneD, NekDouble> uv);
......@@ -94,10 +94,16 @@ public:
return false;
}
cfi::Face *GetCfiPointer()
{
return m_cfiSurface;
}
private:
/// Function which tests the the value of uv used is within the surface
void Test(Array<OneD, NekDouble> uv)
{
boost::ignore_unused(uv);
}
/// CFI object for surface.
cfi::Face *m_cfiSurface;
......
......@@ -37,6 +37,8 @@
#include "CADSurfCFI.h"
#include "CADVertCFI.h"
#include <boost/lexical_cast.hpp>
using namespace std;
namespace Nektar
......@@ -51,27 +53,27 @@ bool CADSystemCFI::LoadCAD()
{
// it is possible to get CFI to lock on to a open gui session
// not sure it ever will with this code
cfiHandel.startServer();
m_cfiHandle.startServer();
if (m_verbose)
{
cout << "cfi loaded in mode: ";
if (cfiHandel.info.mode == cfi::MODE_STANDALONE)
if (m_cfiHandle.info.mode == cfi::MODE_STANDALONE)
{
cout << "standalone" << endl;
}
else if (cfiHandel.info.mode == cfi::MODE_CLIENT)
else if (m_cfiHandle.info.mode == cfi::MODE_CLIENT)
{
cout << "client" << endl;
}
else if (cfiHandel.info.mode == cfi::MODE_SERVER)
else if (m_cfiHandle.info.mode == cfi::MODE_SERVER)
{
cout << "server" << endl;
}
else if (cfiHandel.info.mode == cfi::MODE_BOTH)
else if (m_cfiHandle.info.mode == cfi::MODE_BOTH)
{
cout << "both" << endl;
}
else if (cfiHandel.info.mode == cfi::MODE_PLUGIN)
else if (m_cfiHandle.info.mode == cfi::MODE_PLUGIN)
{
cout << "plugin" << endl;
}
......@@ -80,16 +82,21 @@ bool CADSystemCFI::LoadCAD()
cout << "unknown" << endl;
}
cout << "\tVersion " << cfiHandel.info.version << endl
<< "\tfixno " << cfiHandel.info.fixno << endl
<< "\tubid " << cfiHandel.info.ubid << endl;
cout << "\tVersion " << m_cfiHandle.info.version << endl
<< "\tfixno " << m_cfiHandle.info.fixno << endl
<< "\tubid " << m_cfiHandle.info.ubid << endl;
}
if (m_config.count("UseCFIMesh"))
{
m_useCFIMesh = boost::lexical_cast<bool>(m_config["UseCFIMesh"]);
}
model = cfiHandel.openModelFile(m_name.c_str());
m_model = m_cfiHandle.openModelFile(m_name.c_str());
if (model->getEntityTotal(cfi::TYPE_BODY, cfi::SUBTYPE_ALL) != 1)
if (m_model->getEntityTotal(cfi::TYPE_BODY, cfi::SUBTYPE_ALL) != 1)
{
if (m_cfiMesh)
if (m_useCFIMesh)
{
if (m_verbose)
{
......@@ -101,21 +108,19 @@ bool CADSystemCFI::LoadCAD()
if (m_verbose)
{
cout << "\tHas multibodies and instructions to mesh, this is "
"not "
"possible"
<< endl;
"not possible" << endl;
}
abort();
}
}
vector<cfi::Entity *> *bds =
model->getEntityList(cfi::TYPE_BODY, cfi::SUBTYPE_ALL);
m_model->getEntityList(cfi::TYPE_BODY, cfi::SUBTYPE_ALL);
for (auto &i : *bds)
{
cfi::Body *b = static_cast<cfi::Body *>(i);
bodies.push_back(b);
m_bodies.push_back(b);
}
// cfi doesnt mind stupid units so this scales everything back to meters
......@@ -123,7 +128,7 @@ bool CADSystemCFI::LoadCAD()
// the m_scal object is passed to all cad entities and scales any operation
// before running it.
m_scal = 1.0;
if (model->getUnits() == cfi::UNIT_INCHES)
if (m_model->getUnits() == cfi::UNIT_INCHES)
{
if (m_verbose)
{
......@@ -131,8 +136,8 @@ bool CADSystemCFI::LoadCAD()
}
m_scal = 0.0254;
}
else if (model->getUnits() == cfi::UNIT_MILLIMETERS ||
model->getUnits() == cfi::UNIT_MILLIMETRES</