Nektar++
Public Member Functions | Static Public Member Functions | Private Member Functions | Static Private Member Functions | Private Attributes | Friends | List of all members
Nektar::LibUtilities::SessionReader Class Reference

Reads and parses information from a Nektar++ XML session file. More...

#include <SessionReader.h>

Public Member Functions

 SessionReader (int argc, char *argv[], const std::vector< std::string > &pFilenames, const CommSharedPtr &pComm)
 
 ~SessionReader ()
 Destructor. More...
 
void InitSession (const std::vector< std::string > &filenames=std::vector< std::string >())
 
TiXmlDocument & GetDocument ()
 Provides direct access to the TiXmlDocument object. More...
 
TiXmlElement * GetElement (const std::string &pPath)
 Provides direct access to the TiXmlElement specified. More...
 
bool DefinesElement (const std::string &pPath) const
 Tests if a specified element is defined in the XML document. More...
 
const std::vector< std::string > & GetFilenames () const
 Returns the filename of the loaded XML document. More...
 
const std::string & GetSessionName () const
 Returns the session name of the loaded XML document. More...
 
const std::string GetSessionNameRank () const
 Returns the session name with process rank. More...
 
CommSharedPtr GetComm ()
 Returns the communication object. More...
 
bool GetSharedFilesystem ()
 Returns if file system shared. More...
 
void Finalise ()
 Finalises the session. More...
 
bool DefinesParameter (const std::string &name) const
 Checks if a parameter is specified in the XML document. More...
 
const NekDoubleGetParameter (const std::string &pName) const
 Returns the value of the specified parameter. More...
 
void LoadParameter (const std::string &name, int &var) const
 Load an integer parameter. More...
 
void LoadParameter (const std::string &name, int &var, const int &def) const
 Check for and load an integer parameter. More...
 
void LoadParameter (const std::string &name, NekDouble &var) const
 Load a double precision parameter. More...
 
void LoadParameter (const std::string &name, NekDouble &var, const NekDouble &def) const
 Check for and load a double-precision parameter. More...
 
void SetParameter (const std::string &name, int &var)
 Set an integer parameter. More...
 
void SetParameter (const std::string &name, NekDouble &var)
 Set a double precision parameter. More...
 
bool DefinesSolverInfo (const std::string &name) const
 Checks if a solver info property is specified. More...
 
const std::string & GetSolverInfo (const std::string &pProperty) const
 Returns the value of the specified solver info property. More...
 
void SetSolverInfo (const std::string &pProperty, const std::string &pValue)
 Sets the value of the specified solver info property. More...
 
template<typename T >
const T GetSolverInfoAsEnum (const std::string &pName) const
 Returns the value of the specified solver info property as enum. More...
 
template<typename T >
const T GetValueAsEnum (const std::string &pName, const std::string &vValue) const
 Returns the value of the specified property and value as enum. More...
 
void LoadSolverInfo (const std::string &name, std::string &var, const std::string &def="") const
 Check for and load a solver info property. More...
 
void MatchSolverInfo (const std::string &name, const std::string &trueval, bool &var, const bool &def=false) const
 Check if the value of a solver info property matches. More...
 
bool MatchSolverInfo (const std::string &name, const std::string &trueval) const
 Check if the value of a solver info property matches. More...
 
template<typename T >
bool MatchSolverInfoAsEnum (const std::string &name, const T &trueval) const
 Check if the value of a solver info property matches. More...
 
bool DefinesGlobalSysSolnInfo (const std::string &variable, const std::string &property) const
 
const std::string & GetGlobalSysSolnInfo (const std::string &variable, const std::string &property) const
 
bool DefinesTimeIntScheme () const
 Returns true if the TIMEINTEGRATIONSCHEME section is defined in the session file. More...
 
const TimeIntSchemeGetTimeIntScheme () const
 Returns the time integration scheme structure m_timeIntScheme from the session file. More...
 
std::string GetGeometryType () const
 
bool DefinesGeometricInfo (const std::string &name) const
 Checks if a geometric info property is defined. More...
 
void LoadGeometricInfo (const std::string &name, std::string &var, const std::string &def="") const
 Checks for and load a geometric info string property. More...
 
void LoadGeometricInfo (const std::string &name, bool &var, const bool &def=false) const
 Checks for and loads a geometric info boolean property. More...
 
void LoadGeometricInfo (const std::string &name, NekDouble &var, const NekDouble &def=0.0) const
 Checks for and loads a geometric info double-precision property. More...
 
void MatchGeometricInfo (const std::string &name, const std::string &trueval, bool &var, const bool &def=false) const
 Check if the value of a geometric info string property matches. More...
 
const std::string & GetVariable (const unsigned int &idx) const
 Returns the name of the variable specified by the given index. More...
 
void SetVariable (const unsigned int &idx, std::string newname)
 
std::vector< std::string > GetVariables () const
 Returns the names of all variables. More...
 
bool DefinesFunction (const std::string &name) const
 Checks if a specified function is defined in the XML document. More...
 
bool DefinesFunction (const std::string &name, const std::string &variable, const int pDomain=0) const
 Checks if a specified function has a given variable defined. More...
 
EquationSharedPtr GetFunction (const std::string &name, const std::string &variable, const int pDomain=0) const
 Returns an EquationSharedPtr to a given function variable. More...
 
EquationSharedPtr GetFunction (const std::string &name, const unsigned int &var, const int pDomain=0) const
 Returns an EquationSharedPtr to a given function variable index. More...
 
enum FunctionType GetFunctionType (const std::string &name, const std::string &variable, const int pDomain=0) const
 Returns the type of a given function variable. More...
 
enum FunctionType GetFunctionType (const std::string &pName, const unsigned int &pVar, const int pDomain=0) const
 Returns the type of a given function variable index. More...
 
std::string GetFunctionFilename (const std::string &name, const std::string &variable, const int pDomain=0) const
 Returns the filename to be loaded for a given variable. More...
 
std::string GetFunctionFilename (const std::string &name, const unsigned int &var, const int pDomain=0) const
 Returns the filename to be loaded for a given variable index. More...
 
std::string GetFunctionFilenameVariable (const std::string &name, const std::string &variable, const int pDomain=0) const
 Returns the filename variable to be loaded for a given variable index. More...
 
InterpreterSharedPtr GetInterpreter ()
 Returns the instance of the Interpreter specific to this session. More...
 
bool DefinesTag (const std::string &pName) const
 Checks if a specified tag is defined. More...
 
void SetTag (const std::string &pName, const std::string &pValue)
 Sets a specified tag. More...
 
const std::string & GetTag (const std::string &pName) const
 Returns the value of a specified tag. More...
 
const FilterMapGetFilters () const
 
bool DefinesCmdLineArgument (const std::string &pName) const
 Checks if a specified cmdline argument has been given. More...
 
template<typename T >
GetCmdLineArgument (const std::string &pName) const
 Retrieves a command-line argument value. More...
 
void SubstituteExpressions (std::string &expr)
 Substitutes expressions defined in the XML document. More...
 
void SetUpXmlDoc ()
 
bool GetUpdateOptFile () const
 Get bool to update optimisation file. More...
 
void SetUpdateOptFile (bool flag)
 Set bool to update optimisation file. More...
 

Static Public Member Functions

static SessionReaderSharedPtr CreateInstance (int argc, char *argv[])
 Creates an instance of the SessionReader class. More...
 
static SessionReaderSharedPtr CreateInstance (int argc, char *argv[], std::vector< std::string > &pFilenames, const CommSharedPtr &pComm=CommSharedPtr())
 Creates an instance of the SessionReader class initialised using a separate list of XML documents. More...
 
static std::string RegisterEnumValue (std::string pEnum, std::string pString, int pEnumValue)
 Registers an enumeration value. More...
 
static std::string RegisterDefaultSolverInfo (const std::string &pName, const std::string &pValue)
 Registers the default string value of a solver info property. More...
 
static std::string RegisterCmdLineArgument (const std::string &pName, const std::string &pShortName, const std::string &pDescription)
 Registers a command-line argument with the session reader. More...
 
static std::string RegisterCmdLineFlag (const std::string &pName, const std::string &pShortName, const std::string &pDescription)
 Registers a command-line flag with the session reader. More...
 

Private Member Functions

 SessionReader (int argc, char *argv[])
 Main constructor. More...
 
void TestSharedFilesystem ()
 
std::vector< std::string > ParseCommandLineArguments (int argc, char *argv[])
 Parse the program arguments and fill m_cmdLineOptions. More...
 
std::string ParseSessionName (std::vector< std::string > &filenames)
 Parse the session name. More...
 
void LoadDoc (const std::string &pFilename, TiXmlDocument *pDoc) const
 Loads an xml file into a tinyxml doc and decompresses if needed. More...
 
TiXmlDocument * MergeDoc (const std::vector< std::string > &pFilenames) const
 Creates an XML document from a list of input files. More...
 
void ParseDocument ()
 Loads and parses the specified file. More...
 
void CreateComm (int &argc, char *argv[])
 Loads the given XML document and instantiates an appropriate communication object. More...
 
void PartitionComm ()
 Partitions the comm object based on session parameters. More...
 
void ReadParameters (TiXmlElement *conditions)
 Reads the PARAMETERS section of the XML document. More...
 
void ReadSolverInfo (TiXmlElement *conditions)
 Reads the SOLVERINFO section of the XML document. More...
 
void ReadGlobalSysSolnInfo (TiXmlElement *conditions)
 Reads the GLOBALSYSSOLNINFO section of the XML document. More...
 
void ReadTimeIntScheme (TiXmlElement *conditions)
 Reads the TIMEINTEGRATIONSCHEME section of the XML document. More...
 
void ReadExpressions (TiXmlElement *conditions)
 Reads the EXPRESSIONS section of the XML document. More...
 
void ReadVariables (TiXmlElement *conditions)
 Reads the VARIABLES section of the XML document. More...
 
void ReadFunctions (TiXmlElement *conditions)
 Reads the FUNCTIONS section of the XML document. More...
 
void ReadFilters (TiXmlElement *filters)
 Reads the FILTERS section of the XML document. More...
 
void CmdLineOverride ()
 Enforce parameters from command line arguments. More...
 
void VerifySolverInfo ()
 Check values of solver info options are valid. More...
 
void ParseEquals (const std::string &line, std::string &lhs, std::string &rhs)
 Parse a string in the form lhs = rhs. More...
 

Static Private Member Functions

static EnumMapListGetSolverInfoEnums ()
 String to enumeration map for Solver Info parameters. More...
 
static SolverInfoMapGetSolverInfoDefaults ()
 Default solver info options. More...
 
static GloSysSolnInfoListGetGloSysSolnList ()
 GlobalSysSoln Info map. More...
 
static CmdLineArgMapGetCmdLineArgMap ()
 CmdLine argument map. More...
 

Private Attributes

boost::program_options::variables_map m_cmdLineOptions
 
CommSharedPtr m_comm
 Communication object. More...
 
std::vector< std::string > m_filenames
 Filenames. More...
 
std::string m_sessionName
 Session name of the loaded XML document (filename minus ext). More...
 
TiXmlDocument * m_xmlDoc
 Pointer to the loaded XML document. More...
 
ParameterMap m_parameters
 Parameters. More...
 
SolverInfoMap m_solverInfo
 Solver information properties. More...
 
GeometricInfoMap m_geometricInfo
 Geometric information properties. More...
 
ExpressionMap m_expressions
 Expressions. More...
 
InterpreterSharedPtr m_interpreter
 Interpreter instance. More...
 
FunctionMap m_functions
 Functions. More...
 
VariableList m_variables
 Variables. More...
 
TagMap m_tags
 Custom tags. More...
 
FilterMap m_filters
 Filters map. More...
 
TimeIntScheme m_timeIntScheme
 Time integration scheme information. More...
 
bool m_verbose
 Be verbose. More...
 
bool m_sharedFilesystem
 Running on a shared filesystem. More...
 
bool m_updateOptFile
 Update optimisation file. More...
 

Friends

class MemoryManager< SessionReader >
 Support creation through MemoryManager. More...
 

Detailed Description

Reads and parses information from a Nektar++ XML session file.

This class provides an interface to Nektar++-specific content in a supplied XML document. It also initialises a Nektar++ session including setting up communication for parallel execution and where necessary partitioning the supplied mesh for running across multiple processes.

A session should be initialised at the beginning of a user's application by passing the command-line arguments. This not only allows the SessionReader to extract the name of the XML document to load containing Nektar++ session information, but also supplies the MPI arguments necessary for setting up parallel communication. The SessionReader should be initialised using the CreateInstance function:

static SessionReaderSharedPtr CreateInstance(int argc, char *argv[])
Creates an instance of the SessionReader class.
std::shared_ptr< SessionReader > SessionReaderSharedPtr

The instance vSession can now be passed to other key Nektar++ components during their construction.

Note
At the end of the user application, it is important to call the Finalise routine in order to finalise any MPI communication and correctly free resources.

The SessionReader class provides streamlined, validated access to session parameters, solver information and functions defined within a Nektar++ XML document. The available routines and their usage is documented below.

In the case of solver information properties, the classes to which these parameters are pertinent may register with the SessionReader class the set of valid values for a given property. Such values may also be associated with an enumeration value for more transparent use of the property values in code.

Definition at line 119 of file SessionReader.h.

Constructor & Destructor Documentation

◆ SessionReader() [1/2]

Nektar::LibUtilities::SessionReader::SessionReader ( int  argc,
char *  argv[],
const std::vector< std::string > &  pFilenames,
const CommSharedPtr pComm 
)

Definition at line 217 of file BasicUtils/SessionReader.cpp.

220 {
221  ASSERTL0(pFilenames.size() > 0, "No filenames specified.");
222 
223  ParseCommandLineArguments(argc, argv);
224  m_xmlDoc = 0;
225  m_filenames = pFilenames;
226 
228 
229  // Create communicator
230  if (!pComm.get())
231  {
232  CreateComm(argc, argv);
233  }
234  else
235  {
236  m_comm = pComm;
237  }
238 
240 
241  // If running in parallel change the default global sys solution
242  // type.
243  if (m_comm->GetSize() > 1)
244  {
245  GetSolverInfoDefaults()["GLOBALSYSSOLN"] = "IterativeStaticCond";
246  }
247 
249  m_interpreter->SetRandomSeed((m_comm->GetRank() + 1) *
250  (unsigned int)time(NULL));
251 
252  // Split up the communicator
253  PartitionComm();
254 }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:215
std::vector< std::string > ParseCommandLineArguments(int argc, char *argv[])
Parse the program arguments and fill m_cmdLineOptions.
static SolverInfoMap & GetSolverInfoDefaults()
Default solver info options.
std::string ParseSessionName(std::vector< std::string > &filenames)
Parse the session name.
InterpreterSharedPtr m_interpreter
Interpreter instance.
std::vector< std::string > m_filenames
Filenames.
std::string m_sessionName
Session name of the loaded XML document (filename minus ext).
CommSharedPtr m_comm
Communication object.
void CreateComm(int &argc, char *argv[])
Loads the given XML document and instantiates an appropriate communication object.
void PartitionComm()
Partitions the comm object based on session parameters.
TiXmlDocument * m_xmlDoc
Pointer to the loaded XML document.
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.

References ASSERTL0.

◆ ~SessionReader()

Nektar::LibUtilities::SessionReader::~SessionReader ( )

Destructor.

Definition at line 259 of file BasicUtils/SessionReader.cpp.

260 {
261  if (m_xmlDoc)
262  {
263  delete m_xmlDoc;
264  }
265 }

◆ SessionReader() [2/2]

Nektar::LibUtilities::SessionReader::SessionReader ( int  argc,
char *  argv[] 
)
private

Main constructor.

This constructor parses the command-line arguments given to the user application to set up any MPI communication, read supplied XML session files, and partition meshes where necessary.

Parameters
argcNumber of command-line arguments
argvArray of command-line arguments

Definition at line 185 of file BasicUtils/SessionReader.cpp.

186 {
187  m_xmlDoc = 0;
189 
190  ASSERTL0(m_filenames.size() > 0, "No session file(s) given.");
191 
193 
194  // Create communicator
195  CreateComm(argc, argv);
196 
198 
199  // If running in parallel change the default global sys solution
200  // type.
201  if (m_comm->GetSize() > 1)
202  {
203  GetSolverInfoDefaults()["GLOBALSYSSOLN"] = "IterativeStaticCond";
204  }
205 
207  m_interpreter->SetRandomSeed((m_comm->GetRank() + 1) *
208  (unsigned int)time(NULL));
209 
210  // Split up the communicator
211  PartitionComm();
212 }

References ASSERTL0.

Member Function Documentation

◆ CmdLineOverride()

void Nektar::LibUtilities::SessionReader::CmdLineOverride ( )
private

Enforce parameters from command line arguments.

Definition at line 2442 of file BasicUtils/SessionReader.cpp.

2443 {
2444  // Parse solver info overrides
2445  if (m_cmdLineOptions.count("solverinfo"))
2446  {
2447  std::vector<std::string> solverInfoList =
2448  m_cmdLineOptions["solverinfo"].as<std::vector<std::string>>();
2449 
2450  for (size_t i = 0; i < solverInfoList.size(); ++i)
2451  {
2452  std::string lhs, rhs;
2453 
2454  try
2455  {
2456  ParseEquals(solverInfoList[i], lhs, rhs);
2457  }
2458  catch (...)
2459  {
2460  NEKERROR(ErrorUtil::efatal, "Parse error with command line "
2461  "option: " +
2462  solverInfoList[i]);
2463  }
2464 
2465  std::string lhsUpper = boost::to_upper_copy(lhs);
2466  m_solverInfo[lhsUpper] = rhs;
2467  }
2468  }
2469 
2470  if (m_cmdLineOptions.count("parameter"))
2471  {
2472  std::vector<std::string> parametersList =
2473  m_cmdLineOptions["parameter"].as<std::vector<std::string>>();
2474 
2475  for (size_t i = 0; i < parametersList.size(); ++i)
2476  {
2477  std::string lhs, rhs;
2478 
2479  try
2480  {
2481  ParseEquals(parametersList[i], lhs, rhs);
2482  }
2483  catch (...)
2484  {
2485  NEKERROR(ErrorUtil::efatal, "Parse error with command line "
2486  "option: " +
2487  parametersList[i]);
2488  }
2489 
2490  std::string lhsUpper = boost::to_upper_copy(lhs);
2491 
2492  try
2493  {
2494  m_parameters[lhsUpper] = boost::lexical_cast<NekDouble>(rhs);
2495  }
2496  catch (...)
2497  {
2498  NEKERROR(ErrorUtil::efatal, "Unable to convert string: " + rhs +
2499  "to double value.");
2500  }
2501  }
2502  }
2503 }
#define NEKERROR(type, msg)
Assert Level 0 – Fundamental assert which is used whether in FULLDEBUG, DEBUG or OPT compilation mode...
Definition: ErrorUtil.hpp:209
SolverInfoMap m_solverInfo
Solver information properties.
ParameterMap m_parameters
Parameters.
boost::program_options::variables_map m_cmdLineOptions
void ParseEquals(const std::string &line, std::string &lhs, std::string &rhs)
Parse a string in the form lhs = rhs.

References NEKERROR.

◆ CreateComm()

void Nektar::LibUtilities::SessionReader::CreateComm ( int &  argc,
char *  argv[] 
)
private

Loads the given XML document and instantiates an appropriate communication object.

Definition at line 1561 of file BasicUtils/SessionReader.cpp.

1562 {
1563  if (argc == 0)
1564  {
1565  m_comm = GetCommFactory().CreateInstance("Serial", 0, 0);
1566  }
1567  else
1568  {
1569  string vCommModule("Serial");
1570  if (GetCommFactory().ModuleExists("ParallelMPI"))
1571  {
1572  vCommModule = "ParallelMPI";
1573  }
1574  if (m_cmdLineOptions.count("cwipi") &&
1575  GetCommFactory().ModuleExists("CWIPI"))
1576  {
1577  vCommModule = "CWIPI";
1578  }
1579 
1580  m_comm = GetCommFactory().CreateInstance(vCommModule, argc, argv);
1581  }
1582 }
tBaseSharedPtr CreateInstance(tKey idKey, tParam... args)
Create an instance of the class referred to by idKey.
Definition: NekFactory.hpp:144
CommFactory & GetCommFactory()

References Nektar::LibUtilities::NekFactory< tKey, tBase, tParam >::CreateInstance(), and Nektar::LibUtilities::GetCommFactory().

◆ CreateInstance() [1/2]

static SessionReaderSharedPtr Nektar::LibUtilities::SessionReader::CreateInstance ( int  argc,
char *  argv[] 
)
inlinestatic

Creates an instance of the SessionReader class.

This function should be used by an application to instantiate the session reader. It should be called at the very beginning of the application before any other processing of command-line arguments. After instantiating the class and setting up any parallel communication, it also calls the main initialisation of the object.

Definition at line 135 of file SessionReader.h.

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), and CellMLToNektar.cellml_metadata::p.

Referenced by main(), Nektar::FieldUtils::InputXml::Process(), Nektar::FieldUtils::ProcessDisplacement::Process(), Nektar::FieldUtils::ProcessHomogeneousPlane::Process(), Nektar::FieldUtils::ProcessInterpField::Process(), Nektar::FieldUtils::ProcessInterpPoints::Process(), SessionReader_CreateInstance(), Nektar::SolverUtils::Driver::v_InitObject(), and Nektar::VortexWaveInteraction::VortexWaveInteraction().

◆ CreateInstance() [2/2]

static SessionReaderSharedPtr Nektar::LibUtilities::SessionReader::CreateInstance ( int  argc,
char *  argv[],
std::vector< std::string > &  pFilenames,
const CommSharedPtr pComm = CommSharedPtr() 
)
inlinestatic

Creates an instance of the SessionReader class initialised using a separate list of XML documents.

This function should be used by an application to instantiate the session reader. It may be called after processing of command-line arguments. After instantiating the class and setting up any parallel communication, it also calls the main initialisation of the object.

Definition at line 154 of file SessionReader.h.

157  {
160  argc, argv, pFilenames, pComm);
161  return p;
162  }

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), and CellMLToNektar.cellml_metadata::p.

◆ DefinesCmdLineArgument()

bool Nektar::LibUtilities::SessionReader::DefinesCmdLineArgument ( const std::string &  pName) const

Checks if a specified cmdline argument has been given.

Definition at line 1395 of file BasicUtils/SessionReader.cpp.

1396 {
1397  return (m_cmdLineOptions.find(pName) != m_cmdLineOptions.end());
1398 }

◆ DefinesElement()

bool Nektar::LibUtilities::SessionReader::DefinesElement ( const std::string &  pPath) const

Tests if a specified element is defined in the XML document.

Definition at line 648 of file BasicUtils/SessionReader.cpp.

649 {
650  std::string vPath = boost::to_upper_copy(pPath);
651  std::vector<std::string> st;
652  boost::split(st, vPath, boost::is_any_of("\\/ "));
653  ASSERTL0(st.size() > 0, "No path given in XML element request.");
654 
655  TiXmlElement *vReturn = m_xmlDoc->FirstChildElement(st[0].c_str());
656  ASSERTL0(vReturn,
657  std::string("Cannot find element '") + st[0] + std::string("'."));
658  for (int i = 1; i < st.size(); ++i)
659  {
660  vReturn = vReturn->FirstChildElement(st[i].c_str());
661  if (!vReturn)
662  return false;
663  }
664  return true;
665 }

References ASSERTL0.

◆ DefinesFunction() [1/2]

bool Nektar::LibUtilities::SessionReader::DefinesFunction ( const std::string &  name) const

Checks if a specified function is defined in the XML document.

Definition at line 1143 of file BasicUtils/SessionReader.cpp.

1144 {
1145  std::string vName = boost::to_upper_copy(pName);
1146  return m_functions.find(vName) != m_functions.end();
1147 }
FunctionMap m_functions
Functions.

◆ DefinesFunction() [2/2]

bool Nektar::LibUtilities::SessionReader::DefinesFunction ( const std::string &  name,
const std::string &  variable,
const int  pDomain = 0 
) const

Checks if a specified function has a given variable defined.

Definition at line 1152 of file BasicUtils/SessionReader.cpp.

1155 {
1156  std::string vName = boost::to_upper_copy(pName);
1157 
1158  // Check function exists
1159  auto it1 = m_functions.find(vName);
1160  if (it1 != m_functions.end())
1161  {
1162  pair<std::string, int> key(pVariable, pDomain);
1163  pair<std::string, int> defkey("*", pDomain);
1164  bool varExists = it1->second.find(key) != it1->second.end() ||
1165  it1->second.find(defkey) != it1->second.end();
1166  return varExists;
1167  }
1168  return false;
1169 }

◆ DefinesGeometricInfo()

bool Nektar::LibUtilities::SessionReader::DefinesGeometricInfo ( const std::string &  name) const

Checks if a geometric info property is defined.

Definition at line 1027 of file BasicUtils/SessionReader.cpp.

1028 {
1029  std::string vName = boost::to_upper_copy(pName);
1030  return m_geometricInfo.find(vName) != m_geometricInfo.end();
1031 }
GeometricInfoMap m_geometricInfo
Geometric information properties.

◆ DefinesGlobalSysSolnInfo()

bool Nektar::LibUtilities::SessionReader::DefinesGlobalSysSolnInfo ( const std::string &  variable,
const std::string &  property 
) const

Definition at line 932 of file BasicUtils/SessionReader.cpp.

934 {
935  auto iter = GetGloSysSolnList().find(pVariable);
936  if (iter == GetGloSysSolnList().end())
937  {
938  return false;
939  }
940 
941  std::string vProperty = boost::to_upper_copy(pProperty);
942 
943  auto iter1 = iter->second.find(vProperty);
944  if (iter1 == iter->second.end())
945  {
946  return false;
947  }
948 
949  return true;
950 }
static GloSysSolnInfoList & GetGloSysSolnList()
GlobalSysSoln Info map.

◆ DefinesParameter()

bool Nektar::LibUtilities::SessionReader::DefinesParameter ( const std::string &  name) const

Checks if a parameter is specified in the XML document.

Definition at line 728 of file BasicUtils/SessionReader.cpp.

729 {
730  std::string vName = boost::to_upper_copy(pName);
731  return m_parameters.find(vName) != m_parameters.end();
732 }

Referenced by export_SessionReader().

◆ DefinesSolverInfo()

bool Nektar::LibUtilities::SessionReader::DefinesSolverInfo ( const std::string &  name) const

Checks if a solver info property is specified.

Definition at line 837 of file BasicUtils/SessionReader.cpp.

838 {
839  std::string vName = boost::to_upper_copy(pName);
840  auto infoIter = m_solverInfo.find(vName);
841  return (infoIter != m_solverInfo.end());
842 }

Referenced by GetSolverInfoAsEnum().

◆ DefinesTag()

bool Nektar::LibUtilities::SessionReader::DefinesTag ( const std::string &  pName) const

Checks if a specified tag is defined.

Definition at line 1358 of file BasicUtils/SessionReader.cpp.

1359 {
1360  std::string vName = boost::to_upper_copy(pName);
1361  return m_tags.find(vName) != m_tags.end();
1362 }

◆ DefinesTimeIntScheme()

bool Nektar::LibUtilities::SessionReader::DefinesTimeIntScheme ( ) const

Returns true if the TIMEINTEGRATIONSCHEME section is defined in the session file.

Definition at line 976 of file BasicUtils/SessionReader.cpp.

977 {
978  return m_timeIntScheme.method != "";
979 }
TimeIntScheme m_timeIntScheme
Time integration scheme information.

◆ Finalise()

void Nektar::LibUtilities::SessionReader::Finalise ( )

Finalises the session.

This routine finalises any parallel communication.

Note
This routine should be called at the very end of a users application.

Definition at line 720 of file BasicUtils/SessionReader.cpp.

721 {
722  m_comm->Finalise();
723 }

Referenced by export_SessionReader().

◆ GetCmdLineArgMap()

CmdLineArgMap & Nektar::LibUtilities::SessionReader::GetCmdLineArgMap ( )
staticprivate

CmdLine argument map.

Lists the possible command-line argument which can be specified for this executable.

This list is populated through the RegisterCmdLineArgument static member function which is called statically from various classes to register command-line arguments they need.

Definition at line 171 of file BasicUtils/SessionReader.cpp.

172 {
173  static CmdLineArgMap cmdLineArguments;
174  return cmdLineArguments;
175 }
std::map< std::string, CmdLineArg > CmdLineArgMap
Definition: SessionReader.h:74

Referenced by RegisterCmdLineArgument(), and RegisterCmdLineFlag().

◆ GetCmdLineArgument()

template<typename T >
T Nektar::LibUtilities::SessionReader::GetCmdLineArgument ( const std::string &  pName) const
inline

Retrieves a command-line argument value.

Definition at line 361 of file SessionReader.h.

362  {
363  return m_cmdLineOptions.find(pName)->second.as<T>();
364  }

References m_cmdLineOptions.

◆ GetComm()

CommSharedPtr Nektar::LibUtilities::SessionReader::GetComm ( )

Returns the communication object.

Definition at line 704 of file BasicUtils/SessionReader.cpp.

705 {
706  return m_comm;
707 }

Referenced by export_SessionReader().

◆ GetDocument()

TiXmlDocument & Nektar::LibUtilities::SessionReader::GetDocument ( )

Provides direct access to the TiXmlDocument object.

Definition at line 598 of file BasicUtils/SessionReader.cpp.

599 {
600  ASSERTL1(m_xmlDoc, "XML Document not defined.");
601  return *m_xmlDoc;
602 }
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode....
Definition: ErrorUtil.hpp:249

References ASSERTL1.

◆ GetElement()

TiXmlElement * Nektar::LibUtilities::SessionReader::GetElement ( const std::string &  pPath)

Provides direct access to the TiXmlElement specified.

The single parameter specifies a path to the requested element in a similar format to the filesystem path. Given the following XML:

<NEKTAR>
<CONDITIONS>
<PARAMETERS>
...
</PARAMETERS>
</CONDITIONS>
</NEKTAR>

the PARAMETERS element would be retrieved by requesting the path:

Nektar/Conditions/Parameters
The above copyright notice and this permission notice shall be included.
Definition: CoupledSolver.h:1
Note
Paths are case-insensitive.
Parameters
pPathPath to requested element.
Returns
Direct pointer to requested XML Element.

Definition at line 626 of file BasicUtils/SessionReader.cpp.

627 {
628  std::string vPath = boost::to_upper_copy(pPath);
629  std::vector<std::string> st;
630  boost::split(st, vPath, boost::is_any_of("\\/ "));
631  ASSERTL0(st.size() > 0, "No path given in XML element request.");
632 
633  TiXmlElement *vReturn = m_xmlDoc->FirstChildElement(st[0].c_str());
634  ASSERTL0(vReturn,
635  std::string("Cannot find element '") + st[0] + std::string("'."));
636  for (int i = 1; i < st.size(); ++i)
637  {
638  vReturn = vReturn->FirstChildElement(st[i].c_str());
639  ASSERTL0(vReturn, std::string("Cannot find element '") + st[i] +
640  std::string("'."));
641  }
642  return vReturn;
643 }

References ASSERTL0.

◆ GetFilenames()

const std::vector< std::string > & Nektar::LibUtilities::SessionReader::GetFilenames ( ) const

Returns the filename of the loaded XML document.

Definition at line 670 of file BasicUtils/SessionReader.cpp.

671 {
672  return m_filenames;
673 }

◆ GetFilters()

const FilterMap & Nektar::LibUtilities::SessionReader::GetFilters ( ) const

Definition at line 1387 of file BasicUtils/SessionReader.cpp.

1388 {
1389  return m_filters;
1390 }

◆ GetFunction() [1/2]

EquationSharedPtr Nektar::LibUtilities::SessionReader::GetFunction ( const std::string &  name,
const std::string &  variable,
const int  pDomain = 0 
) const

Returns an EquationSharedPtr to a given function variable.

Definition at line 1174 of file BasicUtils/SessionReader.cpp.

1177 {
1178  std::string vName = boost::to_upper_copy(pName);
1179  auto it1 = m_functions.find(vName);
1180 
1181  ASSERTL0(it1 != m_functions.end(),
1182  std::string("No such function '") + pName +
1183  std::string("' has been defined in the session file."));
1184 
1185  // Check for specific and wildcard definitions
1186  pair<std::string, int> key(pVariable, pDomain);
1187  pair<std::string, int> defkey("*", pDomain);
1188 
1189  auto it2 = it1->second.find(key);
1190  auto it3 = it1->second.find(defkey);
1191  bool specific = it2 != it1->second.end();
1192  bool wildcard = it3 != it1->second.end();
1193 
1194  // Check function is defined somewhere
1195  ASSERTL0(specific || wildcard,
1196  "No such variable " + pVariable + " in domain " +
1197  boost::lexical_cast<string>(pDomain) +
1198  " defined for function " + pName + " in session file.");
1199 
1200  // If not specific, must be wildcard
1201  if (!specific)
1202  {
1203  it2 = it3;
1204  }
1205 
1206  ASSERTL0((it2->second.m_type == eFunctionTypeExpression),
1207  std::string("Function is defined by a file."));
1208  return it2->second.m_expression;
1209 }

References ASSERTL0, and Nektar::LibUtilities::eFunctionTypeExpression.

◆ GetFunction() [2/2]

EquationSharedPtr Nektar::LibUtilities::SessionReader::GetFunction ( const std::string &  name,
const unsigned int &  var,
const int  pDomain = 0 
) const

Returns an EquationSharedPtr to a given function variable index.

Definition at line 1214 of file BasicUtils/SessionReader.cpp.

1217 {
1218  ASSERTL0(pVar < m_variables.size(), "Variable index out of range.");
1219  return GetFunction(pName, m_variables[pVar], pDomain);
1220 }
VariableList m_variables
Variables.
EquationSharedPtr GetFunction(const std::string &name, const std::string &variable, const int pDomain=0) const
Returns an EquationSharedPtr to a given function variable.

References ASSERTL0.

◆ GetFunctionFilename() [1/2]

std::string Nektar::LibUtilities::SessionReader::GetFunctionFilename ( const std::string &  name,
const std::string &  variable,
const int  pDomain = 0 
) const

Returns the filename to be loaded for a given variable.

Definition at line 1273 of file BasicUtils/SessionReader.cpp.

1276 {
1277  std::string vName = boost::to_upper_copy(pName);
1278  auto it1 = m_functions.find(vName);
1279 
1280  ASSERTL0(it1 != m_functions.end(),
1281  std::string("Function '") + pName + std::string("' not found."));
1282 
1283  // Check for specific and wildcard definitions
1284  pair<std::string, int> key(pVariable, pDomain);
1285  pair<std::string, int> defkey("*", pDomain);
1286 
1287  auto it2 = it1->second.find(key);
1288  auto it3 = it1->second.find(defkey);
1289  bool specific = it2 != it1->second.end();
1290  bool wildcard = it3 != it1->second.end();
1291 
1292  // Check function is defined somewhere
1293  ASSERTL0(specific || wildcard,
1294  "No such variable " + pVariable + " in domain " +
1295  boost::lexical_cast<string>(pDomain) +
1296  " defined for function " + pName + " in session file.");
1297 
1298  // If not specific, must be wildcard
1299  if (!specific)
1300  {
1301  it2 = it3;
1302  }
1303 
1304  return it2->second.m_filename;
1305 }

References ASSERTL0.

◆ GetFunctionFilename() [2/2]

std::string Nektar::LibUtilities::SessionReader::GetFunctionFilename ( const std::string &  name,
const unsigned int &  var,
const int  pDomain = 0 
) const

Returns the filename to be loaded for a given variable index.

Definition at line 1310 of file BasicUtils/SessionReader.cpp.

1313 {
1314  ASSERTL0(pVar < m_variables.size(), "Variable index out of range.");
1315  return GetFunctionFilename(pName, m_variables[pVar], pDomain);
1316 }
std::string GetFunctionFilename(const std::string &name, const std::string &variable, const int pDomain=0) const
Returns the filename to be loaded for a given variable.

References ASSERTL0.

◆ GetFunctionFilenameVariable()

std::string Nektar::LibUtilities::SessionReader::GetFunctionFilenameVariable ( const std::string &  name,
const std::string &  variable,
const int  pDomain = 0 
) const

Returns the filename variable to be loaded for a given variable index.

Definition at line 1321 of file BasicUtils/SessionReader.cpp.

1324 {
1325  std::string vName = boost::to_upper_copy(pName);
1326  auto it1 = m_functions.find(vName);
1327 
1328  ASSERTL0(it1 != m_functions.end(),
1329  std::string("Function '") + pName + std::string("' not found."));
1330 
1331  // Check for specific and wildcard definitions
1332  pair<std::string, int> key(pVariable, pDomain);
1333  pair<std::string, int> defkey("*", pDomain);
1334 
1335  auto it2 = it1->second.find(key);
1336  auto it3 = it1->second.find(defkey);
1337  bool specific = it2 != it1->second.end();
1338  bool wildcard = it3 != it1->second.end();
1339 
1340  // Check function is defined somewhere
1341  ASSERTL0(specific || wildcard,
1342  "No such variable " + pVariable + " in domain " +
1343  boost::lexical_cast<string>(pDomain) +
1344  " defined for function " + pName + " in session file.");
1345 
1346  // If not specific, must be wildcard
1347  if (!specific)
1348  {
1349  it2 = it3;
1350  }
1351 
1352  return it2->second.m_fileVariable;
1353 }

References ASSERTL0.

◆ GetFunctionType() [1/2]

enum FunctionType Nektar::LibUtilities::SessionReader::GetFunctionType ( const std::string &  name,
const std::string &  variable,
const int  pDomain = 0 
) const

Returns the type of a given function variable.

Definition at line 1214 of file BasicUtils/SessionReader.cpp.

1228 {
1229  std::string vName = boost::to_upper_copy(pName);
1230  auto it1 = m_functions.find(vName);
1231 
1232  ASSERTL0(it1 != m_functions.end(),
1233  std::string("Function '") + pName + std::string("' not found."));
1234 
1235  // Check for specific and wildcard definitions
1236  pair<std::string, int> key(pVariable, pDomain);
1237  pair<std::string, int> defkey("*", pDomain);
1238 
1239  auto it2 = it1->second.find(key);
1240  auto it3 = it1->second.find(defkey);
1241  bool specific = it2 != it1->second.end();
1242  bool wildcard = it3 != it1->second.end();
1243 
1244  // Check function is defined somewhere
1245  ASSERTL0(specific || wildcard,
1246  "No such variable " + pVariable + " in domain " +
1247  boost::lexical_cast<string>(pDomain) +
1248  " defined for function " + pName + " in session file.");
1249 
1250  // If not specific, must be wildcard
1251  if (!specific)
1252  {
1253  it2 = it3;
1254  }
1255 
1256  return it2->second.m_type;
1257 }

◆ GetFunctionType() [2/2]

enum FunctionType Nektar::LibUtilities::SessionReader::GetFunctionType ( const std::string &  pName,
const unsigned int &  pVar,
const int  pDomain = 0 
) const

Returns the type of a given function variable index.

Definition at line 1214 of file BasicUtils/SessionReader.cpp.

1265 {
1266  ASSERTL0(pVar < m_variables.size(), "Variable index out of range.");
1267  return GetFunctionType(pName, m_variables[pVar], pDomain);
1268 }
enum FunctionType GetFunctionType(const std::string &name, const std::string &variable, const int pDomain=0) const
Returns the type of a given function variable.

◆ GetGeometryType()

std::string Nektar::LibUtilities::SessionReader::GetGeometryType ( ) const

Definition at line 990 of file BasicUtils/SessionReader.cpp.

991 {
992  TiXmlElement *xmlGeom =
993  m_xmlDoc->FirstChildElement("NEKTAR")->FirstChildElement("GEOMETRY");
994  ASSERTL1(xmlGeom, "Failed to find a GEOMETRY section in m_xmlDoc");
995 
996  TiXmlAttribute *attr = xmlGeom->FirstAttribute();
997  while (attr)
998  {
999  std::string attrName(attr->Name());
1000  if (attrName == "HDF5FILE")
1001  {
1002  // there is a file pointer, therefore is HDF5
1003  return "HDF5";
1004  }
1005  // Get the next attribute.
1006  attr = attr->Next();
1007  }
1008 
1009  // Check the VERTEX block. If this is compressed, assume the file is
1010  // compressed, otherwise assume uncompressed.
1011  TiXmlElement *element = xmlGeom->FirstChildElement("VERTEX");
1012  string IsCompressed;
1013  element->QueryStringAttribute("COMPRESSED", &IsCompressed);
1014 
1015  if (IsCompressed.size() > 0)
1016  {
1017  return "XmlCompressed";
1018  }
1019 
1020  // no file pointer or compressed, just standard xml
1021  return "Xml";
1022 }

References ASSERTL1.

◆ GetGlobalSysSolnInfo()

const std::string & Nektar::LibUtilities::SessionReader::GetGlobalSysSolnInfo ( const std::string &  variable,
const std::string &  property 
) const

Definition at line 955 of file BasicUtils/SessionReader.cpp.

957 {
958  auto iter = GetGloSysSolnList().find(pVariable);
959  ASSERTL0(iter != GetGloSysSolnList().end(),
960  "Failed to find variable in GlobalSysSolnInfoList");
961 
962  std::string vProperty = boost::to_upper_copy(pProperty);
963  auto iter1 = iter->second.find(vProperty);
964 
965  ASSERTL0(iter1 != iter->second.end(),
966  "Failed to find property: " + vProperty +
967  " in GlobalSysSolnInfoList");
968 
969  return iter1->second;
970 }

References ASSERTL0.

◆ GetGloSysSolnList()

GloSysSolnInfoList & Nektar::LibUtilities::SessionReader::GetGloSysSolnList ( )
staticprivate

GlobalSysSoln Info map.

List of values for GlobalSysSoln parameters to be used to override details given in SolverInfo

This list is populated by ReadGlobalSysSoln if the GLOBALSYSSOLNINFO section is defined in the input file. This List allows for details to define for the Global Sys solver for each variable.

Definition at line 157 of file BasicUtils/SessionReader.cpp.

158 {
159  static GloSysSolnInfoList gloSysSolnInfoList;
160  return gloSysSolnInfoList;
161 }
std::map< std::string, GloSysInfoMap > GloSysSolnInfoList
Definition: SessionReader.h:80

◆ GetInterpreter()

InterpreterSharedPtr Nektar::LibUtilities::SessionReader::GetInterpreter ( )

Returns the instance of the Interpreter specific to this session.

Definition at line 2529 of file BasicUtils/SessionReader.cpp.

2530 {
2531  return m_interpreter;
2532 }

◆ GetParameter()

const NekDouble & Nektar::LibUtilities::SessionReader::GetParameter ( const std::string &  pName) const

Returns the value of the specified parameter.

If the parameter is not defined, termination occurs. Therefore, the parameters existence should be tested for using DefinesParameter before calling this function.

Parameters
pNameThe name of a floating-point parameter.
Returns
The value of the floating-point parameter.

Definition at line 742 of file BasicUtils/SessionReader.cpp.

743 {
744  std::string vName = boost::to_upper_copy(pName);
745  auto paramIter = m_parameters.find(vName);
746 
747  ASSERTL0(paramIter != m_parameters.end(),
748  "Unable to find requested parameter: " + pName);
749 
750  return paramIter->second;
751 }

References ASSERTL0.

Referenced by export_SessionReader().

◆ GetSessionName()

const std::string & Nektar::LibUtilities::SessionReader::GetSessionName ( ) const

Returns the session name of the loaded XML document.

Definition at line 678 of file BasicUtils/SessionReader.cpp.

679 {
680  return m_sessionName;
681 }

Referenced by export_SessionReader().

◆ GetSessionNameRank()

const std::string Nektar::LibUtilities::SessionReader::GetSessionNameRank ( ) const

Returns the session name with process rank.

Output is of the form [sessionName]_P[idx] where idx is the rank of the process.

Definition at line 687 of file BasicUtils/SessionReader.cpp.

688 {
689  std::string dirname = m_sessionName + "_xml";
690  fs::path pdirname(dirname);
691 
692  std::string vFilename =
693  "P" + boost::lexical_cast<std::string>(m_comm->GetRowComm()->GetRank());
694  fs::path pFilename(vFilename);
695 
696  fs::path fullpath = pdirname / pFilename;
697 
698  return PortablePath(fullpath);
699 }
std::string PortablePath(const boost::filesystem::path &path)
create portable path on different platforms for boost::filesystem path
Definition: FileSystem.cpp:41

References Nektar::LibUtilities::PortablePath().

◆ GetSharedFilesystem()

bool Nektar::LibUtilities::SessionReader::GetSharedFilesystem ( )

Returns if file system shared.

Definition at line 709 of file BasicUtils/SessionReader.cpp.

710 {
711  return m_sharedFilesystem;
712 }
bool m_sharedFilesystem
Running on a shared filesystem.

◆ GetSolverInfo()

const std::string & Nektar::LibUtilities::SessionReader::GetSolverInfo ( const std::string &  pProperty) const

Returns the value of the specified solver info property.

Definition at line 847 of file BasicUtils/SessionReader.cpp.

849 {
850  std::string vProperty = boost::to_upper_copy(pProperty);
851  auto iter = m_solverInfo.find(vProperty);
852 
853  ASSERTL1(iter != m_solverInfo.end(),
854  "Unable to find requested property: " + pProperty);
855 
856  return iter->second;
857 }

References ASSERTL1.

Referenced by GetSolverInfoAsEnum().

◆ GetSolverInfoAsEnum()

template<typename T >
const T Nektar::LibUtilities::SessionReader::GetSolverInfoAsEnum ( const std::string &  pName) const
inline

Returns the value of the specified solver info property as enum.

Definition at line 503 of file SessionReader.h.

505 {
506  std::string vName = boost::to_upper_copy(pName);
508  "Solver info '" + pName + "' not defined.");
509 
510  std::string vValue = GetSolverInfo(vName);
511  auto x = GetSolverInfoEnums().find(vName);
512  ASSERTL0(x != GetSolverInfoEnums().end(),
513  "Enum for SolverInfo property '" + pName + "' not found.");
514 
515  auto y = x->second.find(vValue);
516  ASSERTL0(y != x->second.end(),
517  "Value of SolverInfo property '" + pName + "' is invalid.");
518 
519  return T(y->second);
520 }
bool DefinesSolverInfo(const std::string &name) const
Checks if a solver info property is specified.
static EnumMapList & GetSolverInfoEnums()
String to enumeration map for Solver Info parameters.
const std::string & GetSolverInfo(const std::string &pProperty) const
Returns the value of the specified solver info property.

References ASSERTL0, DefinesSolverInfo(), GetSolverInfo(), and GetSolverInfoEnums().

◆ GetSolverInfoDefaults()

SolverInfoMap & Nektar::LibUtilities::SessionReader::GetSolverInfoDefaults ( )
staticprivate

Default solver info options.

List of default values for solver information parameters to be used in the case of them not being provided.

This list is populated through the RegisterDefaultSolverInfo static member variable which is called statically from various classes to register the default value for a given parameter.

Definition at line 142 of file BasicUtils/SessionReader.cpp.

143 {
144  static SolverInfoMap solverInfoMap;
145  return solverInfoMap;
146 }
std::map< std::string, std::string > SolverInfoMap
Definition: SessionReader.h:58

Referenced by RegisterDefaultSolverInfo().

◆ GetSolverInfoEnums()

EnumMapList & Nektar::LibUtilities::SessionReader::GetSolverInfoEnums ( )
staticprivate

String to enumeration map for Solver Info parameters.

This map of maps stores the list of valid string values for a number of solver information parameters. The top level map connects different parameter names to their list of possible values. The list of possible values is also a map, mapping a valid string to a corresponding enum value.

This list is populated through the RegisterEnumValue static member function which is called statically from various classes to register the valid values for solver info parameters associated with them. The map is therefore fully populated before the SessionReader class is instantiated and a file is read in and parsed.

Definition at line 128 of file BasicUtils/SessionReader.cpp.

129 {
130  static EnumMapList solverInfoEnums;
131  return solverInfoEnums;
132 }
std::map< std::string, EnumMap > EnumMapList
Definition: SessionReader.h:77

Referenced by GetSolverInfoAsEnum(), GetValueAsEnum(), and RegisterEnumValue().

◆ GetTag()

const std::string & Nektar::LibUtilities::SessionReader::GetTag ( const std::string &  pName) const

Returns the value of a specified tag.

Definition at line 1376 of file BasicUtils/SessionReader.cpp.

1377 {
1378  std::string vName = boost::to_upper_copy(pName);
1379  auto vTagIterator = m_tags.find(vName);
1380  ASSERTL0(vTagIterator != m_tags.end(), "Requested tag does not exist.");
1381  return vTagIterator->second;
1382 }

References ASSERTL0.

◆ GetTimeIntScheme()

const TimeIntScheme & Nektar::LibUtilities::SessionReader::GetTimeIntScheme ( ) const

Returns the time integration scheme structure m_timeIntScheme from the session file.

Definition at line 985 of file BasicUtils/SessionReader.cpp.

986 {
987  return m_timeIntScheme;
988 }

◆ GetUpdateOptFile()

bool Nektar::LibUtilities::SessionReader::GetUpdateOptFile ( ) const
inline

Get bool to update optimisation file.

Definition at line 380 of file SessionReader.h.

381  {
382  return m_updateOptFile;
383  }
bool m_updateOptFile
Update optimisation file.

References m_updateOptFile.

◆ GetValueAsEnum()

template<typename T >
const T Nektar::LibUtilities::SessionReader::GetValueAsEnum ( const std::string &  pName,
const std::string &  vValue 
) const
inline

Returns the value of the specified property and value as enum.

Definition at line 526 of file SessionReader.h.

528 {
529  std::string vName = boost::to_upper_copy(pName);
530 
531  auto x = GetSolverInfoEnums().find(vName);
532  ASSERTL0(x != GetSolverInfoEnums().end(),
533  "Enum for property '" + pName + "' not found.");
534 
535  auto y = x->second.find(pValue);
536  ASSERTL0(y != x->second.end(),
537  "Value of property '" + pValue + "' is invalid.");
538  return T(y->second);
539 }

References ASSERTL0, and GetSolverInfoEnums().

◆ GetVariable()

const std::string & Nektar::LibUtilities::SessionReader::GetVariable ( const unsigned int &  idx) const

Returns the name of the variable specified by the given index.

Definition at line 1117 of file BasicUtils/SessionReader.cpp.

1118 {
1119  ASSERTL0(idx < m_variables.size(), "Variable index out of range.");
1120  return m_variables[idx];
1121 }

References ASSERTL0.

Referenced by export_SessionReader().

◆ GetVariables()

std::vector< std::string > Nektar::LibUtilities::SessionReader::GetVariables ( ) const

Returns the names of all variables.

Definition at line 1135 of file BasicUtils/SessionReader.cpp.

1136 {
1137  return m_variables;
1138 }

◆ InitSession()

void Nektar::LibUtilities::SessionReader::InitSession ( const std::vector< std::string > &  filenames = std::vector<std::string>())

Performs the main initialisation of the object. The XML file provided on the command-line is loaded and any mesh partitioning is done. The resulting process-specific XML file (containing the process's geometry partition) is then reloaded and parsed.

Definition at line 273 of file BasicUtils/SessionReader.cpp.

274 {
275  // Re-load filenames for session if required.
276  if (filenames.size() > 0)
277  {
278  m_filenames = filenames;
279  }
280 
281  // check specified opt file
282  std::string optfile;
283  int exists;
284 
285  if (DefinesCmdLineArgument("useoptfile"))
286  {
287  optfile = m_cmdLineOptions.find("useoptfile")->second.as<std::string>();
288  exists = (bool)boost::filesystem::exists(optfile.c_str());
289  ASSERTL0(exists, "A valid .opt file was not specified "
290  "with the --useoptfile command line option");
291 
292  m_filenames.push_back(optfile);
293 
294  // put opt file at beginning
295  std::rotate(m_filenames.rbegin(), m_filenames.rbegin() + 1,
296  m_filenames.rend());
297  }
298  else // check for writeoptfile
299  {
300  // check for opt file
301  optfile = m_sessionName + ".opt";
302  exists = (bool)boost::filesystem::exists(optfile.c_str());
303  if (exists)
304  {
305  m_filenames.push_back(optfile);
306  // rotate order so opt file can be overwritten by
307  // direct choice in xml file
308  std::rotate(m_filenames.rbegin(), m_filenames.rbegin() + 1,
309  m_filenames.rend());
310  }
311  else
312  {
313  m_updateOptFile = true;
314  }
315  }
316 
317  // Merge document if required.
318  if (m_xmlDoc)
319  {
320  delete m_xmlDoc;
321  }
322 
324 
325  // Parse the XML data in #m_xmlDoc
326  ParseDocument();
327 
328  // Override SOLVERINFO and parameters with any specified on the
329  // command line.
330  CmdLineOverride();
331 
332  // Verify SOLVERINFO values
334 
335  // In verbose mode, print out parameters and solver info sections
336  if (m_verbose && m_comm)
337  {
338  if (m_comm->TreatAsRankZero() && m_parameters.size() > 0)
339  {
340  cout << "Parameters:" << endl;
341  for (auto &x : m_parameters)
342  {
343  cout << "\t" << x.first << " = " << x.second << endl;
344  }
345  cout << endl;
346  }
347 
348  if (m_comm->TreatAsRankZero() && m_solverInfo.size() > 0)
349  {
350  cout << "Solver Info:" << endl;
351  for (auto &x : m_solverInfo)
352  {
353  cout << "\t" << x.first << " = " << x.second << endl;
354  }
355  cout << endl;
356  }
357  }
358 }
TiXmlDocument * MergeDoc(const std::vector< std::string > &pFilenames) const
Creates an XML document from a list of input files.
void CmdLineOverride()
Enforce parameters from command line arguments.
void ParseDocument()
Loads and parses the specified file.
void VerifySolverInfo()
Check values of solver info options are valid.
bool DefinesCmdLineArgument(const std::string &pName) const
Checks if a specified cmdline argument has been given.

References ASSERTL0.

◆ LoadDoc()

void Nektar::LibUtilities::SessionReader::LoadDoc ( const std::string &  pFilename,
TiXmlDocument *  pDoc 
) const
private

Loads an xml file into a tinyxml doc and decompresses if needed.

Definition at line 1414 of file BasicUtils/SessionReader.cpp.

1416 {
1417  if (pFilename.size() > 3 &&
1418  pFilename.substr(pFilename.size() - 3, 3) == ".gz")
1419  {
1420  ifstream file(pFilename.c_str(), ios_base::in | ios_base::binary);
1421  ASSERTL0(file.good(), "Unable to open file: " + pFilename);
1422  stringstream ss;
1423  io::filtering_streambuf<io::input> in;
1424  in.push(io::gzip_decompressor());
1425  in.push(file);
1426  try
1427  {
1428  io::copy(in, ss);
1429  ss >> (*pDoc);
1430  }
1431  catch (io::gzip_error &)
1432  {
1434  "Error: File '" + pFilename + "' is corrupt.");
1435  }
1436  }
1437  else if (pFilename.size() > 4 &&
1438  pFilename.substr(pFilename.size() - 4, 4) == "_xml")
1439  {
1440  fs::path pdirname(pFilename);
1441  boost::format pad("P%1$07d.xml");
1442  pad % m_comm->GetRank();
1443  fs::path pRankFilename(pad.str());
1444  fs::path fullpath = pdirname / pRankFilename;
1445 
1446  ifstream file(PortablePath(fullpath).c_str());
1447  ASSERTL0(file.good(), "Unable to open file: " + fullpath.string());
1448  file >> (*pDoc);
1449  }
1450  else
1451  {
1452  ifstream file(pFilename.c_str());
1453  ASSERTL0(file.good(), "Unable to open file: " + pFilename);
1454  file >> (*pDoc);
1455  }
1456 }
def copy(self)
Definition: pycml.py:2663

References ASSERTL0, CellMLToNektar.pycml::copy(), CellMLToNektar.pycml::format, NEKERROR, and Nektar::LibUtilities::PortablePath().

◆ LoadGeometricInfo() [1/3]

void Nektar::LibUtilities::SessionReader::LoadGeometricInfo ( const std::string &  name,
bool &  var,
const bool &  def = false 
) const

Checks for and loads a geometric info boolean property.

Definition at line 1055 of file BasicUtils/SessionReader.cpp.

1057 {
1058  std::string vName = boost::to_upper_copy(pName);
1059  auto iter = m_geometricInfo.find(vName);
1060  if (iter != m_geometricInfo.end())
1061  {
1062  if (iter->second == "TRUE")
1063  {
1064  pVar = true;
1065  }
1066  else
1067  {
1068  pVar = false;
1069  }
1070  }
1071  else
1072  {
1073  pVar = pDefault;
1074  }
1075 }

◆ LoadGeometricInfo() [2/3]

void Nektar::LibUtilities::SessionReader::LoadGeometricInfo ( const std::string &  name,
NekDouble var,
const NekDouble def = 0.0 
) const

Checks for and loads a geometric info double-precision property.

Definition at line 1080 of file BasicUtils/SessionReader.cpp.

1082 {
1083  std::string vName = boost::to_upper_copy(pName);
1084  auto iter = m_geometricInfo.find(vName);
1085  if (iter != m_geometricInfo.end())
1086  {
1087  pVar = std::atoi(iter->second.c_str());
1088  }
1089  else
1090  {
1091  pVar = pDefault;
1092  }
1093 }

◆ LoadGeometricInfo() [3/3]

void Nektar::LibUtilities::SessionReader::LoadGeometricInfo ( const std::string &  name,
std::string &  var,
const std::string &  def = "" 
) const

Checks for and load a geometric info string property.

Definition at line 1036 of file BasicUtils/SessionReader.cpp.

1039 {
1040  std::string vName = boost::to_upper_copy(pName);
1041  auto iter = m_geometricInfo.find(vName);
1042  if (iter != m_geometricInfo.end())
1043  {
1044  pVar = iter->second;
1045  }
1046  else
1047  {
1048  pVar = pDefault;
1049  }
1050 }

◆ LoadParameter() [1/4]

void Nektar::LibUtilities::SessionReader::LoadParameter ( const std::string &  name,
int &  var 
) const

Load an integer parameter.

Definition at line 756 of file BasicUtils/SessionReader.cpp.

757 {
758  std::string vName = boost::to_upper_copy(pName);
759  auto paramIter = m_parameters.find(vName);
760  ASSERTL0(paramIter != m_parameters.end(),
761  "Required parameter '" + pName + "' not specified in session.");
762  NekDouble param = round(paramIter->second);
763  pVar = checked_cast<int>(param);
764 }
double NekDouble

References ASSERTL0.

◆ LoadParameter() [2/4]

void Nektar::LibUtilities::SessionReader::LoadParameter ( const std::string &  name,
int &  var,
const int &  def 
) const

Check for and load an integer parameter.

Definition at line 769 of file BasicUtils/SessionReader.cpp.

771 {
772  std::string vName = boost::to_upper_copy(pName);
773  auto paramIter = m_parameters.find(vName);
774  if (paramIter != m_parameters.end())
775  {
776  NekDouble param = round(paramIter->second);
777  pVar = checked_cast<int>(param);
778  }
779  else
780  {
781  pVar = pDefault;
782  }
783 }

◆ LoadParameter() [3/4]

void Nektar::LibUtilities::SessionReader::LoadParameter ( const std::string &  name,
NekDouble var 
) const

Load a double precision parameter.

Definition at line 788 of file BasicUtils/SessionReader.cpp.

790 {
791  std::string vName = boost::to_upper_copy(pName);
792  auto paramIter = m_parameters.find(vName);
793  ASSERTL0(paramIter != m_parameters.end(),
794  "Required parameter '" + pName + "' not specified in session.");
795  pVar = paramIter->second;
796 }

References ASSERTL0.

◆ LoadParameter() [4/4]

void Nektar::LibUtilities::SessionReader::LoadParameter ( const std::string &  name,
NekDouble var,
const NekDouble def 
) const

Check for and load a double-precision parameter.

Definition at line 801 of file BasicUtils/SessionReader.cpp.

803 {
804  std::string vName = boost::to_upper_copy(pName);
805  auto paramIter = m_parameters.find(vName);
806  if (paramIter != m_parameters.end())
807  {
808  pVar = paramIter->second;
809  }
810  else
811  {
812  pVar = pDefault;
813  }
814 }

◆ LoadSolverInfo()

void Nektar::LibUtilities::SessionReader::LoadSolverInfo ( const std::string &  name,
std::string &  var,
const std::string &  def = "" 
) const

Check for and load a solver info property.

Definition at line 877 of file BasicUtils/SessionReader.cpp.

879 {
880  std::string vName = boost::to_upper_copy(pName);
881  auto infoIter = m_solverInfo.find(vName);
882  if (infoIter != m_solverInfo.end())
883  {
884  pVar = infoIter->second;
885  }
886  else
887  {
888  pVar = pDefault;
889  }
890 }

◆ MatchGeometricInfo()

void Nektar::LibUtilities::SessionReader::MatchGeometricInfo ( const std::string &  name,
const std::string &  trueval,
bool &  var,
const bool &  def = false 
) const

Check if the value of a geometric info string property matches.

Definition at line 1098 of file BasicUtils/SessionReader.cpp.

1101 {
1102  std::string vName = boost::to_upper_copy(pName);
1103  auto iter = m_geometricInfo.find(vName);
1104  if (iter != m_geometricInfo.end())
1105  {
1106  pVar = boost::iequals(iter->second, pTrueVal);
1107  }
1108  else
1109  {
1110  pVar = pDefault;
1111  }
1112 }

◆ MatchSolverInfo() [1/2]

bool Nektar::LibUtilities::SessionReader::MatchSolverInfo ( const std::string &  name,
const std::string &  trueval 
) const

Check if the value of a solver info property matches.

Definition at line 914 of file BasicUtils/SessionReader.cpp.

916 {
917  if (DefinesSolverInfo(pName))
918  {
919  std::string vName = boost::to_upper_copy(pName);
920  auto iter = m_solverInfo.find(vName);
921  if (iter != m_solverInfo.end())
922  {
923  return boost::iequals(iter->second, pTrueVal);
924  }
925  }
926  return false;
927 }

◆ MatchSolverInfo() [2/2]

void Nektar::LibUtilities::SessionReader::MatchSolverInfo ( const std::string &  name,
const std::string &  trueval,
bool &  var,
const bool &  def = false 
) const

Check if the value of a solver info property matches.

Definition at line 895 of file BasicUtils/SessionReader.cpp.

898 {
899  std::string vName = boost::to_upper_copy(pName);
900  auto infoIter = m_solverInfo.find(vName);
901  if (infoIter != m_solverInfo.end())
902  {
903  pVar = boost::iequals(infoIter->second, pTrueVal);
904  }
905  else
906  {
907  pVar = pDefault;
908  }
909 }

◆ MatchSolverInfoAsEnum()

template<typename T >
bool Nektar::LibUtilities::SessionReader::MatchSolverInfoAsEnum ( const std::string &  name,
const T &  trueval 
) const
inline

Check if the value of a solver info property matches.

Definition at line 493 of file SessionReader.h.

495 {
496  return (GetSolverInfoAsEnum<T>(name) == trueval);
497 }

References CellMLToNektar.pycml::name.

◆ MergeDoc()

TiXmlDocument * Nektar::LibUtilities::SessionReader::MergeDoc ( const std::vector< std::string > &  pFilenames) const
private

Creates an XML document from a list of input files.

Definition at line 1461 of file BasicUtils/SessionReader.cpp.

1463 {
1464  ASSERTL0(pFilenames.size() > 0, "No filenames for merging.");
1465 
1466  // Read the first document
1467  TiXmlDocument *vMainDoc = new TiXmlDocument;
1468  LoadDoc(pFilenames[0], vMainDoc);
1469 
1470  TiXmlHandle vMainHandle(vMainDoc);
1471  TiXmlElement *vMainNektar =
1472  GetChildElementOrThrow(pFilenames[0], "NEKTAR", vMainHandle);
1473 
1474  // Read all subsequent XML documents.
1475  // For each element within the NEKTAR tag, use it to replace the
1476  // version already present in the loaded XML data.
1477  for (int i = 1; i < pFilenames.size(); ++i)
1478  {
1479  if ((pFilenames[i].compare(pFilenames[i].size() - 3, 3, "xml") == 0) ||
1480  (pFilenames[i].compare(pFilenames[i].size() - 6, 6, "xml.gz") ==
1481  0) ||
1482  (pFilenames[i].compare(pFilenames[i].size() - 3, 3, "opt") == 0))
1483  {
1484  TiXmlDocument *vTempDoc = new TiXmlDocument;
1485  LoadDoc(pFilenames[i], vTempDoc);
1486 
1487  TiXmlHandle docHandle(vTempDoc);
1488  TiXmlElement *vTempNektar =
1489  GetChildElementOrThrow(pFilenames[i], "NEKTAR", docHandle);
1490  TiXmlElement *p = vTempNektar->FirstChildElement();
1491 
1492  while (p)
1493  {
1494  TiXmlElement *vMainEntry =
1495  vMainNektar->FirstChildElement(p->Value());
1496 
1497  // First check if the new item is in fact blank
1498  // replace if it is a COLLECTIONS section however.
1499 
1500  if (!p->FirstChild() && vMainEntry &&
1501  !boost::iequals(p->Value(), "COLLECTIONS"))
1502  {
1503  std::string warningmsg =
1504  "File " + pFilenames[i] + " contains " +
1505  "an empty XML element " + std::string(p->Value()) +
1506  " which will be ignored.";
1507  NEKERROR(ErrorUtil::ewarning, warningmsg.c_str());
1508  }
1509  else
1510  {
1511  if (vMainEntry)
1512  {
1513  vMainNektar->RemoveChild(vMainEntry);
1514  }
1515  TiXmlElement *q = new TiXmlElement(*p);
1516  vMainNektar->LinkEndChild(q);
1517  }
1518  p = p->NextSiblingElement();
1519  }
1520 
1521  delete vTempDoc;
1522  }
1523  }
1524  return vMainDoc;
1525 }
void LoadDoc(const std::string &pFilename, TiXmlDocument *pDoc) const
Loads an xml file into a tinyxml doc and decompresses if needed.
TiXmlElement * GetChildElementOrThrow(const std::string &filename, std::string elementName, const TiXmlHandle &docHandle)

References ASSERTL0, Nektar::LibUtilities::GetChildElementOrThrow(), NEKERROR, and CellMLToNektar.cellml_metadata::p.

◆ ParseCommandLineArguments()

std::vector< std::string > Nektar::LibUtilities::SessionReader::ParseCommandLineArguments ( int  argc,
char *  argv[] 
)
private

Parse the program arguments and fill m_cmdLineOptions.

Parses the command-line arguments for known options and filenames.

Definition at line 401 of file BasicUtils/SessionReader.cpp.

403 {
404  // List the publically visible options (listed using --help)
405  po::options_description desc("Allowed options");
406 
407  // clang-format off
408  desc.add_options()
409  ("verbose,v", "be verbose")
410  ("version,V", "print version information")
411  ("help,h", "print this help message")
412  ("solverinfo,I", po::value<vector<std::string> >(),
413  "override a SOLVERINFO property")
414  ("parameter,P", po::value<vector<std::string> >(),
415  "override a parameter")
416  ("npx", po::value<int>(),
417  "number of procs in X-dir")
418  ("npy", po::value<int>(),
419  "number of procs in Y-dir")
420  ("npz", po::value<int>(),
421  "number of procs in Z-dir")
422  ("nsz", po::value<int>(),
423  "number of slices in Z-dir")
424  ("part-only", po::value<int>(),
425  "only partition mesh into N partitions.")
426  ("part-only-overlapping", po::value<int>(),
427  "only partition mesh into N overlapping partitions.")
428  ("part-info", "Output partition information")
429 #ifdef NEKTAR_USE_CWIPI
430  ("cwipi", po::value<std::string>(),
431  "set CWIPI name")
432 #endif
433  ("writeoptfile", "write an optimisation file")
434  ("useoptfile", po::value<std::string>(),
435  "use an optimisation file")
436  ;
437  // clang-format on
438 
439  for (auto &cmdIt : GetCmdLineArgMap())
440  {
441  std::string names = cmdIt.first;
442  if (cmdIt.second.shortName != "")
443  {
444  names += "," + cmdIt.second.shortName;
445  }
446  if (cmdIt.second.isFlag)
447  {
448  desc.add_options()(names.c_str(), cmdIt.second.description.c_str());
449  }
450  else
451  {
452  desc.add_options()(names.c_str(), po::value<std::string>(),
453  cmdIt.second.description.c_str());
454  }
455  }
456 
457  // List hidden options (e.g. session file arguments are not actually
458  // specified using the input-file option by the user).
459  po::options_description hidden("Hidden options");
460 
461  // clang-format off
462  hidden.add_options()
463  ("input-file", po::value< vector<string> >(),
464  "input filename")
465  ;
466  // clang-format on
467 
468  // Combine all options for the parser
469  po::options_description all("All options");
470  all.add(desc).add(hidden);
471 
472  // Session file is a positional option
473  po::positional_options_description p;
474  p.add("input-file", -1);
475 
476  // Parse the command-line options
477  po::parsed_options parsed = po::command_line_parser(argc, argv)
478  .options(all)
479  .positional(p)
480  .allow_unregistered()
481  .run();
482 
483  // Extract known options to map and update
484  po::store(parsed, m_cmdLineOptions);
485  po::notify(m_cmdLineOptions);
486 
487  // Help message
488  if (m_cmdLineOptions.count("help"))
489  {
490  cout << desc;
491  exit(0);
492  }
493 
494  // Version information
495  if (m_cmdLineOptions.count("version"))
496  {
497  cout << "Nektar++ version " << NEKTAR_VERSION;
498 
499  if (NekConstants::kGitSha1 != "GITDIR-NOTFOUND")
500  {
501  string sha1(NekConstants::kGitSha1);
502  string branch(NekConstants::kGitBranch);
503  boost::replace_all(branch, "refs/heads/", "");
504 
505  cout << " (git changeset " << sha1.substr(0, 8) << ", ";
506 
507  if (branch == "")
508  {
509  cout << "detached head";
510  }
511  else
512  {
513  cout << "head " << branch;
514  }
515 
516  cout << ")";
517  }
518 
519  cout << endl;
520  exit(0);
521  }
522 
523  // Enable verbose mode
524  if (m_cmdLineOptions.count("verbose"))
525  {
526  m_verbose = true;
527  }
528  else
529  {
530  m_verbose = false;
531  }
532 
533  // Enable update optimisation file
534  if (m_cmdLineOptions.count("writeoptfile"))
535  {
536  m_updateOptFile = true;
537  }
538  else
539  {
540  m_updateOptFile = false;
541  }
542 
543  // Print a warning for unknown options
544  for (auto &x : parsed.options)
545  {
546  if (x.unregistered)
547  {
548  cout << "Warning: Unknown option: " << x.string_key << endl;
549  }
550  }
551 
552  // Return the vector of filename(s) given as positional options
553  if (m_cmdLineOptions.count("input-file"))
554  {
555  return m_cmdLineOptions["input-file"].as<std::vector<std::string>>();
556  }
557  else
558  {
559  return std::vector<std::string>();
560  }
561 }
#define NEKTAR_VERSION
static CmdLineArgMap & GetCmdLineArgMap()
CmdLine argument map.
const std::string kGitBranch
const std::string kGitSha1

References Nektar::NekConstants::kGitBranch, Nektar::NekConstants::kGitSha1, NEKTAR_VERSION, and CellMLToNektar.cellml_metadata::p.

◆ ParseDocument()

void Nektar::LibUtilities::SessionReader::ParseDocument ( )
private

Loads and parses the specified file.

Definition at line 1530 of file BasicUtils/SessionReader.cpp.

1531 {
1532  // Check we actually have a document loaded.
1533  ASSERTL0(m_xmlDoc, "No XML document loaded.");
1534 
1535  // Look for all data in CONDITIONS block.
1536  TiXmlHandle docHandle(m_xmlDoc);
1537  TiXmlElement *e;
1538  e = docHandle.FirstChildElement("NEKTAR")
1539  .FirstChildElement("CONDITIONS")
1540  .Element();
1541 
1542  // Read the various sections of the CONDITIONS block
1543  ReadParameters(e);
1544  ReadSolverInfo(e);
1546  ReadTimeIntScheme(e);
1547  ReadExpressions(e);
1548  ReadVariables(e);
1549  ReadFunctions(e);
1550 
1551  e = docHandle.FirstChildElement("NEKTAR")
1552  .FirstChildElement("FILTERS")
1553  .Element();
1554 
1555  ReadFilters(e);
1556 }
void ReadSolverInfo(TiXmlElement *conditions)
Reads the SOLVERINFO section of the XML document.
void ReadVariables(TiXmlElement *conditions)
Reads the VARIABLES section of the XML document.
void ReadTimeIntScheme(TiXmlElement *conditions)
Reads the TIMEINTEGRATIONSCHEME section of the XML document.
void ReadFilters(TiXmlElement *filters)
Reads the FILTERS section of the XML document.
void ReadGlobalSysSolnInfo(TiXmlElement *conditions)
Reads the GLOBALSYSSOLNINFO section of the XML document.
void ReadFunctions(TiXmlElement *conditions)
Reads the FUNCTIONS section of the XML document.
void ReadParameters(TiXmlElement *conditions)
Reads the PARAMETERS section of the XML document.
void ReadExpressions(TiXmlElement *conditions)
Reads the EXPRESSIONS section of the XML document.

References ASSERTL0.

◆ ParseEquals()

void Nektar::LibUtilities::SessionReader::ParseEquals ( const std::string &  line,
std::string &  lhs,
std::string &  rhs 
)
private

Parse a string in the form lhs = rhs.

Pull out lhs and rhs and eliminate any spaces.

Definition at line 2416 of file BasicUtils/SessionReader.cpp.

2418 {
2419  /// Pull out lhs and rhs and eliminate any spaces.
2420  size_t beg = line.find_first_not_of(" ");
2421  size_t end = line.find_first_of("=");
2422  // Check for no parameter name
2423  if (beg == end)
2424  throw 1;
2425  // Check for no parameter value
2426  if (end != line.find_last_of("="))
2427  throw 1;
2428  // Check for no equals sign
2429  if (end == std::string::npos)
2430  throw 1;
2431 
2432  lhs = line.substr(line.find_first_not_of(" "), end - beg);
2433  lhs = lhs.substr(0, lhs.find_last_not_of(" ") + 1);
2434  rhs = line.substr(line.find_last_of("=") + 1);
2435  rhs = rhs.substr(rhs.find_first_not_of(" "));
2436  rhs = rhs.substr(0, rhs.find_last_not_of(" ") + 1);
2437 }

◆ ParseSessionName()

std::string Nektar::LibUtilities::SessionReader::ParseSessionName ( std::vector< std::string > &  filenames)
private

Parse the session name.

Definition at line 566 of file BasicUtils/SessionReader.cpp.

567 {
568  ASSERTL0(!filenames.empty(), "At least one filename expected.");
569 
570  std::string retval = "";
571 
572  // First input file defines the session name
573  std::string fname = filenames[0];
574 
575  // If loading a pre-partitioned mesh, remove _xml extension
576  if (fname.size() > 4 && fname.substr(fname.size() - 4, 4) == "_xml")
577  {
578  retval = fname.substr(0, fname.find_last_of("_"));
579  }
580  // otherwise remove the .xml extension
581  else if (fname.size() > 4 && fname.substr(fname.size() - 4, 4) == ".xml")
582  {
583  retval = fname.substr(0, fname.find_last_of("."));
584  }
585  // If compressed .xml.gz, remove both extensions
586  else if (fname.size() > 7 && fname.substr(fname.size() - 7, 7) == ".xml.gz")
587  {
588  retval = fname.substr(0, fname.find_last_of("."));
589  retval = retval.substr(0, retval.find_last_of("."));
590  }
591 
592  return retval;
593 }

References ASSERTL0.

◆ PartitionComm()

void Nektar::LibUtilities::SessionReader::PartitionComm ( )
private

Partitions the comm object based on session parameters.

Splits the processes into a cartesian grid and creates communicators for each row and column of the grid. The grid is defined by the PROC_X parameter which, if specified, gives the number of processes spanned by the Fourier direction. PROC_X must exactly divide the total number of processes or an error is thrown.

Definition at line 1591 of file BasicUtils/SessionReader.cpp.

1592 {
1593  if (m_comm->GetSize() > 1)
1594  {
1595  int nProcZ = 1;
1596  int nProcY = 1;
1597  int nProcX = 1;
1598  int nStripZ = 1;
1599  if (DefinesCmdLineArgument("npx"))
1600  {
1601  nProcX = GetCmdLineArgument<int>("npx");
1602  }
1603  if (DefinesCmdLineArgument("npy"))
1604  {
1605  nProcY = GetCmdLineArgument<int>("npy");
1606  }
1607  if (DefinesCmdLineArgument("npz"))
1608  {
1609  nProcZ = GetCmdLineArgument<int>("npz");
1610  }
1611  if (DefinesCmdLineArgument("nsz"))
1612  {
1613  nStripZ = GetCmdLineArgument<int>("nsz");
1614  }
1615  ASSERTL0(m_comm->GetSize() % (nProcZ * nProcY * nProcX) == 0,
1616  "Cannot exactly partition using PROC_Z value.");
1617  ASSERTL0(nProcZ % nProcY == 0,
1618  "Cannot exactly partition using PROC_Y value.");
1619  ASSERTL0(nProcY % nProcX == 0,
1620  "Cannot exactly partition using PROC_X value.");
1621 
1622  // Number of processes associated with the spectral method
1623  int nProcSm = nProcZ * nProcY * nProcX;
1624 
1625  // Number of processes associated with the spectral element
1626  // method.
1627  int nProcSem = m_comm->GetSize() / nProcSm;
1628 
1629  m_comm->SplitComm(nProcSm, nProcSem);
1630  m_comm->GetColumnComm()->SplitComm(nProcZ / nStripZ, nStripZ);
1631  m_comm->GetColumnComm()->GetColumnComm()->SplitComm((nProcY * nProcX),
1632  nProcZ / nStripZ);
1633  m_comm->GetColumnComm()->GetColumnComm()->GetColumnComm()->SplitComm(
1634  nProcX, nProcY);
1635  }
1636 }

References ASSERTL0.

◆ ReadExpressions()

void Nektar::LibUtilities::SessionReader::ReadExpressions ( TiXmlElement *  conditions)
private

Reads the EXPRESSIONS section of the XML document.

Definition at line 2030 of file BasicUtils/SessionReader.cpp.

2031 {
2032  m_expressions.clear();
2033 
2034  if (!conditions)
2035  {
2036  return;
2037  }
2038 
2039  TiXmlElement *expressionsElement =
2040  conditions->FirstChildElement("EXPRESSIONS");
2041 
2042  if (expressionsElement)
2043  {
2044  TiXmlElement *expr = expressionsElement->FirstChildElement("E");
2045 
2046  while (expr)
2047  {
2048  stringstream tagcontent;
2049  tagcontent << *expr;
2050  ASSERTL0(expr->Attribute("NAME"),
2051  "Missing NAME attribute in expression "
2052  "definition: \n\t'" +
2053  tagcontent.str() + "'");
2054  std::string nameString = expr->Attribute("NAME");
2055  ASSERTL0(!nameString.empty(),
2056  "Expressions must have a non-empty name: \n\t'" +
2057  tagcontent.str() + "'");
2058 
2059  ASSERTL0(expr->Attribute("VALUE"),
2060  "Missing VALUE attribute in expression "
2061  "definition: \n\t'" +
2062  tagcontent.str() + "'");
2063  std::string valString = expr->Attribute("VALUE");
2064  ASSERTL0(!valString.empty(),
2065  "Expressions must have a non-empty value: \n\t'" +
2066  tagcontent.str() + "'");
2067 
2068  auto exprIter = m_expressions.find(nameString);
2069  ASSERTL0(exprIter == m_expressions.end(),
2070  std::string("Expression '") + nameString +
2071  std::string("' already specified."));
2072 
2073  m_expressions[nameString] = valString;
2074  expr = expr->NextSiblingElement("E");
2075  }
2076  }
2077 }
ExpressionMap m_expressions
Expressions.

References ASSERTL0.

◆ ReadFilters()

void Nektar::LibUtilities::SessionReader::ReadFilters ( TiXmlElement *  filters)
private

Reads the FILTERS section of the XML document.

Definition at line 2375 of file BasicUtils/SessionReader.cpp.

2376 {
2377  if (!filters)
2378  {
2379  return;
2380  }
2381 
2382  m_filters.clear();
2383 
2384  TiXmlElement *filter = filters->FirstChildElement("FILTER");
2385  while (filter)
2386  {
2387  ASSERTL0(filter->Attribute("TYPE"),
2388  "Missing attribute 'TYPE' for filter.");
2389  std::string typeStr = filter->Attribute("TYPE");
2390 
2391  std::map<std::string, std::string> vParams;
2392 
2393  TiXmlElement *param = filter->FirstChildElement("PARAM");
2394  while (param)
2395  {
2396  ASSERTL0(param->Attribute("NAME"),
2397  "Missing attribute 'NAME' for parameter in filter " +
2398  typeStr + "'.");
2399  std::string nameStr = param->Attribute("NAME");
2400 
2401  ASSERTL0(param->GetText(), "Empty value string for param.");
2402  std::string valueStr = param->GetText();
2403 
2404  vParams[nameStr] = valueStr;
2405 
2406  param = param->NextSiblingElement("PARAM");
2407  }
2408 
2409  m_filters.push_back(
2410  std::pair<std::string, FilterParams>(typeStr, vParams));
2411 
2412  filter = filter->NextSiblingElement("FILTER");
2413  }
2414 }

References ASSERTL0.

◆ ReadFunctions()

void Nektar::LibUtilities::SessionReader::ReadFunctions ( TiXmlElement *  conditions)
private

Reads the FUNCTIONS section of the XML document.

Definition at line 2161 of file BasicUtils/SessionReader.cpp.

2162 {
2163  m_functions.clear();
2164 
2165  if (!conditions)
2166  {
2167  return;
2168  }
2169 
2170  // Scan through conditions section looking for functions.
2171  TiXmlElement *function = conditions->FirstChildElement("FUNCTION");
2172  while (function)
2173  {
2174  stringstream tagcontent;
2175  tagcontent << *function;
2176 
2177  // Every function must have a NAME attribute
2178  ASSERTL0(function->Attribute("NAME"),
2179  "Functions must have a NAME attribute defined in XML "
2180  "element: \n\t'" +
2181  tagcontent.str() + "'");
2182  std::string functionStr = function->Attribute("NAME");
2183  ASSERTL0(!functionStr.empty(),
2184  "Functions must have a non-empty name in XML "
2185  "element: \n\t'" +
2186  tagcontent.str() + "'");
2187 
2188  // Store function names in uppercase to remain case-insensitive.
2189  boost::to_upper(functionStr);
2190 
2191  // Retrieve first entry (variable, or file)
2192  TiXmlElement *variable = function->FirstChildElement();
2193 
2194  // Create new function structure with default type of none.
2195  FunctionVariableMap functionVarMap;
2196 
2197  // Process all entries in the function block
2198  while (variable)
2199  {
2200  FunctionVariableDefinition funcDef;
2201  std::string conditionType = variable->Value();
2202 
2203  // If no var is specified, assume wildcard
2204  std::string variableStr;
2205  if (!variable->Attribute("VAR"))
2206  {
2207  variableStr = "*";
2208  }
2209  else
2210  {
2211  variableStr = variable->Attribute("VAR");
2212  }
2213 
2214  // Parse list of variables
2215  std::vector<std::string> variableList;
2216  ParseUtils::GenerateVector(variableStr, variableList);
2217 
2218  // If no domain string put to 0
2219  std::string domainStr;
2220  if (!variable->Attribute("DOMAIN"))
2221  {
2222  domainStr = "0";
2223  }
2224  else
2225  {
2226  domainStr = variable->Attribute("DOMAIN");
2227  }
2228 
2229  // If no domain string put to 0
2230  std::string evarsStr = "x y z t";
2231  if (variable->Attribute("EVARS"))
2232  {
2233  evarsStr =
2234  evarsStr + std::string(" ") + variable->Attribute("EVARS");
2235  }
2236 
2237  // Parse list of variables
2238  std::vector<std::string> varSplit;
2239  std::vector<unsigned int> domainList;
2240  ParseUtils::GenerateSeqVector(domainStr, domainList);
2241 
2242  // Expressions are denoted by E
2243  if (conditionType == "E")
2244  {
2245  funcDef.m_type = eFunctionTypeExpression;
2246 
2247  // Expression must have a VALUE.
2248  ASSERTL0(variable->Attribute("VALUE"),
2249  "Attribute VALUE expected for function '" +
2250  functionStr + "'.");
2251  std::string fcnStr = variable->Attribute("VALUE");
2252 
2253  ASSERTL0(!fcnStr.empty(),
2254  (std::string("Expression for var: ") + variableStr +
2255  std::string(" must be specified."))
2256  .c_str());
2257 
2258  SubstituteExpressions(fcnStr);
2259 
2260  // set expression
2261  funcDef.m_expression =
2263  m_interpreter, fcnStr, evarsStr);
2264  }
2265 
2266  // Files are denoted by F
2267  else if (conditionType == "F")
2268  {
2269  if (variable->Attribute("TIMEDEPENDENT") &&
2270  boost::lexical_cast<bool>(
2271  variable->Attribute("TIMEDEPENDENT")))
2272  {
2273  funcDef.m_type = eFunctionTypeTransientFile;
2274  }
2275  else
2276  {
2277  funcDef.m_type = eFunctionTypeFile;
2278  }
2279 
2280  // File must have a FILE.
2281  ASSERTL0(variable->Attribute("FILE"),
2282  "Attribute FILE expected for function '" +
2283  functionStr + "'.");
2284  std::string filenameStr = variable->Attribute("FILE");
2285 
2286  ASSERTL0(!filenameStr.empty(),
2287  "A filename must be specified for the FILE "
2288  "attribute of function '" +
2289  functionStr + "'.");
2290 
2291  std::vector<std::string> fSplit;
2292  boost::split(fSplit, filenameStr, boost::is_any_of(":"));
2293 
2294  ASSERTL0(fSplit.size() == 1 || fSplit.size() == 2,
2295  "Incorrect filename specification in function " +
2296  functionStr +
2297  "'. "
2298  "Specify variables inside file as: "
2299  "filename:var1,var2");
2300 
2301  // set the filename
2302  funcDef.m_filename = fSplit[0];
2303 
2304  if (fSplit.size() == 2)
2305  {
2306  ASSERTL0(variableList[0] != "*",
2307  "Filename variable mapping not valid "
2308  "when using * as a variable inside "
2309  "function '" +
2310  functionStr + "'.");
2311 
2312  boost::split(varSplit, fSplit[1], boost::is_any_of(","));
2313  ASSERTL0(varSplit.size() == variableList.size(),
2314  "Filename variables should contain the "
2315  "same number of variables defined in "
2316  "VAR in function " +
2317  functionStr + "'.");
2318  }
2319  }
2320 
2321  // Nothing else supported so throw an error
2322  else
2323  {
2324  stringstream tagcontent;
2325  tagcontent << *variable;
2326 
2328  "Identifier " + conditionType + " in function " +
2329  std::string(function->Attribute("NAME")) +
2330  " is not recognised in XML element: \n\t'" +
2331  tagcontent.str() + "'");
2332  }
2333 
2334  // Add variables to function
2335  for (unsigned int i = 0; i < variableList.size(); ++i)
2336  {
2337  for (unsigned int j = 0; j < domainList.size(); ++j)
2338  {
2339  // Check it has not already been defined
2340  pair<std::string, int> key(variableList[i], domainList[j]);
2341  auto fcnsIter = functionVarMap.find(key);
2342  ASSERTL0(
2343  fcnsIter == functionVarMap.end(),
2344  "Error setting expression '" + variableList[i] +
2345  " in domain " +
2346  boost::lexical_cast<std::string>(domainList[j]) +
2347  "' in function '" + functionStr +
2348  "'. "
2349  "Expression has already been defined.");
2350 
2351  if (varSplit.size() > 0)
2352  {
2353  FunctionVariableDefinition funcDef2 = funcDef;
2354  funcDef2.m_fileVariable = varSplit[i];
2355  functionVarMap[key] = funcDef2;
2356  }
2357  else
2358  {
2359  functionVarMap[key] = funcDef;
2360  }
2361  }
2362  }
2363 
2364  variable = variable->NextSiblingElement();
2365  }
2366  // Add function definition to map
2367  m_functions[functionStr] = functionVarMap;
2368  function = function->NextSiblingElement("FUNCTION");
2369  }
2370 }
void SubstituteExpressions(std::string &expr)
Substitutes expressions defined in the XML document.
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:129
static bool GenerateSeqVector(const std::string &str, std::vector< unsigned int > &out)
Takes a comma-separated compressed string and converts it to entries in a vector.
Definition: ParseUtils.cpp:103
std::map< std::pair< std::string, int >, FunctionVariableDefinition > FunctionVariableMap

References ASSERTL0, Nektar::LibUtilities::eFunctionTypeExpression, Nektar::LibUtilities::eFunctionTypeFile, Nektar::LibUtilities::eFunctionTypeTransientFile, Nektar::LibUtilities::FunctionVariableDefinition::m_expression, Nektar::LibUtilities::FunctionVariableDefinition::m_filename, Nektar::LibUtilities::FunctionVariableDefinition::m_fileVariable, Nektar::LibUtilities::FunctionVariableDefinition::m_type, and NEKERROR.

◆ ReadGlobalSysSolnInfo()

void Nektar::LibUtilities::SessionReader::ReadGlobalSysSolnInfo ( TiXmlElement *  conditions)
private

Reads the GLOBALSYSSOLNINFO section of the XML document.

Definition at line 1798 of file BasicUtils/SessionReader.cpp.

1799 {
1800  GetGloSysSolnList().clear();
1801 
1802  if (!conditions)
1803  {
1804  return;
1805  }
1806 
1807  TiXmlElement *GlobalSys =
1808  conditions->FirstChildElement("GLOBALSYSSOLNINFO");
1809 
1810  if (!GlobalSys)
1811  {
1812  return;
1813  }
1814 
1815  TiXmlElement *VarInfo = GlobalSys->FirstChildElement("V");
1816 
1817  while (VarInfo)
1818  {
1819  std::stringstream tagcontent;
1820  tagcontent << *VarInfo;
1821  ASSERTL0(VarInfo->Attribute("VAR"),
1822  "Missing VAR attribute in GobalSysSolnInfo XML "
1823  "element: \n\t'" +
1824  tagcontent.str() + "'");
1825 
1826  std::string VarList = VarInfo->Attribute("VAR");
1827  ASSERTL0(!VarList.empty(),
1828  "VAR attribute must be non-empty in XML element:\n\t'" +
1829  tagcontent.str() + "'");
1830 
1831  // generate a list of variables.
1832  std::vector<std::string> varStrings;
1833  bool valid = ParseUtils::GenerateVector(VarList, varStrings);
1834 
1835  ASSERTL0(valid, "Unable to process list of variable in XML "
1836  "element \n\t'" +
1837  tagcontent.str() + "'");
1838 
1839  if (varStrings.size())
1840  {
1841  TiXmlElement *SysSolnInfo = VarInfo->FirstChildElement("I");
1842 
1843  while (SysSolnInfo)
1844  {
1845  tagcontent.clear();
1846  tagcontent << *SysSolnInfo;
1847  // read the property name
1848  ASSERTL0(SysSolnInfo->Attribute("PROPERTY"),
1849  "Missing PROPERTY attribute in "
1850  "GlobalSysSolnInfo for variable(s) '" +
1851  VarList + "' in XML element: \n\t'" +
1852  tagcontent.str() + "'");
1853 
1854  std::string SysSolnProperty =
1855  SysSolnInfo->Attribute("PROPERTY");
1856 
1857  ASSERTL0(!SysSolnProperty.empty(),
1858  "GlobalSysSolnIno properties must have a "
1859  "non-empty name for variable(s) : '" +
1860  VarList + "' in XML element: \n\t'" +
1861  tagcontent.str() + "'");
1862 
1863  // make sure that solver property is capitalised
1864  std::string SysSolnPropertyUpper =
1865  boost::to_upper_copy(SysSolnProperty);
1866 
1867  // read the value
1868  ASSERTL0(SysSolnInfo->Attribute("VALUE"),
1869  "Missing VALUE attribute in GlobalSysSolnInfo "
1870  "for variable(s) '" +
1871  VarList + "' in XML element: \n\t" +
1872  tagcontent.str() + "'");
1873 
1874  std::string SysSolnValue = SysSolnInfo->Attribute("VALUE");
1875  ASSERTL0(!SysSolnValue.empty(),
1876  "GlobalSysSolnInfo properties must have a "
1877  "non-empty value for variable(s) '" +
1878  VarList + "' in XML element: \n\t'" +
1879  tagcontent.str() + "'");
1880 
1881  // Store values under variable map.
1882  for (int i = 0; i < varStrings.size(); ++i)
1883  {
1884  auto x = GetGloSysSolnList().find(varStrings[i]);
1885  if (x == GetGloSysSolnList().end())
1886  {
1887  (GetGloSysSolnList()[varStrings[i]])
1888  [SysSolnPropertyUpper] = SysSolnValue;
1889  }
1890  else
1891  {
1892  x->second[SysSolnPropertyUpper] = SysSolnValue;
1893  }
1894  }
1895 
1896  SysSolnInfo = SysSolnInfo->NextSiblingElement("I");
1897  }
1898  VarInfo = VarInfo->NextSiblingElement("V");
1899  }
1900  }
1901 
1902  if (m_verbose && GetGloSysSolnList().size() > 0 && m_comm)
1903  {
1904  if (m_comm->GetRank() == 0)
1905  {
1906  cout << "GlobalSysSoln Info:" << endl;
1907 
1908  for (auto &x : GetGloSysSolnList())
1909  {
1910  cout << "\t Variable: " << x.first << endl;
1911 
1912  for (auto &y : x.second)
1913  {
1914  cout << "\t\t " << y.first << " = " << y.second << endl;
1915  }
1916  }
1917  cout << endl;
1918  }
1919  }
1920 }

References ASSERTL0.

◆ ReadParameters()

void Nektar::LibUtilities::SessionReader::ReadParameters ( TiXmlElement *  conditions)
private

Reads the PARAMETERS section of the XML document.

Definition at line 1641 of file BasicUtils/SessionReader.cpp.

1642 {
1643  m_parameters.clear();
1644 
1645  if (!conditions)
1646  {
1647  return;
1648  }
1649 
1650  TiXmlElement *parametersElement =
1651  conditions->FirstChildElement("PARAMETERS");
1652 
1653  // See if we have parameters defined. They are optional so we go on
1654  // if not.
1655  if (parametersElement)
1656  {
1657  TiXmlElement *parameter = parametersElement->FirstChildElement("P");
1658 
1659  ParameterMap caseSensitiveParameters;
1660 
1661  // Multiple nodes will only occur if there is a comment in
1662  // between definitions.
1663  while (parameter)
1664  {
1665  stringstream tagcontent;
1666  tagcontent << *parameter;
1667  TiXmlNode *node = parameter->FirstChild();
1668 
1669  while (node && node->Type() != TiXmlNode::TINYXML_TEXT)
1670  {
1671  node = node->NextSibling();
1672  }
1673 
1674  if (node)
1675  {
1676  // Format is "paramName = value"
1677  std::string line = node->ToText()->Value(), lhs, rhs;
1678 
1679  try
1680  {
1681  ParseEquals(line, lhs, rhs);
1682  }
1683  catch (...)
1684  {
1686  "Syntax error in parameter expression '" + line +
1687  "' in XML element: \n\t'" + tagcontent.str() +
1688  "'");
1689  }
1690 
1691  // We want the list of parameters to have their RHS
1692  // evaluated, so we use the expression evaluator to do
1693  // the dirty work.
1694  if (!lhs.empty() && !rhs.empty())
1695  {
1696  NekDouble value = 0.0;
1697  try
1698  {
1699  LibUtilities::Equation expession(m_interpreter, rhs);
1700  value = expession.Evaluate();
1701  }
1702  catch (const std::runtime_error &)
1703  {
1705  "Error evaluating parameter expression"
1706  " '" +
1707  rhs + "' in XML element: \n\t'" +
1708  tagcontent.str() + "'");
1709  }
1710  m_interpreter->SetParameter(lhs, value);
1711  caseSensitiveParameters[lhs] = value;
1712  boost::to_upper(lhs);
1713  m_parameters[lhs] = value;
1714  }
1715  }
1716 
1717  parameter = parameter->NextSiblingElement();
1718  }
1719  }
1720 }
std::map< std::string, NekDouble > ParameterMap
Definition: SessionReader.h:59

References Nektar::LibUtilities::Equation::Evaluate(), and NEKERROR.

◆ ReadSolverInfo()

void Nektar::LibUtilities::SessionReader::ReadSolverInfo ( TiXmlElement *  conditions)
private

Reads the SOLVERINFO section of the XML document.

Definition at line 1725 of file BasicUtils/SessionReader.cpp.

1726 {
1727  m_solverInfo.clear();
1729 
1730  if (!conditions)
1731  {
1732  return;
1733  }
1734 
1735  TiXmlElement *solverInfoElement =
1736  conditions->FirstChildElement("SOLVERINFO");
1737 
1738  if (solverInfoElement)
1739  {
1740  TiXmlElement *solverInfo = solverInfoElement->FirstChildElement("I");
1741 
1742  while (solverInfo)
1743  {
1744  std::stringstream tagcontent;
1745  tagcontent << *solverInfo;
1746  // read the property name
1747  ASSERTL0(solverInfo->Attribute("PROPERTY"),
1748  "Missing PROPERTY attribute in solver info "
1749  "XML element: \n\t'" +
1750  tagcontent.str() + "'");
1751  std::string solverProperty = solverInfo->Attribute("PROPERTY");
1752  ASSERTL0(!solverProperty.empty(),
1753  "PROPERTY attribute must be non-empty in XML "
1754  "element: \n\t'" +
1755  tagcontent.str() + "'");
1756 
1757  // make sure that solver property is capitalised
1758  std::string solverPropertyUpper =
1759  boost::to_upper_copy(solverProperty);
1760 
1761  // read the value
1762  ASSERTL0(solverInfo->Attribute("VALUE"),
1763  "Missing VALUE attribute in solver info "
1764  "XML element: \n\t'" +
1765  tagcontent.str() + "'");
1766  std::string solverValue = solverInfo->Attribute("VALUE");
1767  ASSERTL0(!solverValue.empty(),
1768  "VALUE attribute must be non-empty in XML "
1769  "element: \n\t'" +
1770  tagcontent.str() + "'");
1771 
1772  // Set Variable
1773  m_solverInfo[solverPropertyUpper] = solverValue;
1774  solverInfo = solverInfo->NextSiblingElement("I");
1775  }
1776  }
1777 
1778  if (m_comm && m_comm->GetRowComm()->GetSize() > 1)
1779  {
1780  ASSERTL0(
1781  m_solverInfo["GLOBALSYSSOLN"] == "IterativeFull" ||
1782  m_solverInfo["GLOBALSYSSOLN"] == "IterativeStaticCond" ||
1783  m_solverInfo["GLOBALSYSSOLN"] ==
1784  "IterativeMultiLevelStaticCond" ||
1785  m_solverInfo["GLOBALSYSSOLN"] == "XxtFull" ||
1786  m_solverInfo["GLOBALSYSSOLN"] == "XxtStaticCond" ||
1787  m_solverInfo["GLOBALSYSSOLN"] == "XxtMultiLevelStaticCond" ||
1788  m_solverInfo["GLOBALSYSSOLN"] == "PETScFull" ||
1789  m_solverInfo["GLOBALSYSSOLN"] == "PETScStaticCond" ||
1790  m_solverInfo["GLOBALSYSSOLN"] == "PETScMultiLevelStaticCond",
1791  "A parallel solver must be used when run in parallel.");
1792  }
1793 }

References ASSERTL0.

◆ ReadTimeIntScheme()

void Nektar::LibUtilities::SessionReader::ReadTimeIntScheme ( TiXmlElement *  conditions)
private

Reads the TIMEINTEGRATIONSCHEME section of the XML document.

Read the time-integration scheme structure, if present.

Definition at line 1925 of file BasicUtils/SessionReader.cpp.

1926 {
1927  if (!conditions)
1928  {
1929  return;
1930  }
1931 
1932  TiXmlElement *timeInt =
1933  conditions->FirstChildElement("TIMEINTEGRATIONSCHEME");
1934 
1935  if (!timeInt)
1936  {
1937  return;
1938  }
1939 
1940  TiXmlElement *method = timeInt->FirstChildElement("METHOD");
1941  TiXmlElement *variant = timeInt->FirstChildElement("VARIANT");
1942  TiXmlElement *order = timeInt->FirstChildElement("ORDER");
1943  TiXmlElement *params = timeInt->FirstChildElement("FREEPARAMETERS");
1944 
1945  // Only the method and order are required.
1946  ASSERTL0(method, "Missing METHOD tag inside "
1947  "TIMEINTEGRATIONSCHEME.");
1948  ASSERTL0(order, "Missing ORDER tag inside "
1949  "TIMEINTEGRATIONSCHEME.");
1950 
1951  m_timeIntScheme.method = method->GetText();
1952 
1953  std::string orderStr = order->GetText();
1954 
1955  // Only the method and order are required.
1956  ASSERTL0(m_timeIntScheme.method.size() > 0,
1957  "Empty text inside METHOD tag in TIMEINTEGRATIONSCHEME.");
1958  ASSERTL0(orderStr.size() > 0,
1959  "Empty text inside ORDER tag in TIMEINTEGRATIONSCHEME.");
1960  try
1961  {
1962  m_timeIntScheme.order = boost::lexical_cast<unsigned int>(orderStr);
1963  }
1964  catch (...)
1965  {
1966  NEKERROR(ErrorUtil::efatal, "In ORDER tag, unable to convert "
1967  "string '" +
1968  orderStr + "' to an unsigned integer.");
1969  }
1970 
1971  if (variant)
1972  {
1973  m_timeIntScheme.variant = variant->GetText();
1974  }
1975 
1976  if (params)
1977  {
1978  std::string paramsStr = params->GetText();
1979  ASSERTL0(paramsStr.size() > 0,
1980  "Empty text inside FREEPARAMETERS tag in "
1981  "TIMEINTEGRATIONSCHEME.");
1982 
1983  std::vector<std::string> pSplit;
1984  boost::split(pSplit, paramsStr, boost::is_any_of(" "));
1985 
1986  m_timeIntScheme.freeParams.resize(pSplit.size());
1987  for (size_t i = 0; i < pSplit.size(); ++i)
1988  {
1989  try
1990  {
1992  boost::lexical_cast<NekDouble>(pSplit[i]);
1993  }
1994  catch (...)
1995  {
1996  NEKERROR(ErrorUtil::efatal, "In FREEPARAMETERS tag, "
1997  "unable to convert string '" +
1998  pSplit[i] +
1999  "' "
2000  "to a floating-point value.");
2001  }
2002  }
2003  }
2004 
2005  if (m_verbose && m_comm)
2006  {
2007  if (m_comm->GetRank() == 0)
2008  {
2009  cout << "Trying to use time integration scheme:" << endl;
2010  cout << "\t Method : " << m_timeIntScheme.method << endl;
2011  cout << "\t Variant: " << m_timeIntScheme.variant << endl;
2012  cout << "\t Order : " << m_timeIntScheme.order << endl;
2013 
2014  if (m_timeIntScheme.freeParams.size() > 0)
2015  {
2016  cout << "\t Params :";
2017  for (auto &x : m_timeIntScheme.freeParams)
2018  {
2019  cout << " " << x;
2020  }
2021  cout << endl;
2022  }
2023  }
2024  }
2025 }
std::vector< NekDouble > freeParams
Definition: SessionReader.h:90

References ASSERTL0, and NEKERROR.

◆ ReadVariables()

void Nektar::LibUtilities::SessionReader::ReadVariables ( TiXmlElement *  conditions)
private

Reads the VARIABLES section of the XML document.

All elements are of the form: "<V ID="#"> name = value </V>", with ? being the element type.

Definition at line 2082 of file BasicUtils/SessionReader.cpp.

2083 {
2084  m_variables.clear();
2085 
2086  if (!conditions)
2087  {
2088  return;
2089  }
2090 
2091  TiXmlElement *variablesElement = conditions->FirstChildElement("VARIABLES");
2092 
2093  // See if we have parameters defined. They are optional so we go on
2094  // if not.
2095  if (variablesElement)
2096  {
2097  TiXmlElement *varElement = variablesElement->FirstChildElement("V");
2098 
2099  // Sequential counter for the composite numbers.
2100  int nextVariableNumber = -1;
2101 
2102  while (varElement)
2103  {
2104  stringstream tagcontent;
2105  tagcontent << *varElement;
2106 
2107  /// All elements are of the form: "<V ID="#"> name = value
2108  /// </V>", with ? being the element type.
2109  nextVariableNumber++;
2110 
2111  int i;
2112  int err = varElement->QueryIntAttribute("ID", &i);
2113  ASSERTL0(err == TIXML_SUCCESS,
2114  "Variables must have a unique ID number attribute "
2115  "in XML element: \n\t'" +
2116  tagcontent.str() + "'");
2117  ASSERTL0(i == nextVariableNumber,
2118  "ID numbers for variables must begin with zero and"
2119  " be sequential in XML element: \n\t'" +
2120  tagcontent.str() + "'");
2121 
2122  TiXmlNode *varChild = varElement->FirstChild();
2123  // This is primarily to skip comments that may be present.
2124  // Comments appear as nodes just like elements. We are
2125  // specifically looking for text in the body of the
2126  // definition.
2127  while (varChild && varChild->Type() != TiXmlNode::TINYXML_TEXT)
2128  {
2129  varChild = varChild->NextSibling();
2130  }
2131 
2132  ASSERTL0(varChild, "Unable to read variable definition body for "
2133  "variable with ID " +
2134  boost::lexical_cast<string>(i) +
2135  " in XML element: \n\t'" + tagcontent.str() +
2136  "'");
2137  std::string variableName = varChild->ToText()->ValueStr();
2138 
2139  std::istringstream variableStrm(variableName);
2140  variableStrm >> variableName;
2141 
2142  ASSERTL0(std::find(m_variables.begin(), m_variables.end(),
2143  variableName) == m_variables.end(),
2144  "Variable with ID " + boost::lexical_cast<string>(i) +
2145  " in XML element \n\t'" + tagcontent.str() +
2146  "'\nhas already been defined.");
2147 
2148  m_variables.push_back(variableName);
2149 
2150  varElement = varElement->NextSiblingElement("V");
2151  }
2152 
2153  ASSERTL0(nextVariableNumber > -1,
2154  "Number of variables must be greater than zero.");
2155  }
2156 }
InputIterator find(InputIterator first, InputIterator last, InputIterator startingpoint, const EqualityComparable &value)
Definition: StdRegions.hpp:327

References ASSERTL0, and Nektar::StdRegions::find().

◆ RegisterCmdLineArgument()

std::string Nektar::LibUtilities::SessionReader::RegisterCmdLineArgument ( const std::string &  pName,
const std::string &  pShortName,
const std::string &  pDescription 
)
inlinestatic

Registers a command-line argument with the session reader.

Definition at line 611 of file SessionReader.h.

614 {
615  ASSERTL0(!pName.empty(), "Empty name for cmdline argument.");
616  CmdLineArg x;
617  x.shortName = pShortName;
618  x.description = pDescription;
619  x.isFlag = false;
620  GetCmdLineArgMap()[pName] = x;
621  return pName;
622 }

References ASSERTL0, Nektar::LibUtilities::CmdLineArg::description, GetCmdLineArgMap(), Nektar::LibUtilities::CmdLineArg::isFlag, and Nektar::LibUtilities::CmdLineArg::shortName.

◆ RegisterCmdLineFlag()

std::string Nektar::LibUtilities::SessionReader::RegisterCmdLineFlag ( const std::string &  pName,
const std::string &  pShortName,
const std::string &  pDescription 
)
inlinestatic

Registers a command-line flag with the session reader.

Definition at line 627 of file SessionReader.h.

630 {
631  ASSERTL0(!pName.empty(), "Empty name for cmdline argument.");
632  CmdLineArg x;
633  x.shortName = pShortName;
634  x.description = pDescription;
635  x.isFlag = true;
636  GetCmdLineArgMap()[pName] = x;
637  return pName;
638 }

References ASSERTL0, Nektar::LibUtilities::CmdLineArg::description, GetCmdLineArgMap(), Nektar::LibUtilities::CmdLineArg::isFlag, and Nektar::LibUtilities::CmdLineArg::shortName.

◆ RegisterDefaultSolverInfo()

std::string Nektar::LibUtilities::SessionReader::RegisterDefaultSolverInfo ( const std::string &  pName,
const std::string &  pValue 
)
inlinestatic

Registers the default string value of a solver info property.

A default value for a given solver info property may be registered using this function. The property will take this value until it is overwritten by a value specified in the XML document, or specified as a command-line argument. Usage has the form:

std::string GlobalLinSys::def
"GlobalSysSoln","DirectMultiLevelStaticCond");
static std::string RegisterDefaultSolverInfo(const std::string &pName, const std::string &pValue)
Registers the default string value of a solver info property.
Parameters
pNameThe name of the property.
pValueThe default value of the property.
Returns
The default value of the property provided by #pValue.

Definition at line 600 of file SessionReader.h.

602 {
603  std::string vName = boost::to_upper_copy(pName);
604  GetSolverInfoDefaults()[vName] = pValue;
605  return pValue;
606 }

References GetSolverInfoDefaults().

◆ RegisterEnumValue()

std::string Nektar::LibUtilities::SessionReader::RegisterEnumValue ( std::string  pEnum,
std::string  pString,
int  pEnumValue 
)
inlinestatic

Registers an enumeration value.

A set of valid values for a given solver info property may be registered using this function. It must be called statically during the initialisation of a static variable. For example:

std::string GlobalLinSys::lookupIds[2] = {
"GlobalSysSoln",
"DirectFull",
"GlobalSysSoln",
"DirectStaticCond",
}
static std::string RegisterEnumValue(std::string pEnum, std::string pString, int pEnumValue)
Registers an enumeration value.
static std::string lookupIds[]
Definition: GlobalLinSys.h:156
Parameters
pEnumThe name of the property.
pStringA valid value for the property.
pEnumValueAn enumeration value corresponding to this value.
Returns
The value for the property provided by #pString.

Definition at line 566 of file SessionReader.h.

569 {
570  std::string vEnum = boost::to_upper_copy(pEnum);
571  auto x = GetSolverInfoEnums().find(vEnum);
572 
573  if (x == GetSolverInfoEnums().end())
574  {
575  GetSolverInfoEnums()[vEnum] = EnumMap();
576  x = GetSolverInfoEnums().find(vEnum);
577  }
578 
579  x->second[pString] = pEnumValue;
580  return pString;
581 }
std::map< std::string, int > EnumMap
Definition: SessionReader.h:76

References GetSolverInfoEnums().

◆ SetParameter() [1/2]

void Nektar::LibUtilities::SessionReader::SetParameter ( const std::string &  name,
int &  var 
)

Set an integer parameter.

Definition at line 819 of file BasicUtils/SessionReader.cpp.

820 {
821  std::string vName = boost::to_upper_copy(pName);
822  m_parameters[vName] = pVar;
823 }

◆ SetParameter() [2/2]

void Nektar::LibUtilities::SessionReader::SetParameter ( const std::string &  name,
NekDouble var 
)

Set a double precision parameter.

Definition at line 828 of file BasicUtils/SessionReader.cpp.

829 {
830  std::string vName = boost::to_upper_copy(pName);
831  m_parameters[vName] = pVar;
832 }

◆ SetSolverInfo()

void Nektar::LibUtilities::SessionReader::SetSolverInfo ( const std::string &  pProperty,
const std::string &  pValue 
)

Sets the value of the specified solver info property.

Definition at line 862 of file BasicUtils/SessionReader.cpp.

864 {
865  std::string vProperty = boost::to_upper_copy(pProperty);
866  auto iter = m_solverInfo.find(vProperty);
867 
868  ASSERTL1(iter != m_solverInfo.end(),
869  "Unable to find requested property: " + pProperty);
870 
871  iter->second = pValue;
872 }

References ASSERTL1.

◆ SetTag()

void Nektar::LibUtilities::SessionReader::SetTag ( const std::string &  pName,
const std::string &  pValue 
)

Sets a specified tag.

Definition at line 1367 of file BasicUtils/SessionReader.cpp.

1368 {
1369  std::string vName = boost::to_upper_copy(pName);
1370  m_tags[vName] = pValue;
1371 }

◆ SetUpdateOptFile()

void Nektar::LibUtilities::SessionReader::SetUpdateOptFile ( bool  flag)
inline

Set bool to update optimisation file.

Definition at line 386 of file SessionReader.h.

387  {
388  m_updateOptFile = flag;
389  }

References m_updateOptFile.

◆ SetUpXmlDoc()

void Nektar::LibUtilities::SessionReader::SetUpXmlDoc ( void  )

Definition at line 2524 of file BasicUtils/SessionReader.cpp.

2525 {
2527 }

◆ SetVariable()

void Nektar::LibUtilities::SessionReader::SetVariable ( const unsigned int &  idx,
std::string  newname 
)

Definition at line 1126 of file BasicUtils/SessionReader.cpp.

1127 {
1128  ASSERTL0(idx < m_variables.size(), "Variable index out of range.");
1129  m_variables[idx] = newname;
1130 }

References ASSERTL0.

◆ SubstituteExpressions()

void Nektar::LibUtilities::SessionReader::SubstituteExpressions ( std::string &  expr)

Substitutes expressions defined in the XML document.

Definition at line 1403 of file BasicUtils/SessionReader.cpp.

1404 {
1405  for (auto &exprIter : m_expressions)
1406  {
1407  boost::replace_all(pExpr, exprIter.first, exprIter.second);
1408  }
1409 }

◆ TestSharedFilesystem()

void Nektar::LibUtilities::SessionReader::TestSharedFilesystem ( )
private

Definition at line 360 of file BasicUtils/SessionReader.cpp.

361 {
362  m_sharedFilesystem = false;
363 
364  if (m_comm->GetSize() > 1)
365  {
366  if (m_comm->GetRank() == 0)
367  {
368  std::ofstream testfile("shared-fs-testfile");
369  testfile << "" << std::endl;
370  ASSERTL1(!testfile.fail(), "Test file creation failed");
371  testfile.close();
372  }
373  m_comm->Block();
374 
375  int exists = (bool)boost::filesystem::exists("shared-fs-testfile");
376  m_comm->AllReduce(exists, LibUtilities::ReduceSum);
377 
378  m_sharedFilesystem = (exists == m_comm->GetSize());
379 
380  if ((m_sharedFilesystem && m_comm->GetRank() == 0) ||
382  {
383  std::remove("shared-fs-testfile");
384  }
385  }
386  else
387  {
388  m_sharedFilesystem = false;
389  }
390 
391  if (m_verbose && m_comm->GetRank() == 0 && m_sharedFilesystem)
392  {
393  cout << "Shared filesystem detected" << endl;
394  }
395 }

References ASSERTL1, and Nektar::LibUtilities::ReduceSum.

◆ VerifySolverInfo()

void Nektar::LibUtilities::SessionReader::VerifySolverInfo ( )
private

Check values of solver info options are valid.

Definition at line 2505 of file BasicUtils/SessionReader.cpp.

2506 {
2507  for (auto &x : m_solverInfo)
2508  {
2509  std::string solverProperty = x.first;
2510  std::string solverValue = x.second;
2511 
2512  auto propIt = GetSolverInfoEnums().find(solverProperty);
2513  if (propIt != GetSolverInfoEnums().end())
2514  {
2515  auto valIt = propIt->second.find(solverValue);
2516  ASSERTL0(valIt != propIt->second.end(), "Value '" + solverValue +
2517  "' is not valid for "
2518  "property '" +
2519  solverProperty + "'");
2520  }
2521  }
2522 }

References ASSERTL0.

Friends And Related Function Documentation

◆ MemoryManager< SessionReader >

friend class MemoryManager< SessionReader >
friend

Support creation through MemoryManager.

Definition at line 108 of file SessionReader.h.

Member Data Documentation

◆ m_cmdLineOptions

boost::program_options::variables_map Nektar::LibUtilities::SessionReader::m_cmdLineOptions
private

Definition at line 392 of file SessionReader.h.

Referenced by GetCmdLineArgument().

◆ m_comm

CommSharedPtr Nektar::LibUtilities::SessionReader::m_comm
private

Communication object.

Definition at line 395 of file SessionReader.h.

◆ m_expressions

ExpressionMap Nektar::LibUtilities::SessionReader::m_expressions
private

Expressions.

Definition at line 409 of file SessionReader.h.

◆ m_filenames

std::vector<std::string> Nektar::LibUtilities::SessionReader::m_filenames
private

Filenames.

Definition at line 397 of file SessionReader.h.

◆ m_filters

FilterMap Nektar::LibUtilities::SessionReader::m_filters
private

Filters map.

Definition at line 419 of file SessionReader.h.

◆ m_functions

FunctionMap Nektar::LibUtilities::SessionReader::m_functions
private

Functions.

Definition at line 413 of file SessionReader.h.

◆ m_geometricInfo

GeometricInfoMap Nektar::LibUtilities::SessionReader::m_geometricInfo
private

Geometric information properties.

Definition at line 407 of file SessionReader.h.

◆ m_interpreter

InterpreterSharedPtr Nektar::LibUtilities::SessionReader::m_interpreter
private

Interpreter instance.

Definition at line 411 of file SessionReader.h.

◆ m_parameters

ParameterMap Nektar::LibUtilities::SessionReader::m_parameters
private

Parameters.

Definition at line 403 of file SessionReader.h.

◆ m_sessionName

std::string Nektar::LibUtilities::SessionReader::m_sessionName
private

Session name of the loaded XML document (filename minus ext).

Definition at line 399 of file SessionReader.h.

◆ m_sharedFilesystem

bool Nektar::LibUtilities::SessionReader::m_sharedFilesystem
private

Running on a shared filesystem.

Definition at line 425 of file SessionReader.h.

◆ m_solverInfo

SolverInfoMap Nektar::LibUtilities::SessionReader::m_solverInfo
private

Solver information properties.

Definition at line 405 of file SessionReader.h.

◆ m_tags

TagMap Nektar::LibUtilities::SessionReader::m_tags
private

Custom tags.

Definition at line 417 of file SessionReader.h.

◆ m_timeIntScheme

TimeIntScheme Nektar::LibUtilities::SessionReader::m_timeIntScheme
private

Time integration scheme information.

Definition at line 421 of file SessionReader.h.

◆ m_updateOptFile

bool Nektar::LibUtilities::SessionReader::m_updateOptFile
private

Update optimisation file.

Definition at line 427 of file SessionReader.h.

Referenced by GetUpdateOptFile(), and SetUpdateOptFile().

◆ m_variables

VariableList Nektar::LibUtilities::SessionReader::m_variables
private

Variables.

Definition at line 415 of file SessionReader.h.

◆ m_verbose

bool Nektar::LibUtilities::SessionReader::m_verbose
private

Be verbose.

Definition at line 423 of file SessionReader.h.

◆ m_xmlDoc

TiXmlDocument* Nektar::LibUtilities::SessionReader::m_xmlDoc
private

Pointer to the loaded XML document.

Definition at line 401 of file SessionReader.h.