Nektar++
Functions
Python/ExpList.cpp File Reference
#include <MultiRegions/ExpList.h>
#include <LocalRegions/MatrixKey.h>
#include <LibUtilities/Python/NekPyConfig.hpp>
#include <fstream>
#include <sstream>
#include <string>

Go to the source code of this file.

Functions

ExpansionSharedPtr ExpList_GetExp (ExpListSharedPtr exp, int i)
 
void ExpList_WriteVTK (ExpListSharedPtr exp, std::string filename)
 
Array< OneD, NekDoubleExpList_FwdTrans (ExpListSharedPtr exp, const Array< OneD, const NekDouble > &in)
 
Array< OneD, NekDoubleExpList_BwdTrans (ExpListSharedPtr exp, const Array< OneD, const NekDouble > &in)
 
Array< OneD, NekDoubleExpList_IProductWRTBase (ExpListSharedPtr exp, const Array< OneD, const NekDouble > &in)
 
Array< OneD, NekDoubleExpList_MultiplyByInvMassMatrix (ExpListSharedPtr exp, const Array< OneD, const NekDouble > &in)
 
Array< OneD, NekDoubleExpList_HelmSolve (ExpListSharedPtr exp, const Array< OneD, const NekDouble > &in, const py::object constFactorMap, const py::object varCoeffMap)
 
NekDouble ExpList_L2 (ExpListSharedPtr exp, const Array< OneD, const NekDouble > &in)
 
NekDouble ExpList_L2_Error (ExpListSharedPtr exp, const Array< OneD, const NekDouble > &in, const Array< OneD, const NekDouble > &err)
 
NekDouble ExpList_Linf (ExpListSharedPtr exp, const Array< OneD, const NekDouble > &in)
 
NekDouble ExpList_Linf_Error (ExpListSharedPtr exp, const Array< OneD, const NekDouble > &in, const Array< OneD, const NekDouble > &err)
 
py::tuple ExpList_GetCoords (ExpListSharedPtr exp)
 
void ExpList_SetPhysArray (ExpListSharedPtr exp, Array< OneD, NekDouble > inarray)
 
void ExpList_SetPhys (ExpListSharedPtr exp, const Array< OneD, const NekDouble > &inarray)
 
const Array< OneD, const NekDoubleExpList_GetPhys (ExpListSharedPtr exp)
 
NekDouble ExpList_PhysIntegral (ExpListSharedPtr exp)
 
std::string ExpList_GetPhysAddress (ExpListSharedPtr exp)
 
void ExpList_ResetManagers (ExpListSharedPtr exp)
 
void export_ExpList ()
 

Function Documentation

◆ ExpList_BwdTrans()

Array<OneD, NekDouble> ExpList_BwdTrans ( ExpListSharedPtr  exp,
const Array< OneD, const NekDouble > &  in 
)

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

Referenced by export_ExpList().

78 {
79  Array<OneD, NekDouble> out(exp->GetNpoints());
80  exp->BwdTrans(in, out);
81  return out;
82 }

◆ ExpList_FwdTrans()

Array<OneD, NekDouble> ExpList_FwdTrans ( ExpListSharedPtr  exp,
const Array< OneD, const NekDouble > &  in 
)

Definition at line 66 of file Python/ExpList.cpp.

Referenced by export_ExpList().

69 {
70  Array<OneD, NekDouble> out(exp->GetNcoeffs());
71  exp->FwdTrans(in, out);
72  return out;
73 }

◆ ExpList_GetCoords()

py::tuple ExpList_GetCoords ( ExpListSharedPtr  exp)

Definition at line 158 of file Python/ExpList.cpp.

Referenced by export_ExpList().

159 {
160  int nPhys = exp->GetNpoints();
161  int coordim = exp->GetCoordim(0);
162 
163  std::vector<Array<OneD, NekDouble> > coords(coordim);
164  for (int i = 0; i < coordim; ++i)
165  {
166  coords[i] = Array<OneD, NekDouble>(nPhys);
167  }
168 
169  switch (coordim)
170  {
171  case 1:
172  exp->GetCoords(coords[0]);
173  return py::make_tuple(coords[0]);
174  break;
175  case 2:
176  exp->GetCoords(coords[0], coords[1]);
177  return py::make_tuple(coords[0], coords[1]);
178  break;
179  case 3:
180  exp->GetCoords(coords[0], coords[1], coords[2]);
181  return py::make_tuple(coords[0], coords[1], coords[2]);
182  break;
183  }
184 
185  return py::tuple();
186 }

◆ ExpList_GetExp()

ExpansionSharedPtr ExpList_GetExp ( ExpListSharedPtr  exp,
int  i 
)

Definition at line 49 of file Python/ExpList.cpp.

Referenced by export_ExpList().

50 {
51  return exp->GetExp(i);
52 }

◆ ExpList_GetPhys()

const Array<OneD, const NekDouble> ExpList_GetPhys ( ExpListSharedPtr  exp)

Definition at line 202 of file Python/ExpList.cpp.

Referenced by export_ExpList().

203 {
204  return exp->GetPhys();
205 }

◆ ExpList_GetPhysAddress()

std::string ExpList_GetPhysAddress ( ExpListSharedPtr  exp)

Definition at line 212 of file Python/ExpList.cpp.

Referenced by export_ExpList().

213 {
214  std::stringstream ss;
215  ss << static_cast<const void*>(&(exp->GetPhys()[0]));
216  return ss.str();
217 }

◆ ExpList_HelmSolve()

Array<OneD, NekDouble> ExpList_HelmSolve ( ExpListSharedPtr  exp,
const Array< OneD, const NekDouble > &  in,
const py::object  constFactorMap,
const py::object  varCoeffMap 
)

Definition at line 102 of file Python/ExpList.cpp.

References Nektar::StdRegions::NullConstFactorMap, and Nektar::StdRegions::NullVarCoeffMap.

Referenced by export_ExpList().

107 {
108  Array<OneD, NekDouble> out(exp->GetNcoeffs(), 0.0);
109 
112 
113  if (!constFactorMap.is_none())
114  {
115  facMap = py::extract<StdRegions::ConstFactorMap>(constFactorMap);
116  }
117  if (!varCoeffMap.is_none())
118  {
119  coeffMap = py::extract<StdRegions::VarCoeffMap>(varCoeffMap);
120  }
121 
122  FlagList notUsed;
123 
124  exp->HelmSolve(in, out, notUsed, facMap, coeffMap);
125  return out;
126 }
std::map< ConstFactorType, NekDouble > ConstFactorMap
Definition: StdRegions.hpp:294
std::map< StdRegions::VarCoeffType, Array< OneD, NekDouble > > VarCoeffMap
Definition: StdRegions.hpp:264
Defines a list of flags.
static VarCoeffMap NullVarCoeffMap
Definition: StdRegions.hpp:265
static ConstFactorMap NullConstFactorMap
Definition: StdRegions.hpp:295

◆ ExpList_IProductWRTBase()

Array<OneD, NekDouble> ExpList_IProductWRTBase ( ExpListSharedPtr  exp,
const Array< OneD, const NekDouble > &  in 
)

Definition at line 84 of file Python/ExpList.cpp.

Referenced by export_ExpList().

87 {
88  Array<OneD, NekDouble> out(exp->GetNcoeffs());
89  exp->IProductWRTBase(in, out);
90  return out;
91 }

◆ ExpList_L2()

NekDouble ExpList_L2 ( ExpListSharedPtr  exp,
const Array< OneD, const NekDouble > &  in 
)

Definition at line 128 of file Python/ExpList.cpp.

Referenced by export_ExpList().

131 {
132  return exp->L2(in);
133 }

◆ ExpList_L2_Error()

NekDouble ExpList_L2_Error ( ExpListSharedPtr  exp,
const Array< OneD, const NekDouble > &  in,
const Array< OneD, const NekDouble > &  err 
)

Definition at line 135 of file Python/ExpList.cpp.

Referenced by export_ExpList().

139 {
140  return exp->L2(in, err);
141 }

◆ ExpList_Linf()

NekDouble ExpList_Linf ( ExpListSharedPtr  exp,
const Array< OneD, const NekDouble > &  in 
)

Definition at line 143 of file Python/ExpList.cpp.

Referenced by export_ExpList().

146 {
147  return exp->Linf(in);
148 }

◆ ExpList_Linf_Error()

NekDouble ExpList_Linf_Error ( ExpListSharedPtr  exp,
const Array< OneD, const NekDouble > &  in,
const Array< OneD, const NekDouble > &  err 
)

Definition at line 150 of file Python/ExpList.cpp.

Referenced by export_ExpList().

154 {
155  return exp->Linf(in, err);
156 }

◆ ExpList_MultiplyByInvMassMatrix()

Array<OneD, NekDouble> ExpList_MultiplyByInvMassMatrix ( ExpListSharedPtr  exp,
const Array< OneD, const NekDouble > &  in 
)

Definition at line 93 of file Python/ExpList.cpp.

References Nektar::MultiRegions::eLocal.

Referenced by export_ExpList().

96 {
97  Array<OneD, NekDouble> out(exp->GetNcoeffs(), 0.0);
98  exp->MultiplyByInvMassMatrix(in, out, eLocal);
99  return out;
100 }
Local coefficients.

◆ ExpList_PhysIntegral()

NekDouble ExpList_PhysIntegral ( ExpListSharedPtr  exp)

Definition at line 207 of file Python/ExpList.cpp.

Referenced by export_ExpList().

208 {
209  return exp->PhysIntegral();
210 }

◆ ExpList_ResetManagers()

void ExpList_ResetManagers ( ExpListSharedPtr  exp)

Definition at line 219 of file Python/ExpList.cpp.

Referenced by export_ExpList().

220 {
221  exp->ClearGlobalLinSysManager();
224  LocalRegions::MatrixKey::opLess>::ClearManager();
226  LocalRegions::MatrixKey, DNekScalBlkMat,
227  LocalRegions::MatrixKey::opLess>::ClearManager();
228 }
NekMatrix< NekMatrix< NekMatrix< NekDouble, StandardMatrixTag >, ScaledMatrixTag >, BlockMatrixTag > DNekScalBlkMat
Definition: NekTypeDefs.hpp:65
Used to lookup the create function in NekManager.
Definition: MatrixKey.h:67
NekMatrix< NekMatrix< NekDouble, StandardMatrixTag >, ScaledMatrixTag > DNekScalMat

◆ ExpList_SetPhys()

void ExpList_SetPhys ( ExpListSharedPtr  exp,
const Array< OneD, const NekDouble > &  inarray 
)

Definition at line 195 of file Python/ExpList.cpp.

Referenced by export_ExpList().

198 {
199  exp->SetPhys(inarray);
200 }

◆ ExpList_SetPhysArray()

void ExpList_SetPhysArray ( ExpListSharedPtr  exp,
Array< OneD, NekDouble inarray 
)

Definition at line 188 of file Python/ExpList.cpp.

Referenced by export_ExpList().

191 {
192  exp->SetPhysArray(inarray);
193 }

◆ ExpList_WriteVTK()

void ExpList_WriteVTK ( ExpListSharedPtr  exp,
std::string  filename 
)

Definition at line 54 of file Python/ExpList.cpp.

Referenced by export_ExpList().

55 {
56  std::ofstream out(filename.c_str());
57  exp->WriteVtkHeader(out);
58  for (int i = 0; i < exp->GetExpSize(); ++i)
59  {
60  exp->WriteVtkPieceHeader(out, i);
61  exp->WriteVtkPieceFooter(out, i);
62  }
63  exp->WriteVtkFooter(out);
64 }

◆ export_ExpList()

void export_ExpList ( )

Definition at line 230 of file Python/ExpList.cpp.

References ExpList_BwdTrans(), ExpList_FwdTrans(), ExpList_GetCoords(), ExpList_GetExp(), ExpList_GetPhys(), ExpList_GetPhysAddress(), ExpList_HelmSolve(), ExpList_IProductWRTBase(), ExpList_L2(), ExpList_L2_Error(), ExpList_Linf(), ExpList_Linf_Error(), ExpList_MultiplyByInvMassMatrix(), ExpList_PhysIntegral(), ExpList_ResetManagers(), ExpList_SetPhys(), ExpList_SetPhysArray(), ExpList_WriteVTK(), Nektar::MultiRegions::ExpList::GetCoeff_Offset(), Nektar::MultiRegions::ExpList::GetExpSize(), Nektar::MultiRegions::ExpList::GetNcoeffs(), Nektar::MultiRegions::ExpList::GetNpoints(), Nektar::MultiRegions::ExpList::GetPhys_Offset(), Nektar::MultiRegions::ExpList::GetPhysState(), Nektar::MultiRegions::ExpList::PhysEvaluate(), and Nektar::MultiRegions::ExpList::SetPhysState().

Referenced by BOOST_PYTHON_MODULE().

231 {
232  int (ExpList::*GetNcoeffs)() const = &ExpList::GetNcoeffs;
233 
234  py::class_<ExpList,
235  std::shared_ptr<ExpList>,
236  boost::noncopyable>(
237  "ExpList", py::no_init)
238 
239  // Expansions
240  .def("GetExp", ExpList_GetExp)
241  .def("GetExpSize", &ExpList::GetExpSize)
242 
243  // Query points and offset information
244  .def("GetNpoints", &ExpList::GetNpoints)
245  .def("GetNcoeffs", GetNcoeffs)
246  .def("GetCoords", &ExpList_GetCoords)
247  .def("GetPhys_Offset", &ExpList::GetPhys_Offset)
248  .def("GetCoeff_Offset", &ExpList::GetCoeff_Offset)
249 
250  // Evaluations
251  .def("PhysEvaluate", &ExpList::PhysEvaluate)
252 
253  // Operators
254  .def("FwdTrans", &ExpList_FwdTrans)
255  .def("BwdTrans", &ExpList_BwdTrans)
256  .def("IProductWRTBase", &ExpList_IProductWRTBase)
257  .def("MultiplyByInvMassMatrix", &ExpList_MultiplyByInvMassMatrix)
258  .def("HelmSolve", &ExpList_HelmSolve, (
259  py::arg("in"),
260  py::arg("constFactorMap") = py::object(),
261  py::arg("varCoeffMap") = py::object()
262  ))
263 
264  // Error norms
265  .def("L2", &ExpList_L2)
266  .def("L2", &ExpList_L2_Error)
267  .def("Linf", &ExpList_Linf)
268  .def("Linf", &ExpList_Linf_Error)
269 
270  // Storage setups
271  .def("SetPhysArray", &ExpList_SetPhysArray)
272  .def("SetPhys", &ExpList_SetPhys)
273  .def("GetPhys", &ExpList_GetPhys)
274  .def("SetPhysState", &ExpList::SetPhysState)
275  .def("GetPhysState", &ExpList::GetPhysState)
276  .def("PhysIntegral", &ExpList_PhysIntegral)
277  .def("GetPhysAddress", &ExpList_GetPhysAddress)
278 
279  // Misc functions
280  .def("WriteVTK", &ExpList_WriteVTK)
281  .def("ResetManagers", &ExpList_ResetManagers)
282  ;
283 }
Array< OneD, NekDouble > ExpList_MultiplyByInvMassMatrix(ExpListSharedPtr exp, const Array< OneD, const NekDouble > &in)
Array< OneD, NekDouble > ExpList_FwdTrans(ExpListSharedPtr exp, const Array< OneD, const NekDouble > &in)
Array< OneD, NekDouble > ExpList_HelmSolve(ExpListSharedPtr exp, const Array< OneD, const NekDouble > &in, const py::object constFactorMap, const py::object varCoeffMap)
NekDouble ExpList_L2(ExpListSharedPtr exp, const Array< OneD, const NekDouble > &in)
NekDouble ExpList_Linf_Error(ExpListSharedPtr exp, const Array< OneD, const NekDouble > &in, const Array< OneD, const NekDouble > &err)
ExpansionSharedPtr ExpList_GetExp(ExpListSharedPtr exp, int i)
void ExpList_ResetManagers(ExpListSharedPtr exp)
std::string ExpList_GetPhysAddress(ExpListSharedPtr exp)
void ExpList_SetPhys(ExpListSharedPtr exp, const Array< OneD, const NekDouble > &inarray)
NekDouble ExpList_L2_Error(ExpListSharedPtr exp, const Array< OneD, const NekDouble > &in, const Array< OneD, const NekDouble > &err)
Base class for all multi-elemental spectral/hp expansions.
Definition: ExpList.h:103
NekDouble ExpList_PhysIntegral(ExpListSharedPtr exp)
const Array< OneD, const NekDouble > ExpList_GetPhys(ExpListSharedPtr exp)
py::tuple ExpList_GetCoords(ExpListSharedPtr exp)
Array< OneD, NekDouble > ExpList_IProductWRTBase(ExpListSharedPtr exp, const Array< OneD, const NekDouble > &in)
Array< OneD, NekDouble > ExpList_BwdTrans(ExpListSharedPtr exp, const Array< OneD, const NekDouble > &in)
void ExpList_WriteVTK(ExpListSharedPtr exp, std::string filename)
void ExpList_SetPhysArray(ExpListSharedPtr exp, Array< OneD, NekDouble > inarray)
NekDouble ExpList_Linf(ExpListSharedPtr exp, const Array< OneD, const NekDouble > &in)