Nektar++
Loading...
Searching...
No Matches
Functions
StdExpansion.cpp File Reference
#include <StdRegions/StdExpansion.h>
#include <LibUtilities/Python/BasicUtils/SharedArray.hpp>
#include <LibUtilities/Python/LinearAlgebra/NekMatrix.hpp>
#include <LibUtilities/Python/NekPyConfig.hpp>

Go to the source code of this file.

Functions

Array< OneD, NekDoubleStdExpansion_FwdTrans (StdExpansionSharedPtr exp, const Array< OneD, const NekDouble > &in)
 
Array< OneD, NekDoubleStdExpansion_BwdTrans (StdExpansionSharedPtr exp, const Array< OneD, const NekDouble > &in)
 
Array< OneD, NekDoubleStdExpansion_IProductWRTBase (StdExpansionSharedPtr exp, const Array< OneD, const NekDouble > &in)
 
NekDouble StdExpansion_PhysEvaluate (StdExpansionSharedPtr exp, const Array< OneD, const NekDouble > &coords, const Array< OneD, const NekDouble > &physvals)
 
NekDouble StdExpansion_L2 (StdExpansionSharedPtr exp, const Array< OneD, const NekDouble > &in)
 
NekDouble StdExpansion_L2_Error (StdExpansionSharedPtr exp, const Array< OneD, const NekDouble > &in, const Array< OneD, const NekDouble > &err)
 
py::tuple StdExpansion_GetCoords (StdExpansionSharedPtr exp)
 
py::tuple StdExpansion_PhysDeriv (StdExpansionSharedPtr exp, const Array< OneD, const NekDouble > &inarray)
 
py::array StdExpansion_GenMatrix (std::shared_ptr< StdExpansion > &stdExp, const StdMatrixKey &mkey)
 
py::array StdExpansion_GetStdMatrix (std::shared_ptr< StdExpansion > &stdExp, const StdMatrixKey &mkey)
 
void export_StdExpansion (py::module &m)
 

Function Documentation

◆ export_StdExpansion()

void export_StdExpansion ( py::module &  m)

Definition at line 170 of file Python/StdExpansion.cpp.

171{
172 py::class_<StdExpansion, std::shared_ptr<StdExpansion>>(m, "StdExpansion")
173
174 .def("GetNcoeffs", &StdExpansion::GetNcoeffs)
175 .def("GetTotPoints", &StdExpansion::GetTotPoints)
176 .def("GetBasisType", &StdExpansion::GetBasisType)
177 .def("GetPointsType", &StdExpansion::GetPointsType)
178 .def("GetNverts", &StdExpansion::GetNverts)
179 .def("GetNtraces", &StdExpansion::GetNtraces)
180 .def("DetShapeType", &StdExpansion::DetShapeType)
181 .def("GetShapeDimension", &StdExpansion::GetShapeDimension)
182 .def("Integral", &StdExpansion::Integral)
183
184 .def("GetBasis", &StdExpansion::GetBasis)
185
186 .def("GenMatrix", &StdExpansion_GenMatrix)
187 .def("GetStdMatrix", &StdExpansion_GetStdMatrix)
188
189 .def("FwdTrans", &StdExpansion_FwdTrans)
190 .def("BwdTrans", &StdExpansion_BwdTrans)
191 .def("IProductWRTBase", &StdExpansion_IProductWRTBase)
192
193 .def("PhysEvaluate", &StdExpansion_PhysEvaluate)
194 .def("L2", &StdExpansion_L2)
195 .def("L2", &StdExpansion_L2_Error)
196
197 .def("GetCoords", &StdExpansion_GetCoords)
198
199 .def("PhysDeriv", &StdExpansion_PhysDeriv);
200}
NekDouble StdExpansion_L2_Error(StdExpansionSharedPtr exp, const Array< OneD, const NekDouble > &in, const Array< OneD, const NekDouble > &err)
py::tuple StdExpansion_PhysDeriv(StdExpansionSharedPtr exp, const Array< OneD, const NekDouble > &inarray)
py::array StdExpansion_GetStdMatrix(std::shared_ptr< StdExpansion > &stdExp, const StdMatrixKey &mkey)
Array< OneD, NekDouble > StdExpansion_FwdTrans(StdExpansionSharedPtr exp, const Array< OneD, const NekDouble > &in)
NekDouble StdExpansion_L2(StdExpansionSharedPtr exp, const Array< OneD, const NekDouble > &in)
py::tuple StdExpansion_GetCoords(StdExpansionSharedPtr exp)
Array< OneD, NekDouble > StdExpansion_BwdTrans(StdExpansionSharedPtr exp, const Array< OneD, const NekDouble > &in)
NekDouble StdExpansion_PhysEvaluate(StdExpansionSharedPtr exp, const Array< OneD, const NekDouble > &coords, const Array< OneD, const NekDouble > &physvals)
py::array StdExpansion_GenMatrix(std::shared_ptr< StdExpansion > &stdExp, const StdMatrixKey &mkey)
Array< OneD, NekDouble > StdExpansion_IProductWRTBase(StdExpansionSharedPtr exp, const Array< OneD, const NekDouble > &in)
const LibUtilities::BasisSharedPtr & GetBasis(int dir) const
This function gets the shared point to basis in the dir direction.
int GetNcoeffs(void) const
This function returns the total number of coefficients used in the expansion.
int GetTotPoints() const
This function returns the total number of quadrature points used in the element.
LibUtilities::BasisType GetBasisType(const int dir) const
This function returns the type of basis used in the dir direction.
int GetNtraces() const
Returns the number of trace elements connected to this element.
int GetNverts() const
This function returns the number of vertices of the expansion domain.
LibUtilities::ShapeType DetShapeType() const
This function returns the shape of the expansion domain.
LibUtilities::PointsType GetPointsType(const int dir) const
This function returns the type of quadrature points used in the dir direction.
NekDouble Integral(const Array< OneD, const NekDouble > &inarray)
This function integrates the specified function over the domain.

References Nektar::StdRegions::StdExpansion::DetShapeType(), Nektar::StdRegions::StdExpansion::GetBasis(), Nektar::StdRegions::StdExpansion::GetBasisType(), Nektar::StdRegions::StdExpansion::GetNcoeffs(), Nektar::StdRegions::StdExpansion::GetNtraces(), Nektar::StdRegions::StdExpansion::GetNverts(), Nektar::StdRegions::StdExpansion::GetPointsType(), Nektar::StdRegions::StdExpansion::GetShapeDimension(), Nektar::StdRegions::StdExpansion::GetTotPoints(), Nektar::StdRegions::StdExpansion::Integral(), StdExpansion_BwdTrans(), StdExpansion_FwdTrans(), StdExpansion_GenMatrix(), StdExpansion_GetCoords(), StdExpansion_GetStdMatrix(), StdExpansion_IProductWRTBase(), StdExpansion_L2(), StdExpansion_L2_Error(), StdExpansion_PhysDeriv(), and StdExpansion_PhysEvaluate().

Referenced by PYBIND11_MODULE().

◆ StdExpansion_BwdTrans()

Array< OneD, NekDouble > StdExpansion_BwdTrans ( StdExpansionSharedPtr  exp,
const Array< OneD, const NekDouble > &  in 
)

Definition at line 52 of file Python/StdExpansion.cpp.

54{
55 Array<OneD, NekDouble> out(exp->GetTotPoints());
56 exp->BwdTrans(in, out);
57 return out;
58}

Referenced by export_StdExpansion().

◆ StdExpansion_FwdTrans()

Array< OneD, NekDouble > StdExpansion_FwdTrans ( StdExpansionSharedPtr  exp,
const Array< OneD, const NekDouble > &  in 
)

Definition at line 44 of file Python/StdExpansion.cpp.

46{
47 Array<OneD, NekDouble> out(exp->GetNcoeffs());
48 exp->FwdTrans(in, out);
49 return out;
50}

Referenced by export_StdExpansion().

◆ StdExpansion_GenMatrix()

py::array StdExpansion_GenMatrix ( std::shared_ptr< StdExpansion > &  stdExp,
const StdMatrixKey mkey 
)

Definition at line 149 of file Python/StdExpansion.cpp.

151{
152 auto mat = stdExp->GenMatrix(mkey);
153 // A bit wasteful but at the moment we can't easily keep hold of the
154 // NekMatrix returned as a std::shared_ptr, because pybind11 doesn't like
155 // std::shared_ptr around types that have custom casters (see issue #787).
156 return py::array({mat->GetRows(), mat->GetColumns()},
157 {sizeof(NekDouble), mat->GetRows() * sizeof(NekDouble)},
158 mat->GetRawPtr());
159}

Referenced by export_StdExpansion().

◆ StdExpansion_GetCoords()

py::tuple StdExpansion_GetCoords ( StdExpansionSharedPtr  exp)

Definition at line 88 of file Python/StdExpansion.cpp.

89{
90 int nPhys = exp->GetTotPoints();
91 int coordim = exp->GetCoordim();
92
93 std::vector<Array<OneD, NekDouble>> coords(coordim);
94 for (int i = 0; i < coordim; ++i)
95 {
96 coords[i] = Array<OneD, NekDouble>(nPhys);
97 }
98
99 switch (coordim)
100 {
101 case 1:
102 exp->GetCoords(coords[0]);
103 return py::make_tuple(coords[0]);
104 break;
105 case 2:
106 exp->GetCoords(coords[0], coords[1]);
107 return py::make_tuple(coords[0], coords[1]);
108 break;
109 case 3:
110 exp->GetCoords(coords[0], coords[1], coords[2]);
111 return py::make_tuple(coords[0], coords[1], coords[2]);
112 break;
113 }
114
115 return py::tuple();
116}

Referenced by export_StdExpansion().

◆ StdExpansion_GetStdMatrix()

py::array StdExpansion_GetStdMatrix ( std::shared_ptr< StdExpansion > &  stdExp,
const StdMatrixKey mkey 
)

Definition at line 161 of file Python/StdExpansion.cpp.

163{
164 auto mat = stdExp->GetStdMatrix(mkey);
165 return py::array({mat->GetRows(), mat->GetColumns()},
166 {sizeof(NekDouble), mat->GetRows() * sizeof(NekDouble)},
167 mat->GetRawPtr());
168}

Referenced by export_StdExpansion().

◆ StdExpansion_IProductWRTBase()

Array< OneD, NekDouble > StdExpansion_IProductWRTBase ( StdExpansionSharedPtr  exp,
const Array< OneD, const NekDouble > &  in 
)

Definition at line 60 of file Python/StdExpansion.cpp.

62{
63 Array<OneD, NekDouble> out(exp->GetNcoeffs());
64 exp->IProductWRTBase(in, out);
65 return out;
66}

Referenced by export_StdExpansion().

◆ StdExpansion_L2()

NekDouble StdExpansion_L2 ( StdExpansionSharedPtr  exp,
const Array< OneD, const NekDouble > &  in 
)

Definition at line 75 of file Python/StdExpansion.cpp.

77{
78 return exp->L2(in);
79}

Referenced by export_StdExpansion().

◆ StdExpansion_L2_Error()

NekDouble StdExpansion_L2_Error ( StdExpansionSharedPtr  exp,
const Array< OneD, const NekDouble > &  in,
const Array< OneD, const NekDouble > &  err 
)

Definition at line 81 of file Python/StdExpansion.cpp.

84{
85 return exp->L2(in, err);
86}

Referenced by export_StdExpansion().

◆ StdExpansion_PhysDeriv()

py::tuple StdExpansion_PhysDeriv ( StdExpansionSharedPtr  exp,
const Array< OneD, const NekDouble > &  inarray 
)

Definition at line 118 of file Python/StdExpansion.cpp.

120{
121 int nPhys = exp->GetTotPoints();
122 int coordim = exp->GetCoordim();
123
124 std::vector<Array<OneD, NekDouble>> derivs(coordim);
125 for (int i = 0; i < coordim; ++i)
126 {
127 derivs[i] = Array<OneD, NekDouble>(nPhys);
128 }
129
130 switch (coordim)
131 {
132 case 1:
133 exp->PhysDeriv(inarray, derivs[0]);
134 return py::make_tuple(derivs[0]);
135 break;
136 case 2:
137 exp->PhysDeriv(inarray, derivs[0], derivs[1]);
138 return py::make_tuple(derivs[0], derivs[1]);
139 break;
140 case 3:
141 exp->PhysDeriv(inarray, derivs[0], derivs[1], derivs[2]);
142 return py::make_tuple(derivs[0], derivs[1], derivs[2]);
143 break;
144 }
145
146 return py::tuple();
147}

Referenced by export_StdExpansion().

◆ StdExpansion_PhysEvaluate()

NekDouble StdExpansion_PhysEvaluate ( StdExpansionSharedPtr  exp,
const Array< OneD, const NekDouble > &  coords,
const Array< OneD, const NekDouble > &  physvals 
)

Definition at line 68 of file Python/StdExpansion.cpp.

71{
72 return exp->PhysEvaluate(coords, physvals);
73}

Referenced by export_StdExpansion().