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