Nektar++
Functions
Python/ExpList.cpp File Reference
#include <LibUtilities/Python/BasicUtils/SharedArray.hpp>
#include <LibUtilities/Python/NekPyConfig.hpp>
#include <LocalRegions/MatrixKey.h>
#include <MultiRegions/ExpList.h>
#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 StdRegions::ConstFactorMap &constFactorMap, const StdRegions::VarCoeffMap &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)
 
void ExpList_SetCoeffsArray (ExpListSharedPtr exp, Array< OneD, NekDouble > inarray)
 
const Array< OneD, const NekDoubleExpList_GetCoeffs (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 (py::module &m)
 

Function Documentation

◆ ExpList_BwdTrans()

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

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

82{
83 Array<OneD, NekDouble> out(exp->GetNpoints());
84 exp->BwdTrans(in, out);
85 return out;
86}

Referenced by export_ExpList().

◆ ExpList_FwdTrans()

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

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

74{
75 Array<OneD, NekDouble> out(exp->GetNcoeffs());
76 exp->FwdTrans(in, out);
77 return out;
78}

Referenced by export_ExpList().

◆ ExpList_GetCoeffs()

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

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

193{
194 return exp->GetCoeffs();
195}

Referenced by export_ExpList().

◆ ExpList_GetCoords()

py::tuple ExpList_GetCoords ( ExpListSharedPtr  exp)

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

141{
142 size_t nPhys = exp->GetNpoints();
143 size_t coordim = exp->GetCoordim(0);
144
145 std::vector<Array<OneD, NekDouble>> coords(coordim);
146 for (size_t i = 0; i < coordim; ++i)
147 {
148 coords[i] = Array<OneD, NekDouble>(nPhys);
149 }
150
151 switch (coordim)
152 {
153 case 1:
154 exp->GetCoords(coords[0]);
155 return py::make_tuple(coords[0]);
156 break;
157 case 2:
158 exp->GetCoords(coords[0], coords[1]);
159 return py::make_tuple(coords[0], coords[1]);
160 break;
161 case 3:
162 exp->GetCoords(coords[0], coords[1], coords[2]);
163 return py::make_tuple(coords[0], coords[1], coords[2]);
164 break;
165 }
166
167 return py::tuple();
168}

Referenced by export_ExpList().

◆ ExpList_GetExp()

ExpansionSharedPtr ExpList_GetExp ( ExpListSharedPtr  exp,
int  i 
)

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

55{
56 return exp->GetExp(i);
57}

Referenced by export_ExpList().

◆ ExpList_GetPhys()

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

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

182{
183 return exp->GetPhys();
184}

Referenced by export_ExpList().

◆ ExpList_GetPhysAddress()

std::string ExpList_GetPhysAddress ( ExpListSharedPtr  exp)

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

203{
204 std::stringstream ss;
205 ss << static_cast<const void *>(&(exp->GetPhys()[0]));
206 return ss.str();
207}

Referenced by export_ExpList().

◆ ExpList_HelmSolve()

Array< OneD, NekDouble > ExpList_HelmSolve ( ExpListSharedPtr  exp,
const Array< OneD, const NekDouble > &  in,
const StdRegions::ConstFactorMap constFactorMap,
const StdRegions::VarCoeffMap varCoeffMap 
)

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

108{
109 Array<OneD, NekDouble> out(exp->GetNcoeffs(), 0.0);
110 exp->HelmSolve(in, out, constFactorMap, varCoeffMap);
111 return out;
112}

Referenced by export_ExpList().

◆ ExpList_Integral()

NekDouble ExpList_Integral ( ExpListSharedPtr  exp)

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

198{
199 return exp->Integral();
200}

Referenced by export_ExpList().

◆ ExpList_IProductWRTBase()

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

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

90{
91 Array<OneD, NekDouble> out(exp->GetNcoeffs());
92 exp->IProductWRTBase(in, out);
93 return out;
94}

Referenced by export_ExpList().

◆ ExpList_L2()

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

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

116{
117 return exp->L2(in);
118}

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 120 of file Python/ExpList.cpp.

123{
124 return exp->L2(in, err);
125}

Referenced by export_ExpList().

◆ ExpList_Linf()

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

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

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

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 133 of file Python/ExpList.cpp.

136{
137 return exp->Linf(in, err);
138}

Referenced by export_ExpList().

◆ ExpList_LoadField()

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

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

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

References 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 96 of file Python/ExpList.cpp.

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

Referenced by export_ExpList().

◆ ExpList_ResetManagers()

void ExpList_ResetManagers ( ExpListSharedPtr  exp)

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

210{
211 exp->ClearGlobalLinSysManager();
213 LocalRegions::MatrixKey::opLess>::ClearManager();
215 LocalRegions::MatrixKey::opLess>::ClearManager();
216}
NekMatrix< NekMatrix< NekMatrix< NekDouble, StandardMatrixTag >, ScaledMatrixTag >, BlockMatrixTag > DNekScalBlkMat
Definition: NekTypeDefs.hpp:68
NekMatrix< NekMatrix< NekDouble, StandardMatrixTag >, ScaledMatrixTag > DNekScalMat
Used to lookup the create function in NekManager.
Definition: MatrixKey.h:69

Referenced by export_ExpList().

◆ ExpList_SetCoeffsArray()

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

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

188{
189 exp->SetCoeffsArray(inarray);
190}

Referenced by export_ExpList().

◆ ExpList_SetPhys()

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

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

177{
178 exp->SetPhys(inarray);
179}

Referenced by export_ExpList().

◆ ExpList_SetPhysArray()

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

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

171{
172 exp->SetPhysArray(inarray);
173}

Referenced by export_ExpList().

◆ ExpList_WriteVTK()

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

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

60{
61 std::ofstream out(filename.c_str());
62 exp->WriteVtkHeader(out);
63 size_t nExp = exp->GetExpSize();
64 for (size_t i = 0; i < nExp; ++i)
65 {
66 exp->WriteVtkPieceHeader(out, i);
67 exp->WriteVtkPieceFooter(out, i);
68 }
69 exp->WriteVtkFooter(out);
70}

Referenced by export_ExpList().

◆ export_ExpList()

void export_ExpList ( py::module &  m)

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

268{
269 int (ExpList::*GetNcoeffs)() const = &ExpList::GetNcoeffs;
270
271 py::class_<ExpList, std::shared_ptr<ExpList>>(m, "ExpList")
272 .def(py::init<const LibUtilities::SessionReaderSharedPtr &,
274
275 // Query points and offset information
276 .def("GetExp", &ExpList_GetExp)
277 .def("GetExpSize", &ExpList::GetExpSize)
278 .def("GetNpoints", &ExpList::GetNpoints)
279 .def("GetNcoeffs", GetNcoeffs)
280 .def("GetCoords", &ExpList_GetCoords)
281 .def("GetPhys_Offset", &ExpList::GetPhys_Offset)
282 .def("GetCoeff_Offset", &ExpList::GetCoeff_Offset)
283
284 // Evaluations
285 .def("PhysEvaluate", &ExpList::PhysEvaluate)
286 .def("LoadField", &ExpList_LoadField)
287
288 // Operators
289 .def("FwdTrans", &ExpList_FwdTrans)
290 .def("BwdTrans", &ExpList_BwdTrans)
291 .def("IProductWRTBase", &ExpList_IProductWRTBase)
292 .def("MultiplyByInvMassMatrix", &ExpList_MultiplyByInvMassMatrix)
293 .def("HelmSolve", &ExpList_HelmSolve, py::arg("in"),
294 py::arg("constFactorMap") = StdRegions::NullConstFactorMap,
295 py::arg("varCoeffMap") = StdRegions::NullVarCoeffMap)
296
297 // Error norms
298 .def("L2", &ExpList_L2)
299 .def("L2", &ExpList_L2_Error)
300 .def("Linf", &ExpList_Linf)
301 .def("Linf", &ExpList_Linf_Error)
302
303 // Storage setups
304 .def("SetPhysArray", &ExpList_SetPhysArray)
305 .def("SetPhys", &ExpList_SetPhys)
306 .def("GetPhys", &ExpList_GetPhys)
307 .def("SetCoeffsArray", &ExpList_SetCoeffsArray)
308 .def("GetCoeffs", &ExpList_GetCoeffs)
309 .def("SetPhysState", &ExpList::SetPhysState)
310 .def("GetPhysState", &ExpList::GetPhysState)
311 .def("Integral", &ExpList_Integral)
312 .def("GetPhysAddress", &ExpList_GetPhysAddress)
313
314 .def_property("phys", &ExpList_GetPhys, &ExpList_SetPhys)
315 .def_property("coeffs", &ExpList_GetCoeffs, &ExpList_SetCoeffsArray)
316
317 // Misc functions
318 .def("WriteVTK", &ExpList_WriteVTK)
319 .def("ResetManagers", &ExpList_ResetManagers);
320}
const Array< OneD, const NekDouble > ExpList_GetCoeffs(ExpListSharedPtr exp)
NekDouble ExpList_Linf_Error(ExpListSharedPtr exp, const Array< OneD, const NekDouble > &in, const Array< OneD, const NekDouble > &err)
Array< OneD, NekDouble > ExpList_MultiplyByInvMassMatrix(ExpListSharedPtr exp, const Array< OneD, const NekDouble > &in)
void ExpList_ResetManagers(ExpListSharedPtr exp)
NekDouble ExpList_Linf(ExpListSharedPtr exp, const Array< OneD, const NekDouble > &in)
const Array< OneD, const NekDouble > ExpList_GetPhys(ExpListSharedPtr exp)
NekDouble ExpList_Integral(ExpListSharedPtr exp)
Array< OneD, NekDouble > ExpList_IProductWRTBase(ExpListSharedPtr exp, const Array< OneD, const NekDouble > &in)
ExpansionSharedPtr ExpList_GetExp(ExpListSharedPtr exp, int i)
void ExpList_SetPhysArray(ExpListSharedPtr exp, Array< OneD, NekDouble > inarray)
Array< OneD, NekDouble > ExpList_HelmSolve(ExpListSharedPtr exp, const Array< OneD, const NekDouble > &in, const StdRegions::ConstFactorMap &constFactorMap, const StdRegions::VarCoeffMap &varCoeffMap)
Array< OneD, NekDouble > ExpList_BwdTrans(ExpListSharedPtr exp, const Array< OneD, const NekDouble > &in)
void ExpList_LoadField(ExpListSharedPtr exp, std::string filename, std::string varName)
void ExpList_WriteVTK(ExpListSharedPtr exp, std::string filename)
void ExpList_SetPhys(ExpListSharedPtr exp, const Array< OneD, const NekDouble > &inarray)
py::tuple ExpList_GetCoords(ExpListSharedPtr exp)
void ExpList_SetCoeffsArray(ExpListSharedPtr exp, Array< OneD, NekDouble > inarray)
NekDouble ExpList_L2(ExpListSharedPtr exp, const Array< OneD, const NekDouble > &in)
std::string ExpList_GetPhysAddress(ExpListSharedPtr exp)
Array< OneD, NekDouble > ExpList_FwdTrans(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)
Base class for all multi-elemental spectral/hp expansions.
Definition: ExpList.h:99
std::shared_ptr< SessionReader > SessionReaderSharedPtr
std::shared_ptr< MeshGraph > MeshGraphSharedPtr
Definition: MeshGraph.h:174
static ConstFactorMap NullConstFactorMap
Definition: StdRegions.hpp:431
static VarCoeffMap NullVarCoeffMap
Definition: StdRegions.hpp:376

References ExpList_BwdTrans(), ExpList_FwdTrans(), ExpList_GetCoeffs(), 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_SetCoeffsArray(), ExpList_SetPhys(), ExpList_SetPhysArray(), ExpList_WriteVTK(), Nektar::StdRegions::NullConstFactorMap, and Nektar::StdRegions::NullVarCoeffMap.

Referenced by PYBIND11_MODULE().