diff --git a/CHANGELOG.md b/CHANGELOG.md
index 4caa742938733ba61e13228e792cba44dd6bb453..0cf50ff856f0ed30ec4051b0e599fcd1675bc657 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -10,6 +10,8 @@ v5.0.0
- Fix ThridpartyCCM options (!802)
- Fix Windows CRLF tokens in GEO reader and improve comment handling (!805)
- Use chrono in Timer (!807)
+- Fix caching of FUNCTION tags that read from file and provide the same
+ functionality in FUNCTIONs defined for forcings (!759)
**NekMesh**:
- Add feature to read basic 2D geo files as CAD (!731)
diff --git a/library/FieldUtils/Interpolator.h b/library/FieldUtils/Interpolator.h
index 3f9a8348c074aa999a7e5e3aa96f7f53798970dc..7176e2afd7ec8928d817ec3a64906cf8963b8ffd 100644
--- a/library/FieldUtils/Interpolator.h
+++ b/library/FieldUtils/Interpolator.h
@@ -147,7 +147,7 @@ public:
/// Returns the output field
FIELD_UTILS_EXPORT LibUtilities::PtsFieldSharedPtr GetOutField() const;
- /// Print statics of the interpolation weights
+ /// Returns if the weights have already been computed
FIELD_UTILS_EXPORT void PrintStatistics();
/// sets a callback funtion which gets called every time the interpolation
diff --git a/library/SolverUtils/CMakeLists.txt b/library/SolverUtils/CMakeLists.txt
index 677a988393e3cd53e62918445f8b871d6fe86c12..55b4a8de3ba37449def02d597e94a5c04d887a20 100644
--- a/library/SolverUtils/CMakeLists.txt
+++ b/library/SolverUtils/CMakeLists.txt
@@ -1,5 +1,6 @@
SET(SOLVER_UTILS_SOURCES
Core/Misc.cpp
+ Core/SessionFunction.cpp
AdvectionSystem.cpp
Advection/Advection.cpp
Advection/Advection3DHomogeneous1D.cpp
@@ -44,6 +45,7 @@ SET(SOLVER_UTILS_SOURCES
SET(SOLVER_UTILS_HEADERS
Core/Misc.h
+ Core/SessionFunction.h
AdvectionSystem.h
Advection/Advection.h
Advection/AdvectionFR.h
diff --git a/library/SolverUtils/Core/SessionFunction.cpp b/library/SolverUtils/Core/SessionFunction.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..5c4140eb9f2eb1b77bf327c8a82ec0f46e45c9f1
--- /dev/null
+++ b/library/SolverUtils/Core/SessionFunction.cpp
@@ -0,0 +1,500 @@
+////////////////////////////////////////////////////////////////////////////////
+//
+// File: SessionFunction.cpp
+//
+// For more information, please see: http://www.nektar.info/
+//
+// The MIT License
+//
+// Copyright (c) 2017 Kilian Lackhove
+// Copyright (c) 2006 Division of Applied Mathematics, Brown University (USA),
+// Department of Aeronautics, Imperial College London (UK), and Scientific
+// Computing and Imaging Institute, University of Utah (USA).
+//
+// License for the specific language governing rights and limitations under
+// Permission is hereby granted, free of charge, to any person obtaining a
+// copy of this software and associated documentation files (the "Software"),
+// to deal in the Software without restriction, including without limitation
+// the rights to use, copy, modify, merge, publish, distribute, sublicense,
+// and/or sell copies of the Software, and to permit persons to whom the
+// Software is furnished to do so, subject to the following conditions:
+//
+// The above copyright notice and this permission notice shall be included
+// in all copies or substantial portions of the Software.
+//
+// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
+// OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+// FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
+// THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+// LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
+// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
+// DEALINGS IN THE SOFTWARE.
+//
+// Description: Session Function
+//
+////////////////////////////////////////////////////////////////////////////////
+
+#include <SolverUtils/Core/SessionFunction.h>
+
+#include <LibUtilities/BasicUtils/VmathArray.hpp>
+
+#include <boost/format.hpp>
+#include <boost/function.hpp>
+
+using namespace std;
+
+namespace Nektar
+{
+namespace SolverUtils
+{
+
+/**
+ * Representation of a FUNCTION defined in the session xml file.
+ *
+ * @param session The session where the function was defined.
+ * @param field The field the function is defined on.
+ * @param functionName The name of the function.
+ * @param toCache Store the evaluated function for later use.
+ */
+SessionFunction::SessionFunction(const LibUtilities::SessionReaderSharedPtr &session,
+ const MultiRegions::ExpListSharedPtr &field,
+ std::string functionName,
+ bool toCache)
+ : m_session(session), m_field(field), m_name(functionName), m_toCache(toCache)
+{
+ ASSERTL0(m_session->DefinesFunction(m_name),
+ "Function '" + m_name + "' does not exist.");
+}
+
+/**
+ * Evaluates a function defined in the xml session file at each quadrature point.
+ *
+ * @param pArray The array into which to write the values.
+ * @param pTime The time at which to evaluate the function.
+ * @param domain The domain to evaluate the function in.
+ */
+void SessionFunction::Evaluate(Array<OneD, Array<OneD, NekDouble> > &pArray,
+ const NekDouble pTime,
+ const int domain)
+{
+ std::vector<std::string> vFieldNames = m_session->GetVariables();
+
+ for (int i = 0; i < vFieldNames.size(); i++)
+ {
+ Evaluate(vFieldNames[i], pArray[i], pTime, domain);
+ }
+}
+
+/**
+ * Evaluates a function defined in the xml session file at each quadrature point.
+ *
+ * @param pFieldNames The names of the fields to evaluate the function for.
+ * @param pArray The array into which to write the values.
+ * @param pTime The time at which to evaluate the function.
+ * @param domain The domain to evaluate the function in.
+ */
+void SessionFunction::Evaluate(std::vector<std::string> pFieldNames,
+ Array<OneD, Array<OneD, NekDouble> > &pArray,
+ const NekDouble &pTime,
+ const int domain)
+{
+ ASSERTL1(pFieldNames.size() == pArray.num_elements(),
+ "Function '" + m_name +
+ "' variable list size mismatch with array storage.");
+
+ for (int i = 0; i < pFieldNames.size(); i++)
+ {
+ Evaluate(pFieldNames[i], pArray[i], pTime, domain);
+ }
+}
+
+/**
+ * Evaluates a function defined in the xml session file at each quadrature point.
+ *
+ * @param pFieldNames The names of the fields to evaluate the function for.
+ * @param pFields The fields into which to write the values.
+ * @param pTime The time at which to evaluate the function.
+ * @param domain The domain to evaluate the function in.
+ */
+void SessionFunction::Evaluate(
+ std::vector<std::string> pFieldNames,
+ Array<OneD, MultiRegions::ExpListSharedPtr> &pFields,
+ const NekDouble &pTime,
+ const int domain)
+{
+ ASSERTL0(pFieldNames.size() == pFields.num_elements(),
+ "Field list / name list size mismatch.");
+
+ for (int i = 0; i < pFieldNames.size(); i++)
+ {
+ Evaluate(pFieldNames[i], pFields[i]->UpdatePhys(), pTime, domain);
+ pFields[i]->FwdTrans_IterPerExp(pFields[i]->GetPhys(),
+ pFields[i]->UpdateCoeffs());
+ }
+}
+
+/**
+ * Evaluates a function defined in the xml session file at each quadrature point.
+ *
+ * @param pFieldName The name of the field to evaluate the function for.
+ * @param pArray The array into which to write the values.
+ * @param pTime The time at which to evaluate the function.
+ * @param domain The domain to evaluate the function in.
+ */
+void SessionFunction::Evaluate(std::string pFieldName,
+ Array<OneD, NekDouble> &pArray,
+ const NekDouble &pTime,
+ const int domain)
+{
+ LibUtilities::FunctionType vType =
+ m_session->GetFunctionType(m_name, pFieldName, domain);
+
+ unsigned int nq = m_field->GetNpoints();
+
+ std::pair<std::string, int> key(pFieldName, domain);
+ // sorry
+ if ((m_arrays.find(key) != m_arrays.end()) &&
+ (vType == LibUtilities::eFunctionTypeFile ||
+ ((m_lastCached.find(key) != m_lastCached.end()) &&
+ (pTime - m_lastCached[key] < NekConstants::kNekZeroTol))))
+ {
+ // found cached field
+ if (pArray.num_elements() < nq)
+ {
+ pArray = Array<OneD, NekDouble>(nq);
+ }
+ Vmath::Vcopy(nq, m_arrays[key], 1, pArray, 1);
+
+ return;
+ }
+
+ if (vType == LibUtilities::eFunctionTypeExpression)
+ {
+ EvaluateExp(pFieldName, pArray, pTime, domain);
+ }
+ else if (vType == LibUtilities::eFunctionTypeFile ||
+ vType == LibUtilities::eFunctionTypeTransientFile)
+ {
+ std::string filename =
+ m_session->GetFunctionFilename(m_name, pFieldName, domain);
+
+ if (boost::filesystem::path(filename).extension() == ".pts")
+ {
+ EvaluatePts(pFieldName, pArray, pTime, domain);
+ }
+ else
+ {
+ EvaluateFld(pFieldName, pArray, pTime, domain);
+ }
+ }
+
+ if (m_toCache)
+ {
+ m_arrays[key] = Array<OneD, NekDouble>(nq);
+ Vmath::Vcopy(nq, pArray, 1, m_arrays[key], 1);
+
+ m_lastCached[key] = pTime;
+ }
+}
+
+/**
+ * @brief Provide a description of a function for a given field name.
+ *
+ * @param pFieldName Field name.
+ * @param domain The domain to evaluate the function in.
+ */
+std::string SessionFunction::Describe(std::string pFieldName, const int domain)
+{
+ std::string retVal;
+
+ LibUtilities::FunctionType vType =
+ m_session->GetFunctionType(m_name, pFieldName, domain);
+ if (vType == LibUtilities::eFunctionTypeExpression)
+ {
+ LibUtilities::EquationSharedPtr ffunc =
+ m_session->GetFunction(m_name, pFieldName, domain);
+ retVal = ffunc->GetExpression();
+ }
+ else if (vType == LibUtilities::eFunctionTypeFile ||
+ LibUtilities::eFunctionTypeTransientFile)
+ {
+ std::string filename =
+ m_session->GetFunctionFilename(m_name, pFieldName, domain);
+ retVal = "from file " + filename;
+ }
+
+ return retVal;
+}
+
+/**
+ * Evaluates a function from expression
+ *
+ * @param pFieldName The name of the field to evaluate the function for.
+ * @param pArray The array into which to write the values.
+ * @param pTime The time at which to evaluate the function.
+ * @param domain The domain to evaluate the function in.
+ */
+void SessionFunction::EvaluateExp(string pFieldName,
+ Array<OneD, NekDouble> &pArray,
+ const NekDouble &pTime,
+ const int domain)
+{
+ unsigned int nq = m_field->GetNpoints();
+ if (pArray.num_elements() < nq)
+ {
+ pArray = Array<OneD, NekDouble>(nq);
+ }
+
+ Array<OneD, NekDouble> x0(nq);
+ Array<OneD, NekDouble> x1(nq);
+ Array<OneD, NekDouble> x2(nq);
+
+ // Get the coordinates (assuming all fields have the same
+ // discretisation)
+ m_field->GetCoords(x0, x1, x2);
+ LibUtilities::EquationSharedPtr ffunc =
+ m_session->GetFunction(m_name, pFieldName, domain);
+
+ ffunc->Evaluate(x0, x1, x2, pTime, pArray);
+}
+
+/**
+ * Evaluates a function from fld file
+ *
+ * @param pFieldName The name of the field to evaluate the function for.
+ * @param pArray The array into which to write the values.
+ * @param pTime The time at which to evaluate the function.
+ * @param domain The domain to evaluate the function in.
+ */
+void SessionFunction::EvaluateFld(string pFieldName,
+ Array<OneD, NekDouble> &pArray,
+ const NekDouble &pTime,
+ const int domain)
+{
+ unsigned int nq = m_field->GetNpoints();
+ if (pArray.num_elements() < nq)
+ {
+ pArray = Array<OneD, NekDouble>(nq);
+ }
+
+ std::string filename =
+ m_session->GetFunctionFilename(m_name, pFieldName, domain);
+ std::string fileVar =
+ m_session->GetFunctionFilenameVariable(m_name, pFieldName, domain);
+
+ if (fileVar.length() == 0)
+ {
+ fileVar = pFieldName;
+ }
+
+ // In case of eFunctionTypeTransientFile, generate filename from
+ // format string
+ LibUtilities::FunctionType vType =
+ m_session->GetFunctionType(m_name, pFieldName, domain);
+ if (vType == LibUtilities::eFunctionTypeTransientFile)
+ {
+ try
+ {
+#if (defined _WIN32 && _MSC_VER < 1900)
+ // We need this to make sure boost::format has always
+ // two digits in the exponents of Scientific notation.
+ unsigned int old_exponent_format;
+ old_exponent_format = _set_output_format(_TWO_DIGIT_EXPONENT);
+ filename = boost::str(boost::format(filename) % pTime);
+ _set_output_format(old_exponent_format);
+#else
+ filename = boost::str(boost::format(filename) % pTime);
+#endif
+ }
+ catch (...)
+ {
+ ASSERTL0(false,
+ "Invalid Filename in function \"" + m_name +
+ "\", variable \"" + fileVar + "\"")
+ }
+ }
+
+ // Define list of global element ids
+ int numexp = m_field->GetExpSize();
+ Array<OneD, int> ElementGIDs(numexp);
+ for (int i = 0; i < numexp; ++i)
+ {
+ ElementGIDs[i] = m_field->GetExp(i)->GetGeom()->GetGlobalID();
+ }
+
+ std::vector<LibUtilities::FieldDefinitionsSharedPtr> FieldDef;
+ std::vector<std::vector<NekDouble> > FieldData;
+ LibUtilities::FieldIOSharedPtr fldIO =
+ LibUtilities::FieldIO::CreateForFile(m_session, filename);
+ fldIO->Import(filename,
+ FieldDef,
+ FieldData,
+ LibUtilities::NullFieldMetaDataMap,
+ ElementGIDs);
+
+ int idx = -1;
+ Array<OneD, NekDouble> vCoeffs(m_field->GetNcoeffs(), 0.0);
+ // Loop over all the expansions
+ for (int i = 0; i < FieldDef.size(); ++i)
+ {
+ // Find the index of the required field in the
+ // expansion segment
+ for (int j = 0; j < FieldDef[i]->m_fields.size(); ++j)
+ {
+ if (FieldDef[i]->m_fields[j] == fileVar)
+ {
+ idx = j;
+ }
+ }
+
+ if (idx >= 0)
+ {
+ m_field->ExtractDataToCoeffs(
+ FieldDef[i], FieldData[i], FieldDef[i]->m_fields[idx], vCoeffs);
+ }
+ else
+ {
+ cout << "Field " + fileVar + " not found." << endl;
+ }
+ }
+
+ m_field->BwdTrans_IterPerExp(vCoeffs, pArray);
+}
+
+/**
+ * Evaluates a function from pts file
+ *
+ * @param pFieldName The name of the field to evaluate the function for.
+ * @param pArray The array into which to write the values.
+ * @param pTime The time at which to evaluate the function.
+ * @param domain The domain to evaluate the function in.
+ */
+void SessionFunction::EvaluatePts(string pFieldName,
+ Array<OneD, NekDouble> &pArray,
+ const NekDouble &pTime,
+ const int domain)
+{
+ unsigned int nq = m_field->GetNpoints();
+ if (pArray.num_elements() < nq)
+ {
+ pArray = Array<OneD, NekDouble>(nq);
+ }
+
+ std::string filename =
+ m_session->GetFunctionFilename(m_name, pFieldName, domain);
+ std::string fileVar =
+ m_session->GetFunctionFilenameVariable(m_name, pFieldName, domain);
+
+ if (fileVar.length() == 0)
+ {
+ fileVar = pFieldName;
+ }
+
+ // In case of eFunctionTypeTransientFile, generate filename from
+ // format string
+ LibUtilities::FunctionType vType =
+ m_session->GetFunctionType(m_name, pFieldName, domain);
+ if (vType == LibUtilities::eFunctionTypeTransientFile)
+ {
+ try
+ {
+#if (defined _WIN32 && _MSC_VER < 1900)
+ // We need this to make sure boost::format has always
+ // two digits in the exponents of Scientific notation.
+ unsigned int old_exponent_format;
+ old_exponent_format = _set_output_format(_TWO_DIGIT_EXPONENT);
+ filename = boost::str(boost::format(filename) % pTime);
+ _set_output_format(old_exponent_format);
+#else
+ filename = boost::str(boost::format(filename) % pTime);
+#endif
+ }
+ catch (...)
+ {
+ ASSERTL0(false,
+ "Invalid Filename in function \"" + m_name +
+ "\", variable \"" + fileVar + "\"")
+ }
+ }
+
+ LibUtilities::PtsFieldSharedPtr outPts;
+ // check if we already loaded this file. For transient files,
+ // funcFilename != filename so we can make sure we only keep the
+ // latest pts field per funcFilename.
+ std::string funcFilename =
+ m_session->GetFunctionFilename(m_name, pFieldName, domain);
+
+ LibUtilities::PtsFieldSharedPtr inPts;
+ LibUtilities::PtsIO ptsIO(m_session->GetComm());
+ ptsIO.Import(filename, inPts);
+
+ Array<OneD, Array<OneD, NekDouble> > pts(inPts->GetDim() +
+ inPts->GetNFields());
+ for (int i = 0; i < inPts->GetDim() + inPts->GetNFields(); ++i)
+ {
+ pts[i] = Array<OneD, NekDouble>(nq);
+ }
+ if (inPts->GetDim() == 1)
+ {
+ m_field->GetCoords(pts[0]);
+ }
+ else if (inPts->GetDim() == 2)
+ {
+ m_field->GetCoords(pts[0], pts[1]);
+ }
+ else if (inPts->GetDim() == 3)
+ {
+ m_field->GetCoords(pts[0], pts[1], pts[2]);
+ }
+ outPts = MemoryManager<LibUtilities::PtsField>::AllocateSharedPtr(
+ inPts->GetDim(), inPts->GetFieldNames(), pts);
+
+ FieldUtils::Interpolator interp;
+ if (m_interpolators.find(funcFilename) != m_interpolators.end())
+ {
+ interp = m_interpolators[funcFilename];
+ }
+ else
+ {
+ interp = FieldUtils::Interpolator(Nektar::FieldUtils::eShepard);
+ if (m_session->GetComm()->GetRank() == 0)
+ {
+ interp.SetProgressCallback(&SessionFunction::PrintProgressbar,
+ this);
+ }
+ interp.CalcWeights(inPts, outPts);
+ if (m_session->GetComm()->GetRank() == 0)
+ {
+ cout << endl;
+ if (m_session->DefinesCmdLineArgument("verbose"))
+ {
+ interp.PrintStatistics();
+ }
+ }
+ }
+
+ if (m_toCache)
+ {
+ m_interpolators[funcFilename] = interp;
+ }
+
+ // TODO: only interpolate the field we actually want
+ interp.Interpolate(inPts, outPts);
+
+ int fieldInd;
+ vector<string> fieldNames = outPts->GetFieldNames();
+ for (fieldInd = 0; fieldInd < fieldNames.size(); ++fieldInd)
+ {
+ if (outPts->GetFieldName(fieldInd) == pFieldName)
+ {
+ break;
+ }
+ }
+ ASSERTL0(fieldInd != fieldNames.size(), "field not found");
+
+ pArray = outPts->GetPts(fieldInd + outPts->GetDim());
+}
+
+// end of namespaces
+}
+}
diff --git a/library/SolverUtils/Core/SessionFunction.h b/library/SolverUtils/Core/SessionFunction.h
new file mode 100644
index 0000000000000000000000000000000000000000..246a923e843b6c92a64960d1375293df0bab77b8
--- /dev/null
+++ b/library/SolverUtils/Core/SessionFunction.h
@@ -0,0 +1,152 @@
+///////////////////////////////////////////////////////////////////////////////
+//
+// File SessionFunction.h
+//
+// For more information, please see: http://www.nektar.info
+//
+// The MIT License
+//
+// Copyright (c) 2017 Kilian Lackhove
+// Copyright (c) 2006 Division of Applied Mathematics, Brown University (USA),
+// Department of Aeronautics, Imperial College London (UK), and Scientific
+// Computing and Imaging Institute, University of Utah (USA).
+//
+// License for the specific language governing rights and limitations under
+// Permission is hereby granted, free of charge, to any person obtaining a
+// copy of this software and associated documentation files (the "Software"),
+// to deal in the Software without restriction, including without limitation
+// the rights to use, copy, modify, merge, publish, distribute, sublicense,
+// and/or sell copies of the Software, and to permit persons to whom the
+// Software is furnished to do so, subject to the following conditions:
+//
+// The above copyright notice and this permission notice shall be included
+// in all copies or substantial portions of the Software.
+//
+// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
+// OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+// FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
+// THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+// LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
+// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
+// DEALINGS IN THE SOFTWARE.
+//
+// Description: Session Function
+//
+///////////////////////////////////////////////////////////////////////////////
+
+#ifndef NEKTAR_SOLVERUTILS_SESSIONFUNCTION_H
+#define NEKTAR_SOLVERUTILS_SESSIONFUNCTION_H
+
+#include <FieldUtils/Interpolator.h>
+#include <LibUtilities/BasicUtils/FieldIO.h>
+#include <LibUtilities/BasicUtils/FileSystem.h>
+#include <LibUtilities/BasicUtils/NekFactory.hpp>
+#include <LibUtilities/BasicUtils/Progressbar.hpp>
+#include <LibUtilities/BasicUtils/PtsField.h>
+#include <LibUtilities/BasicUtils/PtsIO.h>
+#include <LibUtilities/BasicUtils/SharedArray.hpp>
+#include <MultiRegions/ExpList.h>
+#include <SolverUtils/SolverUtilsDeclspec.h>
+
+namespace Nektar
+{
+namespace SolverUtils
+{
+
+class SessionFunction
+{
+public:
+ /// Representation of a FUNCTION defined in the session xml file.
+ SOLVER_UTILS_EXPORT SessionFunction(
+ const LibUtilities::SessionReaderSharedPtr &session,
+ const MultiRegions::ExpListSharedPtr &field,
+ std::string functionName,
+ bool toCache = false);
+
+ /// Evaluates a function defined in the xml session file at each quadrature point.
+ SOLVER_UTILS_EXPORT void Evaluate(
+ Array<OneD, Array<OneD, NekDouble> > &pArray,
+ const NekDouble pTime = 0.0,
+ const int domain = 0);
+
+ /// Evaluates a function defined in the xml session file at each quadrature point.
+ SOLVER_UTILS_EXPORT void Evaluate(
+ std::vector<std::string> pFieldNames,
+ Array<OneD, Array<OneD, NekDouble> > &pArray,
+ const NekDouble &pTime = 0.0,
+ const int domain = 0);
+
+ /// Evaluates a function defined in the xml session file at each quadrature point.
+ SOLVER_UTILS_EXPORT void Evaluate(
+ std::vector<std::string> pFieldNames,
+ Array<OneD, MultiRegions::ExpListSharedPtr> &pFields,
+ const NekDouble &pTime = 0.0,
+ const int domain = 0);
+
+ // Evaluates a function defined in the xml session file at each quadrature point.
+ SOLVER_UTILS_EXPORT void Evaluate(std::string pFieldName,
+ Array<OneD, NekDouble> &pArray,
+ const NekDouble &pTime = 0.0,
+ const int domain = 0);
+
+ /// Provide a description of a function for a given field name.
+ SOLVER_UTILS_EXPORT std::string Describe(std::string pFieldName,
+ const int domain = 0);
+
+ SOLVER_UTILS_EXPORT const LibUtilities::SessionReaderSharedPtr & GetSession()
+ {
+ return m_session;
+ }
+
+ SOLVER_UTILS_EXPORT const MultiRegions::ExpListSharedPtr & GetExpansion()
+ {
+ return m_field;
+ }
+
+private:
+ /// The session reader
+ LibUtilities::SessionReaderSharedPtr m_session;
+ /// The expansion we want to evaluate this function for
+ MultiRegions::ExpListSharedPtr m_field;
+ // Name of this function
+ std::string m_name;
+ /// Store resulting arrays (and interpolators)
+ bool m_toCache;
+ /// Last time the cache for this variable & domain combo was updated
+ std::map<std::pair<std::string, int>, NekDouble> m_lastCached;
+ /// Interpolator for pts file input for a variable & domain combination
+ std::map<std::string, FieldUtils::Interpolator> m_interpolators;
+ /// Cached result arrays
+ std::map<std::pair<std::string, int>, Array<OneD, NekDouble> > m_arrays;
+
+ // Evaluates a function from expression
+ SOLVER_UTILS_EXPORT void EvaluateExp(std::string pFieldName,
+ Array<OneD, NekDouble> &pArray,
+ const NekDouble &pTime = 0.0,
+ const int domain = 0);
+
+ // Evaluates a function from fld file
+ SOLVER_UTILS_EXPORT void EvaluateFld(std::string pFieldName,
+ Array<OneD, NekDouble> &pArray,
+ const NekDouble &pTime = 0.0,
+ const int domain = 0);
+
+ /// Evaluates a function from pts file
+ SOLVER_UTILS_EXPORT void EvaluatePts(std::string pFieldName,
+ Array<OneD, NekDouble> &pArray,
+ const NekDouble &pTime = 0.0,
+ const int domain = 0);
+
+ SOLVER_UTILS_EXPORT void PrintProgressbar(const int position,
+ const int goal) const
+ {
+ LibUtilities::PrintProgressbar(position, goal, "Interpolating");
+ }
+};
+
+typedef boost::shared_ptr<SessionFunction> SessionFunctionSharedPtr;
+static SessionFunctionSharedPtr NullSessionFunction;
+}
+}
+
+#endif
diff --git a/library/SolverUtils/EquationSystem.cpp b/library/SolverUtils/EquationSystem.cpp
index 866026f52632bcb926611246a2987c3084cb4dc7..cec6ab95bd732af67ddf1003e120509d6d9a15d4 100644
--- a/library/SolverUtils/EquationSystem.cpp
+++ b/library/SolverUtils/EquationSystem.cpp
@@ -676,469 +676,35 @@ namespace Nektar
DNekScalBlkMat, LocalRegions::MatrixKey::opLess>::ClearManager();
}
- /**
- * Evaluates a physical function at each quadrature point in the domain.
- *
- * @param pArray The array into which to write the values.
- * @param pEqn The equation to evaluate.
- */
- void EquationSystem::EvaluateFunction(
- Array<OneD, Array<OneD, NekDouble> >& pArray,
- std::string pFunctionName,
- const NekDouble pTime,
- const int domain)
- {
- ASSERTL0(m_session->DefinesFunction(pFunctionName),
- "Function '" + pFunctionName + "' does not exist.");
-
- std::vector<std::string> vFieldNames = m_session->GetVariables();
-
- for(int i = 0 ; i < vFieldNames.size(); i++)
- {
- EvaluateFunction(vFieldNames[i], pArray[i], pFunctionName,
- pTime, domain);
- }
- }
-
- /**
- * Populates a forcing function for each of the dependent variables
- * using the expression provided by the BoundaryConditions object.
- * @param force Array of fields to assign forcing.
- */
- void EquationSystem::EvaluateFunction(
- std::vector<std::string> pFieldNames,
- Array<OneD, Array<OneD, NekDouble> > &pFields,
- const std::string &pFunctionName,
- const NekDouble &pTime,
- const int domain)
- {
- ASSERTL1(pFieldNames.size() == pFields.num_elements(),
- "Function '" + pFunctionName +
- "' variable list size mismatch with array storage.");
- ASSERTL0(m_session->DefinesFunction(pFunctionName),
- "Function '" + pFunctionName + "' does not exist.");
-
- for (int i = 0; i < pFieldNames.size(); i++)
- {
- EvaluateFunction(pFieldNames[i], pFields[i], pFunctionName, pTime,
- domain);
- }
- }
-
- /**
- * Populates a function for each of the dependent variables using
- * the expression or filenames provided by the SessionReader object.
- * @param force Array of fields to assign forcing.
- */
- void EquationSystem::EvaluateFunction(
- std::vector<std::string> pFieldNames,
- Array<OneD, MultiRegions::ExpListSharedPtr> &pFields,
- const std::string &pFunctionName,
- const NekDouble &pTime,
- const int domain)
- {
- ASSERTL0(m_session->DefinesFunction(pFunctionName),
- "Function '" + pFunctionName + "' does not exist.");
- ASSERTL0(pFieldNames.size() == pFields.num_elements(),
- "Field list / name list size mismatch.");
-
- for (int i = 0; i < pFieldNames.size(); i++)
- {
- EvaluateFunction(pFieldNames[i], pFields[i]->UpdatePhys(),
- pFunctionName, pTime, domain);
- pFields[i]->FwdTrans_IterPerExp(pFields[i]->GetPhys(),
- pFields[i]->UpdateCoeffs());
- }
- }
-
- void EquationSystem::EvaluateFunction(std::string pFieldName,
- Array<OneD, NekDouble> &pArray,
- const std::string &pFunctionName,
- const NekDouble &pTime,
- const int domain)
- {
- ASSERTL0(m_session->DefinesFunction(pFunctionName),
- "Function '" + pFunctionName + "' does not exist.");
-
- LibUtilities::FunctionType vType;
- vType = m_session->GetFunctionType(pFunctionName, pFieldName, domain);
- if (vType == LibUtilities::eFunctionTypeExpression)
- {
- EvaluateFunctionExp(pFieldName, pArray, pFunctionName, pTime, domain);
- }
- else if (vType == LibUtilities::eFunctionTypeFile ||
- vType == LibUtilities::eFunctionTypeTransientFile)
- {
- std::string filename =
- m_session->GetFunctionFilename(pFunctionName, pFieldName, domain);
-
- if (boost::filesystem::path(filename).extension() != ".pts")
- {
- EvaluateFunctionFld(pFieldName, pArray, pFunctionName, pTime,
- domain);
- }
- else
- {
- EvaluateFunctionPts(pFieldName, pArray, pFunctionName, pTime,
- domain);
- }
- }
- }
-
- void EquationSystem::EvaluateFunctionExp(string pFieldName,
- Array<OneD, NekDouble> &pArray,
- const string &pFunctionName,
- const NekDouble &pTime,
- const int domain)
- {
- ASSERTL0(m_session->DefinesFunction(pFunctionName),
- "Function '" + pFunctionName + "' does not exist.");
-
- unsigned int nq = m_fields[0]->GetNpoints();
- if (pArray.num_elements() < nq)
- {
- pArray = Array<OneD, NekDouble>(nq);
- }
-
- LibUtilities::FunctionType vType;
- vType = m_session->GetFunctionType(pFunctionName, pFieldName, domain);
- ASSERTL0(vType == LibUtilities::eFunctionTypeExpression,
- "vType not eFunctionTypeExpression");
-
- Array<OneD, NekDouble> x0(nq);
- Array<OneD, NekDouble> x1(nq);
- Array<OneD, NekDouble> x2(nq);
-
- // Get the coordinates (assuming all fields have the same
- // discretisation)
- m_fields[0]->GetCoords(x0, x1, x2);
- LibUtilities::EquationSharedPtr ffunc =
- m_session->GetFunction(pFunctionName, pFieldName, domain);
-
- ffunc->Evaluate(x0, x1, x2, pTime, pArray);
- }
-
- void EquationSystem::EvaluateFunctionFld(string pFieldName,
- Array<OneD, NekDouble> &pArray,
- const string &pFunctionName,
- const NekDouble &pTime,
- const int domain)
- {
-
- ASSERTL0(m_session->DefinesFunction(pFunctionName),
- "Function '" + pFunctionName + "' does not exist.");
-
- unsigned int nq = m_fields[0]->GetNpoints();
- if (pArray.num_elements() < nq)
- {
- pArray = Array<OneD, NekDouble>(nq);
- }
-
- LibUtilities::FunctionType vType;
- vType = m_session->GetFunctionType(pFunctionName, pFieldName, domain);
- ASSERTL0(vType == LibUtilities::eFunctionTypeFile ||
- vType == LibUtilities::eFunctionTypeTransientFile,
- "vType not eFunctionTypeFile or eFunctionTypeTransientFile");
-
- std::string filename =
- m_session->GetFunctionFilename(pFunctionName, pFieldName, domain);
- std::string fileVar = m_session->GetFunctionFilenameVariable(
- pFunctionName, pFieldName, domain);
-
- if (fileVar.length() == 0)
- {
- fileVar = pFieldName;
- }
-
- // In case of eFunctionTypeTransientFile, generate filename from
- // format string
- if (vType == LibUtilities::eFunctionTypeTransientFile)
- {
- try
- {
-#if (defined _WIN32 && _MSC_VER < 1900)
- // We need this to make sure boost::format has always
- // two digits in the exponents of Scientific notation.
- unsigned int old_exponent_format;
- old_exponent_format = _set_output_format(_TWO_DIGIT_EXPONENT);
- filename = boost::str(boost::format(filename) % pTime);
- _set_output_format(old_exponent_format);
-#else
- filename = boost::str(boost::format(filename) % pTime);
-#endif
- }
- catch (...)
- {
- ASSERTL0(false, "Invalid Filename in function \"" + pFunctionName +
- "\", variable \"" + fileVar + "\"")
- }
- }
-
- std::vector<LibUtilities::FieldDefinitionsSharedPtr> FieldDef;
- std::vector<std::vector<NekDouble> > FieldData;
-
- // Define list of global element ids
- int numexp = m_fields[0]->GetExpSize();
- Array<OneD, int> ElementGIDs(numexp);
- for (int i = 0; i < numexp; ++i)
- {
- ElementGIDs[i] = m_fields[0]->GetExp(i)->GetGeom()->GetGlobalID();
- }
-
- // check if we already loaded this file. For transient files,
- // funcFilename != filename so we can make sure we only keep the
- // latest field per funcFilename.
- std::string funcFilename =
- m_session->GetFunctionFilename(pFunctionName, pFieldName, domain);
- if (m_loadedFldFields.find(funcFilename) != m_loadedFldFields.end())
- {
- // found
- if (m_loadedFldFields[funcFilename].first == filename)
- {
- // found
- FieldDef = m_loadedFldFields[funcFilename].second.fieldDef;
- FieldData = m_loadedFldFields[funcFilename].second.fieldData;
- }
- else
- {
- LibUtilities::FieldIOSharedPtr pts_fld =
- LibUtilities::FieldIO::CreateForFile(m_session, filename);
- pts_fld->Import(filename, FieldDef, FieldData,
- LibUtilities::NullFieldMetaDataMap,
- ElementGIDs);
- }
- }
- else
- {
- LibUtilities::FieldIOSharedPtr pts_fld =
- LibUtilities::FieldIO::CreateForFile(m_session, filename);
- pts_fld->Import(filename, FieldDef, FieldData,
- LibUtilities::NullFieldMetaDataMap,
- ElementGIDs);
- }
- // Now we overwrite the field we had in
- // m_loadedFldFields[funcFilename] before, making sure we only keep
- // one field per funcFilename in memory
- m_loadedFldFields[funcFilename].first = filename;
- m_loadedFldFields[funcFilename].second.fieldDef = FieldDef;
- m_loadedFldFields[funcFilename].second.fieldData = FieldData;
-
- int idx = -1;
- Array<OneD, NekDouble> vCoeffs(m_fields[0]->GetNcoeffs(), 0.0);
- // Loop over all the expansions
- for (int i = 0; i < FieldDef.size(); ++i)
- {
- // Find the index of the required field in the
- // expansion segment
- for (int j = 0; j < FieldDef[i]->m_fields.size(); ++j)
- {
- if (FieldDef[i]->m_fields[j] == fileVar)
- {
- idx = j;
- }
- }
-
- if (idx >= 0)
- {
- m_fields[0]->ExtractDataToCoeffs(
- FieldDef[i], FieldData[i], FieldDef[i]->m_fields[idx], vCoeffs);
- }
- else
- {
- cout << "Field " + fileVar + " not found." << endl;
- }
- }
-
- m_fields[0]->BwdTrans_IterPerExp(vCoeffs, pArray);
- }
-
- void EquationSystem::EvaluateFunctionPts(string pFieldName,
- Array<OneD, NekDouble> &pArray,
- const string &pFunctionName,
- const NekDouble &pTime,
- const int domain)
+ SessionFunctionSharedPtr EquationSystem::GetFunction(
+ std::string name,
+ const MultiRegions::ExpListSharedPtr &field,
+ bool cache)
{
-
- ASSERTL0(m_session->DefinesFunction(pFunctionName),
- "Function '" + pFunctionName + "' does not exist.");
-
- unsigned int nq = m_fields[0]->GetNpoints();
- if (pArray.num_elements() < nq)
- {
- pArray = Array<OneD, NekDouble>(nq);
- }
-
- LibUtilities::FunctionType vType;
- vType = m_session->GetFunctionType(pFunctionName, pFieldName, domain);
- ASSERTL0(vType == LibUtilities::eFunctionTypeFile ||
- vType == LibUtilities::eFunctionTypeTransientFile,
- "vType not eFunctionTypeFile or eFunctionTypeTransientFile");
-
- std::string filename =
- m_session->GetFunctionFilename(pFunctionName, pFieldName, domain);
- std::string fileVar = m_session->GetFunctionFilenameVariable(
- pFunctionName, pFieldName, domain);
-
- if (fileVar.length() == 0)
+ MultiRegions::ExpListSharedPtr vField = field;
+ if (!field)
{
- fileVar = pFieldName;
+ vField = m_fields[0];
}
- // In case of eFunctionTypeTransientFile, generate filename from
- // format string
- if (vType == LibUtilities::eFunctionTypeTransientFile)
+ if (cache)
{
- try
+ if ((m_sessionFunctions.find(name) == m_sessionFunctions.end())
+ || (m_sessionFunctions[name]->GetSession() != m_session)
+ || (m_sessionFunctions[name]->GetExpansion() != vField)
+ )
{
-#if (defined _WIN32 && _MSC_VER < 1900)
- // We need this to make sure boost::format has always
- // two digits in the exponents of Scientific notation.
- unsigned int old_exponent_format;
- old_exponent_format = _set_output_format(_TWO_DIGIT_EXPONENT);
- filename = boost::str(boost::format(filename) % pTime);
- _set_output_format(old_exponent_format);
-#else
- filename = boost::str(boost::format(filename) % pTime);
-#endif
+ m_sessionFunctions[name] =
+ MemoryManager<SessionFunction>::AllocateSharedPtr(
+ m_session, vField, name, cache);
}
- catch (...)
- {
- ASSERTL0(false, "Invalid Filename in function \"" + pFunctionName +
- "\", variable \"" + fileVar + "\"")
- }
- }
- LibUtilities::PtsFieldSharedPtr outPts;
- // check if we already loaded this file. For transient files,
- // funcFilename != filename so we can make sure we only keep the
- // latest pts field per funcFilename.
- std::string funcFilename =
- m_session->GetFunctionFilename(pFunctionName, pFieldName, domain);
-
- if (m_loadedPtsFields.find(funcFilename) != m_loadedPtsFields.end())
- {
- // found
- if (m_loadedPtsFields[funcFilename].first == filename)
- {
- // found
- outPts = m_loadedPtsFields[funcFilename].second;
- }
- else
- {
- LoadPts(funcFilename, filename, outPts);
- }
+ return m_sessionFunctions[name];
}
else
{
- LoadPts(funcFilename, filename, outPts);
- }
- // Now we overwrite the field we had in
- // m_loadedPtsFields[funcFilename] before, making sure we only keep
- // one field per funcFilename in memory
- m_loadedPtsFields[funcFilename].first = filename;
- m_loadedPtsFields[funcFilename].second = outPts;
-
- int fieldInd;
- vector<string> fieldNames = outPts->GetFieldNames();
- for (fieldInd = 0; fieldInd < fieldNames.size(); ++fieldInd)
- {
- if (outPts->GetFieldName(fieldInd) == pFieldName)
- {
- break;
- }
- }
- ASSERTL0(fieldInd != fieldNames.size(), "field not found");
-
- pArray = outPts->GetPts(fieldInd + outPts->GetDim());
- }
-
- void EquationSystem::LoadPts(std::string funcFilename,
- std::string filename,
- LibUtilities::PtsFieldSharedPtr &outPts)
- {
- unsigned int nq = m_fields[0]->GetNpoints();
-
- LibUtilities::PtsFieldSharedPtr inPts;
- LibUtilities::PtsIO ptsIO(m_session->GetComm());
- ptsIO.Import(filename, inPts);
-
- Array<OneD, Array<OneD, NekDouble> > pts(inPts->GetDim() +
- inPts->GetNFields());
- for (int i = 0; i < inPts->GetDim() + inPts->GetNFields(); ++i)
- {
- pts[i] = Array<OneD, NekDouble>(nq);
- }
- if (inPts->GetDim() == 1)
- {
- m_fields[0]->GetCoords(pts[0]);
- }
- else if (inPts->GetDim() == 2)
- {
- m_fields[0]->GetCoords(pts[0], pts[1]);
- }
- else if (inPts->GetDim() == 3)
- {
- m_fields[0]->GetCoords(pts[0], pts[1], pts[2]);
- }
- outPts = MemoryManager<LibUtilities::PtsField>::AllocateSharedPtr(
- inPts->GetDim(), inPts->GetFieldNames(), pts);
-
- // check if we already have an interolator for this funcFilename
- if (m_interpolators.find(funcFilename) == m_interpolators.end())
- {
- m_interpolators[funcFilename] =
- FieldUtils::Interpolator(Nektar::FieldUtils::eShepard);
- if (m_comm->GetRank() == 0)
- {
- m_interpolators[funcFilename].SetProgressCallback(
- &EquationSystem::PrintProgressbar, this);
- }
- m_interpolators[funcFilename].CalcWeights(inPts, outPts);
- if (m_comm->GetRank() == 0)
- {
- cout << endl;
- if (GetSession()->DefinesCmdLineArgument("verbose"))
- {
- m_interpolators[funcFilename].PrintStatistics();
- }
- }
- }
- m_interpolators[funcFilename].Interpolate(inPts, outPts);
- }
-
-
- /**
- * @brief Provide a description of a function for a given field name.
- *
- * @param pFieldName Field name.
- * @param pFunctionName Function name.
- */
- std::string EquationSystem::DescribeFunction(
- std::string pFieldName,
- const std::string &pFunctionName,
- const int domain)
- {
- ASSERTL0(m_session->DefinesFunction(pFunctionName),
- "Function '" + pFunctionName + "' does not exist.");
-
- std::string retVal;
- LibUtilities::FunctionType vType;
-
- vType = m_session->GetFunctionType(pFunctionName, pFieldName);
- if (vType == LibUtilities::eFunctionTypeExpression)
- {
- LibUtilities::EquationSharedPtr ffunc
- = m_session->GetFunction(pFunctionName, pFieldName,domain);
- retVal = ffunc->GetExpression();
+ return SessionFunctionSharedPtr(new SessionFunction(m_session, vField, name, cache));
}
- else if (vType == LibUtilities::eFunctionTypeFile)
- {
- std::string filename
- = m_session->GetFunctionFilename(pFunctionName, pFieldName,domain);
- retVal = "from file " + filename;
- }
-
- return retVal;
}
/**
@@ -1188,8 +754,7 @@ namespace Nektar
Array<OneD, NekDouble>
exactsoln(m_fields[field]->GetNpoints());
- EvaluateFunction(m_session->GetVariable(field), exactsoln,
- "ExactSolution", m_time);
+ GetFunction("ExactSolution")->Evaluate(m_session->GetVariable(field), exactsoln, m_time);
L2error = m_fields[field]->L2(m_fields[field]->GetPhys(), exactsoln);
}
@@ -1247,8 +812,7 @@ namespace Nektar
Array<OneD, NekDouble>
exactsoln(m_fields[field]->GetNpoints());
- EvaluateFunction(m_session->GetVariable(field), exactsoln,
- "ExactSolution", m_time);
+ GetFunction("ExactSolution")->Evaluate(m_session->GetVariable(field), exactsoln, m_time);
Linferror = m_fields[field]->Linf(m_fields[field]->GetPhys(), exactsoln);
}
@@ -1356,8 +920,7 @@ namespace Nektar
if (m_session->DefinesFunction("InitialConditions"))
{
- EvaluateFunction(m_session->GetVariables(), m_fields,
- "InitialConditions", m_time, domain);
+ GetFunction("InitialConditions")->Evaluate(m_session->GetVariables(), m_fields, m_time, domain);
if (m_session->GetComm()->GetRank() == 0)
{
@@ -1366,7 +929,7 @@ namespace Nektar
{
std::string varName = m_session->GetVariable(i);
cout << " - Field " << varName << ": "
- << DescribeFunction(varName, "InitialConditions",domain)
+ << GetFunction("InitialConditions")->Describe(varName, domain)
<< endl;
}
}
@@ -1406,8 +969,7 @@ namespace Nektar
Vmath::Zero(outfield.num_elements(), outfield, 1);
if (m_session->DefinesFunction("ExactSolution"))
{
- EvaluateFunction(m_session->GetVariable(field), outfield,
- "ExactSolution", time);
+ GetFunction("ExactSolution")->Evaluate(m_session->GetVariable(field), outfield, time);
}
}
@@ -1431,7 +993,7 @@ namespace Nektar
vel.push_back("Vz");
vel.resize(m_spacedim);
SetUpBaseFields(m_graph);
- EvaluateFunction(vel, base, "BaseFlow");
+ GetFunction("BaseFlow")->Evaluate(vel, base);
}
void EquationSystem::SetUpBaseFields(
diff --git a/library/SolverUtils/EquationSystem.h b/library/SolverUtils/EquationSystem.h
index 6e8777469cd8ccfbc31f8adc071868fc0fa70f5d..fc0063a52cbb19e720d33919665e6a4fd2cec12e 100644
--- a/library/SolverUtils/EquationSystem.h
+++ b/library/SolverUtils/EquationSystem.h
@@ -48,6 +48,7 @@
#include <MultiRegions/ExpList.h>
#include <SolverUtils/SolverUtilsDeclspec.h>
#include <SolverUtils/Core/Misc.h>
+#include <SolverUtils/Core/SessionFunction.h>
namespace Nektar
{
@@ -68,11 +69,6 @@ class Interpolator;
> EquationSystemFactory;
SOLVER_UTILS_EXPORT EquationSystemFactory& GetEquationSystemFactory();
- struct loadedFldField {
- std::vector<LibUtilities::FieldDefinitionsSharedPtr> fieldDef;
- std::vector<std::vector<NekDouble> > fieldData;
- } ;
-
/// A base class for describing how to solve specific equations.
class EquationSystem : public boost::enable_shared_from_this<EquationSystem>
{
@@ -143,43 +139,12 @@ class Interpolator;
/// Set parameter m_lambda
SOLVER_UTILS_EXPORT inline void SetLambda(NekDouble lambda);
-
- /// Evaluates a function as specified in the session file.
- SOLVER_UTILS_EXPORT void EvaluateFunction(
- Array<OneD, Array<OneD, NekDouble> >& pArray,
- std::string pFunctionName,
- const NekDouble pTime = 0.0,
- const int domain = 0);
-
- /// Populate given fields with the function from session.
- SOLVER_UTILS_EXPORT void EvaluateFunction(
- std::vector<std::string> pFieldNames,
- Array<OneD, Array<OneD, NekDouble> > &pFields,
- const std::string& pName,
- const NekDouble& pTime = 0.0,
- const int domain = 0);
-
- /// Populate given fields with the function from session.
- SOLVER_UTILS_EXPORT void EvaluateFunction(
- std::vector<std::string> pFieldNames,
- Array<OneD, MultiRegions::ExpListSharedPtr> &pFields,
- const std::string& pName,
- const NekDouble& pTime = 0.0,
- const int domain = 0);
-
- // Populate an array with a function variable from session.
- SOLVER_UTILS_EXPORT void EvaluateFunction(
- std::string pFieldName,
- Array<OneD, NekDouble>& pArray,
- const std::string& pFunctionName,
- const NekDouble& pTime = 0.0,
- const int domain = 0);
-
- // Describe a function.
- SOLVER_UTILS_EXPORT std::string DescribeFunction(
- std::string pFieldName,
- const std::string &pFunctionName,
- const int domain);
+
+ /// Get a SessionFunction by name
+ SOLVER_UTILS_EXPORT SessionFunctionSharedPtr GetFunction(
+ std::string name,
+ const MultiRegions::ExpListSharedPtr &field = MultiRegions::NullExpListSharedPtr,
+ bool cache = false);
/// Perform initialisation of the base flow.
SOLVER_UTILS_EXPORT void InitialiseBaseFlow(
@@ -458,14 +423,10 @@ class Interpolator;
LibUtilities::CommSharedPtr m_comm;
/// The session reader
LibUtilities::SessionReaderSharedPtr m_session;
+ /// Map of known SessionFunctions
+ std::map<std::string, SolverUtils::SessionFunctionSharedPtr> m_sessionFunctions;
/// Field input/output
LibUtilities::FieldIOSharedPtr m_fld;
- /// Map of interpolator objects
- std::map<std::string, FieldUtils::Interpolator > m_interpolators;
- /// pts fields we already read from disk: {funcFilename: (filename, ptsfield)}
- std::map<std::string, std::pair<std::string, LibUtilities::PtsFieldSharedPtr> > m_loadedPtsFields;
- // fld fiels already loaded from disk: {funcFilename: (filename, loadedFldField)}
- std::map<std::string, std::pair<std::string, loadedFldField> > m_loadedFldFields;
/// Array holding all dependent variables.
Array<OneD, MultiRegions::ExpListSharedPtr> m_fields;
/// Base fields.
@@ -605,34 +566,6 @@ class Interpolator;
unsigned int field,
Array<OneD, NekDouble> &outfield,
const NekDouble time);
-
- // Populate an array with a function variable from session.
- SOLVER_UTILS_EXPORT void EvaluateFunctionExp(
- std::string pFieldName,
- Array<OneD, NekDouble>& pArray,
- const std::string& pFunctionName,
- const NekDouble& pTime = 0.0,
- const int domain = 0);
-
- // Populate an array with a function variable from session.
- SOLVER_UTILS_EXPORT void EvaluateFunctionFld(
- std::string pFieldName,
- Array<OneD, NekDouble>& pArray,
- const std::string& pFunctionName,
- const NekDouble& pTime = 0.0,
- const int domain = 0);
-
- SOLVER_UTILS_EXPORT void EvaluateFunctionPts(
- std::string pFieldName,
- Array<OneD, NekDouble>& pArray,
- const std::string& pFunctionName,
- const NekDouble& pTime = 0.0,
- const int domain = 0);
-
- SOLVER_UTILS_EXPORT void LoadPts(
- std::string funcFilename,
- std::string filename,
- Nektar::LibUtilities::PtsFieldSharedPtr &outPts);
//Initialise m_base in order to store the base flow from a file
SOLVER_UTILS_EXPORT void SetUpBaseFields(SpatialDomains::MeshGraphSharedPtr &mesh);
diff --git a/library/SolverUtils/Forcing/Forcing.cpp b/library/SolverUtils/Forcing/Forcing.cpp
index 1fd32714ecc61d5e86eda74be21eec587cc313b9..dc2c0da7911717cc7f2eada892cbc0eaa03be2a4 100644
--- a/library/SolverUtils/Forcing/Forcing.cpp
+++ b/library/SolverUtils/Forcing/Forcing.cpp
@@ -153,79 +153,29 @@ namespace Nektar
pEqn->Evaluate(x0, x0, x0, pTime, pArray);
}
-
-
- void Forcing::EvaluateFunction(
- Array<OneD, MultiRegions::ExpListSharedPtr> pFields,
- LibUtilities::SessionReaderSharedPtr pSession,
- std::string pFieldName,
- Array<OneD, NekDouble>& pArray,
- const std::string& pFunctionName,
- NekDouble pTime)
+ SessionFunctionSharedPtr Forcing::GetFunction(
+ const Array<OneD, MultiRegions::ExpListSharedPtr> &pFields,
+ const LibUtilities::SessionReaderSharedPtr &pSession,
+ std::string pName,
+ bool pCache)
{
- ASSERTL0(pSession->DefinesFunction(pFunctionName),
- "Function '" + pFunctionName + "' does not exist.");
-
- unsigned int nq = pFields[0]->GetNpoints();
- if (pArray.num_elements() != nq)
+ if (pCache)
{
- pArray = Array<OneD, NekDouble> (nq);
- }
+ if ((m_sessionFunctions.find(pName) == m_sessionFunctions.end())
+ || (m_sessionFunctions[pName]->GetSession() != pSession)
+ || (m_sessionFunctions[pName]->GetExpansion() != pFields[0])
+ )
+ {
+ m_sessionFunctions[pName] =
+ MemoryManager<SessionFunction>::AllocateSharedPtr(
+ pSession, pFields[0], pName, pCache);
+ }
- LibUtilities::FunctionType vType;
- vType = pSession->GetFunctionType(pFunctionName, pFieldName);
- if (vType == LibUtilities::eFunctionTypeExpression)
- {
- Array<OneD, NekDouble> x0(nq);
- Array<OneD, NekDouble> x1(nq);
- Array<OneD, NekDouble> x2(nq);
-
- pFields[0]->GetCoords(x0, x1, x2);
- LibUtilities::EquationSharedPtr ffunc =
- pSession->GetFunction(pFunctionName, pFieldName);
-
- ffunc->Evaluate(x0, x1, x2, pTime, pArray);
+ return m_sessionFunctions[pName];
}
- else if (vType == LibUtilities::eFunctionTypeFile)
+ else
{
- std::string filename = pSession->GetFunctionFilename(
- pFunctionName,
- pFieldName);
-
- std::vector<LibUtilities::FieldDefinitionsSharedPtr> FieldDef;
- std::vector<std::vector<NekDouble> > FieldData;
- Array<OneD, NekDouble> vCoeffs(pFields[0]->GetNcoeffs());
- Vmath::Zero(vCoeffs.num_elements(), vCoeffs, 1);
-
- LibUtilities::FieldIOSharedPtr fld =
- LibUtilities::FieldIO::CreateForFile(m_session, filename);
- fld->Import(filename, FieldDef, FieldData);
-
- int idx = -1;
- for (int i = 0; i < FieldDef.size(); ++i)
- {
- for (int j = 0; j < FieldDef[i]->m_fields.size(); ++j)
- {
- if (FieldDef[i]->m_fields[j] == pFieldName)
- {
- idx = j;
- }
- }
-
- if (idx >= 0)
- {
- pFields[0]->ExtractDataToCoeffs(
- FieldDef[i],
- FieldData[i],
- FieldDef[i]->m_fields[idx],
- vCoeffs);
- }
- else
- {
- cout << "Field " + pFieldName + " not found." << endl;
- }
- }
- pFields[0]->BwdTrans_IterPerExp(vCoeffs, pArray);
+ return SessionFunctionSharedPtr(new SessionFunction(pSession, pFields[0], pName, pCache));
}
}
diff --git a/library/SolverUtils/Forcing/Forcing.h b/library/SolverUtils/Forcing/Forcing.h
index 905dee1effd5e652819cc90d2cddb4683e1e475c..3b0c726fa4c2a902a6183f63725b4671503f66d0 100644
--- a/library/SolverUtils/Forcing/Forcing.h
+++ b/library/SolverUtils/Forcing/Forcing.h
@@ -41,6 +41,7 @@
#include <LibUtilities/BasicUtils/NekFactory.hpp>
#include <LibUtilities/BasicUtils/SharedArray.hpp>
#include <MultiRegions/ExpList.h>
+#include <SolverUtils/Core/SessionFunction.h>
#include <SolverUtils/SolverUtilsDeclspec.h>
namespace Nektar
@@ -102,6 +103,8 @@ namespace SolverUtils
Array<OneD, Array<OneD, NekDouble> > m_Forcing;
/// Number of variables
int m_NumVariable;
+ /// Map of known SessionFunctions
+ std::map<std::string, SolverUtils::SessionFunctionSharedPtr> m_sessionFunctions;
/// Constructor
SOLVER_UTILS_EXPORT Forcing(
@@ -118,13 +121,12 @@ namespace SolverUtils
Array<OneD, Array<OneD, NekDouble> > &outarray,
const NekDouble &time)=0;
- SOLVER_UTILS_EXPORT void EvaluateFunction(
- Array<OneD, MultiRegions::ExpListSharedPtr> pFields,
- LibUtilities::SessionReaderSharedPtr pSession,
- std::string pFieldName,
- Array<OneD, NekDouble>& pArray,
- const std::string& pFunctionName,
- NekDouble pTime = NekDouble(0));
+ /// Get a SessionFunction by name
+ SOLVER_UTILS_EXPORT SessionFunctionSharedPtr GetFunction(
+ const Array<OneD, MultiRegions::ExpListSharedPtr> &pFields,
+ const LibUtilities::SessionReaderSharedPtr &pSession,
+ std::string pName,
+ bool pCache = false);
SOLVER_UTILS_EXPORT void EvaluateTimeFunction(
LibUtilities::SessionReaderSharedPtr pSession,
diff --git a/library/SolverUtils/Forcing/ForcingAbsorption.cpp b/library/SolverUtils/Forcing/ForcingAbsorption.cpp
index 6b632e9882e24e11c32fd399a4cb88fc71f01695..dc62e044ed753d3c5e21c1ee6655f941f15eb128 100644
--- a/library/SolverUtils/Forcing/ForcingAbsorption.cpp
+++ b/library/SolverUtils/Forcing/ForcingAbsorption.cpp
@@ -81,8 +81,7 @@ namespace SolverUtils
"Variable '" + s_FieldStr + "' not defined.");
m_Absorption[i] = Array<OneD, NekDouble> (npts, 0.0);
m_Forcing[i] = Array<OneD, NekDouble> (npts, 0.0);
- EvaluateFunction(pFields, m_session, s_FieldStr,
- m_Absorption[i], funcName);
+ GetFunction(pFields, m_session, funcName)->Evaluate(s_FieldStr, m_Absorption[i]);
}
funcNameElmt = pForce->FirstChildElement("REFFLOW");
@@ -99,8 +98,7 @@ namespace SolverUtils
ASSERTL0(m_session->DefinesFunction(funcName, s_FieldStr),
"Variable '" + s_FieldStr + "' not defined.");
m_Refflow[i] = Array<OneD, NekDouble> (npts, 0.0);
- EvaluateFunction(pFields, m_session, s_FieldStr,
- m_Refflow[i], funcName);
+ GetFunction(pFields, m_session, funcName)->Evaluate(s_FieldStr, m_Refflow[i]);
}
m_hasRefFlow = true;
}
diff --git a/library/SolverUtils/Forcing/ForcingBody.cpp b/library/SolverUtils/Forcing/ForcingBody.cpp
index de4463a683c05542eb7f2de61caf10d520b7a792..657d1dcb078d7e1769c899978a3747b1bdeea472 100644
--- a/library/SolverUtils/Forcing/ForcingBody.cpp
+++ b/library/SolverUtils/Forcing/ForcingBody.cpp
@@ -129,8 +129,7 @@ namespace SolverUtils
std::string s_FieldStr = m_session->GetVariable(i);
ASSERTL0(m_session->DefinesFunction(m_funcName, s_FieldStr),
"Variable '" + s_FieldStr + "' not defined.");
- EvaluateFunction(pFields, m_session, s_FieldStr,
- m_Forcing[i], m_funcName, time);
+ GetFunction(pFields, m_session, m_funcName, true)->Evaluate(s_FieldStr, m_Forcing[i], time);
}
// If singleMode or halfMode, transform the forcing term to be in
diff --git a/solvers/ADRSolver/EquationSystems/CFLtester.cpp b/solvers/ADRSolver/EquationSystems/CFLtester.cpp
index e66f4ca49a32775f192ca30ccb06de0fe065026e..0426e238714762740ca75a1f5266b2445e5837ac 100644
--- a/solvers/ADRSolver/EquationSystems/CFLtester.cpp
+++ b/solvers/ADRSolver/EquationSystems/CFLtester.cpp
@@ -62,7 +62,7 @@ namespace Nektar
vel.push_back("Vz");
vel.resize(m_spacedim);
- EvaluateFunction(vel, m_velocity, "AdvectionVelocity");
+ GetFunction( "AdvectionVelocity")->Evaluate(vel, m_velocity);
if (m_explicitAdvection)
{
diff --git a/solvers/ADRSolver/EquationSystems/EigenValuesAdvection.cpp b/solvers/ADRSolver/EquationSystems/EigenValuesAdvection.cpp
index 379e67d1725dbac06a00b35a7cd0ab87afa81f46..5d036d20e7d9a59e76684d964868bbcd8c20c7ac 100644
--- a/solvers/ADRSolver/EquationSystems/EigenValuesAdvection.cpp
+++ b/solvers/ADRSolver/EquationSystems/EigenValuesAdvection.cpp
@@ -61,7 +61,7 @@ namespace Nektar
vel.push_back("Vz");
vel.resize(m_spacedim);
- EvaluateFunction(vel, m_velocity, "AdvectionVelocity");
+ GetFunction( "AdvectionVelocity")->Evaluate(vel, m_velocity);
}
void EigenValuesAdvection::v_DoInitialise()
diff --git a/solvers/ADRSolver/EquationSystems/Poisson.cpp b/solvers/ADRSolver/EquationSystems/Poisson.cpp
index 579909c24b7f5d14d1d7dff31f2ff07b33620a35..9a7e8448f54f871d67a974e5ba27dd2f99e9a3ba 100644
--- a/solvers/ADRSolver/EquationSystems/Poisson.cpp
+++ b/solvers/ADRSolver/EquationSystems/Poisson.cpp
@@ -52,7 +52,7 @@ namespace Nektar
{
Laplace::v_InitObject();
- EvaluateFunction(m_session->GetVariables(), m_fields, "Forcing");
+ GetFunction("Forcing")->Evaluate(m_session->GetVariables(), m_fields);
}
Poisson::~Poisson()
diff --git a/solvers/ADRSolver/EquationSystems/Projection.cpp b/solvers/ADRSolver/EquationSystems/Projection.cpp
index 78d86fccead0f59f77e3c892b7ea200aff0ed5ff..0c0ca8e56129acac445ff6a88679a79e0e3b4be2 100644
--- a/solvers/ADRSolver/EquationSystems/Projection.cpp
+++ b/solvers/ADRSolver/EquationSystems/Projection.cpp
@@ -51,7 +51,7 @@ void Projection::v_InitObject()
{
EquationSystem::v_InitObject();
- EvaluateFunction(m_session->GetVariables(), m_fields, "Forcing");
+ GetFunction("Forcing")->Evaluate(m_session->GetVariables(), m_fields);
}
Projection::~Projection()
diff --git a/solvers/ADRSolver/EquationSystems/SteadyAdvectionDiffusion.cpp b/solvers/ADRSolver/EquationSystems/SteadyAdvectionDiffusion.cpp
index c603b166f6df2b144787fe74cbc4ec4adbd13b64..b886cfa30a2a025ddfd63740086c480c58d9ce8f 100644
--- a/solvers/ADRSolver/EquationSystems/SteadyAdvectionDiffusion.cpp
+++ b/solvers/ADRSolver/EquationSystems/SteadyAdvectionDiffusion.cpp
@@ -77,7 +77,7 @@ namespace Nektar
void SteadyAdvectionDiffusion::v_DoInitialise()
{
// set initial forcing from session file
- EvaluateFunction(m_session->GetVariables(), m_fields, "Forcing");
+ GetFunction("Forcing")->Evaluate(m_session->GetVariables(), m_fields);
}
void SteadyAdvectionDiffusion::v_DoSolve()
diff --git a/solvers/ADRSolver/EquationSystems/UnsteadyAdvection.cpp b/solvers/ADRSolver/EquationSystems/UnsteadyAdvection.cpp
index 6f9feed1efc252d1f850742908965aa1d2276658..7dc5aad8c0834618a544cc2dadf7923bb47d8628 100644
--- a/solvers/ADRSolver/EquationSystems/UnsteadyAdvection.cpp
+++ b/solvers/ADRSolver/EquationSystems/UnsteadyAdvection.cpp
@@ -74,7 +74,7 @@ namespace Nektar
// Store in the global variable m_velocity the advection velocities
m_velocity = Array<OneD, Array<OneD, NekDouble> >(m_spacedim);
- EvaluateFunction(vel, m_velocity, "AdvectionVelocity");
+ GetFunction( "AdvectionVelocity")->Evaluate(vel, m_velocity);
// Type of advection class to be used
switch(m_projectionType)
diff --git a/solvers/ADRSolver/EquationSystems/UnsteadyAdvectionDiffusion.cpp b/solvers/ADRSolver/EquationSystems/UnsteadyAdvectionDiffusion.cpp
index ac0afcce65a5d62cb2dac661c8606cb95ebfa0f0..7eeac15b60b7fa891c98004a77dc7ef2e53b55d9 100644
--- a/solvers/ADRSolver/EquationSystems/UnsteadyAdvectionDiffusion.cpp
+++ b/solvers/ADRSolver/EquationSystems/UnsteadyAdvectionDiffusion.cpp
@@ -77,7 +77,7 @@ namespace Nektar
vel.push_back("Vz");
vel.resize(m_spacedim);
- EvaluateFunction(vel, m_velocity, "AdvectionVelocity");
+ GetFunction( "AdvectionVelocity")->Evaluate(vel, m_velocity);
m_session->MatchSolverInfo(
"SpectralVanishingViscosity", "True", m_useSpecVanVisc, false);
diff --git a/solvers/APESolver/EquationSystems/APE.cpp b/solvers/APESolver/EquationSystems/APE.cpp
index 92661ea3cf994155fb30e1e62044a22dda9e4ab2..16da1cab5b4b0b3eb9149dabb402ecd7ed0643bc 100755
--- a/solvers/APESolver/EquationSystems/APE.cpp
+++ b/solvers/APESolver/EquationSystems/APE.cpp
@@ -137,7 +137,7 @@ void APE::v_InitObject()
{
m_bf[i] = Array<OneD, NekDouble>(GetTotPoints());
}
- EvaluateFunction(m_bfNames, m_bf, "Baseflow", m_time);
+ GetFunction("Baseflow", m_bfField, true)->Evaluate(m_bfNames, m_bf, m_time);
m_forcing = SolverUtils::Forcing::Load(m_session, m_fields, m_spacedim + 1);
@@ -301,7 +301,7 @@ void APE::GetFluxVector(
*/
bool APE::v_PreIntegrate(int step)
{
- EvaluateFunction(m_bfNames, m_bf, "Baseflow", m_time);
+ GetFunction("Baseflow", m_bfField, true)->Evaluate(m_bfNames, m_bf, m_time);
Array<OneD, NekDouble> tmpC(GetNcoeffs());
std::vector<SolverUtils::ForcingSharedPtr>::const_iterator x;
diff --git a/solvers/CardiacEPSolver/EquationSystems/BidomainRoth.cpp b/solvers/CardiacEPSolver/EquationSystems/BidomainRoth.cpp
index 1897d4e9cd365a114f1fc39c181893f8da6d9fcc..39ef519890064090721bab0f4e889d8e831f2116 100644
--- a/solvers/CardiacEPSolver/EquationSystems/BidomainRoth.cpp
+++ b/solvers/CardiacEPSolver/EquationSystems/BidomainRoth.cpp
@@ -158,8 +158,8 @@ void BidomainRoth::v_InitObject()
aniso_var[j]),
"Function 'AnisotropicConductivity' not correctly "
"defined.");
- EvaluateFunction(aniso_var[j], vTemp_j,
- "ExtracellularAnisotropicConductivity");
+
+ GetFunction("ExtracellularAnisotropicConductivity")->Evaluate(aniso_var[j], vTemp_j);
// Loop through rows of D
for (int i = 0; i < j + 1; ++i)
@@ -170,8 +170,7 @@ void BidomainRoth::v_InitObject()
"Function 'ExtracellularAnisotropicConductivity' not "
"correctly defined.");
- EvaluateFunction(aniso_var[i], vTemp_i,
- "ExtracellularAnisotropicConductivity");
+ GetFunction("ExtracellularAnisotropicConductivity")->Evaluate(aniso_var[i], vTemp_i);
Vmath::Vmul(nq, vTemp_i, 1, vTemp_j, 1,
m_vardiffe[varCoeffEnum[k]], 1);
@@ -232,8 +231,7 @@ void BidomainRoth::v_InitObject()
"Function 'IntracellularAnisotropicConductivity' not "
"correctly defined.");
- EvaluateFunction(aniso_var[j], vTemp_j,
- "IntracellularAnisotropicConductivity");
+ GetFunction("IntracellularAnisotropicConductivity")->Evaluate(aniso_var[j], vTemp_j);
// Loop through rows of D
for (int i = 0; i < j + 1; ++i)
@@ -243,8 +241,7 @@ void BidomainRoth::v_InitObject()
aniso_var[i]),
"Function 'IntracellularAnisotropicConductivity' not "
"correctly defined.");
- EvaluateFunction(aniso_var[i], vTemp_i,
- "IntracellularAnisotropicConductivity");
+ GetFunction("IntracellularAnisotropicConductivity")->Evaluate(aniso_var[i], vTemp_i);
Vmath::Vmul(nq, vTemp_i, 1, vTemp_j, 1,
m_vardiffi[varCoeffEnum[k]], 1);
diff --git a/solvers/CardiacEPSolver/EquationSystems/Monodomain.cpp b/solvers/CardiacEPSolver/EquationSystems/Monodomain.cpp
index 58dd3ccb210d2c0348dd69bf7143fdb316b41ece..2eec9c40f9e2eb6e06d77b24581d0caad32417b2 100644
--- a/solvers/CardiacEPSolver/EquationSystems/Monodomain.cpp
+++ b/solvers/CardiacEPSolver/EquationSystems/Monodomain.cpp
@@ -168,8 +168,7 @@ namespace Nektar
aniso_var[j]),
"Function 'AnisotropicConductivity' not correctly "
"defined.");
- EvaluateFunction(aniso_var[j], vTemp_j,
- "AnisotropicConductivity");
+ GetFunction("AnisotropicConductivity")->Evaluate(aniso_var[j], vTemp_j);
// Loop through rows of D
for (int i = 0; i < j + 1; ++i)
@@ -178,8 +177,7 @@ namespace Nektar
"AnisotropicConductivity",aniso_var[i]),
"Function 'AnisotropicConductivity' not correctly "
"defined.");
- EvaluateFunction(aniso_var[i], vTemp_i,
- "AnisotropicConductivity");
+ GetFunction("AnisotropicConductivity")->Evaluate(aniso_var[i], vTemp_i);
Vmath::Vmul(nq, vTemp_i, 1, vTemp_j, 1,
m_vardiff[varCoeffEnum[k]], 1);
@@ -223,7 +221,7 @@ namespace Nektar
const std::string varName = "intensity";
Array<OneD, NekDouble> vTemp;
- EvaluateFunction(varName, vTemp, "IsotropicConductivity");
+ GetFunction( "IsotropicConductivity")->Evaluate(varName, vTemp);
// If the d_min and d_max parameters are defined, then we need to
// rescale the isotropic conductivity to convert from the source
diff --git a/solvers/CompressibleFlowSolver/Forcing/ForcingQuasi1D.cpp b/solvers/CompressibleFlowSolver/Forcing/ForcingQuasi1D.cpp
index 6f8e6db25833f8a021a288de3dd4e0f18306f706..92794660ad05da1e0a90c9d3aa7c75fdad7bc29b 100644
--- a/solvers/CompressibleFlowSolver/Forcing/ForcingQuasi1D.cpp
+++ b/solvers/CompressibleFlowSolver/Forcing/ForcingQuasi1D.cpp
@@ -80,8 +80,8 @@ void ForcingQuasi1D::v_InitObject(
std::string sFieldStr = m_session->GetVariable(0);
ASSERTL0(m_session->DefinesFunction(funcName, sFieldStr),
"Variable '" + sFieldStr + "' not defined.");
- EvaluateFunction(pFields, m_session, sFieldStr,
- m_geomFactor, funcName, 0.0);
+ GetFunction(pFields, m_session, funcName, true)
+ ->Evaluate(sFieldStr, m_geomFactor, 0.0);
pFields[0]->PhysDeriv(MultiRegions::DirCartesianMap[0], m_geomFactor, tmp);
Vmath::Vdiv(pFields[0]->GetTotPoints(), tmp, 1,
m_geomFactor, 1,
diff --git a/solvers/IncNavierStokesSolver/EquationSystems/CoupledLinearNS.cpp b/solvers/IncNavierStokesSolver/EquationSystems/CoupledLinearNS.cpp
index 50c241722f3a4c255eb843df98b8952e10d91d1c..0ac98a87eb0a9ff629ae6c1ff898dd5e9282ba9f 100644
--- a/solvers/IncNavierStokesSolver/EquationSystems/CoupledLinearNS.cpp
+++ b/solvers/IncNavierStokesSolver/EquationSystems/CoupledLinearNS.cpp
@@ -1258,7 +1258,7 @@ namespace Nektar
{
fieldStr.push_back(m_boundaryConditions->GetVariable(m_velocity[i]));
}
- EvaluateFunction(fieldStr,AdvField,"AdvectionVelocity");
+ GetFunction("AdvectionVelocity")->Evaluate(fieldStr, AdvField);
SetUpCoupledMatrix(0.0,AdvField,false);
}
@@ -1290,7 +1290,7 @@ namespace Nektar
{
fieldStr.push_back(m_boundaryConditions->GetVariable(m_velocity[i]));
}
- EvaluateFunction(fieldStr, Restart, "Restart");
+ GetFunction( "Restart")->Evaluate(fieldStr, Restart);
for(int i = 0; i < m_velocity.num_elements(); ++i)
{
@@ -1338,7 +1338,7 @@ namespace Nektar
{
fieldStr.push_back(m_boundaryConditions->GetVariable(m_velocity[i]));
}
- EvaluateFunction(fieldStr,AdvField,"AdvectionVelocity");
+ GetFunction("AdvectionVelocity")->Evaluate(fieldStr, AdvField);
SetUpCoupledMatrix(m_lambda,AdvField,true);
}
@@ -1548,7 +1548,7 @@ namespace Nektar
{
fieldStr.push_back(m_boundaryConditions->GetVariable(m_velocity[i]));
}
- EvaluateFunction(fieldStr, m_ForcingTerm, "ForcingTerm");
+ GetFunction( "ForcingTerm")->Evaluate(fieldStr, m_ForcingTerm);
for(int i = 0; i < m_velocity.num_elements(); ++i)
{
m_fields[m_velocity[i]]->FwdTrans_IterPerExp(m_ForcingTerm[i], m_ForcingTerm_Coeffs[i]);
diff --git a/solvers/IncNavierStokesSolver/Forcing/ForcingMovingBody.cpp b/solvers/IncNavierStokesSolver/Forcing/ForcingMovingBody.cpp
index ac734c19edc8500f139b31b16a1094c7bb3dcf4b..d8ab9ed5b114c654da8cfb7c1e96307b0800ad18 100644
--- a/solvers/IncNavierStokesSolver/Forcing/ForcingMovingBody.cpp
+++ b/solvers/IncNavierStokesSolver/Forcing/ForcingMovingBody.cpp
@@ -155,10 +155,8 @@ void ForcingMovingBody::v_Apply(
}
else
{
- EvaluateFunction(pFields, m_session, m_motion[0], m_zta[j],
- m_funcName[j], time);
- EvaluateFunction(pFields, m_session, m_motion[1], m_eta[j],
- m_funcName[j], time);
+ GetFunction(pFields, m_session, m_funcName[j], true)->Evaluate(m_motion[0], m_zta[j], time);
+ GetFunction(pFields, m_session, m_funcName[j], true)->Evaluate(m_motion[1], m_eta[j], time);
cnt = cnt + 2;
}
}
diff --git a/solvers/LinearElasticSolver/EquationSystems/LinearElasticSystem.cpp b/solvers/LinearElasticSolver/EquationSystems/LinearElasticSystem.cpp
index 66001cbba9e467e58c2286c192fec81c584edfac..042ebe4211583eca5fd41218339ac9e5745d614a 100644
--- a/solvers/LinearElasticSolver/EquationSystems/LinearElasticSystem.cpp
+++ b/solvers/LinearElasticSolver/EquationSystems/LinearElasticSystem.cpp
@@ -453,7 +453,7 @@ void LinearElasticSystem::v_DoSolve()
// Evaluate the forcing function from the XML file.
Array<OneD, Array<OneD, NekDouble> > forcing(nVel);
- EvaluateFunction(forcing, "Forcing");
+ GetFunction("Forcing")->Evaluate(forcing);
// Add temperature term
string tempEval;
diff --git a/solvers/PulseWaveSolver/EquationSystems/PulseWaveSystem.cpp b/solvers/PulseWaveSolver/EquationSystems/PulseWaveSystem.cpp
index b641aba351b265bc4e746353cc0825e49e6479ec..86acd3bb134292b1c45e4095070e4ec325244131 100755
--- a/solvers/PulseWaveSolver/EquationSystems/PulseWaveSystem.cpp
+++ b/solvers/PulseWaveSolver/EquationSystems/PulseWaveSystem.cpp
@@ -220,10 +220,10 @@ namespace Nektar
m_fields[0] = m_vessels[2*omega];
m_beta[omega] = Array<OneD, NekDouble>(nq);
- EvaluateFunction("beta", m_beta[omega], "MaterialProperties", m_time, omega);
+ GetFunction("MaterialProperties")->Evaluate("beta", m_beta[omega], m_time, omega);
m_A_0[omega] = Array<OneD, NekDouble>(nq);
- EvaluateFunction("A_0", m_A_0[omega], "A_0", m_time, omega);
+ GetFunction("A_0")->Evaluate("A_0", m_A_0[omega], m_time, omega);
int nqTrace = GetTraceTotPoints();
@@ -1187,8 +1187,7 @@ namespace Nektar
LibUtilities::EquationSharedPtr vEqu
= m_session->GetFunction("ExactSolution",field,omega);
- EvaluateFunction(m_session->GetVariable(field),exactsoln,"ExactSolution",
- m_time);
+ GetFunction("ExactSolution")->Evaluate(m_session->GetVariable(field), exactsoln, m_time);
L2error_dom = m_vessels[vesselid]->L2(
m_vessels[vesselid]->GetPhys(),
@@ -1266,8 +1265,7 @@ namespace Nektar
{
Array<OneD, NekDouble> exactsoln(m_vessels[vesselid]->GetNpoints());
- EvaluateFunction(m_session->GetVariable(field),exactsoln,"ExactSolution",
- m_time);
+ GetFunction("ExactSolution")->Evaluate(m_session->GetVariable(field), exactsoln, m_time);
LinferrorDom = m_vessels[vesselid]->Linf(
m_vessels[vesselid]->GetPhys(),
diff --git a/solvers/ShallowWaterSolver/EquationSystems/ShallowWaterSystem.cpp b/solvers/ShallowWaterSolver/EquationSystems/ShallowWaterSystem.cpp
index d412e8f1f0f5d6176b8307657a275f0abd23b74b..a4d8af8caf66647b616e0c4a6aeceb38a6d73c80 100644
--- a/solvers/ShallowWaterSolver/EquationSystems/ShallowWaterSystem.cpp
+++ b/solvers/ShallowWaterSolver/EquationSystems/ShallowWaterSystem.cpp
@@ -166,13 +166,13 @@ namespace Nektar
void ShallowWaterSystem::EvaluateWaterDepth(void)
{
- EvaluateFunction("d",m_depth,"WaterDepth");
+ GetFunction("WaterDepth")->Evaluate("d", m_depth);
}
void ShallowWaterSystem::EvaluateCoriolis(void)
{
- EvaluateFunction("f",m_coriolis,"Coriolis");
+ GetFunction("Coriolis")->Evaluate("f", m_coriolis);
}
void ShallowWaterSystem::CopyBoundaryTrace(const Array<OneD, NekDouble>&Fwd, Array<OneD, NekDouble>&Bwd)
diff --git a/solvers/VortexWaveInteraction/VortexWaveInteraction.cpp b/solvers/VortexWaveInteraction/VortexWaveInteraction.cpp
index 010074f590a93cd572aa179a6c8727c57748b55c..d72e783bb56f293ca6dd133191768c078deca878 100644
--- a/solvers/VortexWaveInteraction/VortexWaveInteraction.cpp
+++ b/solvers/VortexWaveInteraction/VortexWaveInteraction.cpp
@@ -393,7 +393,7 @@ namespace Nektar
std::vector<std::string> vFieldNames = m_sessionRoll->GetVariables();
vFieldNames.erase(vFieldNames.end()-1);
- m_solverRoll->EvaluateFunction(vFieldNames, m_vwiForcingObj->UpdateForces(), "BodyForce");
+ m_solverRoll->GetFunction("BodyForce")->Evaluate(vFieldNames, m_vwiForcingObj->UpdateForces());
// Scale forcing
for(int i = 0; i < m_vwiForcingObj->UpdateForces().num_elements(); ++i)