From b1df5e3b65eeb1f4d80c6854482e6dddaa2d2078 Mon Sep 17 00:00:00 2001 From: Giacomo Longo Date: Tue, 5 Dec 2023 16:05:55 +0100 Subject: [PATCH] Allow to set CVODE parameters for ME FMUs --- src/OMSimulatorLib/Flags.cpp | 28 ++++++++++++++++++++++++++++ src/OMSimulatorLib/Flags.h | 16 ++++++++++++++++ src/OMSimulatorLib/SystemSC.cpp | 16 ++++++++-------- 3 files changed, 52 insertions(+), 8 deletions(-) diff --git a/src/OMSimulatorLib/Flags.cpp b/src/OMSimulatorLib/Flags.cpp index 984624cfa..e4b8ee7c5 100644 --- a/src/OMSimulatorLib/Flags.cpp +++ b/src/OMSimulatorLib/Flags.cpp @@ -65,6 +65,10 @@ void oms::Flags::setDefaults() { addParametersToCSV = false; algLoopSolver = oms_alg_solver_kinsol; + cvodeMaxErrTestFails = 100; + cvodeMaxNLSFails = 100; + cvodeMaxNLSIterations = 5; + cvodeMaxSteps = 1000; defaultModeIsCS = false; deleteTempFiles = true; directionalDerivatives = true; @@ -203,6 +207,30 @@ oms_status_enu_t oms::Flags::ClearAllOptions(const std::string& value) return oms_status_ok; } +oms_status_enu_t oms::Flags::CVODEMaxErrTestFails(const std::string& value) +{ + GetInstance().cvodeMaxErrTestFails = atoi(value.c_str()); + return oms_status_ok; +} + +oms_status_enu_t oms::Flags::CVODEMaxNLSFailures(const std::string& value) +{ + GetInstance().cvodeMaxNLSFails = atoi(value.c_str()); + return oms_status_ok; +} + +oms_status_enu_t oms::Flags::CVODEMaxNLSIterations(const std::string& value) +{ + GetInstance().cvodeMaxNLSIterations = atoi(value.c_str()); + return oms_status_ok; +} + +oms_status_enu_t oms::Flags::CVODEMaxSteps(const std::string& value) +{ + GetInstance().cvodeMaxSteps = atoi(value.c_str()); + return oms_status_ok; +} + oms_status_enu_t oms::Flags::DeleteTempFiles(const std::string& value) { GetInstance().deleteTempFiles = (value == "true"); diff --git a/src/OMSimulatorLib/Flags.h b/src/OMSimulatorLib/Flags.h index 9a643f965..0e62690d9 100644 --- a/src/OMSimulatorLib/Flags.h +++ b/src/OMSimulatorLib/Flags.h @@ -56,6 +56,10 @@ namespace oms static oms_status_enu_t SetCommandLineOption(const std::string& cmd); static bool AddParametersToCSV() {return GetInstance().addParametersToCSV;} + static int CVODEMaxErrTestFails() {return GetInstance().cvodeMaxErrTestFails;} + static int CVODEMaxNLSFailures() {return GetInstance().cvodeMaxNLSFails;} + static int CVODEMaxNLSIterations() {return GetInstance().cvodeMaxNLSIterations;} + static int CVODEMaxSteps() {return GetInstance().cvodeMaxSteps;} static bool DefaultModeIsCS() {return GetInstance().defaultModeIsCS;} static bool DeleteTempFiles() {return GetInstance().deleteTempFiles;} static bool DirectionalDerivatives() {return GetInstance().directionalDerivatives;} @@ -89,6 +93,10 @@ namespace oms private: bool addParametersToCSV; + int cvodeMaxErrTestFails; + int cvodeMaxNLSFails; + int cvodeMaxNLSIterations; + int cvodeMaxSteps; bool defaultModeIsCS; bool deleteTempFiles; bool directionalDerivatives; @@ -148,6 +156,10 @@ namespace oms {"--addParametersToCSV", "", "Export parameters to .csv file (true, [false])", re_default, Flags::AddParametersToCSV, false}, {"--algLoopSolver", "", "Specifies the alg. loop solver method (fixedpoint, [kinsol]) used for algebraic loops spanning over multiple components.", re_default, Flags::AlgLoopSolver, false}, {"--clearAllOptions", "", "Reset all flags to default values", re_void, Flags::ClearAllOptions, false}, + {"--CVODEMaxErrTestFails", "", "Maximum number of error test failures for CVODE", re_number, Flags::CVODEMaxErrTestFails, false}, + {"--CVODEMaxNLSFailures", "", "Maximum number of nonlinear convergence failures for CVODE", re_number, Flags::CVODEMaxNLSFailures, false}, + {"--CVODEMaxNLSIterations", "", "Maximum number of nonlinear solver iterations for CVODE", re_number, Flags::CVODEMaxNLSIterations, false}, + {"--CVODEMaxSteps", "", "Maximum number of steps for CVODE", re_number, Flags::CVODEMaxSteps, false}, {"--deleteTempFiles", "", "Deletes temp files as soon as they are no longer needed ([true], false)", re_bool, Flags::DeleteTempFiles, false}, {"--directionalDerivatives", "", "Specifies whether directional derivatives should be used to calculate the Jacobian for alg. loops or if a numerical approximation should be used instead ([true], false)", re_bool, Flags::DirectionalDerivatives, false}, {"--dumpAlgLoops", "", "Dump information for alg loops (true, [false])", re_bool, Flags::DumpAlgLoops, false}, @@ -186,6 +198,10 @@ namespace oms static oms_status_enu_t AddParametersToCSV(const std::string& value); static oms_status_enu_t AlgLoopSolver(const std::string& value); static oms_status_enu_t ClearAllOptions(const std::string& value); + static oms_status_enu_t CVODEMaxErrTestFails(const std::string& value); + static oms_status_enu_t CVODEMaxNLSFailures(const std::string& value); + static oms_status_enu_t CVODEMaxNLSIterations(const std::string& value); + static oms_status_enu_t CVODEMaxSteps(const std::string& value); static oms_status_enu_t DeleteTempFiles(const std::string& value); static oms_status_enu_t DirectionalDerivatives(const std::string& value); static oms_status_enu_t DumpAlgLoops(const std::string& value); diff --git a/src/OMSimulatorLib/SystemSC.cpp b/src/OMSimulatorLib/SystemSC.cpp index aaa260973..0b5f11226 100644 --- a/src/OMSimulatorLib/SystemSC.cpp +++ b/src/OMSimulatorLib/SystemSC.cpp @@ -348,21 +348,21 @@ oms_status_enu_t oms::SystemSC::initialize() if (flag < 0) logError("SUNDIALS_ERROR: CVodeSetMaxStep() failed with flag = " + std::to_string(flag)); // further settings from cpp runtime - flag = CVodeSetInitStep(solverData.cvode.mem, initialStepSize); // INITIAL STEPSIZE + flag = CVodeSetInitStep(solverData.cvode.mem, initialStepSize); // INITIAL STEPSIZE if (flag < 0) logError("SUNDIALS_ERROR: CVodeSetInitStep() failed with flag = " + std::to_string(flag)); - flag = CVodeSetMaxOrd(solverData.cvode.mem, 5); // MAXIMUM ORDER + flag = CVodeSetMaxOrd(solverData.cvode.mem, 5); // MAXIMUM ORDER if (flag < 0) logError("SUNDIALS_ERROR: CVodeSetMaxOrd() failed with flag = " + std::to_string(flag)); - flag = CVodeSetMaxConvFails(solverData.cvode.mem, 100); // MAXIMUM NUMBER OF NONLINEAR CONVERGENCE FAILURES + flag = CVodeSetMaxConvFails(solverData.cvode.mem, Flags::CVODEMaxNLSFailures()); // MAXIMUM NUMBER OF NONLINEAR CONVERGENCE FAILURES if (flag < 0) logError("SUNDIALS_ERROR: CVodeSetMaxConvFails() failed with flag = " + std::to_string(flag)); - flag = CVodeSetStabLimDet(solverData.cvode.mem, true); // STABILITY DETECTION + flag = CVodeSetStabLimDet(solverData.cvode.mem, true); // STABILITY DETECTION if (flag < 0) logError("SUNDIALS_ERROR: CVodeSetStabLimDet() failed with flag = " + std::to_string(flag)); - flag = CVodeSetMinStep(solverData.cvode.mem, minimumStepSize); // MINIMUM STEPSIZE + flag = CVodeSetMinStep(solverData.cvode.mem, minimumStepSize); // MINIMUM STEPSIZE if (flag < 0) logError("SUNDIALS_ERROR: CVodeSetMinStep() failed with flag = " + std::to_string(flag)); - flag = CVodeSetMaxNonlinIters(solverData.cvode.mem, 5); // MAXIMUM NUMBER OF ITERATIONS + flag = CVodeSetMaxNonlinIters(solverData.cvode.mem, Flags::CVODEMaxNLSIterations()); // MAXIMUM NUMBER OF ITERATIONS if (flag < 0) logError("SUNDIALS_ERROR: CVodeSetMaxNonlinIters() failed with flag = " + std::to_string(flag)); - flag = CVodeSetMaxErrTestFails(solverData.cvode.mem, 100); // MAXIMUM NUMBER OF ERROR TEST FAILURES + flag = CVodeSetMaxErrTestFails(solverData.cvode.mem, Flags::CVODEMaxErrTestFails()); // MAXIMUM NUMBER OF ERROR TEST FAILURES if (flag < 0) logError("SUNDIALS_ERROR: CVodeSetMaxErrTestFails() failed with flag = " + std::to_string(flag)); - flag = CVodeSetMaxNumSteps(solverData.cvode.mem, 1000); // MAXIMUM NUMBER OF STEPS + flag = CVodeSetMaxNumSteps(solverData.cvode.mem, Flags::CVODEMaxSteps()); // MAXIMUM NUMBER OF STEPS if (flag < 0) logError("SUNDIALS_ERROR: CVodeSetMaxNumSteps() failed with flag = " + std::to_string(flag)); }