Nektar++
FilterFieldConvert.cpp
Go to the documentation of this file.
1///////////////////////////////////////////////////////////////////////////////
2//
3// File: FilterFieldConvert.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: Base clase for filters performing operations on samples
32// of the field.
33//
34///////////////////////////////////////////////////////////////////////////////
35
37#include <boost/algorithm/string/classification.hpp>
38#include <boost/algorithm/string/predicate.hpp>
39#include <boost/algorithm/string/split.hpp>
40#include <boost/program_options.hpp>
41
42namespace Nektar::SolverUtils
43{
47
50 const std::shared_ptr<EquationSystem> &pEquation, const ParamMap &pParams)
51 : Filter(pSession, pEquation)
52{
53 m_dt = m_session->GetParameter("TimeStep");
54
55 // OutputFile
56 std::string ext = ".fld";
57 m_outputFile = Filter::SetupOutput(ext, pParams);
58
59 // Time after which we need to write checkfiles
60 auto it = pParams.find("OutputStartTime");
61 if (it == pParams.end())
62 {
64 }
65 else
66 {
67 LibUtilities::Equation equ(m_session->GetInterpreter(), it->second);
69 }
70
71 // Restart file
72 it = pParams.find("RestartFile");
73 if (it == pParams.end())
74 {
75 m_restartFile = "";
76 }
77 else
78 {
79 ASSERTL0(it->second.length() > 0, "Missing parameter 'RestartFile'.");
80 if (it->second.find_last_of('.') != std::string::npos)
81 {
82 m_restartFile = it->second;
83 }
84 else
85 {
86 std::stringstream outname;
87 outname << it->second << ".fld";
88 m_restartFile = outname.str();
89 }
90 }
91
92 // OutputFrequency
93 it = pParams.find("OutputFrequency");
94 if (it == pParams.end())
95 {
96 m_outputFrequency = m_session->GetParameter("NumSteps");
97 }
98 else
99 {
100 LibUtilities::Equation equ(m_session->GetInterpreter(), it->second);
101 m_outputFrequency = round(equ.Evaluate());
102 }
103
104 // The base class can use SampleFrequency = OutputFrequency
105 // (Derived classes need to override this if needed)
107
108 // Phase sampling option
109 it = pParams.find("PhaseAverage");
110 if (it == pParams.end())
111 {
112 m_phaseSample = false;
113 }
114 else
115 {
116 std::string sOption = it->second.c_str();
117 m_phaseSample = (boost::iequals(sOption, "true")) ||
118 (boost::iequals(sOption, "yes"));
119 }
120
121 if (m_phaseSample)
122 {
123 auto itPeriod = pParams.find("PhaseAveragePeriod");
124 auto itPhase = pParams.find("PhaseAveragePhase");
125
126 // Error if only one of the required params for PhaseAverage is present
127 ASSERTL0(
128 (itPeriod != pParams.end() && itPhase != pParams.end()),
129 "The phase sampling feature requires both 'PhaseAveragePeriod' and "
130 "'PhaseAveragePhase' to be set.");
131
132 LibUtilities::Equation equPeriod(m_session->GetInterpreter(),
133 itPeriod->second);
134 m_phaseSamplePeriod = equPeriod.Evaluate();
135
136 LibUtilities::Equation equPhase(m_session->GetInterpreter(),
137 itPhase->second);
138 m_phaseSamplePhase = equPhase.Evaluate();
139
140 // Check that phase and period are within required limits
142 "PhaseAveragePeriod must be greater than 0.");
144 "PhaseAveragePhase must be between 0 and 1.");
145
146 // Load sampling frequency, overriding the previous value
147 it = pParams.find("SampleFrequency");
148 if (it == pParams.end())
149 {
151 }
152 else
153 {
154 LibUtilities::Equation equ(m_session->GetInterpreter(), it->second);
155 m_sampleFrequency = round(equ.Evaluate());
156 }
157
158 // Compute tolerance within which sampling occurs.
160
161 // Display worst case scenario sampling tolerance for exact phase, if
162 // verbose option is active
163 if (m_session->GetComm()->GetRank() == 0 &&
164 m_session->DefinesCmdLineArgument("verbose"))
165 {
166 std::cout << "Phase sampling activated with period "
167 << m_phaseSamplePeriod << " and phase "
168 << m_phaseSamplePhase << "." << std::endl
169 << "Sampling within a tolerance of "
170 << std::setprecision(6) << m_phaseTolerance << "."
171 << std::endl;
172 }
173 }
174
175 m_numSamples = 0;
176 m_index = 0;
177 m_outputIndex = 0;
178
179 //
180 // FieldConvert modules
181 //
182 m_f = std::shared_ptr<Field>(new Field());
183 std::vector<std::string> modcmds;
184 // Process modules
185 std::stringstream moduleStream;
186 it = pParams.find("Modules");
187 if (it != pParams.end())
188 {
189 moduleStream.str(it->second);
190 }
191 while (!moduleStream.fail())
192 {
193 std::string sMod;
194 moduleStream >> sMod;
195 if (!moduleStream.fail())
196 {
197 modcmds.push_back(sMod);
198 }
199 }
200 // Output module
201 modcmds.push_back(m_outputFile);
202 // Create modules
203 CreateModules(modcmds);
204 // Strip options from m_outputFile
205 std::vector<std::string> tmp;
206 boost::split(tmp, m_outputFile, boost::is_any_of(":"));
207 std::string outName = tmp[0];
208 ext = fs::path(outName).extension().string();
209 m_outputFile = Filter::SetupOutput(ext, outName);
210
211 // Prevent checking before overwriting
212 it = pParams.find("options");
213 if (it != pParams.end())
214 {
215 int argc = 0;
216 std::string strargv;
217 std::vector<char *> argv;
218
219 strargv = "dummy " + it->second;
220 strargv.push_back((char)0);
221 bool flag = true;
222 for (size_t i = 0; strargv[i]; ++i)
223 {
224 if (strargv[i] != ' ' && flag)
225 {
226 argv.push_back(&strargv[i]);
227 ++argc;
228 flag = false;
229 }
230 if (strargv[i] == ' ')
231 {
232 flag = true;
233 strargv[i] = 0;
234 }
235 }
236
237 po::options_description desc("Available options");
238
239 // clang-format off
240 desc.add_options()
241 ("help,h",
242 "Produce this help message.")
243 ("modules-list,l",
244 "Print the list of available modules.")
245 ("output-points,n", po::value<int>(),
246 "Output at n equipspaced points along the "
247 "collapsed coordinates (for .dat, .vtu).")
248 ("output-points-hom-z", po::value<int>(),
249 "Number of planes in the z-direction for output of "
250 "Homogeneous 1D expansion(for .dat, .vtu).")
251 ("error,e",
252 "Write error of fields for regression checking")
253 ("force-output,f",
254 "Force the output to be written without any checks")
255 ("range,r", po::value<std::string>(),
256 "Define output range i.e. (-r xmin,xmax,ymin,ymax,zmin,zmax) "
257 "in which any vertex is contained.")
258 ("no-equispaced",
259 "Do not use equispaced output.")
260 ("nparts", po::value<int>(),
261 "Define nparts if running serial problem to mimic "
262 "parallel run with many partitions.")
263 ("npz", po::value<int>(),
264 "Used to define number of partitions in z for Homogeneous1D "
265 "expansions for parallel runs.")
266 ("npt", po::value<int>(),
267 "Used to define number of partitions in time for Parareal runs. ")
268 ("onlyshape", po::value<std::string>(),
269 "Only use element with defined shape type i.e. -onlyshape "
270 " Tetrahedron")
271 ("part-only", po::value<int>(),
272 "Partition into specified npart partitions and exit")
273 ("part-only-overlapping", po::value<int>(),
274 "Partition into specified npart overlapping partitions and exit")
275 ("modules-opt,p", po::value<std::string>(),
276 "Print options for a module.")
277 ("module,m", po::value<std::vector<std::string> >(),
278 "Specify modules which are to be used.")
279 ("use-session-variables",
280 "Use variables defined in session for output")
281 ("use-session-expansion",
282 "Use expansion defined in session.")
283 ("verbose,v",
284 "Enable verbose mode.");
285 // clang-format on
286
287 po::options_description hidden("Hidden options");
288
289 // clang-format off
290 hidden.add_options()
291 ("input-file", po::value<std::vector<std::string> >(),
292 "Input filename");
293 // clang-format on
294
295 po::options_description cmdline_options;
296 cmdline_options.add(hidden).add(desc);
297
298 po::options_description visible("Allowed options");
299 visible.add(desc);
300
301 po::positional_options_description p;
302 p.add("input-file", -1);
303
304 try
305 {
306 po::store(po::command_line_parser(argc, &(argv[0]))
307 .options(cmdline_options)
308 .positional(p)
309 .run(),
310 m_vm);
311 po::notify(m_vm);
312 }
313 catch (const std::exception &e)
314 {
315 std::cerr << e.what() << std::endl;
316 std::cerr << desc;
317 }
318 }
319 m_vm.insert(std::make_pair("force-output", po::variable_value()));
320}
321
323{
324}
325
328 const NekDouble &time)
329{
330 v_FillVariablesName(pFields);
331
332 // m_variables need to be filled by a derived class
333 m_outFields.resize(m_variables.size());
334 int nfield;
335
336 for (int n = 0; n < m_variables.size(); ++n)
337 {
338 // if n >= pFields.size() assum we have used n=0 field
339 nfield = (n < pFields.size()) ? n : 0;
340
341 m_outFields[n] =
342 Array<OneD, NekDouble>(pFields[nfield]->GetNcoeffs(), 0.0);
343 }
344
345 m_fieldMetaData["InitialTime"] = boost::lexical_cast<std::string>(time);
346
347 // Load restart file if necessary
348 if (m_restartFile != "")
349 {
350 // Load file
351 std::vector<LibUtilities::FieldDefinitionsSharedPtr> fieldDef;
352 std::vector<std::vector<NekDouble>> fieldData;
353 LibUtilities::FieldMetaDataMap fieldMetaData;
356 fld->Import(m_restartFile, fieldDef, fieldData, fieldMetaData);
357
358 // Extract fields to output
359 int nfield = -1, k;
360 for (int j = 0; j < m_variables.size(); ++j)
361 {
362 // see if m_variables is part of pFields definition and if
363 // so use that field for extract
364 for (k = 0; k < pFields.size(); ++k)
365 {
366 if (pFields[k]->GetSession()->GetVariable(k) == m_variables[j])
367 {
368 nfield = k;
369 break;
370 }
371 }
372 if (nfield == -1)
373 {
374 nfield = 0;
375 }
376
377 for (int i = 0; i < fieldData.size(); ++i)
378 {
379 pFields[nfield]->ExtractDataToCoeffs(
380 fieldDef[i], fieldData[i], m_variables[j], m_outFields[j]);
381 }
382 }
383
384 // Load information for numSamples
385 if (fieldMetaData.count("NumberOfFieldDumps"))
386 {
387 m_numSamples = atoi(fieldMetaData["NumberOfFieldDumps"].c_str());
388 }
389 else
390 {
391 m_numSamples = 1;
392 }
393
394 if (fieldMetaData.count("InitialTime"))
395 {
396 m_fieldMetaData["InitialTime"] = fieldMetaData["InitialTime"];
397 }
398
399 // Load information for outputIndex
400 if (fieldMetaData.count("FilterFileNum"))
401 {
402 m_outputIndex = atoi(fieldMetaData["FilterFileNum"].c_str());
403 }
404
405 // Divide by scale
406 NekDouble scale = v_GetScale();
407 for (int n = 0; n < m_outFields.size(); ++n)
408 {
409 Vmath::Smul(m_outFields[n].size(), 1.0 / scale, m_outFields[n], 1,
410 m_outFields[n], 1);
411 }
412 }
413}
414
417{
418 int nfield = pFields.size();
419 m_variables.resize(pFields.size());
420 for (int n = 0; n < nfield; ++n)
421 {
422 m_variables[n] = pFields[n]->GetSession()->GetVariable(n);
423 }
424
425 // Need to create a dummy coeffs vector to get extra variables names...
426 std::vector<Array<OneD, NekDouble>> coeffs(nfield);
427 for (int n = 0; n < nfield; ++n)
428 {
429 coeffs[n] = pFields[n]->GetCoeffs();
430 }
431 // Get extra variables
432 auto equ = m_equ.lock();
433 ASSERTL0(equ, "Weak pointer expired");
434 equ->ExtraFldOutput(coeffs, m_variables);
435}
436
439 const NekDouble &time)
440{
441 m_index++;
442
443 if (m_index % m_sampleFrequency > 0 ||
444 (time - m_outputStartTime) < -1.0e-07)
445 {
446 return;
447 }
448
449 // Append extra fields
450 int nfield = pFields.size();
451 std::vector<Array<OneD, NekDouble>> coeffs(nfield);
452 for (int n = 0; n < nfield; ++n)
453 {
454 coeffs[n] = pFields[n]->GetCoeffs();
455 }
456 std::vector<std::string> variables = m_variables;
457 auto equ = m_equ.lock();
458 ASSERTL0(equ, "Weak pointer expired");
459 equ->ExtraFldOutput(coeffs, variables);
460
461 if (m_phaseSample)
462 {
463 // The sample is added to the filter only if the current time
464 // corresponds to the correct phase. Introducing M as number of
465 // cycles and N nondimensional phase (between 0 and 1):
466 // t = M * m_phaseSamplePeriod + N * m_phaseSamplePeriod
467 int currentCycle = floor(time / m_phaseSamplePeriod);
468 NekDouble currentPhase = time / m_phaseSamplePeriod - currentCycle;
469
470 // Evaluate phase relative to the requested value.
471 NekDouble relativePhase = fabs(m_phaseSamplePhase - currentPhase);
472
473 // Check if relative phase is within required tolerance and sample.
474 // Care must be taken to handle the cases at phase 0 as the sample might
475 // have to be taken at the very end of the previous cycle instead.
476 if (relativePhase < m_phaseTolerance ||
477 fabs(relativePhase - 1) < m_phaseTolerance)
478 {
479 m_numSamples++;
480 v_ProcessSample(pFields, coeffs, time);
481 if (m_session->GetComm()->GetRank() == 0 &&
482 m_session->DefinesCmdLineArgument("verbose"))
483 {
484 std::cout << "Sample: " << std::setw(8) << std::left
485 << m_numSamples << "Phase: " << std::setw(8)
486 << std::left << currentPhase << std::endl;
487 }
488 }
489 }
490 else
491 {
492 m_numSamples++;
493 v_ProcessSample(pFields, coeffs, time);
494 }
495
496 if (m_index % m_outputFrequency == 0)
497 {
498 m_fieldMetaData["FinalTime"] = boost::lexical_cast<std::string>(time);
499 v_PrepareOutput(pFields, time);
500 m_fieldMetaData["FilterFileNum"] = std::to_string(++m_outputIndex);
501 OutputField(pFields, m_outputIndex);
502 }
503}
504
507 const NekDouble &time)
508{
509 m_fieldMetaData["FinalTime"] = boost::lexical_cast<std::string>(time);
510 v_PrepareOutput(pFields, time);
511 OutputField(pFields);
512}
513
516 &pFields,
517 std::vector<Array<OneD, NekDouble>> &fieldcoeffs,
518 [[maybe_unused]] const NekDouble &time)
519{
520 for (int n = 0; n < m_outFields.size(); ++n)
521 {
522 Vmath::Vcopy(m_outFields[n].size(), fieldcoeffs[n], 1, m_outFields[n],
523 1);
524 }
525}
526
529{
530 NekDouble scale = v_GetScale();
531 for (int n = 0; n < m_outFields.size(); ++n)
532 {
533 Vmath::Smul(m_outFields[n].size(), scale, m_outFields[n], 1,
534 m_outFields[n], 1);
535 }
536
537 CreateFields(pFields);
538
539 // Determine new file name
540 std::stringstream tmpOutname;
541 std::string outname;
542 int dot = m_outputFile.find_last_of('.');
543 std::string name = m_outputFile.substr(0, dot);
544 std::string ext = m_outputFile.substr(dot, m_outputFile.length() - dot);
545 std::string suffix = v_GetFileSuffix();
546
547 if (dump == -1) // final dump
548 {
549 tmpOutname << name << suffix << ext;
550 }
551 else
552 {
553 tmpOutname << name << "_" << dump << suffix << ext;
554 }
555 outname = Filter::SetupOutput(ext, tmpOutname.str());
556 m_modules[m_modules.size() - 1]->RegisterConfig("outfile", outname);
557
558 // Run field process.
559 for (int n = 0; n < SIZE_ModulePriority; ++n)
560 {
561 ModulePriority priority = static_cast<ModulePriority>(n);
562 for (int i = 0; i < m_modules.size(); ++i)
563 {
564 if (m_modules[i]->GetModulePriority() == priority)
565 {
566 m_modules[i]->Process(m_vm);
567 }
568 }
569 }
570
571 // Empty m_f to save memory
572 m_f->ClearField();
573
574 if (dump != -1) // not final dump so rescale
575 {
576 for (int n = 0; n < m_outFields.size(); ++n)
577 {
578 Vmath::Smul(m_outFields[n].size(), 1.0 / scale, m_outFields[n], 1,
579 m_outFields[n], 1);
580 }
581 }
582}
583
585{
586 return true;
587}
588
589void FilterFieldConvert::CreateModules(std::vector<std::string> &modcmds)
590{
591 for (int i = 0; i < modcmds.size(); ++i)
592 {
593 // First split each command by the colon separator.
594 std::vector<std::string> tmp1;
595 ModuleKey module;
596 int offset = 1;
597
598 boost::split(tmp1, modcmds[i], boost::is_any_of(":"));
599
600 if (i == modcmds.size() - 1)
601 {
602 module.first = eOutputModule;
603
604 // If no colon detected, automatically detect mesh type from
605 // file extension. Otherwise override and use tmp1[1] as the
606 // module to load. This also allows us to pass options to
607 // input/output modules. So, for example, to override
608 // filename.xml to be read as vtk, you use:
609 //
610 // filename.xml:vtk:opt1=arg1:opt2=arg2
611 if (tmp1.size() == 1)
612 {
613 int dot = tmp1[0].find_last_of('.') + 1;
614 std::string ext = tmp1[0].substr(dot, tmp1[0].length() - dot);
615 module.second = ext;
616 tmp1.push_back(std::string("outfile=") + tmp1[0]);
617 }
618 else
619 {
620 module.second = tmp1[1];
621 tmp1.push_back(std::string("outfile=") + tmp1[0]);
622 offset++;
623 }
624 }
625 else
626 {
627 module.first = eProcessModule;
628 module.second = tmp1[0];
629 }
630
631 // Create modules
632 ModuleSharedPtr mod;
633 mod = GetModuleFactory().CreateInstance(module, m_f);
634 m_modules.push_back(mod);
635
636 // Set options for this module.
637 for (int j = offset; j < tmp1.size(); ++j)
638 {
639 std::vector<std::string> tmp2;
640 boost::split(tmp2, tmp1[j], boost::is_any_of("="));
641
642 if (tmp2.size() == 1)
643 {
644 mod->RegisterConfig(tmp2[0]);
645 }
646 else if (tmp2.size() == 2)
647 {
648 mod->RegisterConfig(tmp2[0], tmp2[1]);
649 }
650 else
651 {
652 std::cerr << "ERROR: Invalid module configuration: format is "
653 << "either :arg or :arg=val" << std::endl;
654 abort();
655 }
656 }
657
658 // Ensure configuration options have been set.
659 mod->SetDefaults();
660 }
661
662 // Include equispaced output if needed
663 Array<OneD, int> modulesCount(SIZE_ModulePriority, 0);
664 for (int i = 0; i < m_modules.size(); ++i)
665 {
666 ++modulesCount[m_modules[i]->GetModulePriority()];
667 }
668 if (modulesCount[eModifyPts] != 0 && modulesCount[eCreatePts] == 0 &&
669 modulesCount[eConvertExpToPts] == 0)
670 {
671 ModuleKey module;
672 ModuleSharedPtr mod;
673 module.first = eProcessModule;
674 module.second = std::string("equispacedoutput");
675 mod = GetModuleFactory().CreateInstance(module, m_f);
676 m_modules.insert(m_modules.end() - 1, mod);
677 mod->SetDefaults();
678 }
679
680 // Check if modules provided are compatible
682}
683
686{
687 // Fill some parameters of m_f
688 m_f->m_session = m_session;
689 m_f->m_graph = pFields[0]->GetGraph();
690 m_f->m_comm = m_f->m_session->GetComm();
691 m_f->m_fieldMetaDataMap = m_fieldMetaData;
692 m_f->m_fieldPts = LibUtilities::NullPtsField;
693 // Create m_f->m_exp
694 m_f->m_numHomogeneousDir = 0;
695 if (pFields[0]->GetExpType() == MultiRegions::e3DH1D)
696 {
697 m_f->m_numHomogeneousDir = 1;
698 }
699 else if (pFields[0]->GetExpType() == MultiRegions::e3DH2D)
700 {
701 m_f->m_numHomogeneousDir = 2;
702 }
703
704 m_f->m_exp.resize(m_variables.size());
705 m_f->m_exp[0] = pFields[0];
706 int nfield;
707 for (int n = 0; n < m_variables.size(); ++n)
708 {
709 // if n >= pFields.size() assume we have used n=0 field
710 nfield = (n < pFields.size()) ? n : 0;
711
712 m_f->m_exp[n] =
713 m_f->AppendExpList(m_f->m_numHomogeneousDir, m_variables[0]);
714 m_f->m_exp[n]->SetWaveSpace(false);
715
716 ASSERTL1(pFields[nfield]->GetNcoeffs() == m_outFields[n].size(),
717 "pFields[nfield] does not have the "
718 "same number of coefficients as m_outFields[n]");
719
720 m_f->m_exp[n]->ExtractCoeffsToCoeffs(pFields[nfield], m_outFields[n],
721 m_f->m_exp[n]->UpdateCoeffs());
722
723 m_f->m_exp[n]->BwdTrans(m_f->m_exp[n]->GetCoeffs(),
724 m_f->m_exp[n]->UpdatePhys());
725 }
726 m_f->m_variables = m_variables;
727}
728
729// This function checks validity conditions for the list of modules provided
730void FilterFieldConvert::CheckModules(std::vector<ModuleSharedPtr> &modules)
731{
732 // Count number of modules by priority
733 Array<OneD, int> modulesCount(SIZE_ModulePriority, 0);
734 for (int i = 0; i < modules.size(); ++i)
735 {
736 ++modulesCount[modules[i]->GetModulePriority()];
737 }
738
739 // FilterFieldConvert already starts with m_exp, so anything before
740 // eModifyExp is not valid, and also eCreatePts
741 if (modulesCount[eCreateGraph] != 0 ||
742 modulesCount[eCreateFieldData] != 0 ||
743 modulesCount[eModifyFieldData] != 0 || modulesCount[eCreateExp] != 0 ||
744 modulesCount[eFillExp] != 0 || modulesCount[eCreatePts] != 0)
745 {
746 std::stringstream ss;
747 ss << "Module(s): ";
748 for (int i = 0; i < modules.size(); ++i)
749 {
750 if (modules[i]->GetModulePriority() == eCreateGraph ||
751 modules[i]->GetModulePriority() == eCreateFieldData ||
752 modules[i]->GetModulePriority() == eModifyFieldData ||
753 modules[i]->GetModulePriority() == eCreateExp ||
754 modules[i]->GetModulePriority() == eFillExp ||
755 modules[i]->GetModulePriority() == eCreatePts)
756 {
757 ss << modules[i]->GetModuleName() << " ";
758 }
759 }
760 ss << "not compatible with FilterFieldConvert.";
761 ASSERTL0(false, ss.str());
762 }
763
764 // Modules of type eConvertExpToPts are not compatible with eBndExtraction
765 if (modulesCount[eConvertExpToPts] != 0 &&
766 modulesCount[eBndExtraction] != 0)
767 {
768 std::stringstream ss;
769 ss << "Module(s): ";
770 for (int i = 0; i < modules.size(); ++i)
771 {
772 if (modules[i]->GetModulePriority() == eBndExtraction)
773 {
774 ss << modules[i]->GetModuleName() << " ";
775 }
776 }
777 ss << "is not compatible with module(s): ";
778 for (int i = 0; i < modules.size(); ++i)
779 {
780 if (modules[i]->GetModulePriority() == eConvertExpToPts)
781 {
782 ss << modules[i]->GetModuleName() << " ";
783 }
784 }
785 ss << ".";
786 ASSERTL0(false, ss.str());
787 }
788}
789} // namespace Nektar::SolverUtils
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:208
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode....
Definition: ErrorUtil.hpp:242
static std::shared_ptr< FieldIO > CreateForFile(const LibUtilities::SessionReaderSharedPtr session, const std::string &filename)
Construct a FieldIO object for a given input filename.
Definition: FieldIO.cpp:223
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, std::string pDesc="")
Register a class with the factory.
tBaseSharedPtr CreateInstance(tKey idKey, tParam... args)
Create an instance of the class referred to by idKey.
SOLVER_UTILS_EXPORT bool v_IsTimeDependent() override
virtual SOLVER_UTILS_EXPORT void v_ProcessSample(const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields, std::vector< Array< OneD, NekDouble > > &fieldcoeffs, const NekDouble &time)
static std::string className
Name of the class.
void CreateFields(const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields)
void CheckModules(std::vector< ModuleSharedPtr > &modules)
SOLVER_UTILS_EXPORT void v_Finalise(const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields, const NekDouble &time) override
SOLVER_UTILS_EXPORT void v_Update(const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields, const NekDouble &time) override
virtual SOLVER_UTILS_EXPORT void v_FillVariablesName(const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields)
virtual SOLVER_UTILS_EXPORT void v_PrepareOutput(const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields, const NekDouble &time)
SOLVER_UTILS_EXPORT ~FilterFieldConvert() override
void CreateModules(std::vector< std::string > &modcmds)
static FilterSharedPtr create(const LibUtilities::SessionReaderSharedPtr &pSession, const std::shared_ptr< EquationSystem > &pEquation, const std::map< std::string, std::string > &pParams)
Creates an instance of this class.
SOLVER_UTILS_EXPORT FilterFieldConvert(const LibUtilities::SessionReaderSharedPtr &pSession, const std::shared_ptr< EquationSystem > &pEquation, const ParamMap &pParams)
void OutputField(const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields, int dump=-1)
LibUtilities::FieldMetaDataMap m_fieldMetaData
std::vector< Array< OneD, NekDouble > > m_outFields
virtual SOLVER_UTILS_EXPORT NekDouble v_GetScale()
virtual SOLVER_UTILS_EXPORT std::string v_GetFileSuffix()
std::vector< ModuleSharedPtr > m_modules
SOLVER_UTILS_EXPORT void v_Initialise(const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields, const NekDouble &time) override
SOLVER_UTILS_EXPORT std::string SetupOutput(const std::string ext, const ParamMap &pParams)
Definition: Filter.h:139
LibUtilities::SessionReaderSharedPtr m_session
Definition: Filter.h:93
const std::weak_ptr< EquationSystem > m_equ
Definition: Filter.h:94
std::map< std::string, std::string > ParamMap
Definition: Filter.h:66
std::pair< ModuleType, std::string > ModuleKey
Definition: Module.h:180
std::shared_ptr< Module > ModuleSharedPtr
Definition: Module.h:329
ModuleFactory & GetModuleFactory()
Definition: Module.cpp:47
std::shared_ptr< FieldIO > FieldIOSharedPtr
Definition: FieldIO.h:322
std::map< std::string, std::string > FieldMetaDataMap
Definition: FieldIO.h:50
std::shared_ptr< SessionReader > SessionReaderSharedPtr
static PtsFieldSharedPtr NullPtsField
Definition: PtsField.h:185
FilterFactory & GetFilterFactory()
double NekDouble
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 Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.hpp:825