Nektar++
ProcessInterpPoints.cpp
Go to the documentation of this file.
1////////////////////////////////////////////////////////////////////////////////
2//
3// File: ProcessInterpPoints.cpp
4//
5// For more information, please see: http://www.nektar.info/
6//
7// The MIT License
8//
9// Copyright (c) 2006 Division of Applied Mathematics, Brown University (USA),
10// Department of Aeronautics, Imperial College London (UK), and Scientific
11// Computing and Imaging Institute, University of Utah (USA).
12//
13// Permission is hereby granted, free of charge, to any person obtaining a
14// copy of this software and associated documentation files (the "Software"),
15// to deal in the Software without restriction, including without limitation
16// the rights to use, copy, modify, merge, publish, distribute, sublicense,
17// and/or sell copies of the Software, and to permit persons to whom the
18// Software is furnished to do so, subject to the following conditions:
19//
20// The above copyright notice and this permission notice shall be included
21// in all copies or substantial portions of the Software.
22//
23// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
24// OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
25// FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
26// THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
27// LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
28// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
29// DEALINGS IN THE SOFTWARE.
30//
31// Description: Interpolate field to a series of specified points.
32//
33////////////////////////////////////////////////////////////////////////////////
34
35#include <iostream>
36#include <string>
37
38#include <boost/geometry.hpp>
39
47
48#include "ProcessInterpPoints.h"
49
50using namespace std;
51
52namespace bg = boost::geometry;
53namespace bgi = boost::geometry::index;
54
55namespace Nektar::FieldUtils
56{
57
61 "Interpolates a field to a set of points. Requires fromfld, fromxml "
62 "to be defined, and a topts, line, plane or block of target points ");
63
65{
66 m_config["fromxml"] = ConfigOption(
67 false, "NotSet", "Xml file from which to interpolate field");
68
69 m_config["fromfld"] = ConfigOption(
70 false, "NotSet", "Fld file from which to interpolate field");
71
72 m_config["topts"] =
73 ConfigOption(false, "NotSet", "Pts file to which interpolate field");
74 m_config["line"] = ConfigOption(false, "NotSet",
75 "Specify a line of N points using "
76 "line=N,x0,y0,z0,z1,y1,z1");
77 m_config["plane"] =
78 ConfigOption(false, "NotSet",
79 "Specify a plane of N1 x N2 points using "
80 "plane=N1,N2,x0,y0,z0,z1,y1,z1,x2,y2,z2,x3,y3,z3");
81 m_config["box"] =
82 ConfigOption(false, "NotSet",
83 "Specify a rectangular box of N1 x N2 x N3 points "
84 "using a box of points limited by box="
85 "N1,N2,N3,xmin,xmax,ymin,ymax,zmin,zmax");
86
87 m_config["clamptolowervalue"] =
88 ConfigOption(false, "-10000000", "Lower bound for interpolation value");
89 m_config["clamptouppervalue"] =
90 ConfigOption(false, "10000000", "Upper bound for interpolation value");
91 m_config["defaultvalue"] =
92 ConfigOption(false, "0", "Default value if point is outside domain");
93
94 m_config["cp"] =
95 ConfigOption(false, "NotSet",
96 "Parameters p0 and q to determine pressure coefficients");
97 m_config["realmodetoimag"] =
98 ConfigOption(false, "NotSet", "Take fields as sin mode");
99 m_config["distTolerance"] =
100 ConfigOption(false, "1e-5", "The maximum acceptable distance.");
101}
102
104{
105}
106
107void ProcessInterpPoints::v_Process(po::variables_map &vm)
108{
109 m_f->SetUpExp(vm);
110
111 CreateFieldPts(vm);
112
113 FieldSharedPtr fromField = std::shared_ptr<Field>(new Field());
114 std::vector<std::string> files;
115 ParseUtils::GenerateVector(m_config["fromxml"].as<string>(), files);
116
117 // set up session file for from field
118 char *argv[] = {const_cast<char *>("FieldConvert"), nullptr};
119 fromField->m_session = LibUtilities::SessionReader::CreateInstance(
120 1, argv, files,
121 LibUtilities::GetCommFactory().CreateInstance("Serial", 0, 0));
122
123 // Set up range based on min and max of local parallel partition
126
127 int coordim = m_f->m_fieldPts->GetDim();
128 int npts = m_f->m_fieldPts->GetNpoints();
129 std::vector<std::string> fieldNames = m_f->m_fieldPts->GetFieldNames();
130 for (auto &it : fieldNames)
131 {
132 m_f->m_fieldPts->RemoveField(it);
133 }
134
136 m_f->m_fieldPts->GetPts(pts);
137
138 rng->m_checkShape = false;
139 rng->m_zmin = -1;
140 rng->m_zmax = 1;
141 rng->m_ymin = -1;
142 rng->m_ymax = 1;
143 switch (coordim)
144 {
145 case 3:
146 rng->m_doZrange = true;
147 rng->m_zmin = Vmath::Vmin(npts, pts[2], 1);
148 rng->m_zmax = Vmath::Vmax(npts, pts[2], 1);
149 if (rng->m_zmax == rng->m_zmin)
150 {
151 rng->m_zmin -= 1;
152 rng->m_zmax += 1;
153 }
154 /* Falls through. */
155 case 2:
156 rng->m_doYrange = true;
157 rng->m_ymin = Vmath::Vmin(npts, pts[1], 1);
158 rng->m_ymax = Vmath::Vmax(npts, pts[1], 1);
159 /* Falls through. */
160 case 1:
161 rng->m_doXrange = true;
162 rng->m_xmin = Vmath::Vmin(npts, pts[0], 1);
163 rng->m_xmax = Vmath::Vmax(npts, pts[0], 1);
164 break;
165 default:
166 NEKERROR(ErrorUtil::efatal, "Too many values specified in range");
167 }
168
169 // setup rng parameters.
170 fromField->m_graph =
171 SpatialDomains::MeshGraphIO::Read(fromField->m_session, rng);
172
173 // Read in local from field partitions
174 const SpatialDomains::ExpansionInfoMap &expansions =
175 fromField->m_graph->GetExpansionInfo();
176 Array<OneD, int> ElementGIDs(expansions.size());
177
178 int i = 0;
179 for (auto &expIt : expansions)
180 {
181 ElementGIDs[i++] = expIt.second->m_geomShPtr->GetGlobalID();
182 }
183 // check to see that we do have some element in the domain since
184 // possibly all points could be outside of the domain
185 ASSERTL0(i > 0, "No elements are set. Are the interpolated points "
186 "within the domain given by the xml files?");
187 string fromfld = m_config["fromfld"].as<string>();
188 m_f->FieldIOForFile(fromfld)->Import(
189 fromfld, fromField->m_fielddef, fromField->m_data,
191 int NumHomogeneousDir = fromField->m_fielddef[0]->m_numHomogeneousDir;
192 for (i = 0; i < fromField->m_fielddef.size(); ++i)
193 {
194 int d1 = fromField->m_fielddef[i]->m_basis.size();
195 d1 -= 1;
196 if (d1 >= 0 && (fromField->m_fielddef[i]->m_basis[d1] ==
198 fromField->m_fielddef[i]->m_basis[d1] ==
200 {
201 fromField->m_fielddef[i]->m_homogeneousZIDs[0] += 2;
202 fromField->m_fielddef[i]->m_numModes[d1] = 4;
203 fromField->m_fielddef[i]->m_basis[d1] = LibUtilities::eFourier;
204 }
205 }
206
207 //----------------------------------------------
208 // Set up Expansion information to use mode order from field
209 fromField->m_graph->SetExpansionInfo(fromField->m_fielddef);
210 int nfields = fromField->m_fielddef[0]->m_fields.size();
211 fromField->m_exp.resize(nfields);
212 fromField->m_exp[0] = fromField->SetUpFirstExpList(NumHomogeneousDir, true);
213 m_f->m_exp.resize(nfields);
214
215 // declare auxiliary fields.
216 for (i = 1; i < nfields; ++i)
217 {
218 fromField->m_exp[i] = fromField->AppendExpList(NumHomogeneousDir);
219 }
220
221 // load field into expansion in fromfield.
222 set<int> sinmode;
223 if (m_config["realmodetoimag"].as<string>().compare("NotSet"))
224 {
225 ParseUtils::GenerateVariableSet(m_config["realmodetoimag"].as<string>(),
226 m_f->m_variables, sinmode);
227 }
228 for (int j = 0; j < nfields; ++j)
229 {
230 for (i = 0; i < fromField->m_fielddef.size(); i++)
231 {
232 fromField->m_exp[j]->ExtractDataToCoeffs(
233 fromField->m_fielddef[i], fromField->m_data[i],
234 fromField->m_fielddef[0]->m_fields[j],
235 fromField->m_exp[j]->UpdateCoeffs());
236 }
237 if (NumHomogeneousDir == 1)
238 {
239 fromField->m_exp[j]->SetWaveSpace(true);
240 if (sinmode.count(j))
241 {
242 int Ncoeff = fromField->m_exp[j]->GetPlane(2)->GetNcoeffs();
244 Ncoeff, -1., fromField->m_exp[j]->GetPlane(2)->GetCoeffs(),
245 1, fromField->m_exp[j]->GetPlane(3)->UpdateCoeffs(), 1);
246 Vmath::Zero(Ncoeff,
247 fromField->m_exp[j]->GetPlane(2)->UpdateCoeffs(),
248 1);
249 }
250 }
251 fromField->m_exp[j]->BwdTrans(fromField->m_exp[j]->GetCoeffs(),
252 fromField->m_exp[j]->UpdatePhys());
253
254 Array<OneD, NekDouble> newPts(m_f->m_fieldPts->GetNpoints());
255 m_f->m_fieldPts->AddField(newPts,
256 fromField->m_fielddef[0]->m_fields[j]);
257 m_f->m_variables.push_back(fromField->m_fielddef[0]->m_fields[j]);
258 }
259
260 NekDouble clamp_low = m_config["clamptolowervalue"].as<NekDouble>();
261 NekDouble clamp_up = m_config["clamptouppervalue"].as<NekDouble>();
262 NekDouble def_value = m_config["defaultvalue"].as<NekDouble>();
263 NekDouble tolerance = m_config["distTolerance"].as<NekDouble>();
264
265 InterpolateFieldToPts(fromField->m_exp, m_f->m_fieldPts, clamp_low,
266 clamp_up, def_value, tolerance);
267
268 if (!boost::iequals(m_config["cp"].as<string>(), "NotSet"))
269 {
270 calcCp0();
271 }
272}
273
274void ProcessInterpPoints::CreateFieldPts([[maybe_unused]] po::variables_map &vm)
275{
276 int rank = m_f->m_comm->GetSpaceComm()->GetRank();
277 int nprocs = m_f->m_comm->GetSpaceComm()->GetSize();
278 // Check for command line point specification
279
280 if (m_config["topts"].as<string>().compare("NotSet") != 0)
281 {
282 string inFile = m_config["topts"].as<string>();
283
284 if (fs::path(inFile).extension() == ".pts")
285 {
288 m_f->m_comm);
289
290 ptsIO->Import(inFile, m_f->m_fieldPts);
291 }
292 else if (fs::path(inFile).extension() == ".csv")
293 {
296 m_f->m_comm);
297
298 csvIO->Import(inFile, m_f->m_fieldPts);
299 }
300 else
301 {
302 ASSERTL0(false, "unknown topts file type");
303 }
304 }
305 else if (m_config["line"].as<string>().compare("NotSet") != 0)
306 {
307 vector<NekDouble> values;
308 ASSERTL0(
309 ParseUtils::GenerateVector(m_config["line"].as<string>(), values),
310 "Failed to interpret line string");
311
312 ASSERTL0(values.size() > 2, "line string should contain 2*Dim+1 values "
313 "N,x0,y0,z0,x1,y1,z1");
314
315 double tmp;
316 ASSERTL0(std::modf(values[0], &tmp) == 0.0, "N is not an integer");
317 ASSERTL0(values[0] > 1, "N is not a valid number");
318
319 int dim = (values.size() - 1) / 2;
320 int npts = values[0];
321
322 // Information for partitioning
323 int ptsPerProc = npts / nprocs;
324 int extraPts = (rank < nprocs - 1) ? 0 : npts % nprocs;
325 int locPts = ptsPerProc + extraPts;
326 int start = rank * ptsPerProc;
327 int end = start + locPts;
328
330 Array<OneD, NekDouble> delta(dim);
331 for (int i = 0; i < dim; ++i)
332 {
333 pts[i] = Array<OneD, NekDouble>(locPts);
334 delta[i] = (values[dim + i + 1] - values[i + 1]) / (npts - 1);
335 }
336
337 for (int i = 0, cntLoc = 0; i < npts; ++i)
338 {
339 if (i >= start && i < end)
340 {
341 for (int n = 0; n < dim; ++n)
342 {
343 pts[n][cntLoc] = values[n + 1] + i * delta[n];
344 }
345 ++cntLoc;
346 }
347 }
348
349 vector<size_t> ppe;
350 ppe.push_back(npts);
351 m_f->m_fieldPts =
353 m_f->m_fieldPts->SetPtsType(LibUtilities::ePtsLine);
354 m_f->m_fieldPts->SetPointsPerEdge(ppe);
355 }
356 else if (m_config["plane"].as<string>().compare("NotSet") != 0)
357 {
358 vector<NekDouble> values;
359 ASSERTL0(
360 ParseUtils::GenerateVector(m_config["plane"].as<string>(), values),
361 "Failed to interpret plane string");
362
363 ASSERTL0(values.size() > 9,
364 "plane string should contain 4 Dim+2 values "
365 "N1,N2,x0,y0,z0,x1,y1,z1,x2,y2,z2,x3,y3,z3");
366
367 double tmp;
368 ASSERTL0(std::modf(values[0], &tmp) == 0.0, "N1 is not an integer");
369 ASSERTL0(std::modf(values[1], &tmp) == 0.0, "N2 is not an integer");
370
371 ASSERTL0(values[0] > 1, "N1 is not a valid number");
372 ASSERTL0(values[1] > 1, "N2 is not a valid number");
373
374 int dim = (values.size() - 2) / 4;
375
376 Array<OneD, int> npts(2);
377 npts[0] = values[0];
378 npts[1] = values[1];
379
380 int totpts = npts[0] * npts[1];
381
382 // Information for partitioning
383 int ptsPerProc = totpts / nprocs;
384 int extraPts = (rank < nprocs - 1) ? 0 : totpts % nprocs;
385 int locPts = ptsPerProc + extraPts;
386 int start = rank * ptsPerProc;
387 int end = start + locPts;
388
390 Array<OneD, NekDouble> delta1(dim);
391 Array<OneD, NekDouble> delta2(dim);
392 for (int i = 0; i < dim; ++i)
393 {
394 pts[i] = Array<OneD, NekDouble>(locPts);
395 delta1[i] = (values[2 + 1 * dim + i] - values[2 + 0 * dim + i]) /
396 (npts[0] - 1);
397 delta2[i] = (values[2 + 2 * dim + i] - values[2 + 3 * dim + i]) /
398 (npts[0] - 1);
399 }
400
401 for (int j = 0, cnt = 0, cntLoc = 0; j < npts[1]; ++j)
402 {
403 for (int i = 0; i < npts[0]; ++i, ++cnt)
404 {
405 if (cnt >= start && cnt < end)
406 {
407 for (int n = 0; n < dim; ++n)
408 {
409 pts[n][cntLoc] =
410 (values[2 + n] + i * delta1[n]) *
411 (1.0 - j / ((NekDouble)(npts[1] - 1))) +
412 (values[2 + 3 * dim + n] + i * delta2[n]) *
413 (j / ((NekDouble)(npts[1] - 1)));
414 }
415 ++cntLoc;
416 }
417 }
418 }
419
420 vector<size_t> ppe;
421 ppe.push_back(npts[0]);
422 ppe.push_back(npts[1]);
423 m_f->m_fieldPts =
425 m_f->m_fieldPts->SetPtsType(LibUtilities::ePtsPlane);
426 m_f->m_fieldPts->SetPointsPerEdge(ppe);
427 }
428 else if (m_config["box"].as<string>().compare("NotSet") != 0)
429 {
430 vector<NekDouble> values;
431 ASSERTL0(
432 ParseUtils::GenerateVector(m_config["box"].as<string>(), values),
433 "Failed to interpret box string");
434
435 ASSERTL0(values.size() == 9, "box string should contain 9 values "
436 "N1,N2,N3,xmin,xmax,ymin,ymax,zmin,zmax");
437
438 int dim = 3;
439 Array<OneD, int> npts(3);
440 npts[0] = values[0];
441 npts[1] = values[1];
442 npts[2] = values[2];
443
444 int totpts = npts[0] * npts[1] * npts[2];
445
447 Array<OneD, NekDouble> delta(dim);
448
449 // Information for partitioning
450 int ptsPerProc = totpts / nprocs;
451 int extraPts = (rank < nprocs - 1) ? 0 : totpts % nprocs;
452 int locPts = ptsPerProc + extraPts;
453 int start = rank * ptsPerProc;
454 int end = start + locPts;
455
456 for (int i = 0; i < dim; ++i)
457 {
458 pts[i] = Array<OneD, NekDouble>(locPts);
459 delta[i] = (values[4 + 2 * i] - values[3 + 2 * i]) / (npts[i] - 1);
460 }
461
462 for (int k = 0, cnt = 0, cntLoc = 0; k < npts[2]; ++k)
463 {
464 for (int j = 0; j < npts[1]; ++j)
465 {
466 for (int i = 0; i < npts[0]; ++i, ++cnt)
467 {
468 if (cnt >= start && cnt < end)
469 {
470 pts[0][cntLoc] = values[3] + i * delta[0];
471 pts[1][cntLoc] = values[5] + j * delta[1];
472 pts[2][cntLoc] = values[7] + k * delta[2];
473 ++cntLoc;
474 }
475 }
476 }
477 }
478
479 vector<size_t> ppe;
480 ppe.push_back(npts[0]);
481 ppe.push_back(npts[1]);
482 ppe.push_back(npts[2]);
483 m_f->m_fieldPts =
485 m_f->m_fieldPts->SetPtsType(LibUtilities::ePtsBox);
486 m_f->m_fieldPts->SetPointsPerEdge(ppe);
487 vector<NekDouble> boxdim;
488 boxdim.assign(&values[3], &values[3] + 6);
489 m_f->m_fieldPts->SetBoxSize(boxdim);
490 }
491 else
492 {
493 ASSERTL0(false, "Missing target points for ProcessInterpPoints.");
494 }
495}
496
498 vector<MultiRegions::ExpListSharedPtr> &field0,
500 NekDouble clamp_up, NekDouble def_value, NekDouble tolerance)
501{
502 ASSERTL0(pts->GetNFields() == field0.size(), "ptField has too few fields");
503
504 int nfields = field0.size();
505
507 if (m_f->m_comm->GetRank() == 0)
508 {
510 this);
511 }
512 interp.Interpolate(field0, pts, def_value, tolerance);
513
514 if (m_f->m_comm->GetRank() == 0)
515 {
516 cout << endl;
517 }
518
519 for (int f = 0; f < nfields; ++f)
520 {
521 for (int i = 0; i < pts->GetNpoints(); ++i)
522 {
523 if (pts->GetPointVal(f, i) > clamp_up)
524 {
525 pts->SetPointVal(f, i, clamp_up);
526 }
527 else if (pts->GetPointVal(f, i) < clamp_low)
528 {
529 pts->SetPointVal(f, i, clamp_low);
530 }
531 }
532 }
533}
534
536{
537 LibUtilities::PtsFieldSharedPtr pts = m_f->m_fieldPts;
538 int dim = pts->GetDim();
539 int nq1 = pts->GetNpoints();
540 int r, f;
541 int pfield = -1;
542 NekDouble p0, qinv;
543 vector<int> velid;
544
545 vector<NekDouble> values;
546 ASSERTL0(ParseUtils::GenerateVector(m_config["cp"].as<string>(), values),
547 "Failed to interpret cp string");
548
549 ASSERTL0(values.size() == 2, "cp string should contain 2 values "
550 "p0 and q (=1/2 rho u^2)");
551
552 p0 = values[0];
553 qinv = 1.0 / values[1];
554
555 for (int i = 0; i < pts->GetNFields(); ++i)
556 {
557 if (boost::iequals(pts->GetFieldName(i), "p"))
558 {
559 pfield = i;
560 }
561
562 if (boost::iequals(pts->GetFieldName(i), "u") ||
563 boost::iequals(pts->GetFieldName(i), "v") ||
564 boost::iequals(pts->GetFieldName(i), "w"))
565 {
566 velid.push_back(i);
567 }
568 }
569
570 if (pfield != -1)
571 {
572 if (!velid.size())
573 {
574 WARNINGL0(false, "Did not find velocity components for Cp0");
575 }
576 }
577 else
578 {
579 WARNINGL0(false, "Failed to find 'p' field to determine cp0");
580 }
581
582 // Allocate data storage
584
585 for (f = 0; f < 2; ++f)
586 {
587 data[f] = Array<OneD, NekDouble>(nq1, 0.0);
588 }
589
590 for (r = 0; r < nq1; r++)
591 {
592 if (pfield != -1) // calculate cp
593 {
594 data[0][r] = qinv * (pts->GetPointVal(dim + pfield, r) - p0);
595
596 if (velid.size()) // calculate cp0
597 {
598 NekDouble q = 0;
599 for (int i = 0; i < velid.size(); ++i)
600 {
601 q += 0.5 * pts->GetPointVal(dim + velid[i], r) *
602 pts->GetPointVal(dim + velid[i], r);
603 }
604 data[1][r] =
605 qinv * (pts->GetPointVal(dim + pfield, r) + q - p0);
606 }
607 }
608 }
609
610 if (pfield != -1)
611 {
612 pts->AddField(data[0], "Cp");
613 m_f->m_variables.push_back("Cp");
614 if (velid.size())
615 {
616 pts->AddField(data[1], "Cp0");
617 m_f->m_variables.push_back("Cp0");
618 }
619 }
620}
621
623 const int goal) const
624{
625 LibUtilities::PrintProgressbar(position, goal, "Interpolating");
626}
627} // namespace Nektar::FieldUtils
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:208
#define NEKERROR(type, msg)
Assert Level 0 – Fundamental assert which is used whether in FULLDEBUG, DEBUG or OPT compilation mode...
Definition: ErrorUtil.hpp:202
#define WARNINGL0(condition, msg)
Definition: ErrorUtil.hpp:215
A class that contains algorithms for interpolation between pts fields, expansions and different meshe...
FIELD_UTILS_EXPORT void Interpolate(const T expInField, T &expOutField, NekDouble def_value=0., NekDouble tolerance=NekConstants::kFindDistanceMin)
Interpolate from an expansion to an expansion.
FieldSharedPtr m_f
Field object.
Definition: Module.h:239
std::map< std::string, ConfigOption > m_config
List of configuration values.
Definition: Module.h:272
void InterpolateFieldToPts(std::vector< MultiRegions::ExpListSharedPtr > &field0, LibUtilities::PtsFieldSharedPtr &pts, NekDouble clamp_low, NekDouble clamp_up, NekDouble def_value, NekDouble tolerance)
void PrintProgressbar(const int position, const int goal) const
static std::shared_ptr< Module > create(FieldSharedPtr f)
Creates an instance of this class.
void v_Process(po::variables_map &vm) override
Write mesh to output file.
Abstract base class for processing modules.
Definition: Module.h:301
void SetProgressCallback(FuncPointerT func, ObjectPointerT obj)
sets a callback funtion which gets called every time the interpolation progresses
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, std::string pDesc="")
Register a class with the factory.
static SessionReaderSharedPtr CreateInstance(int argc, char *argv[])
Creates an instance of the SessionReader class.
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.
static bool GenerateVariableSet(const std::string &str, const std::vector< std::string > &variables, std::set< int > &out)
Generate a set of variable locations.
Definition: ParseUtils.cpp:168
static bool GenerateVector(const std::string &str, std::vector< T > &out)
Takes a comma-separated string and converts it to entries in a vector.
Definition: ParseUtils.cpp:130
static MeshGraphSharedPtr Read(const LibUtilities::SessionReaderSharedPtr pSession, LibUtilities::DomainRangeShPtr rng=LibUtilities::NullDomainRangeShPtr, bool fillGraph=true, SpatialDomains::MeshGraphSharedPtr partitionedGraph=nullptr)
Definition: MeshGraphIO.cpp:51
std::shared_ptr< Field > FieldSharedPtr
Definition: Field.hpp:1026
std::pair< ModuleType, std::string > ModuleKey
Definition: Module.h:180
ModuleFactory & GetModuleFactory()
Definition: Module.cpp:47
int PrintProgressbar(const int position, const int goal, const std::string message, int lastprogress=-1)
Prints a progressbar.
Definition: Progressbar.hpp:65
static FieldMetaDataMap NullFieldMetaDataMap
Definition: FieldIO.h:51
std::shared_ptr< PtsField > PtsFieldSharedPtr
Definition: PtsField.h:184
std::shared_ptr< DomainRange > DomainRangeShPtr
Definition: DomainRange.h:64
std::shared_ptr< CsvIO > CsvIOSharedPtr
Definition: CsvIO.h:74
std::shared_ptr< PtsIO > PtsIOSharedPtr
Definition: PtsIO.h:90
CommFactory & GetCommFactory()
@ eFourierHalfModeIm
Fourier Modified expansions with just the imaginary part of the first mode .
Definition: BasisType.h:68
@ eFourierHalfModeRe
Fourier Modified expansions with just the real part of the first mode .
Definition: BasisType.h:66
@ eFourier
Fourier Expansion .
Definition: BasisType.h:55
std::map< int, ExpansionInfoShPtr > ExpansionInfoMap
Definition: MeshGraph.h:141
std::vector< double > q(NPUPPER *NPUPPER)
double NekDouble
T Vmin(int n, const T *x, const int incx)
Return the minimum element in x - called vmin to avoid conflict with min.
Definition: Vmath.hpp:725
void Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*x.
Definition: Vmath.hpp:100
void Zero(int n, T *x, const int incx)
Zero vector.
Definition: Vmath.hpp:273
T Vmax(int n, const T *x, const int incx)
Return the maximum element in x – called vmax to avoid conflict with max.
Definition: Vmath.hpp:644
STL namespace.
Represents a command-line configuration option.
Definition: Module.h:129