Nektar++
Loading...
Searching...
No Matches
Typedefs | Functions
GeomElements.cpp File Reference
#include <LibUtilities/Python/NekPyConfig.hpp>
#include <SpatialDomains/Python/SpatialDomains.h>
#include <SpatialDomains/HexGeom.h>
#include <SpatialDomains/PointGeom.h>
#include <SpatialDomains/PrismGeom.h>
#include <SpatialDomains/PyrGeom.h>
#include <SpatialDomains/QuadGeom.h>
#include <SpatialDomains/SegGeom.h>
#include <SpatialDomains/TetGeom.h>
#include <SpatialDomains/TriGeom.h>

Go to the source code of this file.

Typedefs

using NekError = Nektar::ErrorUtil::NekError
 

Functions

template<class T , class S >
unique_ptr_objpool< T > Geometry_Init (int id, py::list &facets)
 
template<class T , class S >
unique_ptr_objpool< T > Geometry_Init_Curved (int id, py::list &facets, CurveSharedPtr curve)
 
template<class T , class S >
void export_Geom_2d (py::module &m, const char *name)
 
template<class T , class S >
void export_Geom_3d (py::module &m, const char *name)
 
unique_ptr_objpool< SegGeomSegGeom_Init (int id, int coordim, py::list &points, CurveSharedPtr curve)
 
unique_ptr_objpool< PointGeomPointGeom_Init (int coordim, int id, NekDouble x, NekDouble y, NekDouble z)
 
py::tuple PointGeom_GetCoordinates (const PointGeom &self)
 
void export_GeomElements (py::module &m)
 

Typedef Documentation

◆ NekError

Definition at line 50 of file GeomElements.cpp.

Function Documentation

◆ export_Geom_2d()

template<class T , class S >
void export_Geom_2d ( py::module &  m,
const char *  name 
)

Definition at line 89 of file GeomElements.cpp.

90{
91 py::classh<T, Geometry2D>(m, name)
92 .def(py::init<>(&Geometry_Init<T, S>), py::arg("id"),
93 py::arg("segments") = py::list())
94 .def(py::init<>(&Geometry_Init_Curved<T, S>), py::arg("id"),
95 py::arg("segments"), py::arg("curve"));
96}

◆ export_Geom_3d()

template<class T , class S >
void export_Geom_3d ( py::module &  m,
const char *  name 
)

Definition at line 98 of file GeomElements.cpp.

99{
100 py::classh<T, Geometry3D>(m, name).def(py::init<>(&Geometry_Init<T, S>),
101 py::arg("id"),
102 py::arg("segments") = py::list());
103}

◆ export_GeomElements()

void export_GeomElements ( py::module &  m)

Definition at line 142 of file GeomElements.cpp.

143{
144 // Geometry dimensioned base classes
145 py::classh<Geometry1D, Geometry>(m, "Geometry1D");
146 py::classh<Geometry2D, Geometry>(m, "Geometry2D")
147 .def("GetCurve", &Geometry2D::GetCurve,
148 py::return_value_policy::reference);
149 py::classh<Geometry3D, Geometry>(m, "Geometry3D");
150
151 // Point geometries
152 py::classh<PointGeom, Geometry>(m, "PointGeom")
153 .def(py::init<>(&PointGeom_Init), py::arg("coordim"), py::arg("id"),
154 py::arg("x"), py::arg("y"), py::arg("z"))
155 .def("GetCoordinates", &PointGeom_GetCoordinates);
156
157 // Segment geometries
158 py::classh<SegGeom, Geometry>(m, "SegGeom")
159 .def(py::init<>(&SegGeom_Init), py::arg("id"), py::arg("coordim"),
160 py::arg("points") = py::list(),
161 py::arg("curve") = CurveSharedPtr())
162 .def("GetCurve", &SegGeom::GetCurve);
163
164 export_Geom_2d<TriGeom, SegGeom>(m, "TriGeom");
165 export_Geom_2d<QuadGeom, SegGeom>(m, "QuadGeom");
166 export_Geom_3d<TetGeom, TriGeom>(m, "TetGeom");
167 export_Geom_3d<PrismGeom, Geometry2D>(m, "PrismGeom");
168 export_Geom_3d<PyrGeom, Geometry2D>(m, "PyrGeom");
169 export_Geom_3d<HexGeom, QuadGeom>(m, "HexGeom");
170}
py::tuple PointGeom_GetCoordinates(const PointGeom &self)
unique_ptr_objpool< SegGeom > SegGeom_Init(int id, int coordim, py::list &points, CurveSharedPtr curve)
unique_ptr_objpool< PointGeom > PointGeom_Init(int coordim, int id, NekDouble x, NekDouble y, NekDouble z)
CurveSharedPtr GetCurve()
Definition SegGeom.h:80
std::shared_ptr< Curve > CurveSharedPtr
Definition Curve.hpp:60

References Nektar::SpatialDomains::Geometry2D::GetCurve(), Nektar::SpatialDomains::SegGeom::GetCurve(), PointGeom_GetCoordinates(), PointGeom_Init(), and SegGeom_Init().

Referenced by PYBIND11_MODULE().

◆ Geometry_Init()

template<class T , class S >
unique_ptr_objpool< T > Geometry_Init ( int  id,
py::list &  facets 
)

Definition at line 53 of file GeomElements.cpp.

54{
55 std::array<S *, T::kNfacets> geomArr;
56
57 if (py::len(facets) != T::kNfacets)
58 {
59 throw new NekError("Incorrect number of facets for this geometry");
60 }
61
62 for (int i = 0; i < T::kNfacets; ++i)
63 {
64 geomArr[i] = py::cast<S *>(facets[i]);
65 }
66
67 return ObjPoolManager<T>::AllocateUniquePtr(id, geomArr);
68}
Nektar::ErrorUtil::NekError NekError
static std::unique_ptr< DataType, UniquePtrDeleter > AllocateUniquePtr(const Args &...args)

References Nektar::ObjPoolManager< DataType >::AllocateUniquePtr().

◆ Geometry_Init_Curved()

template<class T , class S >
unique_ptr_objpool< T > Geometry_Init_Curved ( int  id,
py::list &  facets,
CurveSharedPtr  curve 
)

Definition at line 71 of file GeomElements.cpp.

73{
74 std::array<S *, T::kNfacets> geomArr;
75
76 if (py::len(facets) != T::kNfacets)
77 {
78 throw new NekError("Incorrect number of facets for this geometry");
79 }
80
81 for (int i = 0; i < T::kNfacets; ++i)
82 {
83 geomArr[i] = py::cast<S *>(facets[i]);
84 }
85
86 return ObjPoolManager<T>::AllocateUniquePtr(id, geomArr, curve);
87}

References Nektar::ObjPoolManager< DataType >::AllocateUniquePtr().

◆ PointGeom_GetCoordinates()

py::tuple PointGeom_GetCoordinates ( const PointGeom self)

Definition at line 137 of file GeomElements.cpp.

138{
139 return py::make_tuple(self.x(), self.y(), self.z());
140}
boost::call_traits< DataType >::const_reference x() const
Definition NekPoint.hpp:160
boost::call_traits< DataType >::const_reference z() const
Definition NekPoint.hpp:172
boost::call_traits< DataType >::const_reference y() const
Definition NekPoint.hpp:166

References Nektar::NekPoint< data_type >::x(), Nektar::NekPoint< data_type >::y(), and Nektar::NekPoint< data_type >::z().

Referenced by export_GeomElements().

◆ PointGeom_Init()

unique_ptr_objpool< PointGeom > PointGeom_Init ( int  coordim,
int  id,
NekDouble  x,
NekDouble  y,
NekDouble  z 
)

Definition at line 131 of file GeomElements.cpp.

133{
134 return ObjPoolManager<PointGeom>::AllocateUniquePtr(coordim, id, x, y, z);
135}

References Nektar::ObjPoolManager< DataType >::AllocateUniquePtr().

Referenced by export_GeomElements().

◆ SegGeom_Init()

unique_ptr_objpool< SegGeom > SegGeom_Init ( int  id,
int  coordim,
py::list &  points,
CurveSharedPtr  curve 
)

Definition at line 105 of file GeomElements.cpp.

107{
108 std::array<PointGeom *, 2> geomArr;
109
110 if (py::len(points) != 2)
111 {
112 throw new NekError("Incorrect number of facets for this geometry");
113 }
114
115 for (int i = 0; i < 2; ++i)
116 {
117 geomArr[i] = py::cast<PointGeom *>(points[i]);
118 }
119
120 if (!curve)
121 {
122 return ObjPoolManager<SegGeom>::AllocateUniquePtr(id, coordim, geomArr);
123 }
124 else
125 {
126 return ObjPoolManager<SegGeom>::AllocateUniquePtr(id, coordim, geomArr,
127 curve);
128 }
129}

References Nektar::ObjPoolManager< DataType >::AllocateUniquePtr().

Referenced by export_GeomElements().