Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
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>

Inheritance diagram for Nektar::LibUtilities::SessionReader:
Inheritance graph
[legend]
Collaboration diagram for Nektar::LibUtilities::SessionReader:
Collaboration graph
[legend]

Public Member Functions

 SessionReader (int argc, char *argv[], const std::vector< std::string > &pFilenames, const CommSharedPtr &pComm)
 
 ~SessionReader ()
 Destructor. More...
 
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...
 
CommSharedPtrGetComm ()
 Returns the communication object. 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 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...
 
AnalyticExpressionEvaluatorGetExpressionEvaluator ()
 Returns the instance of AnalyticExpressionEvaluator 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...
 
CompositeOrdering GetCompositeOrdering () const
 
BndRegionOrdering GetBndRegionOrdering () const
 
void SetUpXmlDoc ()
 

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 InitSession ()
 
SessionReaderSharedPtr GetSharedThisPtr ()
 Returns a shared pointer to the current object. More...
 
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 PartitionMesh ()
 Partitions the mesh when running in parallel. 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 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...
 
AnalyticExpressionEvaluator m_exprEvaluator
 Analytic expression evaluator instance. More...
 
FunctionMap m_functions
 Functions. More...
 
VariableList m_variables
 Variables. More...
 
TagMap m_tags
 Custom tags. More...
 
FilterMap m_filters
 Filters map. More...
 
bool m_verbose
 Be verbose. More...
 
CompositeOrdering m_compOrder
 Map of original composite ordering for parallel periodic bcs. More...
 
BndRegionOrdering m_bndRegOrder
 Map of original boundary region ordering for parallel periodic bcs. 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:

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 123 of file SessionReader.h.

Constructor & Destructor Documentation

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

Definition at line 210 of file SessionReader.cpp.

References ASSERTL0.

215  {
216  ASSERTL0(pFilenames.size() > 0, "No filenames specified.");
217 
218  ParseCommandLineArguments(argc, argv);
219  m_xmlDoc = 0;
220  m_filenames = pFilenames;
221 
223 
224  // Create communicator
225  if (!pComm.get())
226  {
227  CreateComm(argc, argv);
228  }
229  else
230  {
231  m_comm = pComm;
232 
233  if (m_comm->GetSize() > 1)
234  {
235  GetSolverInfoDefaults()["GLOBALSYSSOLN"] =
236  "IterativeStaticCond";
237  }
238  }
239 
240  // If running in parallel change the default global sys solution
241  // type.
242  if (m_comm->GetSize() > 1)
243  {
244  GetSolverInfoDefaults()["GLOBALSYSSOLN"] =
245  "IterativeStaticCond";
246  }
247  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
std::vector< std::string > ParseCommandLineArguments(int argc, char *argv[])
Parse the program arguments and fill m_cmdLineOptions.
CommSharedPtr m_comm
Communication object.
std::vector< std::string > m_filenames
Filenames.
void CreateComm(int &argc, char *argv[])
Loads the given XML document and instantiates an appropriate communication object.
std::string ParseSessionName(std::vector< std::string > &filenames)
Parse the session name.
static SolverInfoMap & GetSolverInfoDefaults()
Default solver info options.
TiXmlDocument * m_xmlDoc
Pointer to the loaded XML document.
std::string m_sessionName
Session name of the loaded XML document (filename minus ext).
Nektar::LibUtilities::SessionReader::~SessionReader ( )

Destructor.

Definition at line 253 of file SessionReader.cpp.

254  {
255  delete m_xmlDoc;
256  }
TiXmlDocument * m_xmlDoc
Pointer to the loaded XML document.
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 SessionReader.cpp.

References ASSERTL0.

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 
197  // If running in parallel change the default global sys solution
198  // type.
199  if (m_comm->GetSize() > 1)
200  {
201  GetSolverInfoDefaults()["GLOBALSYSSOLN"] =
202  "IterativeStaticCond";
203  }
204  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
std::vector< std::string > ParseCommandLineArguments(int argc, char *argv[])
Parse the program arguments and fill m_cmdLineOptions.
CommSharedPtr m_comm
Communication object.
std::vector< std::string > m_filenames
Filenames.
void CreateComm(int &argc, char *argv[])
Loads the given XML document and instantiates an appropriate communication object.
std::string ParseSessionName(std::vector< std::string > &filenames)
Parse the session name.
static SolverInfoMap & GetSolverInfoDefaults()
Default solver info options.
TiXmlDocument * m_xmlDoc
Pointer to the loaded XML document.
std::string m_sessionName
Session name of the loaded XML document (filename minus ext).

Member Function Documentation

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

Enforce parameters from command line arguments.

Definition at line 2588 of file SessionReader.cpp.

References ASSERTL0, and Nektar::lhs.

2589  {
2590  // Parse solver info overrides
2591  if (m_cmdLineOptions.count("solverinfo"))
2592  {
2593  std::vector<std::string> solverInfoList =
2594  m_cmdLineOptions["solverinfo"].as<
2595  std::vector<std::string> >();
2596 
2597  for (int i = 0; i < solverInfoList.size(); ++i)
2598  {
2599  std::string lhs, rhs;
2600 
2601  try
2602  {
2603  ParseEquals(solverInfoList[i], lhs, rhs);
2604  }
2605  catch (...)
2606  {
2607  ASSERTL0(false, "Parse error with command line "
2608  "option: "+solverInfoList[i]);
2609  }
2610 
2611  std::string lhsUpper = boost::to_upper_copy(lhs);
2612  m_solverInfo[lhsUpper] = rhs;
2613  }
2614  }
2615 
2616  if (m_cmdLineOptions.count("parameter"))
2617  {
2618  std::vector<std::string> parametersList =
2619  m_cmdLineOptions["parameter"].as<
2620  std::vector<std::string> >();
2621 
2622  for (int i = 0; i < parametersList.size(); ++i)
2623  {
2624  std::string lhs, rhs;
2625 
2626  try
2627  {
2628  ParseEquals(parametersList[i], lhs, rhs);
2629  }
2630  catch (...)
2631  {
2632  ASSERTL0(false, "Parse error with command line "
2633  "option: "+parametersList[i]);
2634  }
2635 
2636  std::string lhsUpper = boost::to_upper_copy(lhs);
2637 
2638  try
2639  {
2640  m_parameters[lhsUpper] =
2641  boost::lexical_cast<NekDouble>(rhs);
2642  }
2643  catch (...)
2644  {
2645  ASSERTL0(false, "Unable to convert string: "+rhs+
2646  "to double value.");
2647  }
2648  }
2649  }
2650  }
ParameterMap m_parameters
Parameters.
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
SolverInfoMap m_solverInfo
Solver information properties.
StandardMatrixTag & lhs
void ParseEquals(const std::string &line, std::string &lhs, std::string &rhs)
Parse a string in the form lhs = rhs.
boost::program_options::variables_map m_cmdLineOptions
double NekDouble
void Nektar::LibUtilities::SessionReader::CreateComm ( int &  argc,
char *  argv[] 
)
private

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

Definition at line 1528 of file SessionReader.cpp.

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

1531  {
1532  if (argc == 0)
1533  {
1534  m_comm = GetCommFactory().CreateInstance("Serial", 0, 0);
1535  }
1536  else
1537  {
1538  string vCommModule("Serial");
1539  if (GetCommFactory().ModuleExists("ParallelMPI"))
1540  {
1541  vCommModule = "ParallelMPI";
1542  }
1543 
1544  m_comm = GetCommFactory().CreateInstance(vCommModule,argc,argv);
1545  }
1546  }
tBaseSharedPtr CreateInstance(tKey idKey BOOST_PP_COMMA_IF(MAX_PARAM) BOOST_PP_ENUM_BINARY_PARAMS(MAX_PARAM, tParam, x))
Create an instance of the class referred to by idKey.
Definition: NekFactory.hpp:162
CommSharedPtr m_comm
Communication object.
CommFactory & GetCommFactory()
Definition: Comm.cpp:64
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 140 of file SessionReader.h.

Referenced by main(), Nektar::Utilities::InputCAD::Process(), Nektar::Utilities::InputNekpp::Process(), Nektar::Utilities::InputXml::Process(), Nektar::Utilities::ProcessDisplacement::Process(), Nektar::Utilities::ProcessInterpField::Process(), Nektar::Utilities::ProcessInterpPoints::Process(), Nektar::Utilities::OutputNekpp::Process(), Nektar::SolverUtils::Driver::v_InitObject(), and Nektar::VortexWaveInteraction::VortexWaveInteraction().

142  {
143  SessionReaderSharedPtr p = MemoryManager<
144  LibUtilities::SessionReader>::AllocateSharedPtr(argc, argv);
145  p->InitSession();
146  return p;
147  }
boost::shared_ptr< SessionReader > SessionReaderSharedPtr
Definition: MeshPartition.h:51
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 159 of file SessionReader.h.

164  {
165  SessionReaderSharedPtr p = MemoryManager<
166  LibUtilities::SessionReader>
167  ::AllocateSharedPtr(argc, argv, pFilenames, pComm);
168  p->InitSession();
169  return p;
170  }
boost::shared_ptr< SessionReader > SessionReaderSharedPtr
Definition: MeshPartition.h:51
bool Nektar::LibUtilities::SessionReader::DefinesCmdLineArgument ( const std::string &  pName) const

Checks if a specified cmdline argument has been given.

Definition at line 1350 of file SessionReader.cpp.

1352  {
1353  return (m_cmdLineOptions.find(pName) != m_cmdLineOptions.end());
1354  }
boost::program_options::variables_map m_cmdLineOptions
bool Nektar::LibUtilities::SessionReader::DefinesElement ( const std::string &  pPath) const

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

Definition at line 560 of file SessionReader.cpp.

References ASSERTL0.

561  {
562  std::string vPath = boost::to_upper_copy(pPath);
563  std::vector<std::string> st;
564  boost::split(st, vPath, boost::is_any_of("\\/ "));
565  ASSERTL0(st.size() > 0, "No path given in XML element request.");
566 
567  TiXmlElement* vReturn = m_xmlDoc->FirstChildElement(st[0].c_str());
568  ASSERTL0(vReturn, std::string("Cannot find element '")
569  + st[0] + std::string("'."));
570  for (int i = 1; i < st.size(); ++i)
571  {
572  vReturn = vReturn->FirstChildElement(st[i].c_str());
573  if (!vReturn) return false;
574  }
575  return true;
576  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
TiXmlDocument * m_xmlDoc
Pointer to the loaded XML document.
bool Nektar::LibUtilities::SessionReader::DefinesFunction ( const std::string &  name) const

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

Definition at line 1046 of file SessionReader.cpp.

1047  {
1048  FunctionMap::const_iterator it1;
1049  std::string vName = boost::to_upper_copy(pName);
1050 
1051  if ((it1 = m_functions.find(vName)) != m_functions.end())
1052  {
1053  return true;
1054  }
1055  return false;
1056  }
FunctionMap m_functions
Functions.
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 1062 of file SessionReader.cpp.

1066  {
1067  FunctionMap::const_iterator it1;
1068  FunctionVariableMap::const_iterator it2;
1069  std::string vName = boost::to_upper_copy(pName);
1070 
1071  // Check function exists
1072  if ((it1 = m_functions.find(vName)) != m_functions.end())
1073  {
1074  pair<std::string, int> key(pVariable,pDomain);
1075  pair<std::string, int> defkey("*",pDomain);
1076  bool varExists =
1077  (it2 = it1->second.find(key)) != it1->second.end() ||
1078  (it2 = it1->second.find(defkey)) != it1->second.end();
1079  return varExists;
1080  }
1081  return false;
1082  }
FunctionMap m_functions
Functions.
bool Nektar::LibUtilities::SessionReader::DefinesGeometricInfo ( const std::string &  name) const

Checks if a geometric info property is defined.

Definition at line 911 of file SessionReader.cpp.

912  {
913  std::string vName = boost::to_upper_copy(pName);
914  GeometricInfoMap::const_iterator iter = m_geometricInfo.find(vName);
915  return (iter != m_geometricInfo.end());
916  }
GeometricInfoMap m_geometricInfo
Geometric information properties.
bool Nektar::LibUtilities::SessionReader::DefinesGlobalSysSolnInfo ( const std::string &  variable,
const std::string &  property 
) const

Definition at line 865 of file SessionReader.cpp.

867  {
868 
869  GloSysSolnInfoList::const_iterator iter =
870  GetGloSysSolnList().find(pVariable);
871  if(iter == GetGloSysSolnList().end())
872  {
873  return false;
874  }
875 
876  std::string vProperty = boost::to_upper_copy(pProperty);
877 
878  GloSysInfoMap::const_iterator iter1 = iter->second.find(vProperty);
879  if(iter1 == iter->second.end())
880  {
881  return false;
882  }
883 
884  return true;
885  }
static GloSysSolnInfoList & GetGloSysSolnList()
GlobalSysSoln Info map.
bool Nektar::LibUtilities::SessionReader::DefinesParameter ( const std::string &  name) const

Checks if a parameter is specified in the XML document.

Definition at line 638 of file SessionReader.cpp.

639  {
640  std::string vName = boost::to_upper_copy(pName);
641  ParameterMap::const_iterator paramIter = m_parameters.find(vName);
642  return (paramIter != m_parameters.end());
643  }
ParameterMap m_parameters
Parameters.
bool Nektar::LibUtilities::SessionReader::DefinesSolverInfo ( const std::string &  name) const

Checks if a solver info property is specified.

Definition at line 761 of file SessionReader.cpp.

Referenced by GetSolverInfoAsEnum().

762  {
763  std::string vName = boost::to_upper_copy(pName);
764  SolverInfoMap::const_iterator infoIter = m_solverInfo.find(vName);
765  return (infoIter != m_solverInfo.end());
766  }
SolverInfoMap m_solverInfo
Solver information properties.
bool Nektar::LibUtilities::SessionReader::DefinesTag ( const std::string &  pName) const

Checks if a specified tag is defined.

Definition at line 1305 of file SessionReader.cpp.

1306  {
1307  std::string vName = boost::to_upper_copy(pName);
1308  TagMap::const_iterator vTagIterator = m_tags.find(vName);
1309  return (vTagIterator != m_tags.end());
1310  }
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 629 of file SessionReader.cpp.

630  {
631  m_comm->Finalise();
632  }
CommSharedPtr m_comm
Communication object.
BndRegionOrdering Nektar::LibUtilities::SessionReader::GetBndRegionOrdering ( ) const

Definition at line 1375 of file SessionReader.cpp.

1376  {
1377  return m_bndRegOrder;
1378  }
BndRegionOrdering m_bndRegOrder
Map of original boundary region ordering for parallel periodic bcs.
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 170 of file SessionReader.cpp.

Referenced by RegisterCmdLineArgument(), and RegisterCmdLineFlag().

171  {
172  static CmdLineArgMap cmdLineArguments;
173  return cmdLineArguments;
174  }
std::map< std::string, CmdLineArg > CmdLineArgMap
Definition: SessionReader.h:75
template<typename T >
T Nektar::LibUtilities::SessionReader::GetCmdLineArgument ( const std::string &  pName) const
inline

Retrieves a command-line argument value.

Definition at line 404 of file SessionReader.h.

References m_cmdLineOptions.

406  {
407  return m_cmdLineOptions.find(pName)->second.as<T>();
408  }
boost::program_options::variables_map m_cmdLineOptions
CommSharedPtr & Nektar::LibUtilities::SessionReader::GetComm ( )

Returns the communication object.

Definition at line 617 of file SessionReader.cpp.

618  {
619  return m_comm;
620  }
CommSharedPtr m_comm
Communication object.
CompositeOrdering Nektar::LibUtilities::SessionReader::GetCompositeOrdering ( ) const

Definition at line 1370 of file SessionReader.cpp.

1371  {
1372  return m_compOrder;
1373  }
CompositeOrdering m_compOrder
Map of original composite ordering for parallel periodic bcs.
TiXmlDocument & Nektar::LibUtilities::SessionReader::GetDocument ( )

Provides direct access to the TiXmlDocument object.

Definition at line 508 of file SessionReader.cpp.

References ASSERTL1.

509  {
510  ASSERTL1(m_xmlDoc, "XML Document not defined.");
511  return *m_xmlDoc;
512  }
TiXmlDocument * m_xmlDoc
Pointer to the loaded XML document.
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:191
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
Note
Paths are case-insensitive.
Parameters
pPathPath to requested element.
Returns
Direct pointer to requested XML Element.

Definition at line 537 of file SessionReader.cpp.

References ASSERTL0.

538  {
539  std::string vPath = boost::to_upper_copy(pPath);
540  std::vector<std::string> st;
541  boost::split(st, vPath, boost::is_any_of("\\/ "));
542  ASSERTL0(st.size() > 0, "No path given in XML element request.");
543 
544  TiXmlElement* vReturn = m_xmlDoc->FirstChildElement(st[0].c_str());
545  ASSERTL0(vReturn, std::string("Cannot find element '")
546  + st[0] + std::string("'."));
547  for (int i = 1; i < st.size(); ++i)
548  {
549  vReturn = vReturn->FirstChildElement(st[i].c_str());
550  ASSERTL0(vReturn, std::string("Cannot find element '")
551  + st[i] + std::string("'."));
552  }
553  return vReturn;
554  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
TiXmlDocument * m_xmlDoc
Pointer to the loaded XML document.
AnalyticExpressionEvaluator & Nektar::LibUtilities::SessionReader::GetExpressionEvaluator ( )

Returns the instance of AnalyticExpressionEvaluator specific to this session.

Definition at line 1296 of file SessionReader.cpp.

1297  {
1298  return m_exprEvaluator;
1299  }
AnalyticExpressionEvaluator m_exprEvaluator
Analytic expression evaluator instance.
const std::vector< std::string > & Nektar::LibUtilities::SessionReader::GetFilenames ( ) const

Returns the filename of the loaded XML document.

Definition at line 582 of file SessionReader.cpp.

583  {
584  return m_filenames;
585  }
std::vector< std::string > m_filenames
Filenames.
const FilterMap & Nektar::LibUtilities::SessionReader::GetFilters ( ) const

Definition at line 1341 of file SessionReader.cpp.

1342  {
1343  return m_filters;
1344  }
FilterMap m_filters
Filters map.
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 1088 of file SessionReader.cpp.

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

1092  {
1093  FunctionMap::const_iterator it1;
1094  FunctionVariableMap::const_iterator it2, it3;
1095  std::string vName = boost::to_upper_copy(pName);
1096 
1097  ASSERTL0((it1 = m_functions.find(vName)) != m_functions.end(),
1098  std::string("No such function '") + pName
1099  + std::string("' has been defined in the session file."));
1100 
1101  // Check for specific and wildcard definitions
1102  pair<std::string,int> key(pVariable,pDomain);
1103  pair<std::string,int> defkey("*",pDomain);
1104  bool specific = (it2 = it1->second.find(key)) !=
1105  it1->second.end();
1106  bool wildcard = (it3 = it1->second.find(defkey)) !=
1107  it1->second.end();
1108 
1109  // Check function is defined somewhere
1110  ASSERTL0(specific || wildcard,
1111  "No such variable " + pVariable
1112  + " in domain " + boost::lexical_cast<string>(pDomain)
1113  + " defined for function " + pName
1114  + " in session file.");
1115 
1116  // If not specific, must be wildcard
1117  if (!specific)
1118  {
1119  it2 = it3;
1120  }
1121 
1122  ASSERTL0((it2->second.m_type == eFunctionTypeExpression),
1123  std::string("Function is defined by a file."));
1124  return it2->second.m_expression;
1125  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
FunctionMap m_functions
Functions.
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 1131 of file SessionReader.cpp.

References ASSERTL0.

1135  {
1136  ASSERTL0(pVar < m_variables.size(), "Variable index out of range.");
1137  return GetFunction(pName, m_variables[pVar],pDomain);
1138  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
EquationSharedPtr GetFunction(const std::string &name, const std::string &variable, const int pDomain=0) const
Returns an EquationSharedPtr to a given function variable.
VariableList m_variables
Variables.
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 1199 of file SessionReader.cpp.

References ASSERTL0.

1203  {
1204  FunctionMap::const_iterator it1;
1205  FunctionVariableMap::const_iterator it2, it3;
1206  std::string vName = boost::to_upper_copy(pName);
1207 
1208  it1 = m_functions.find(vName);
1209  ASSERTL0 (it1 != m_functions.end(),
1210  std::string("Function '") + pName
1211  + std::string("' not found."));
1212 
1213  // Check for specific and wildcard definitions
1214  pair<std::string,int> key(pVariable,pDomain);
1215  pair<std::string,int> defkey("*",pDomain);
1216  bool specific = (it2 = it1->second.find(key)) !=
1217  it1->second.end();
1218  bool wildcard = (it3 = it1->second.find(defkey)) !=
1219  it1->second.end();
1220 
1221  // Check function is defined somewhere
1222  ASSERTL0(specific || wildcard,
1223  "No such variable " + pVariable
1224  + " in domain " + boost::lexical_cast<string>(pDomain)
1225  + " defined for function " + pName
1226  + " in session file.");
1227 
1228  // If not specific, must be wildcard
1229  if (!specific)
1230  {
1231  it2 = it3;
1232  }
1233 
1234  return it2->second.m_filename;
1235  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
FunctionMap m_functions
Functions.
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 1241 of file SessionReader.cpp.

References ASSERTL0.

1245  {
1246  ASSERTL0(pVar < m_variables.size(), "Variable index out of range.");
1247  return GetFunctionFilename(pName, m_variables[pVar],pDomain);
1248  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
VariableList m_variables
Variables.
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.
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 1254 of file SessionReader.cpp.

References ASSERTL0.

1258  {
1259  FunctionMap::const_iterator it1;
1260  FunctionVariableMap::const_iterator it2, it3;
1261  std::string vName = boost::to_upper_copy(pName);
1262 
1263  it1 = m_functions.find(vName);
1264  ASSERTL0 (it1 != m_functions.end(),
1265  std::string("Function '") + pName
1266  + std::string("' not found."));
1267 
1268  // Check for specific and wildcard definitions
1269  pair<std::string,int> key(pVariable,pDomain);
1270  pair<std::string,int> defkey("*",pDomain);
1271  bool specific = (it2 = it1->second.find(key)) !=
1272  it1->second.end();
1273  bool wildcard = (it3 = it1->second.find(defkey)) !=
1274  it1->second.end();
1275 
1276  // Check function is defined somewhere
1277  ASSERTL0(specific || wildcard,
1278  "No such variable " + pVariable
1279  + " in domain " + boost::lexical_cast<string>(pDomain)
1280  + " defined for function " + pName
1281  + " in session file.");
1282 
1283  // If not specific, must be wildcard
1284  if (!specific)
1285  {
1286  it2 = it3;
1287  }
1288 
1289  return it2->second.m_fileVariable;
1290  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
FunctionMap m_functions
Functions.
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 1144 of file SessionReader.cpp.

References ASSERTL0.

1148  {
1149  FunctionMap::const_iterator it1;
1150  FunctionVariableMap::const_iterator it2, it3;
1151  std::string vName = boost::to_upper_copy(pName);
1152 
1153  it1 = m_functions.find(vName);
1154  ASSERTL0 (it1 != m_functions.end(),
1155  std::string("Function '") + pName
1156  + std::string("' not found."));
1157 
1158  // Check for specific and wildcard definitions
1159  pair<std::string,int> key(pVariable,pDomain);
1160  pair<std::string,int> defkey("*",pDomain);
1161  bool specific = (it2 = it1->second.find(key)) !=
1162  it1->second.end();
1163  bool wildcard = (it3 = it1->second.find(defkey)) !=
1164  it1->second.end();
1165 
1166  // Check function is defined somewhere
1167  ASSERTL0(specific || wildcard,
1168  "No such variable " + pVariable
1169  + " in domain " + boost::lexical_cast<string>(pDomain)
1170  + " defined for function " + pName
1171  + " in session file.");
1172 
1173  // If not specific, must be wildcard
1174  if (!specific)
1175  {
1176  it2 = it3;
1177  }
1178 
1179  return it2->second.m_type;
1180  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
FunctionMap m_functions
Functions.
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 1186 of file SessionReader.cpp.

References ASSERTL0.

1190  {
1191  ASSERTL0(pVar < m_variables.size(), "Variable index out of range.");
1192  return GetFunctionType(pName, m_variables[pVar],pDomain);
1193  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
VariableList m_variables
Variables.
enum FunctionType GetFunctionType(const std::string &name, const std::string &variable, const int pDomain=0) const
Returns the type of a given function variable.
const std::string & Nektar::LibUtilities::SessionReader::GetGlobalSysSolnInfo ( const std::string &  variable,
const std::string &  property 
) const

Definition at line 891 of file SessionReader.cpp.

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

892  {
893  GloSysSolnInfoList::const_iterator iter;
894 
895  ASSERTL0( (iter = GetGloSysSolnList().find(pVariable)) !=
896  GetGloSysSolnList().end(),
897  "Failed to find variable in GlobalSysSolnInfoList");
898 
899  std::string vProperty = boost::to_upper_copy(pProperty);
900  GloSysInfoMap::const_iterator iter1;
901 
902  ASSERTL0( (iter1 = iter->second.find(vProperty)) != iter->second.end(),
903  "Failed to find property: " + vProperty + " in GlobalSysSolnInfoList");
904 
905  return iter1->second;
906  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
static GloSysSolnInfoList & GetGloSysSolnList()
GlobalSysSoln Info map.
InputIterator find(InputIterator first, InputIterator last, InputIterator startingpoint, const EqualityComparable &value)
Definition: StdRegions.hpp:315
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 156 of file SessionReader.cpp.

157  {
158  static GloSysSolnInfoList gloSysSolnInfoList;
159  return gloSysSolnInfoList;
160  }
std::map< std::string, GloSysInfoMap > GloSysSolnInfoList
Definition: SessionReader.h:81
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 654 of file SessionReader.cpp.

References ASSERTL0.

656  {
657  std::string vName = boost::to_upper_copy(pName);
658  ParameterMap::const_iterator paramIter = m_parameters.find(vName);
659 
660  ASSERTL0(paramIter != m_parameters.end(),
661  "Unable to find requested parameter: " + pName);
662 
663  return paramIter->second;
664  }
ParameterMap m_parameters
Parameters.
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
const std::string & Nektar::LibUtilities::SessionReader::GetSessionName ( ) const

Returns the session name of the loaded XML document.

Definition at line 591 of file SessionReader.cpp.

592  {
593  return m_sessionName;
594  }
std::string m_sessionName
Session name of the loaded XML document (filename minus ext).
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 601 of file SessionReader.cpp.

References Nektar::LibUtilities::PortablePath().

602  {
603  std::string dirname = m_sessionName + "_xml";
604  fs::path pdirname(dirname);
605 
606  std::string vFilename = "P" + boost::lexical_cast<std::string>(m_comm->GetRowComm()->GetRank());
607  fs::path pFilename(vFilename);
608 
609  fs::path fullpath = pdirname / pFilename;
610 
611  return PortablePath(fullpath);
612  }
CommSharedPtr m_comm
Communication object.
std::string PortablePath(const boost::filesystem::path &path)
create portable path on different platforms for boost::filesystem path
Definition: FileSystem.cpp:41
std::string m_sessionName
Session name of the loaded XML document (filename minus ext).
SessionReaderSharedPtr Nektar::LibUtilities::SessionReader::GetSharedThisPtr ( )
inlineprivate

Returns a shared pointer to the current object.

This allows a member function to pass a shared pointer to itself during a call to another function.

Definition at line 708 of file SessionReader.h.

709  {
710  return shared_from_this();
711  }
const std::string & Nektar::LibUtilities::SessionReader::GetSolverInfo ( const std::string &  pProperty) const

Returns the value of the specified solver info property.

Definition at line 772 of file SessionReader.cpp.

References ASSERTL1.

Referenced by GetSolverInfoAsEnum().

774  {
775  std::string vProperty = boost::to_upper_copy(pProperty);
776  SolverInfoMap::const_iterator iter = m_solverInfo.find(vProperty);
777 
778  ASSERTL1(iter != m_solverInfo.end(),
779  "Unable to find requested property: " + pProperty);
780 
781  return iter->second;
782  }
SolverInfoMap m_solverInfo
Solver information properties.
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:191
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 554 of file SessionReader.h.

References ASSERTL0, DefinesSolverInfo(), Nektar::StdRegions::find(), GetSolverInfo(), GetSolverInfoEnums(), and Nektar::iterator.

556  {
557  std::string vName = boost::to_upper_copy(pName);
559  "Solver info '" + pName + "' not defined.");
560 
561  std::string vValue = GetSolverInfo(vName);
563  ASSERTL0((x = GetSolverInfoEnums().find(vName)) !=
564  GetSolverInfoEnums().end(),
565  "Enum for SolverInfo property '" + pName + "' not found.");
566 
568  ASSERTL0((y = x->second.find(vValue)) != x->second.end(),
569  "Value of SolverInfo property '" + pName +
570  "' is invalid.");
571 
572  return T(y->second);
573  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
const std::string & GetSolverInfo(const std::string &pProperty) const
Returns the value of the specified solver info property.
bool DefinesSolverInfo(const std::string &name) const
Checks if a solver info property is specified.
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
static EnumMapList & GetSolverInfoEnums()
String to enumeration map for Solver Info parameters.
InputIterator find(InputIterator first, InputIterator last, InputIterator startingpoint, const EqualityComparable &value)
Definition: StdRegions.hpp:315
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 140 of file SessionReader.cpp.

Referenced by RegisterDefaultSolverInfo().

141  {
142  static SolverInfoMap solverInfoMap;
143  return solverInfoMap;
144  }
std::map< std::string, std::string > SolverInfoMap
Definition: SessionReader.h:58
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 125 of file SessionReader.cpp.

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

126  {
127  static EnumMapList solverInfoEnums;
128  return solverInfoEnums;
129  }
std::map< std::string, EnumMap > EnumMapList
Definition: SessionReader.h:78
const std::string & Nektar::LibUtilities::SessionReader::GetTag ( const std::string &  pName) const

Returns the value of a specified tag.

Definition at line 1328 of file SessionReader.cpp.

References ASSERTL0.

1329  {
1330  std::string vName = boost::to_upper_copy(pName);
1331  TagMap::const_iterator vTagIterator = m_tags.find(vName);
1332  ASSERTL0(vTagIterator != m_tags.end(),
1333  "Requested tag does not exist.");
1334  return vTagIterator->second;
1335  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
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 581 of file SessionReader.h.

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

584  {
585  std::string vName = boost::to_upper_copy(pName);
586 
588  ASSERTL0((x = GetSolverInfoEnums().find(vName)) !=
589  GetSolverInfoEnums().end(),
590  "Enum for property '" + pName + "' not found.");
591 
593  ASSERTL0((y = x->second.find(pValue)) != x->second.end(),
594  "Value of property '" + pValue + "' is invalid.");
595  return T(y->second);
596  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
static EnumMapList & GetSolverInfoEnums()
String to enumeration map for Solver Info parameters.
InputIterator find(InputIterator first, InputIterator last, InputIterator startingpoint, const EqualityComparable &value)
Definition: StdRegions.hpp:315
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 1014 of file SessionReader.cpp.

References ASSERTL0.

1016  {
1017  ASSERTL0(idx < m_variables.size(), "Variable index out of range.");
1018  return m_variables[idx];
1019  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
VariableList m_variables
Variables.
std::vector< std::string > Nektar::LibUtilities::SessionReader::GetVariables ( ) const

Returns the names of all variables.

Definition at line 1037 of file SessionReader.cpp.

1038  {
1039  return m_variables;
1040  }
VariableList m_variables
Variables.
void Nektar::LibUtilities::SessionReader::InitSession ( )
private

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 265 of file SessionReader.cpp.

References Nektar::iterator.

266  {
267  m_exprEvaluator.SetRandomSeed((m_comm->GetRank() + 1) * time(NULL));
268 
269  // Split up the communicator
270  PartitionComm();
271 
272  // Partition mesh
273  PartitionMesh();
274 
275  // Parse the XML data in #m_xmlDoc
276  ParseDocument();
277 
278  // Override SOLVERINFO and parameters with any specified on the
279  // command line.
280  CmdLineOverride();
281 
282  // Verify SOLVERINFO values
284 
285  // In verbose mode, print out parameters and solver info sections
286  if (m_verbose && m_comm)
287  {
288  if (m_comm->GetRank() == 0 && m_parameters.size() > 0)
289  {
290  cout << "Parameters:" << endl;
292  for (x = m_parameters.begin(); x != m_parameters.end(); ++x)
293  {
294  cout << "\t" << x->first << " = " << x->second << endl;
295  }
296  cout << endl;
297  }
298 
299  if (m_comm->GetRank() == 0 && m_solverInfo.size() > 0)
300  {
301  cout << "Solver Info:" << endl;
303  for (x = m_solverInfo.begin(); x != m_solverInfo.end(); ++x)
304  {
305  cout << "\t" << x->first << " = " << x->second << endl;
306  }
307  cout << endl;
308  }
309  }
310  }
ParameterMap m_parameters
Parameters.
void CmdLineOverride()
Enforce parameters from command line arguments.
SolverInfoMap m_solverInfo
Solver information properties.
AnalyticExpressionEvaluator m_exprEvaluator
Analytic expression evaluator instance.
void PartitionComm()
Partitions the comm object based on session parameters.
CommSharedPtr m_comm
Communication object.
void VerifySolverInfo()
Check values of solver info options are valid.
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
void ParseDocument()
Loads and parses the specified file.
void PartitionMesh()
Partitions the mesh when running in parallel.
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 1383 of file SessionReader.cpp.

References ASSERTL0, and Nektar::LibUtilities::PortablePath().

1386  {
1387  if (pFilename.size() > 3 &&
1388  pFilename.substr(pFilename.size() - 3, 3) == ".gz")
1389  {
1390  ifstream file(pFilename.c_str(),
1391  ios_base::in | ios_base::binary);
1392  ASSERTL0(file.good(), "Unable to open file: " + pFilename);
1393  stringstream ss;
1394  io::filtering_streambuf<io::input> in;
1395  in.push(io::gzip_decompressor());
1396  in.push(file);
1397  try
1398  {
1399  io::copy(in, ss);
1400  ss >> (*pDoc);
1401  }
1402  catch (io::gzip_error& e)
1403  {
1404  ASSERTL0(false,
1405  "Error: File '" + pFilename + "' is corrupt.");
1406  }
1407  }
1408  else if (pFilename.size() > 4 &&
1409  pFilename.substr(pFilename.size() - 4, 4) == "_xml")
1410  {
1411  fs::path pdirname(pFilename);
1412  boost::format pad("P%1$07d.xml");
1413  pad % m_comm->GetRank();
1414  fs::path pRankFilename(pad.str());
1415  fs::path fullpath = pdirname / pRankFilename;
1416 
1417  ifstream file(PortablePath(fullpath).c_str());
1418  ASSERTL0(file.good(), "Unable to open file: " + fullpath.string());
1419  file >> (*pDoc);
1420  }
1421  else
1422  {
1423  ifstream file(pFilename.c_str());
1424  ASSERTL0(file.good(), "Unable to open file: " + pFilename);
1425  file >> (*pDoc);
1426  }
1427  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
CommSharedPtr m_comm
Communication object.
std::string PortablePath(const boost::filesystem::path &path)
create portable path on different platforms for boost::filesystem path
Definition: FileSystem.cpp:41
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 922 of file SessionReader.cpp.

926  {
927  std::string vName = boost::to_upper_copy(pName);
928  GeometricInfoMap::const_iterator iter = m_geometricInfo.find(vName);
929  if(iter != m_geometricInfo.end())
930  {
931  pVar = iter->second;
932  }
933  else
934  {
935  pVar = pDefault;
936  }
937  }
GeometricInfoMap m_geometricInfo
Geometric information properties.
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 943 of file SessionReader.cpp.

947  {
948  std::string vName = boost::to_upper_copy(pName);
949  GeometricInfoMap::const_iterator iter = m_geometricInfo.find(vName);
950  if(iter != m_geometricInfo.end())
951  {
952  if (iter->second == "TRUE")
953  {
954  pVar = true;
955  }
956  else
957  {
958  pVar = false;
959  }
960  }
961  else
962  {
963  pVar = pDefault;
964  }
965  }
GeometricInfoMap m_geometricInfo
Geometric information properties.
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 971 of file SessionReader.cpp.

975  {
976  std::string vName = boost::to_upper_copy(pName);
977  GeometricInfoMap::const_iterator iter = m_geometricInfo.find(vName);
978  if(iter != m_geometricInfo.end())
979  {
980  pVar = std::atoi(iter->second.c_str());
981  }
982  else
983  {
984  pVar = pDefault;
985  }
986  }
GeometricInfoMap m_geometricInfo
Geometric information properties.
void Nektar::LibUtilities::SessionReader::LoadParameter ( const std::string &  name,
int &  var 
) const

Load an integer parameter.

Definition at line 670 of file SessionReader.cpp.

References ASSERTL0.

672  {
673  std::string vName = boost::to_upper_copy(pName);
674  ParameterMap::const_iterator paramIter = m_parameters.find(vName);
675  ASSERTL0(paramIter != m_parameters.end(), "Required parameter '" +
676  pName + "' not specified in session.");
677  pVar = (int)floor(paramIter->second);
678  }
ParameterMap m_parameters
Parameters.
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
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 684 of file SessionReader.cpp.

686  {
687  std::string vName = boost::to_upper_copy(pName);
688  ParameterMap::const_iterator paramIter = m_parameters.find(vName);
689  if(paramIter != m_parameters.end())
690  {
691  pVar = (int)floor(paramIter->second);
692  }
693  else
694  {
695  pVar = pDefault;
696  }
697  }
ParameterMap m_parameters
Parameters.
void Nektar::LibUtilities::SessionReader::LoadParameter ( const std::string &  name,
NekDouble var 
) const

Load a double precision parameter.

Definition at line 703 of file SessionReader.cpp.

References ASSERTL0.

705  {
706  std::string vName = boost::to_upper_copy(pName);
707  ParameterMap::const_iterator paramIter = m_parameters.find(vName);
708  ASSERTL0(paramIter != m_parameters.end(), "Required parameter '" +
709  pName + "' not specified in session.");
710  pVar = paramIter->second;
711  }
ParameterMap m_parameters
Parameters.
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
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 717 of file SessionReader.cpp.

721  {
722  std::string vName = boost::to_upper_copy(pName);
723  ParameterMap::const_iterator paramIter = m_parameters.find(vName);
724  if(paramIter != m_parameters.end())
725  {
726  pVar = paramIter->second;
727  }
728  else
729  {
730  pVar = pDefault;
731  }
732  }
ParameterMap m_parameters
Parameters.
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 802 of file SessionReader.cpp.

806  {
807  std::string vName = boost::to_upper_copy(pName);
808  SolverInfoMap::const_iterator infoIter = m_solverInfo.find(vName);
809  if(infoIter != m_solverInfo.end())
810  {
811  pVar = infoIter->second;
812  }
813  else
814  {
815  pVar = pDefault;
816  }
817  }
SolverInfoMap m_solverInfo
Solver information properties.
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 992 of file SessionReader.cpp.

997  {
998  std::string vName = boost::to_upper_copy(pName);
999  GeometricInfoMap::const_iterator iter = m_geometricInfo.find(vName);
1000  if(iter != m_geometricInfo.end())
1001  {
1002  pVar = boost::iequals(iter->second, pTrueVal);
1003  }
1004  else
1005  {
1006  pVar = pDefault;
1007  }
1008  }
GeometricInfoMap m_geometricInfo
Geometric information properties.
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 823 of file SessionReader.cpp.

828  {
829  std::string vName = boost::to_upper_copy(pName);
830  SolverInfoMap::const_iterator infoIter = m_solverInfo.find(vName);
831  if(infoIter != m_solverInfo.end())
832  {
833  pVar = boost::iequals(infoIter->second, pTrueVal);
834  }
835  else
836  {
837  pVar = pDefault;
838  }
839  }
SolverInfoMap m_solverInfo
Solver information properties.
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 845 of file SessionReader.cpp.

848  {
849  if (DefinesSolverInfo(pName))
850  {
851  std::string vName = boost::to_upper_copy(pName);
852  SolverInfoMap::const_iterator iter = m_solverInfo.find(vName);
853  if(iter != m_solverInfo.end())
854  {
855  return true;
856  }
857  }
858  return false;
859  }
SolverInfoMap m_solverInfo
Solver information properties.
bool DefinesSolverInfo(const std::string &name) const
Checks if a solver info property is specified.
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 543 of file SessionReader.h.

545  {
546  return (GetSolverInfoAsEnum<T>(name) == trueval);
547  }
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 1432 of file SessionReader.cpp.

References ASSERTL0, and WARNINGL0.

1434  {
1435  ASSERTL0(pFilenames.size() > 0, "No filenames for merging.");
1436 
1437  // Read the first document
1438  TiXmlDocument *vMainDoc = new TiXmlDocument;
1439  LoadDoc(pFilenames[0], vMainDoc);
1440 
1441  TiXmlHandle vMainHandle(vMainDoc);
1442  TiXmlElement* vMainNektar =
1443  vMainHandle.FirstChildElement("NEKTAR").Element();
1444 
1445  // Read all subsequent XML documents.
1446  // For each element within the NEKTAR tag, use it to replace the
1447  // version already present in the loaded XML data.
1448  for (int i = 1; i < pFilenames.size(); ++i)
1449  {
1450  if((pFilenames[i].compare(pFilenames[i].size()-3,3,"xml") == 0)
1451  ||(pFilenames[i].compare(pFilenames[i].size()-6,6,"xml.gz") == 0))
1452  {
1453  TiXmlDocument* vTempDoc = new TiXmlDocument;
1454  LoadDoc(pFilenames[i], vTempDoc);
1455 
1456  TiXmlHandle docHandle(vTempDoc);
1457  TiXmlElement* vTempNektar;
1458  vTempNektar = docHandle.FirstChildElement("NEKTAR").Element();
1459  ASSERTL0(vTempNektar, "Unable to find NEKTAR tag in file.");
1460  TiXmlElement* p = vTempNektar->FirstChildElement();
1461 
1462  while (p)
1463  {
1464  TiXmlElement *vMainEntry =
1465  vMainNektar->FirstChildElement(p->Value());
1466 
1467  // First check if the new item is in fact blank
1468  if (!p->FirstChild() && vMainEntry)
1469  {
1470  std::string warningmsg =
1471  "File " + pFilenames[i] + " contains " +
1472  "an empty XML element " +
1473  std::string(p->Value()) +
1474  " which will be ignored.";
1475  WARNINGL0(false, warningmsg.c_str());
1476  }
1477  else
1478  {
1479  if (vMainEntry)
1480  {
1481  vMainNektar->RemoveChild(vMainEntry);
1482  }
1483  TiXmlElement *q = new TiXmlElement(*p);
1484  vMainNektar->LinkEndChild(q);
1485  }
1486  p = p->NextSiblingElement();
1487  }
1488 
1489  delete vTempDoc;
1490  }
1491  }
1492  return vMainDoc;
1493  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
void LoadDoc(const std::string &pFilename, TiXmlDocument *pDoc) const
Loads an xml file into a tinyxml doc and decompresses if needed.
#define WARNINGL0(condition, msg)
Definition: ErrorUtil.hpp:167
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 317 of file SessionReader.cpp.

References Nektar::iterator, Nektar::NekConstants::kGitBranch, Nektar::NekConstants::kGitSha1, and NEKTAR_VERSION.

319  {
320  // List the publically visible options (listed using --help)
321  po::options_description desc("Allowed options");
322  desc.add_options()
323  ("verbose,v", "be verbose")
324  ("version,V", "print version information")
325  ("help,h", "print this help message")
326  ("solverinfo,I", po::value<vector<std::string> >(),
327  "override a SOLVERINFO property")
328  ("parameter,P", po::value<vector<std::string> >(),
329  "override a parameter")
330  ("shared-filesystem,s", "Using shared filesystem.")
331  ("npx", po::value<int>(),
332  "number of procs in X-dir")
333  ("npy", po::value<int>(),
334  "number of procs in Y-dir")
335  ("npz", po::value<int>(),
336  "number of procs in Z-dir")
337  ("nsz", po::value<int>(),
338  "number of slices in Z-dir")
339  ("part-only", po::value<int>(),
340  "only partition mesh into N partitions.")
341  ("part-only-overlapping", po::value<int>(),
342  "only partition mesh into N overlapping partitions.")
343  ("part-info", "Output partition information")
344  ;
345 
346  CmdLineArgMap::const_iterator cmdIt;
347  for (cmdIt = GetCmdLineArgMap().begin();
348  cmdIt != GetCmdLineArgMap().end(); ++cmdIt)
349  {
350  std::string names = cmdIt->first;
351  if (cmdIt->second.shortName != "")
352  {
353  names += "," + cmdIt->second.shortName;
354  }
355  if (cmdIt->second.isFlag)
356  {
357  desc.add_options()
358  (names.c_str(), cmdIt->second.description.c_str())
359  ;
360  }
361  else
362  {
363  desc.add_options()
364  (names.c_str(), po::value<std::string>(),
365  cmdIt->second.description.c_str())
366  ;
367  }
368  }
369 
370  // List hidden options (e.g. session file arguments are not actually
371  // specified using the input-file option by the user).
372  po::options_description hidden("Hidden options");
373  hidden.add_options()
374  ("input-file", po::value< vector<string> >(),
375  "input filename")
376  ;
377 
378  // Combine all options for the parser
379  po::options_description all("All options");
380  all.add(desc).add(hidden);
381 
382  // Session file is a positional option
383  po::positional_options_description p;
384  p.add("input-file", -1);
385 
386  // Parse the command-line options
387  po::parsed_options parsed = po::command_line_parser(argc, argv).
388  options(all).
389  positional(p).
390  allow_unregistered().
391  run();
392 
393  // Extract known options to map and update
394  po::store(parsed, m_cmdLineOptions);
395  po::notify(m_cmdLineOptions);
396 
397  // Help message
398  if (m_cmdLineOptions.count("help"))
399  {
400  cout << desc;
401  exit(0);
402  }
403 
404  // Version information
405  if (m_cmdLineOptions.count("version"))
406  {
407  cout << "Nektar++ version " << NEKTAR_VERSION;
408 
409  if (NekConstants::kGitSha1 != "GITDIR-NOTFOUND")
410  {
411  string sha1(NekConstants::kGitSha1);
412  string branch(NekConstants::kGitBranch);
413  boost::replace_all(branch, "refs/heads/", "");
414 
415  cout << " (git changeset " << sha1.substr(0, 8) << ", ";
416 
417  if (branch == "")
418  {
419  cout << "detached head";
420  }
421  else
422  {
423  cout << "head " << branch;
424  }
425 
426  cout << ")";
427  }
428 
429  cout << endl;
430  exit(0);
431  }
432 
433  // Enable verbose mode
434  if (m_cmdLineOptions.count("verbose"))
435  {
436  m_verbose = true;
437  }
438  else
439  {
440  m_verbose = false;
441  }
442 
443  // Print a warning for unknown options
444  std::vector< po::basic_option<char> >::iterator x;
445  for (x = parsed.options.begin(); x != parsed.options.end(); ++x)
446  {
447  if (x->unregistered)
448  {
449  cout << "Warning: Unknown option: " << x->string_key
450  << endl;
451  }
452  }
453 
454  // Return the vector of filename(s) given as positional options
455  if (m_cmdLineOptions.count("input-file"))
456  {
457  return m_cmdLineOptions["input-file"].as<
458  std::vector<std::string> >();
459  }
460  else
461  {
462  return std::vector<std::string>();
463  }
464  }
#define NEKTAR_VERSION
const std::string kGitSha1
boost::program_options::variables_map m_cmdLineOptions
const std::string kGitBranch
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
static CmdLineArgMap & GetCmdLineArgMap()
CmdLine argument map.
void Nektar::LibUtilities::SessionReader::ParseDocument ( )
private

Loads and parses the specified file.

Definition at line 1499 of file SessionReader.cpp.

References ASSERTL0.

1500  {
1501  // Check we actually have a document loaded.
1502  ASSERTL0(m_xmlDoc, "No XML document loaded.");
1503 
1504  // Look for all data in CONDITIONS block.
1505  TiXmlHandle docHandle(m_xmlDoc);
1506  TiXmlElement* e;
1507  e = docHandle.FirstChildElement("NEKTAR").
1508  FirstChildElement("CONDITIONS").Element();
1509 
1510  // Read the various sections of the CONDITIONS block
1511  ReadParameters (e);
1512  ReadSolverInfo (e);
1514  ReadExpressions (e);
1515  ReadVariables (e);
1516  ReadFunctions (e);
1517 
1518  e = docHandle.FirstChildElement("NEKTAR").
1519  FirstChildElement("FILTERS").Element();
1520 
1521  ReadFilters(e);
1522  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
void ReadExpressions(TiXmlElement *conditions)
Reads the EXPRESSIONS section of the XML document.
void ReadSolverInfo(TiXmlElement *conditions)
Reads the SOLVERINFO 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 ReadVariables(TiXmlElement *conditions)
Reads the VARIABLES section of the XML document.
void ReadFilters(TiXmlElement *filters)
Reads the FILTERS section of the XML document.
TiXmlDocument * m_xmlDoc
Pointer to the loaded XML document.
void ReadParameters(TiXmlElement *conditions)
Reads the PARAMETERS section of the XML document.
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 2562 of file SessionReader.cpp.

2566  {
2567  /// Pull out lhs and rhs and eliminate any spaces.
2568  int beg = line.find_first_not_of(" ");
2569  int end = line.find_first_of("=");
2570  // Check for no parameter name
2571  if (beg == end) throw 1;
2572  // Check for no parameter value
2573  if (end != line.find_last_of("=")) throw 1;
2574  // Check for no equals sign
2575  if (end == std::string::npos) throw 1;
2576 
2577  lhs = line.substr(line.find_first_not_of(" "),
2578  end-beg);
2579  lhs = lhs .substr(0, lhs.find_last_not_of(" ")+1);
2580  rhs = line.substr(line.find_last_of("=")+1);
2581  rhs = rhs .substr(rhs.find_first_not_of(" "));
2582  rhs = rhs .substr(0, rhs.find_last_not_of(" ")+1);
2583  }
StandardMatrixTag & lhs
std::string Nektar::LibUtilities::SessionReader::ParseSessionName ( std::vector< std::string > &  filenames)
private

Parse the session name.

Definition at line 470 of file SessionReader.cpp.

References ASSERTL0.

472  {
473  ASSERTL0(!filenames.empty(),
474  "At least one filename expected.");
475 
476  std::string retval = "";
477 
478  // First input file defines the session name
479  std::string fname = filenames[0];
480 
481  // If loading a pre-partitioned mesh, remove _xml extension
482  if (fname.size() > 4 &&
483  fname.substr(fname.size() - 4, 4) == "_xml")
484  {
485  retval = fname.substr(0, fname.find_last_of("_"));
486  }
487  // otherwise remove the .xml extension
488  else if (fname.size() > 4 &&
489  fname.substr(fname.size() - 4, 4) == ".xml")
490  {
491  retval = fname.substr(0, fname.find_last_of("."));
492  }
493  // If compressed .xml.gz, remove both extensions
494  else if (fname.size() > 7 &&
495  fname.substr(fname.size() - 7, 7) == ".xml.gz")
496  {
497  retval = fname.substr(0, fname.find_last_of("."));
498  retval = retval.substr(0, retval.find_last_of("."));
499  }
500 
501  return retval;
502  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
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 1842 of file SessionReader.cpp.

References ASSERTL0.

1843  {
1844  if (m_comm->GetSize() > 1)
1845  {
1846  int nProcZ = 1;
1847  int nProcY = 1;
1848  int nProcX = 1;
1849  int nStripZ = 1;
1850  if (DefinesCmdLineArgument("npx")) {
1851  nProcX = GetCmdLineArgument<int>("npx");
1852  }
1853  if (DefinesCmdLineArgument("npy")) {
1854  nProcY = GetCmdLineArgument<int>("npy");
1855  }
1856  if (DefinesCmdLineArgument("npz")) {
1857  nProcZ = GetCmdLineArgument<int>("npz");
1858  }
1859  if (DefinesCmdLineArgument("nsz")) {
1860  nStripZ = GetCmdLineArgument<int>("nsz");
1861  }
1862  ASSERTL0(m_comm->GetSize() % (nProcZ*nProcY*nProcX) == 0,
1863  "Cannot exactly partition using PROC_Z value.");
1864  ASSERTL0(nProcZ % nProcY == 0,
1865  "Cannot exactly partition using PROC_Y value.");
1866  ASSERTL0(nProcY % nProcX == 0,
1867  "Cannot exactly partition using PROC_X value.");
1868 
1869  // Number of processes associated with the spectral method
1870  int nProcSm = nProcZ * nProcY * nProcX;
1871 
1872  // Number of processes associated with the spectral element
1873  // method.
1874  int nProcSem = m_comm->GetSize() / nProcSm;
1875 
1876  m_comm->SplitComm(nProcSm,nProcSem);
1877  m_comm->GetColumnComm()->SplitComm(nProcZ/nStripZ,nStripZ);
1878  m_comm->GetColumnComm()->GetColumnComm()->SplitComm(
1879  (nProcY*nProcX),nProcZ/nStripZ);
1880  m_comm->GetColumnComm()->GetColumnComm()->GetColumnComm()
1881  ->SplitComm(nProcX,nProcY);
1882  }
1883  }
bool DefinesCmdLineArgument(const std::string &pName) const
Checks if a specified cmdline argument has been given.
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
CommSharedPtr m_comm
Communication object.
void Nektar::LibUtilities::SessionReader::PartitionMesh ( )
private

Partitions the mesh when running in parallel.

Definition at line 1552 of file SessionReader.cpp.

References ASSERTL0, Nektar::LibUtilities::NekFactory< tKey, tBase, >::CreateInstance(), Gs::Finalise(), Nektar::LibUtilities::GetMeshPartitionFactory(), Nektar::LibUtilities::GetSize(), Nektar::iterator, Nektar::LibUtilities::PortablePath(), and Nektar::LibUtilities::ReduceMax.

1553  {
1554  ASSERTL0(m_comm.get(), "Communication not initialised.");
1555 
1556  // Get row of comm, or the whole comm if not split
1557  CommSharedPtr vCommMesh = m_comm->GetRowComm();
1558  const bool isRoot = (m_comm->GetRank() == 0);
1559 
1560  // Delete any existing loaded mesh
1561  if (m_xmlDoc)
1562  {
1563  delete m_xmlDoc;
1564  }
1565 
1566  // Load file for root process only (since this is always needed)
1567  // and determine if the provided geometry has already been
1568  // partitioned. This will be the case if the user provides the
1569  // directory of mesh partitions as an input. Partitioned geometries
1570  // have the attribute
1571  // PARTITION=X
1572  // where X is the number of the partition (and should match the
1573  // process rank). The result is shared with all other processes.
1574  int isPartitioned = 0;
1575  if (isRoot)
1576  {
1578  if (DefinesElement("Nektar/Geometry"))
1579  {
1580  if (GetElement("Nektar/Geometry")->Attribute("PARTITION"))
1581  {
1582  cout << "Using pre-partitioned mesh." << endl;
1583  isPartitioned = 1;
1584  }
1585  }
1586  }
1587  GetComm()->AllReduce(isPartitioned, LibUtilities::ReduceMax);
1588 
1589  // If the mesh is already partitioned, we are done. Remaining
1590  // processes must load their partitions.
1591  if (isPartitioned) {
1592  if (!isRoot)
1593  {
1595  }
1596  return;
1597  }
1598 
1599  // Default partitioner to use is Metis. Use Scotch as default
1600  // if it is installed. Override default with command-line flags
1601  // if they are set.
1602  string vPartitionerName = "Metis";
1603  if (GetMeshPartitionFactory().ModuleExists("Scotch"))
1604  {
1605  vPartitionerName = "Scotch";
1606  }
1607  if (DefinesCmdLineArgument("use-metis"))
1608  {
1609  vPartitionerName = "Metis";
1610  }
1611  if (DefinesCmdLineArgument("use-scotch"))
1612  {
1613  vPartitionerName = "Scotch";
1614  }
1615 
1616  // Mesh has not been partitioned so do partitioning if required.
1617  // Note in the serial case nothing is done as we have already loaded
1618  // the mesh.
1619  if (DefinesCmdLineArgument("part-only")||
1620  DefinesCmdLineArgument("part-only-overlapping"))
1621  {
1622  // Perform partitioning of the mesh only. For this we insist
1623  // the code is run in serial (parallel execution is pointless).
1624  ASSERTL0(GetComm()->GetSize() == 1,
1625  "The 'part-only' option should be used in serial.");
1626 
1627  // Number of partitions is specified by the parameter.
1628  int nParts;
1630  MeshPartitionSharedPtr vPartitioner =
1632  vPartitionerName, vSession);
1633  if(DefinesCmdLineArgument("part-only"))
1634  {
1635  nParts = GetCmdLineArgument<int>("part-only");
1636  vPartitioner->PartitionMesh(nParts, true);
1637  }
1638  else
1639  {
1640  nParts = GetCmdLineArgument<int>("part-only-overlapping");
1641  vPartitioner->PartitionMesh(nParts, true, true);
1642  }
1643  vPartitioner->WriteAllPartitions(vSession);
1644  vPartitioner->GetCompositeOrdering(m_compOrder);
1645  vPartitioner->GetBndRegionOrdering(m_bndRegOrder);
1646 
1647  if (isRoot && DefinesCmdLineArgument("part-info"))
1648  {
1649  vPartitioner->PrintPartInfo(std::cout);
1650  }
1651 
1652  Finalise();
1653  exit(0);
1654  }
1655  else if (vCommMesh->GetSize() > 1)
1656  {
1658  int nParts = vCommMesh->GetSize();
1659  if (DefinesCmdLineArgument("shared-filesystem"))
1660  {
1661  CommSharedPtr vComm = GetComm();
1662  vector<unsigned int> keys, vals;
1663  int i;
1664 
1665  if (vComm->GetRank() == 0)
1666  {
1668 
1669  MeshPartitionSharedPtr vPartitioner =
1671  vPartitionerName, vSession);
1672  vPartitioner->PartitionMesh(nParts, true);
1673  vPartitioner->WriteAllPartitions(vSession);
1674  vPartitioner->GetCompositeOrdering(m_compOrder);
1675  vPartitioner->GetBndRegionOrdering(m_bndRegOrder);
1676 
1677  // Communicate orderings to the other processors.
1678 
1679  // First send sizes of the orderings and boundary
1680  // regions to allocate storage on the remote end.
1681  keys.resize(2);
1682  keys[0] = m_compOrder.size();
1683  keys[1] = m_bndRegOrder.size();
1684 
1685  for (i = 1; i < vComm->GetSize(); ++i)
1686  {
1687  vComm->Send(i, keys);
1688  }
1689 
1690  // Construct the keys and sizes of values for composite
1691  // ordering
1693  keys.resize(m_compOrder.size());
1694  vals.resize(m_compOrder.size());
1695 
1696  for (cIt = m_compOrder.begin(), i = 0;
1697  cIt != m_compOrder.end(); ++cIt, ++i)
1698  {
1699  keys[i] = cIt->first;
1700  vals[i] = cIt->second.size();
1701  }
1702 
1703  // Send across data.
1704  for (i = 1; i < vComm->GetSize(); ++i)
1705  {
1706  vComm->Send(i, keys);
1707  vComm->Send(i, vals);
1708 
1709  for (cIt = m_compOrder.begin();
1710  cIt != m_compOrder.end(); ++cIt)
1711  {
1712  vComm->Send(i, cIt->second);
1713  }
1714  }
1715 
1716  // Construct the keys and sizes of values for composite
1717  // ordering
1719  keys.resize(m_bndRegOrder.size());
1720  vals.resize(m_bndRegOrder.size());
1721 
1722  for (bIt = m_bndRegOrder.begin(), i = 0;
1723  bIt != m_bndRegOrder.end(); ++bIt, ++i)
1724  {
1725  keys[i] = bIt->first;
1726  vals[i] = bIt->second.size();
1727  }
1728 
1729  // Send across data.
1730  for (i = 1; i < vComm->GetSize(); ++i)
1731  {
1732  vComm->Send(i, keys);
1733  vComm->Send(i, vals);
1734 
1735  for (bIt = m_bndRegOrder.begin();
1736  bIt != m_bndRegOrder.end(); ++bIt)
1737  {
1738  vComm->Send(i, bIt->second);
1739  }
1740  }
1741 
1742  if (DefinesCmdLineArgument("part-info"))
1743  {
1744  vPartitioner->PrintPartInfo(std::cout);
1745  }
1746  }
1747  else
1748  {
1749  keys.resize(2);
1750  vComm->Recv(0, keys);
1751 
1752  int cmpSize = keys[0];
1753  int bndSize = keys[1];
1754 
1755  keys.resize(cmpSize);
1756  vals.resize(cmpSize);
1757  vComm->Recv(0, keys);
1758  vComm->Recv(0, vals);
1759 
1760  for (int i = 0; i < keys.size(); ++i)
1761  {
1762  vector<unsigned int> tmp(vals[i]);
1763  vComm->Recv(0, tmp);
1764  m_compOrder[keys[i]] = tmp;
1765  }
1766 
1767  keys.resize(bndSize);
1768  vals.resize(bndSize);
1769  vComm->Recv(0, keys);
1770  vComm->Recv(0, vals);
1771 
1772  for (int i = 0; i < keys.size(); ++i)
1773  {
1774  vector<unsigned int> tmp(vals[i]);
1775  vComm->Recv(0, tmp);
1776  m_bndRegOrder[keys[i]] = tmp;
1777  }
1778  }
1779  }
1780  else
1781  {
1782  // Need to load mesh on non-root processes.
1783  if (!isRoot)
1784  {
1786  }
1787 
1788  // Partitioner now operates in parallel
1789  // Each process receives partitioning over interconnect
1790  // and writes its own session file to the working directory.
1791  MeshPartitionSharedPtr vPartitioner =
1793  vPartitionerName, vSession);
1794  vPartitioner->PartitionMesh(nParts, false);
1795  vPartitioner->WriteLocalPartition(vSession);
1796  vPartitioner->GetCompositeOrdering(m_compOrder);
1797  vPartitioner->GetBndRegionOrdering(m_bndRegOrder);
1798 
1799  if (DefinesCmdLineArgument("part-info") && isRoot)
1800  {
1801  vPartitioner->PrintPartInfo(std::cout);
1802  }
1803  }
1804  m_comm->Block();
1805 
1806  std::string dirname = GetSessionName() + "_xml";
1807  fs::path pdirname(dirname);
1808  boost::format pad("P%1$07d.xml");
1809  pad % m_comm->GetRowComm()->GetRank();
1810  fs::path pFilename(pad.str());
1811  fs::path fullpath = pdirname / pFilename;
1812 
1813  std::string vFilename = PortablePath(fullpath);
1814 
1815  if (m_xmlDoc)
1816  {
1817  delete m_xmlDoc;
1818  }
1819  m_xmlDoc = new TiXmlDocument(vFilename);
1820 
1821  ASSERTL0(m_xmlDoc, "Failed to create XML document object.");
1822 
1823  bool loadOkay = m_xmlDoc->LoadFile(vFilename);
1824  ASSERTL0(loadOkay, "Unable to load file: " + vFilename +
1825  ". Check XML standards compliance. Error on line: " +
1826  boost::lexical_cast<std::string>(m_xmlDoc->Row()));
1827  }
1828  else
1829  {
1831  }
1832  }
bool DefinesCmdLineArgument(const std::string &pName) const
Checks if a specified cmdline argument has been given.
TiXmlElement * GetElement(const std::string &pPath)
Provides direct access to the TiXmlElement specified.
TiXmlDocument * MergeDoc(const std::vector< std::string > &pFilenames) const
Creates an XML document from a list of input files.
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
CommSharedPtr & GetComm()
Returns the communication object.
tBaseSharedPtr CreateInstance(tKey idKey BOOST_PP_COMMA_IF(MAX_PARAM) BOOST_PP_ENUM_BINARY_PARAMS(MAX_PARAM, tParam, x))
Create an instance of the class referred to by idKey.
Definition: NekFactory.hpp:162
boost::shared_ptr< MeshPartition > MeshPartitionSharedPtr
SessionReaderSharedPtr GetSharedThisPtr()
Returns a shared pointer to the current object.
CommSharedPtr m_comm
Communication object.
boost::shared_ptr< SessionReader > SessionReaderSharedPtr
Definition: MeshPartition.h:51
boost::shared_ptr< Comm > CommSharedPtr
Pointer to a Communicator object.
Definition: Comm.h:53
void Finalise()
Finalises the session.
const std::string & GetSessionName() const
Returns the session name of the loaded XML document.
std::vector< std::string > m_filenames
Filenames.
int GetSize(const Array< OneD, const NekDouble > &x)
Definition: NodalUtil.cpp:111
std::string PortablePath(const boost::filesystem::path &path)
create portable path on different platforms for boost::filesystem path
Definition: FileSystem.cpp:41
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
TiXmlDocument * m_xmlDoc
Pointer to the loaded XML document.
bool DefinesElement(const std::string &pPath) const
Tests if a specified element is defined in the XML document.
CompositeOrdering m_compOrder
Map of original composite ordering for parallel periodic bcs.
MeshPartitionFactory & GetMeshPartitionFactory()
BndRegionOrdering m_bndRegOrder
Map of original boundary region ordering for parallel periodic bcs.
void Nektar::LibUtilities::SessionReader::ReadExpressions ( TiXmlElement *  conditions)
private

Reads the EXPRESSIONS section of the XML document.

Definition at line 2183 of file SessionReader.cpp.

References ASSERTL0, and Nektar::iterator.

2184  {
2185  m_expressions.clear();
2186 
2187  if (!conditions)
2188  {
2189  return;
2190  }
2191 
2192  TiXmlElement *expressionsElement =
2193  conditions->FirstChildElement("EXPRESSIONS");
2194 
2195  if (expressionsElement)
2196  {
2197  TiXmlElement *expr = expressionsElement->FirstChildElement("E");
2198 
2199  while (expr)
2200  {
2201  stringstream tagcontent;
2202  tagcontent << *expr;
2203  ASSERTL0(expr->Attribute("NAME"),
2204  "Missing NAME attribute in expression "
2205  "definition: \n\t'" + tagcontent.str() + "'");
2206  std::string nameString = expr->Attribute("NAME");
2207  ASSERTL0(!nameString.empty(),
2208  "Expressions must have a non-empty name: \n\t'"
2209  + tagcontent.str() + "'");
2210 
2211  ASSERTL0(expr->Attribute("VALUE"),
2212  "Missing VALUE attribute in expression "
2213  "definition: \n\t'" + tagcontent.str() + "'");
2214  std::string valString = expr->Attribute("VALUE");
2215  ASSERTL0(!valString.empty(),
2216  "Expressions must have a non-empty value: \n\t'"
2217  + tagcontent.str() + "'");
2218 
2219  ExpressionMap::iterator exprIter
2220  = m_expressions.find(nameString);
2221  ASSERTL0(exprIter == m_expressions.end(),
2222  std::string("Expression '") + nameString
2223  + std::string("' already specified."));
2224 
2225  m_expressions[nameString] = valString;
2226  expr = expr->NextSiblingElement("E");
2227  }
2228  }
2229  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
ExpressionMap m_expressions
Expressions.
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
void Nektar::LibUtilities::SessionReader::ReadFilters ( TiXmlElement *  filters)
private

Reads the FILTERS section of the XML document.

Definition at line 2521 of file SessionReader.cpp.

References ASSERTL0.

2522  {
2523  if (!filters)
2524  {
2525  return;
2526  }
2527 
2528  m_filters.clear();
2529 
2530  TiXmlElement *filter = filters->FirstChildElement("FILTER");
2531  while (filter)
2532  {
2533  ASSERTL0(filter->Attribute("TYPE"),
2534  "Missing attribute 'TYPE' for filter.");
2535  std::string typeStr = filter->Attribute("TYPE");
2536 
2537  std::map<std::string, std::string> vParams;
2538 
2539  TiXmlElement *param = filter->FirstChildElement("PARAM");
2540  while (param)
2541  {
2542  ASSERTL0(param->Attribute("NAME"),
2543  "Missing attribute 'NAME' for parameter in filter "
2544  + typeStr + "'.");
2545  std::string nameStr = param->Attribute("NAME");
2546 
2547  ASSERTL0(param->GetText(), "Empty value string for param.");
2548  std::string valueStr = param->GetText();
2549 
2550  vParams[nameStr] = valueStr;
2551 
2552  param = param->NextSiblingElement("PARAM");
2553  }
2554 
2555  m_filters.push_back(
2556  std::pair<std::string, FilterParams>(typeStr, vParams));
2557 
2558  filter = filter->NextSiblingElement("FILTER");
2559  }
2560  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
FilterMap m_filters
Filters map.
void Nektar::LibUtilities::SessionReader::ReadFunctions ( TiXmlElement *  conditions)
private

Reads the FUNCTIONS section of the XML document.

Definition at line 2318 of file SessionReader.cpp.

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

2319  {
2320  m_functions.clear();
2321 
2322  if (!conditions)
2323  {
2324  return;
2325  }
2326 
2327  // Scan through conditions section looking for functions.
2328  TiXmlElement *function = conditions->FirstChildElement("FUNCTION");
2329  while (function)
2330  {
2331  stringstream tagcontent;
2332  tagcontent << *function;
2333 
2334  // Every function must have a NAME attribute
2335  ASSERTL0(function->Attribute("NAME"),
2336  "Functions must have a NAME attribute defined in XML "
2337  "element: \n\t'" + tagcontent.str() + "'");
2338  std::string functionStr = function->Attribute("NAME");
2339  ASSERTL0(!functionStr.empty(),
2340  "Functions must have a non-empty name in XML "
2341  "element: \n\t'" + tagcontent.str() + "'");
2342 
2343  // Store function names in uppercase to remain case-insensitive.
2344  boost::to_upper(functionStr);
2345 
2346  // Retrieve first entry (variable, or file)
2347  TiXmlElement *variable = function->FirstChildElement();
2348 
2349  // Create new function structure with default type of none.
2350  FunctionVariableMap functionVarMap;
2351 
2352  // Process all entries in the function block
2353  while (variable)
2354  {
2355  FunctionVariableDefinition funcDef;
2356  std::string conditionType = variable->Value();
2357 
2358  // If no var is specified, assume wildcard
2359  std::string variableStr;
2360  if (!variable->Attribute("VAR"))
2361  {
2362  variableStr = "*";
2363  }
2364  else
2365  {
2366  variableStr = variable->Attribute("VAR");
2367  }
2368 
2369  // Parse list of variables
2370  std::vector<std::string> variableList;
2371  ParseUtils::GenerateOrderedStringVector(variableStr.c_str(),
2372  variableList);
2373 
2374  // If no domain string put to 0
2375  std::string domainStr;
2376  if (!variable->Attribute("DOMAIN"))
2377  {
2378  domainStr = "0";
2379  }
2380  else
2381  {
2382  domainStr = variable->Attribute("DOMAIN");
2383  }
2384 
2385  // Parse list of variables
2386  std::vector<std::string> varSplit;
2387  std::vector<unsigned int> domainList;
2388  ParseUtils::GenerateSeqVector(domainStr.c_str(), domainList);
2389 
2390  // Expressions are denoted by E
2391  if (conditionType == "E")
2392  {
2393  funcDef.m_type = eFunctionTypeExpression;
2394 
2395  // Expression must have a VALUE.
2396  ASSERTL0(variable->Attribute("VALUE"),
2397  "Attribute VALUE expected for function '"
2398  + functionStr + "'.");
2399  std::string fcnStr = variable->Attribute("VALUE");
2400 
2401  ASSERTL0(!fcnStr.empty(),
2402  (std::string("Expression for var: ")
2403  + variableStr
2404  + std::string(" must be specified.")).c_str());
2405 
2406  SubstituteExpressions(fcnStr);
2407 
2408  // set expression
2409  funcDef.m_expression = MemoryManager<Equation>
2411  }
2412 
2413  // Files are denoted by F
2414  else if (conditionType == "F")
2415  {
2416  if (variable->Attribute("TIMEDEPENDENT") &&
2417  boost::lexical_cast<bool>(variable->Attribute("TIMEDEPENDENT")))
2418  {
2419  funcDef.m_type = eFunctionTypeTransientFile;
2420  }
2421  else
2422  {
2423  funcDef.m_type = eFunctionTypeFile;
2424  }
2425 
2426  // File must have a FILE.
2427  ASSERTL0(variable->Attribute("FILE"),
2428  "Attribute FILE expected for function '"
2429  + functionStr + "'.");
2430  std::string filenameStr = variable->Attribute("FILE");
2431 
2432  ASSERTL0(!filenameStr.empty(),
2433  "A filename must be specified for the FILE "
2434  "attribute of function '" + functionStr
2435  + "'.");
2436 
2437  std::vector<std::string> fSplit;
2438  boost::split(fSplit, filenameStr, boost::is_any_of(":"));
2439 
2440  ASSERTL0(fSplit.size() == 1 || fSplit.size() == 2,
2441  "Incorrect filename specification in function "
2442  + functionStr + "'. "
2443  "Specify variables inside file as: "
2444  "filename:var1,var2");
2445 
2446  // set the filename
2447  funcDef.m_filename = fSplit[0];
2448 
2449  if (fSplit.size() == 2)
2450  {
2451  ASSERTL0(variableList[0] != "*",
2452  "Filename variable mapping not valid "
2453  "when using * as a variable inside "
2454  "function '" + functionStr + "'.");
2455 
2456  boost::split(
2457  varSplit, fSplit[1], boost::is_any_of(","));
2458  ASSERTL0(varSplit.size() == variableList.size(),
2459  "Filename variables should contain the "
2460  "same number of variables defined in "
2461  "VAR in function " + functionStr + "'.");
2462  }
2463  }
2464 
2465  // Nothing else supported so throw an error
2466  else
2467  {
2468  stringstream tagcontent;
2469  tagcontent << *variable;
2470 
2471  ASSERTL0(false,
2472  "Identifier " + conditionType + " in function "
2473  + std::string(function->Attribute("NAME"))
2474  + " is not recognised in XML element: \n\t'"
2475  + tagcontent.str() + "'");
2476  }
2477 
2478 
2479 
2480  // Add variables to function
2481  for (unsigned int i = 0; i < variableList.size(); ++i)
2482  {
2483  for(unsigned int j = 0; j < domainList.size(); ++j)
2484  {
2485  // Check it has not already been defined
2486  pair<std::string,int> key(variableList[i],domainList[j]);
2488  = functionVarMap.find(key);
2489  ASSERTL0(fcnsIter == functionVarMap.end(),
2490  "Error setting expression '" + variableList[i]
2491  + " in domain "
2492  + boost::lexical_cast<std::string>(domainList[j])
2493  + "' in function '" + functionStr + "'. "
2494  "Expression has already been defined.");
2495 
2496  if (varSplit.size() > 0)
2497  {
2498  FunctionVariableDefinition funcDef2 = funcDef;
2499  funcDef2.m_fileVariable = varSplit[i];
2500  functionVarMap[key] = funcDef2;
2501  }
2502  else
2503  {
2504  functionVarMap[key] = funcDef;
2505  }
2506  }
2507  }
2508 
2509  variable = variable->NextSiblingElement();
2510  }
2511  // Add function definition to map
2512  m_functions[functionStr] = functionVarMap;
2513  function = function->NextSiblingElement("FUNCTION");
2514  }
2515  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
static bool GenerateOrderedStringVector(const char *const str, std::vector< std::string > &vec)
Definition: ParseUtils.hpp:143
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
std::map< std::pair< std::string, int >, FunctionVariableDefinition > FunctionVariableMap
void SubstituteExpressions(std::string &expr)
Substitutes expressions defined in the XML document.
SessionReaderSharedPtr GetSharedThisPtr()
Returns a shared pointer to the current object.
static bool GenerateSeqVector(const char *const str, std::vector< unsigned int > &vec)
Definition: ParseUtils.hpp:79
FunctionMap m_functions
Functions.
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
void Nektar::LibUtilities::SessionReader::ReadGlobalSysSolnInfo ( TiXmlElement *  conditions)
private

Reads the GLOBALSYSSOLNINFO section of the XML document.

Definition at line 2049 of file SessionReader.cpp.

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

2050  {
2051  GetGloSysSolnList().clear();
2052 
2053  if (!conditions)
2054  {
2055  return;
2056  }
2057 
2058  TiXmlElement *GlobalSys =
2059  conditions->FirstChildElement("GLOBALSYSSOLNINFO");
2060 
2061  if(!GlobalSys)
2062  {
2063  return;
2064  }
2065 
2066  TiXmlElement *VarInfo = GlobalSys->FirstChildElement("V");
2067 
2068  while (VarInfo)
2069  {
2070  std::stringstream tagcontent;
2071  tagcontent << *VarInfo;
2072  ASSERTL0(VarInfo->Attribute("VAR"),
2073  "Missing VAR attribute in GobalSysSolnInfo XML "
2074  "element: \n\t'" + tagcontent.str() + "'");
2075 
2076  std::string VarList = VarInfo->Attribute("VAR");
2077  ASSERTL0(!VarList.empty(),
2078  "VAR attribute must be non-empty in XML element:\n\t'"
2079  + tagcontent.str() + "'");
2080 
2081  // generate a list of variables.
2082  std::vector<std::string> varStrings;
2084  VarList.c_str(),varStrings);
2085 
2086  ASSERTL0(valid,"Unable to process list of variable in XML "
2087  "element \n\t'" + tagcontent.str() + "'");
2088 
2089  if(varStrings.size())
2090  {
2091  TiXmlElement *SysSolnInfo = VarInfo->FirstChildElement("I");
2092 
2093  while (SysSolnInfo)
2094  {
2095  tagcontent.clear();
2096  tagcontent << *SysSolnInfo;
2097  // read the property name
2098  ASSERTL0(SysSolnInfo->Attribute("PROPERTY"),
2099  "Missing PROPERTY attribute in "
2100  "GlobalSysSolnInfo for variable(s) '"
2101  + VarList + "' in XML element: \n\t'"
2102  + tagcontent.str() + "'");
2103 
2104  std::string SysSolnProperty =
2105  SysSolnInfo->Attribute("PROPERTY");
2106 
2107  ASSERTL0(!SysSolnProperty.empty(),
2108  "GlobalSysSolnIno properties must have a "
2109  "non-empty name for variable(s) : '"
2110  + VarList + "' in XML element: \n\t'"
2111  + tagcontent.str() + "'");
2112 
2113  // make sure that solver property is capitalised
2114  std::string SysSolnPropertyUpper =
2115  boost::to_upper_copy(SysSolnProperty);
2116 
2117  // read the value
2118  ASSERTL0(SysSolnInfo->Attribute("VALUE"),
2119  "Missing VALUE attribute in GlobalSysSolnInfo "
2120  "for variable(s) '" + VarList
2121  + "' in XML element: \n\t"
2122  + tagcontent.str() + "'");
2123 
2124  std::string SysSolnValue =
2125  SysSolnInfo->Attribute("VALUE");
2126  ASSERTL0(!SysSolnValue.empty(),
2127  "GlobalSysSolnInfo properties must have a "
2128  "non-empty value for variable(s) '"
2129  + VarList + "' in XML element: \n\t'"
2130  + tagcontent.str() + "'");
2131 
2132  // Store values under variable map.
2133  for(int i = 0; i < varStrings.size(); ++i)
2134  {
2136  if ((x = GetGloSysSolnList().find(varStrings[i])) ==
2137  GetGloSysSolnList().end())
2138  {
2139  (GetGloSysSolnList()[varStrings[i]])[
2140  SysSolnPropertyUpper] = SysSolnValue;
2141  }
2142  else
2143  {
2144  x->second[SysSolnPropertyUpper] = SysSolnValue;
2145  }
2146  }
2147 
2148  SysSolnInfo = SysSolnInfo->NextSiblingElement("I");
2149  }
2150  VarInfo = VarInfo->NextSiblingElement("V");
2151  }
2152  }
2153 
2154  if (m_verbose && GetGloSysSolnList().size() > 0 && m_comm)
2155  {
2156  if(m_comm->GetRank() == 0)
2157  {
2158  cout << "GlobalSysSoln Info:" << endl;
2159 
2161  for (x = GetGloSysSolnList().begin();
2162  x != GetGloSysSolnList().end();
2163  ++x)
2164  {
2165  cout << "\t Variable: " << x->first << endl;
2166 
2168  for (y = x->second.begin(); y != x->second.end(); ++y)
2169  {
2170  cout << "\t\t " << y->first << " = " << y->second
2171  << endl;
2172  }
2173  }
2174  cout << endl;
2175  }
2176  }
2177  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
static bool GenerateOrderedStringVector(const char *const str, std::vector< std::string > &vec)
Definition: ParseUtils.hpp:143
static GloSysSolnInfoList & GetGloSysSolnList()
GlobalSysSoln Info map.
CommSharedPtr m_comm
Communication object.
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
InputIterator find(InputIterator first, InputIterator last, InputIterator startingpoint, const EqualityComparable &value)
Definition: StdRegions.hpp:315
void Nektar::LibUtilities::SessionReader::ReadParameters ( TiXmlElement *  conditions)
private

Reads the PARAMETERS section of the XML document.

Definition at line 1889 of file SessionReader.cpp.

References ASSERTL0, Nektar::LibUtilities::Equation::Evaluate(), and Nektar::lhs.

1890  {
1891  m_parameters.clear();
1892 
1893  if (!conditions)
1894  {
1895  return;
1896  }
1897 
1898  TiXmlElement *parametersElement = conditions->FirstChildElement(
1899  "PARAMETERS");
1900 
1901  // See if we have parameters defined. They are optional so we go on
1902  // if not.
1903  if (parametersElement)
1904  {
1905  TiXmlElement *parameter =
1906  parametersElement->FirstChildElement("P");
1907 
1908  ParameterMap caseSensitiveParameters;
1909 
1910  // Multiple nodes will only occur if there is a comment in
1911  // between definitions.
1912  while (parameter)
1913  {
1914  stringstream tagcontent;
1915  tagcontent << *parameter;
1916  TiXmlNode *node = parameter->FirstChild();
1917 
1918  while (node && node->Type() != TiXmlNode::TINYXML_TEXT)
1919  {
1920  node = node->NextSibling();
1921  }
1922 
1923  if (node)
1924  {
1925  // Format is "paramName = value"
1926  std::string line = node->ToText()->Value(), lhs, rhs;
1927 
1928  try {
1929  ParseEquals(line, lhs, rhs);
1930  }
1931  catch (...)
1932  {
1933  ASSERTL0(false, "Syntax error in parameter "
1934  "expression '" + line
1935  + "' in XML element: \n\t'"
1936  + tagcontent.str() + "'");
1937  }
1938 
1939  // We want the list of parameters to have their RHS
1940  // evaluated, so we use the expression evaluator to do
1941  // the dirty work.
1942  if (!lhs.empty() && !rhs.empty())
1943  {
1944  NekDouble value=0.0;
1945  try
1946  {
1947  LibUtilities::Equation expession(
1948  GetSharedThisPtr(), rhs);
1949  value = expession.Evaluate();
1950  }
1951  catch (const std::runtime_error &)
1952  {
1953  ASSERTL0(false,
1954  "Error evaluating parameter expression"
1955  " '" + rhs + "' in XML element: \n\t'"
1956  + tagcontent.str() + "'");
1957  }
1959  caseSensitiveParameters[lhs] = value;
1960  boost::to_upper(lhs);
1961  m_parameters[lhs] = value;
1962  }
1963  }
1964 
1965  parameter = parameter->NextSiblingElement();
1966  }
1967  }
1968  }
ParameterMap m_parameters
Parameters.
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
void SetParameter(std::string const &name, NekDouble value)
This function behaves in the same way as SetParameters, but it only adds one parameter and it does no...
AnalyticExpressionEvaluator m_exprEvaluator
Analytic expression evaluator instance.
SessionReaderSharedPtr GetSharedThisPtr()
Returns a shared pointer to the current object.
StandardMatrixTag & lhs
void ParseEquals(const std::string &line, std::string &lhs, std::string &rhs)
Parse a string in the form lhs = rhs.
double NekDouble
std::map< std::string, NekDouble > ParameterMap
Definition: SessionReader.h:59
void Nektar::LibUtilities::SessionReader::ReadSolverInfo ( TiXmlElement *  conditions)
private

Reads the SOLVERINFO section of the XML document.

Definition at line 1974 of file SessionReader.cpp.

References ASSERTL0.

1975  {
1976  m_solverInfo.clear();
1978 
1979  if (!conditions)
1980  {
1981  return;
1982  }
1983 
1984  TiXmlElement *solverInfoElement =
1985  conditions->FirstChildElement("SOLVERINFO");
1986 
1987  if (solverInfoElement)
1988  {
1989  TiXmlElement *solverInfo =
1990  solverInfoElement->FirstChildElement("I");
1991 
1992  while (solverInfo)
1993  {
1994  std::stringstream tagcontent;
1995  tagcontent << *solverInfo;
1996  // read the property name
1997  ASSERTL0(solverInfo->Attribute("PROPERTY"),
1998  "Missing PROPERTY attribute in solver info "
1999  "XML element: \n\t'" + tagcontent.str() + "'");
2000  std::string solverProperty =
2001  solverInfo->Attribute("PROPERTY");
2002  ASSERTL0(!solverProperty.empty(),
2003  "PROPERTY attribute must be non-empty in XML "
2004  "element: \n\t'" + tagcontent.str() + "'");
2005 
2006  // make sure that solver property is capitalised
2007  std::string solverPropertyUpper =
2008  boost::to_upper_copy(solverProperty);
2009 
2010  // read the value
2011  ASSERTL0(solverInfo->Attribute("VALUE"),
2012  "Missing VALUE attribute in solver info "
2013  "XML element: \n\t'" + tagcontent.str() + "'");
2014  std::string solverValue = solverInfo->Attribute("VALUE");
2015  ASSERTL0(!solverValue.empty(),
2016  "VALUE attribute must be non-empty in XML "
2017  "element: \n\t'" + tagcontent.str() + "'");
2018 
2019  // Set Variable
2020  m_solverInfo[solverPropertyUpper] = solverValue;
2021  solverInfo = solverInfo->NextSiblingElement("I");
2022  }
2023  }
2024 
2025  if (m_comm && m_comm->GetRowComm()->GetSize() > 1)
2026  {
2027  ASSERTL0(
2028  m_solverInfo["GLOBALSYSSOLN"] == "IterativeFull" ||
2029  m_solverInfo["GLOBALSYSSOLN"] == "IterativeStaticCond" ||
2030  m_solverInfo["GLOBALSYSSOLN"] ==
2031  "IterativeMultiLevelStaticCond" ||
2032  m_solverInfo["GLOBALSYSSOLN"] == "XxtFull" ||
2033  m_solverInfo["GLOBALSYSSOLN"] == "XxtStaticCond" ||
2034  m_solverInfo["GLOBALSYSSOLN"] ==
2035  "XxtMultiLevelStaticCond" ||
2036  m_solverInfo["GLOBALSYSSOLN"] == "PETScFull" ||
2037  m_solverInfo["GLOBALSYSSOLN"] == "PETScStaticCond" ||
2038  m_solverInfo["GLOBALSYSSOLN"] ==
2039  "PETScMultiLevelStaticCond",
2040  "A parallel solver must be used when run in parallel.");
2041  }
2042  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
SolverInfoMap m_solverInfo
Solver information properties.
CommSharedPtr m_comm
Communication object.
static SolverInfoMap & GetSolverInfoDefaults()
Default solver info options.
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 2235 of file SessionReader.cpp.

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

2236  {
2237  m_variables.clear();
2238 
2239  if (!conditions)
2240  {
2241  return;
2242  }
2243 
2244  TiXmlElement *variablesElement =
2245  conditions->FirstChildElement("VARIABLES");
2246 
2247  // See if we have parameters defined. They are optional so we go on
2248  // if not.
2249  if (variablesElement)
2250  {
2251  TiXmlElement *varElement =
2252  variablesElement->FirstChildElement("V");
2253 
2254  // Sequential counter for the composite numbers.
2255  int nextVariableNumber = -1;
2256 
2257  while (varElement)
2258  {
2259  stringstream tagcontent;
2260  tagcontent << *varElement;
2261 
2262  /// All elements are of the form: "<V ID="#"> name = value
2263  /// </V>", with ? being the element type.
2264  nextVariableNumber++;
2265 
2266  int i;
2267  int err = varElement->QueryIntAttribute("ID", &i);
2268  ASSERTL0(err == TIXML_SUCCESS,
2269  "Variables must have a unique ID number attribute "
2270  "in XML element: \n\t'" + tagcontent.str() + "'");
2271  ASSERTL0(i == nextVariableNumber,
2272  "ID numbers for variables must begin with zero and"
2273  " be sequential in XML element: \n\t'"
2274  + tagcontent.str() + "'");
2275 
2276  TiXmlNode* varChild = varElement->FirstChild();
2277  // This is primarily to skip comments that may be present.
2278  // Comments appear as nodes just like elements. We are
2279  // specifically looking for text in the body of the
2280  // definition.
2281  while(varChild && varChild->Type() != TiXmlNode::TINYXML_TEXT)
2282  {
2283  varChild = varChild->NextSibling();
2284  }
2285 
2286  ASSERTL0(varChild,
2287  "Unable to read variable definition body for "
2288  "variable with ID "
2289  + boost::lexical_cast<string>(i)
2290  + " in XML element: \n\t'"
2291  + tagcontent.str() + "'");
2292  std::string variableName = varChild->ToText()->ValueStr();
2293 
2294  std::istringstream variableStrm(variableName);
2295  variableStrm >> variableName;
2296 
2297  ASSERTL0(std::find(m_variables.begin(), m_variables.end(),
2298  variableName) == m_variables.end(),
2299  "Variable with ID "
2300  + boost::lexical_cast<string>(i)
2301  + " in XML element \n\t'" + tagcontent.str()
2302  + "'\nhas already been defined.");
2303 
2304  m_variables.push_back(variableName);
2305 
2306  varElement = varElement->NextSiblingElement("V");
2307  }
2308 
2309  ASSERTL0(nextVariableNumber > -1,
2310  "Number of variables must be greater than zero.");
2311  }
2312  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
VariableList m_variables
Variables.
InputIterator find(InputIterator first, InputIterator last, InputIterator startingpoint, const EqualityComparable &value)
Definition: StdRegions.hpp:315
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 671 of file SessionReader.h.

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

675  {
676  ASSERTL0(!pName.empty(), "Empty name for cmdline argument.");
677  CmdLineArg x;
678  x.shortName = pShortName;
679  x.description = pDescription;
680  x.isFlag = false;
681  GetCmdLineArgMap()[pName] = x;
682  return pName;
683  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
static CmdLineArgMap & GetCmdLineArgMap()
CmdLine argument map.
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 689 of file SessionReader.h.

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

693  {
694  ASSERTL0(!pName.empty(), "Empty name for cmdline argument.");
695  CmdLineArg x;
696  x.shortName = pShortName;
697  x.description = pDescription;
698  x.isFlag = true;
699  GetCmdLineArgMap()[pName] = x;
700  return pName;
701  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
static CmdLineArgMap & GetCmdLineArgMap()
CmdLine argument map.
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");
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 658 of file SessionReader.h.

References GetSolverInfoDefaults().

661  {
662  std::string vName = boost::to_upper_copy(pName);
663  GetSolverInfoDefaults()[vName] = pValue;
664  return pValue;
665  }
static SolverInfoMap & GetSolverInfoDefaults()
Default solver info options.
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",
}
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 625 of file SessionReader.h.

References Nektar::StdRegions::find(), GetSolverInfoEnums(), and Nektar::iterator.

627  {
628  std::string vEnum = boost::to_upper_copy(pEnum);
630  if ((x = GetSolverInfoEnums().find(vEnum)) ==
631  GetSolverInfoEnums().end())
632  {
633  GetSolverInfoEnums()[vEnum] = EnumMap();
634  x = GetSolverInfoEnums().find(vEnum);
635  }
636  x->second[pString] = pEnumValue;
637  return pString;
638  }
std::map< std::string, int > EnumMap
Definition: SessionReader.h:77
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
static EnumMapList & GetSolverInfoEnums()
String to enumeration map for Solver Info parameters.
InputIterator find(InputIterator first, InputIterator last, InputIterator startingpoint, const EqualityComparable &value)
Definition: StdRegions.hpp:315
void Nektar::LibUtilities::SessionReader::SetParameter ( const std::string &  name,
int &  var 
)

Set an integer parameter.

Definition at line 739 of file SessionReader.cpp.

740  {
741  std::string vName = boost::to_upper_copy(pName);
742  m_parameters[vName] = pVar;
743  }
ParameterMap m_parameters
Parameters.
void Nektar::LibUtilities::SessionReader::SetParameter ( const std::string &  name,
NekDouble var 
)

Set a double precision parameter.

Definition at line 749 of file SessionReader.cpp.

751  {
752  std::string vName = boost::to_upper_copy(pName);
753  m_parameters[vName] = pVar;
754  }
ParameterMap m_parameters
Parameters.
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 787 of file SessionReader.cpp.

References ASSERTL1, and Nektar::iterator.

789  {
790  std::string vProperty = boost::to_upper_copy(pProperty);
791  SolverInfoMap::iterator iter = m_solverInfo.find(vProperty);
792 
793  ASSERTL1(iter != m_solverInfo.end(),
794  "Unable to find requested property: " + pProperty);
795 
796  iter->second = pValue;
797  }
SolverInfoMap m_solverInfo
Solver information properties.
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:191
void Nektar::LibUtilities::SessionReader::SetTag ( const std::string &  pName,
const std::string &  pValue 
)

Sets a specified tag.

Definition at line 1316 of file SessionReader.cpp.

1319  {
1320  std::string vName = boost::to_upper_copy(pName);
1321  m_tags[vName] = pValue;
1322  }
void Nektar::LibUtilities::SessionReader::SetUpXmlDoc ( void  )

Definition at line 2674 of file SessionReader.cpp.

2675  {
2677  }
TiXmlDocument * MergeDoc(const std::vector< std::string > &pFilenames) const
Creates an XML document from a list of input files.
std::vector< std::string > m_filenames
Filenames.
TiXmlDocument * m_xmlDoc
Pointer to the loaded XML document.
void Nektar::LibUtilities::SessionReader::SetVariable ( const unsigned int &  idx,
std::string  newname 
)

Definition at line 1026 of file SessionReader.cpp.

References ASSERTL0.

1028  {
1029  ASSERTL0(idx < m_variables.size(), "Variable index out of range.");
1030  m_variables[idx] = newname;
1031  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
VariableList m_variables
Variables.
void Nektar::LibUtilities::SessionReader::SubstituteExpressions ( std::string &  expr)

Substitutes expressions defined in the XML document.

Definition at line 1360 of file SessionReader.cpp.

References Nektar::iterator.

1361  {
1362  ExpressionMap::iterator exprIter;
1363  for (exprIter = m_expressions.begin();
1364  exprIter != m_expressions.end(); ++exprIter)
1365  {
1366  boost::replace_all(pExpr, exprIter->first, exprIter->second);
1367  }
1368  }
ExpressionMap m_expressions
Expressions.
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
void Nektar::LibUtilities::SessionReader::VerifySolverInfo ( )
private

Check values of solver info options are valid.

Definition at line 2652 of file SessionReader.cpp.

References ASSERTL0.

2653  {
2654  SolverInfoMap::const_iterator x;
2655  for (x = m_solverInfo.begin(); x != m_solverInfo.end(); ++x)
2656  {
2657  std::string solverProperty = x->first;
2658  std::string solverValue = x->second;
2659 
2660  EnumMapList::const_iterator propIt =
2661  GetSolverInfoEnums().find(solverProperty);
2662  if (propIt != GetSolverInfoEnums().end())
2663  {
2664  EnumMap::const_iterator valIt =
2665  propIt->second.find(solverValue);
2666  ASSERTL0(valIt != propIt->second.end(),
2667  "Value '" + solverValue + "' is not valid for "
2668  "property '" + solverProperty + "'");
2669  }
2670  }
2671  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
SolverInfoMap m_solverInfo
Solver information properties.
static EnumMapList & GetSolverInfoEnums()
String to enumeration map for Solver Info parameters.

Friends And Related Function Documentation

friend class MemoryManager< SessionReader >
friend

Support creation through MemoryManager.

Definition at line 128 of file SessionReader.h.

Member Data Documentation

BndRegionOrdering Nektar::LibUtilities::SessionReader::m_bndRegOrder
private

Map of original boundary region ordering for parallel periodic bcs.

Definition at line 464 of file SessionReader.h.

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

Definition at line 430 of file SessionReader.h.

Referenced by GetCmdLineArgument().

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

Communication object.

Definition at line 433 of file SessionReader.h.

CompositeOrdering Nektar::LibUtilities::SessionReader::m_compOrder
private

Map of original composite ordering for parallel periodic bcs.

Definition at line 461 of file SessionReader.h.

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

Expressions.

Definition at line 447 of file SessionReader.h.

AnalyticExpressionEvaluator Nektar::LibUtilities::SessionReader::m_exprEvaluator
private

Analytic expression evaluator instance.

Definition at line 449 of file SessionReader.h.

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

Filenames.

Definition at line 435 of file SessionReader.h.

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

Filters map.

Definition at line 457 of file SessionReader.h.

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

Functions.

Definition at line 451 of file SessionReader.h.

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

Geometric information properties.

Definition at line 445 of file SessionReader.h.

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

Parameters.

Definition at line 441 of file SessionReader.h.

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

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

Definition at line 437 of file SessionReader.h.

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

Solver information properties.

Definition at line 443 of file SessionReader.h.

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

Custom tags.

Definition at line 455 of file SessionReader.h.

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

Variables.

Definition at line 453 of file SessionReader.h.

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

Be verbose.

Definition at line 459 of file SessionReader.h.

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

Pointer to the loaded XML document.

Definition at line 439 of file SessionReader.h.