Nektar++
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
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_Integral (ExpListSharedPtr exp)
 
std::string ExpList_GetPhysAddress (ExpListSharedPtr exp)
 
void ExpList_ResetManagers (ExpListSharedPtr exp)
 
void ExpList_LoadField (ExpListSharedPtr exp, std::string filename, std::string varName)
 
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.

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

Referenced by export_ExpList().

◆ ExpList_FwdTrans()

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

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

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

Referenced by export_ExpList().

◆ ExpList_GetCoords()

py::tuple ExpList_GetCoords ( ExpListSharedPtr  exp)

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

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

Referenced by export_ExpList().

◆ ExpList_GetExp()

ExpansionSharedPtr ExpList_GetExp ( ExpListSharedPtr  exp,
int  i 
)

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

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

Referenced by export_ExpList().

◆ ExpList_GetPhys()

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

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

201 {
202  return exp->GetPhys();
203 }

Referenced by export_ExpList().

◆ ExpList_GetPhysAddress()

std::string ExpList_GetPhysAddress ( ExpListSharedPtr  exp)

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

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

Referenced by export_ExpList().

◆ 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.

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  exp->HelmSolve(in, out, facMap, coeffMap);
123  return out;
124 }
std::map< StdRegions::VarCoeffType, Array< OneD, NekDouble > > VarCoeffMap
Definition: StdRegions.hpp:272
std::map< ConstFactorType, NekDouble > ConstFactorMap
Definition: StdRegions.hpp:314
static ConstFactorMap NullConstFactorMap
Definition: StdRegions.hpp:315
static VarCoeffMap NullVarCoeffMap
Definition: StdRegions.hpp:273

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

Referenced by export_ExpList().

◆ ExpList_Integral()

NekDouble ExpList_Integral ( ExpListSharedPtr  exp)

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

206 {
207  return exp->Integral();
208 }

Referenced by export_ExpList().

◆ ExpList_IProductWRTBase()

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

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

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

Referenced by export_ExpList().

◆ ExpList_L2()

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

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

129 {
130  return exp->L2(in);
131 }

Referenced by export_ExpList().

◆ ExpList_L2_Error()

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

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

137 {
138  return exp->L2(in, err);
139 }

Referenced by export_ExpList().

◆ ExpList_Linf()

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

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

144 {
145  return exp->Linf(in);
146 }

Referenced by export_ExpList().

◆ ExpList_Linf_Error()

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

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

152 {
153  return exp->Linf(in, err);
154 }

Referenced by export_ExpList().

◆ ExpList_LoadField()

void ExpList_LoadField ( ExpListSharedPtr  exp,
std::string  filename,
std::string  varName 
)

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

229 {
230  int nExp = exp->GetExpSize();
231  Array<OneD, int> elementGIDs(nExp);
232 
233  // Construct a map from element locations to global IDs
234  for (int i = 0; i < nExp; ++i)
235  {
236  elementGIDs[i] = exp->GetExp(i)->GetGeom()->GetGlobalID();
237  }
238 
239  std::vector<LibUtilities::FieldDefinitionsSharedPtr> def;
240  std::vector<std::vector<NekDouble> > data;
242  LibUtilities::FieldIO::CreateForFile(exp->GetSession(), filename);
243  fldIO->Import(filename, def, data, LibUtilities::NullFieldMetaDataMap,
244  elementGIDs);
245 
246  int idx = -1;
247 
248  Vmath::Zero(exp->GetNcoeffs(), exp->UpdateCoeffs(), 1);
249 
250  // Loop over all the expansions
251  for (int i = 0; i < def.size(); ++i)
252  {
253  // Find the index of the required field in the expansion segment
254  for (int j = 0; j < def[i]->m_fields.size(); ++j)
255  {
256  if (def[i]->m_fields[j] == varName)
257  {
258  idx = j;
259  }
260  }
261 
262  if (idx >= 0)
263  {
264  exp->ExtractDataToCoeffs(
265  def[i], data[i], def[i]->m_fields[idx], exp->UpdateCoeffs());
266  }
267  else
268  {
269  std::cout << "Field " + varName + " not found." << std::endl;
270  }
271  }
272 
273  exp->BwdTrans_IterPerExp(exp->GetCoeffs(), exp->UpdatePhys());
274 }
std::shared_ptr< FieldIO > FieldIOSharedPtr
Definition: FieldIO.h:306
static FieldMetaDataMap NullFieldMetaDataMap
Definition: FieldIO.h:53
void Zero(int n, T *x, const int incx)
Zero vector.
Definition: Vmath.cpp:436

References Nektar::LibUtilities::FieldIO::CreateForFile(), Nektar::LibUtilities::NullFieldMetaDataMap, and Vmath::Zero().

Referenced by export_ExpList().

◆ ExpList_MultiplyByInvMassMatrix()

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

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

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

Referenced by export_ExpList().

◆ ExpList_ResetManagers()

void ExpList_ResetManagers ( ExpListSharedPtr  exp)

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

218 {
219  exp->ClearGlobalLinSysManager();
222  LocalRegions::MatrixKey::opLess>::ClearManager();
225  LocalRegions::MatrixKey::opLess>::ClearManager();
226 }
NekMatrix< NekMatrix< NekDouble, StandardMatrixTag >, ScaledMatrixTag > DNekScalMat
NekMatrix< NekMatrix< NekMatrix< NekDouble, StandardMatrixTag >, ScaledMatrixTag >, BlockMatrixTag > DNekScalBlkMat
Definition: NekTypeDefs.hpp:65
Used to lookup the create function in NekManager.
Definition: MatrixKey.h:68

Referenced by export_ExpList().

◆ ExpList_SetPhys()

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

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

196 {
197  exp->SetPhys(inarray);
198 }

Referenced by export_ExpList().

◆ ExpList_SetPhysArray()

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

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

189 {
190  exp->SetPhysArray(inarray);
191 }

Referenced by export_ExpList().

◆ ExpList_WriteVTK()

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

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

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 }

Referenced by export_ExpList().

◆ export_ExpList()

void export_ExpList ( )

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

277 {
278  int (ExpList::*GetNcoeffs)() const = &ExpList::GetNcoeffs;
279 
280  py::class_<ExpList,
281  std::shared_ptr<ExpList>,
282  boost::noncopyable>(
283  "ExpList", py::init<
286 
287  // Query points and offset information
288  .def("GetExp", &ExpList_GetExp)
289  .def("GetExpSize", &ExpList::GetExpSize)
290  .def("GetNpoints", &ExpList::GetNpoints)
291  .def("GetNcoeffs", GetNcoeffs)
292  .def("GetCoords", &ExpList_GetCoords)
293  .def("GetPhys_Offset", &ExpList::GetPhys_Offset)
294  .def("GetCoeff_Offset", &ExpList::GetCoeff_Offset)
295 
296  // Evaluations
297  .def("PhysEvaluate", &ExpList::PhysEvaluate)
298  .def("LoadField", &ExpList_LoadField)
299 
300  // Operators
301  .def("FwdTrans", &ExpList_FwdTrans)
302  .def("BwdTrans", &ExpList_BwdTrans)
303  .def("IProductWRTBase", &ExpList_IProductWRTBase)
304  .def("MultiplyByInvMassMatrix", &ExpList_MultiplyByInvMassMatrix)
305  .def("HelmSolve", &ExpList_HelmSolve, (
306  py::arg("in"),
307  py::arg("constFactorMap") = py::object(),
308  py::arg("varCoeffMap") = py::object()
309  ))
310 
311  // Error norms
312  .def("L2", &ExpList_L2)
313  .def("L2", &ExpList_L2_Error)
314  .def("Linf", &ExpList_Linf)
315  .def("Linf", &ExpList_Linf_Error)
316 
317  // Storage setups
318  .def("SetPhysArray", &ExpList_SetPhysArray)
319  .def("SetPhys", &ExpList_SetPhys)
320  .def("GetPhys", &ExpList_GetPhys)
321  .def("SetPhysState", &ExpList::SetPhysState)
322  .def("GetPhysState", &ExpList::GetPhysState)
323  .def("Integral", &ExpList_Integral)
324  .def("GetPhysAddress", &ExpList_GetPhysAddress)
325 
326  // Misc functions
327  .def("WriteVTK", &ExpList_WriteVTK)
328  .def("ResetManagers", &ExpList_ResetManagers)
329  ;
330 }
NekDouble ExpList_Linf_Error(ExpListSharedPtr exp, const Array< OneD, const NekDouble > &in, const Array< OneD, const NekDouble > &err)
void ExpList_ResetManagers(ExpListSharedPtr exp)
NekDouble ExpList_Linf(ExpListSharedPtr exp, const Array< OneD, const NekDouble > &in)
Array< OneD, NekDouble > ExpList_IProductWRTBase(ExpListSharedPtr exp, const Array< OneD, const NekDouble > &in)
NekDouble ExpList_Integral(ExpListSharedPtr exp)
ExpansionSharedPtr ExpList_GetExp(ExpListSharedPtr exp, int i)
void ExpList_SetPhysArray(ExpListSharedPtr exp, Array< OneD, NekDouble > inarray)
const Array< OneD, const NekDouble > ExpList_GetPhys(ExpListSharedPtr exp)
void ExpList_LoadField(ExpListSharedPtr exp, std::string filename, std::string varName)
Array< OneD, NekDouble > ExpList_MultiplyByInvMassMatrix(ExpListSharedPtr exp, const Array< OneD, const NekDouble > &in)
void ExpList_WriteVTK(ExpListSharedPtr exp, std::string filename)
Array< OneD, NekDouble > ExpList_FwdTrans(ExpListSharedPtr exp, const Array< OneD, const NekDouble > &in)
void ExpList_SetPhys(ExpListSharedPtr exp, const Array< OneD, const NekDouble > &inarray)
py::tuple ExpList_GetCoords(ExpListSharedPtr exp)
NekDouble ExpList_L2(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)
Array< OneD, NekDouble > ExpList_BwdTrans(ExpListSharedPtr exp, const Array< OneD, const NekDouble > &in)
std::string ExpList_GetPhysAddress(ExpListSharedPtr exp)
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:107
std::shared_ptr< SessionReader > SessionReaderSharedPtr
std::shared_ptr< MeshGraph > MeshGraphSharedPtr
Definition: MeshGraph.h:174

References ExpList_BwdTrans(), ExpList_FwdTrans(), ExpList_GetCoords(), ExpList_GetExp(), ExpList_GetPhys(), ExpList_GetPhysAddress(), ExpList_HelmSolve(), ExpList_Integral(), ExpList_IProductWRTBase(), ExpList_L2(), ExpList_L2_Error(), ExpList_Linf(), ExpList_Linf_Error(), ExpList_LoadField(), ExpList_MultiplyByInvMassMatrix(), 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().