- added rounding of stepSize to hit the numberOfInterval point.

 - added cruntime flag(-noEventEmit) to suppress output at events.


git-svn-id: https://openmodelica.org/svn/OpenModelica/trunk@18102 f25d12d1-65f4-0310-ae8a-bbce733d8d8e

SimulationRuntime/c/Makefile.common

@@ -31,6 +31,7 @@ RUNTIME_HEADERS = $(LIBF2CHEADER) \
 ./openmodelica.h \
 ./simulation_data.h \
 ./simulation/modelinfo.h \
+./simulation/options.h \
 ./simulation/simulation_info_xml.h \
 ./simulation/simulation_runtime.h \
 ./simulation/solver/delay.h \
@@ -57,6 +58,7 @@ RUNTIME_HEADERS = $(LIBF2CHEADER) \
 ./util/real_array.h \
 ./util/ringbuffer.h \
 ./util/rtclock.h \
+./util/simulation_options.h \
 ./util/string_array.h \
 ./util/utility.h \
 ./util/varinfo.h \

SimulationRuntime/c/simulation/solver/dassl.c

@@ -270,7 +270,16 @@ int dasrt_step(DATA* simData, SOLVER_INFO* solverInfo)
 
   /* Calculate time steps until TOUT is reached
    * (DASSL calculates beyond TOUT unless info[6] is set to 1!) */
-  tout = solverInfo->currentTime + solverInfo->currentStepSize;
+
+  /* adjust next time step by rounding. No rounding, 
+   * when time event is activated, since then currectStepSize is exact */
+  if (simData->simulationInfo.sampleActivated){
+    tout = solverInfo->currentTime + solverInfo->currentStepSize;
+  } else {
+    tout = _omc_round(solverInfo->currentTime + solverInfo->currentStepSize, solverInfo->stepPrecision);
+    solverInfo->currentStepSize = tout - solverInfo->currentTime;
+  }
+
   /* Check that tout is not less than timeValue
    * else will dassl get in trouble. If that is the case we skip the current step. */
   if(solverInfo->currentTime - tout >= -1e-13)

SimulationRuntime/c/simulation/solver/events.c

@@ -30,6 +30,7 @@
 
 #include "events.h"
 #include "omc_error.h"
+#include "options.h"
 #include "simulation_data.h"
 #include "simulation_result.h"
 #include "openmodelica.h"         /* for modelica types */
@@ -192,7 +193,9 @@ void handleEvents(DATA* data, LIST* eventLst, double *eventTime, SOLVER_INFO* so
   long i;
   LIST_NODE* it;
 
-  sim_result.emit(&sim_result,data);
+  /* prevent emit if noeventemit flag is used */
+  if (!(omc_flag[FLAG_NOEVENTEMIT]))
+    sim_result.emit(&sim_result,data);
 
   /* time event */
   if(data->simulationInfo.sampleActivated)

SimulationRuntime/c/simulation/solver/model_help.c

@@ -32,6 +32,7 @@
 #include <stdlib.h>
 #include <string.h>
 #include <float.h>
+#include <math.h>
 
 #include "simulation_data.h"
 #include "openmodelica_func.h"
@@ -1177,5 +1178,18 @@ modelica_real _event_div_real(modelica_real x1, modelica_real x2, modelica_integ
 }
 
 
-
+/*! \fn _omc_round
+ *
+ *  \param [in]  [x]
+ *  \param [in]  [precision]
+ *
+ * Returns a round value of x with a certain presision.
+ */
+modelica_real _omc_round(modelica_real x, modelica_integer precision)
+{
+  ASSERT(precision<20, "_omc_round fails with presision > 19");
+  double p = (double) pow(10, precision);
+  double trunc = round(x * p);
+  return (modelica_real) (trunc / p);
+}
 

SimulationRuntime/c/simulation/solver/model_help.h

@@ -138,6 +138,7 @@ modelica_real _event_div_real(modelica_real x1, modelica_real x2, modelica_integ
 modelica_integer _event_mod_integer(modelica_integer x1, modelica_integer x2, modelica_integer index, DATA *data);
 modelica_integer _event_rem_integer(modelica_integer x1, modelica_integer x2, modelica_integer index, DATA *data);
 
+modelica_real _omc_round(modelica_real x, modelica_integer prec);
 /* functions used for relation which
  * are not used as zero-crossings
  */

SimulationRuntime/c/simulation/solver/perform_simulation.c

@@ -40,6 +40,7 @@
 #include "openmodelica_func.h"
 
 #include "omc_error.h"
+#include "options.h"
 #include <math.h>
 #include <string.h>
 #include <errno.h>
@@ -274,7 +275,13 @@ int performSimulation(DATA* data, SOLVER_INFO* solverInfo)
           fmt = NULL;
         }
       }
-      sim_result.emit(&sim_result,data);
+      /* prevent emit if noeventemit flag is used, if it's an event */
+      if ((omc_flag[FLAG_NOEVENTEMIT])){
+        if (solverInfo->didEventStep == 0)
+          sim_result.emit(&sim_result,data);
+      } else {
+        sim_result.emit(&sim_result,data);
+      }
 
       printAllVarsDebug(data, 0, LOG_DEBUG);  /* ??? */
 

SimulationRuntime/c/simulation/solver/solver_main.c

@@ -53,6 +53,7 @@
  */
 #include "rtclock.h"
 #include "omc_error.h"
+#include "options.h"
 #include <math.h>
 #include <string.h>
 #include <errno.h>
@@ -146,6 +147,7 @@ int initializeSolverData(DATA* data, SOLVER_INFO* solverInfo)
   solverInfo->didEventStep = 0;
   solverInfo->stateEvents = 0;
   solverInfo->sampleEvents = 0;
+  solverInfo->stepPrecision = 15;
 
   if(solverInfo->solverMethod == 2)
   {
@@ -421,7 +423,11 @@ int finishSimulation(DATA* data, SOLVER_INFO* solverInfo, const char* outputVari
   {
     data->simulationInfo.terminal = 1;
     updateDiscreteSystem(data);
-    sim_result.emit(&sim_result,data);
+
+    /* prevent emit if noeventemit flag is used */
+    if (!(omc_flag[FLAG_NOEVENTEMIT]))
+      sim_result.emit(&sim_result,data);
+
     data->simulationInfo.terminal = 0;
   }
   communicateStatus("Finished", 1);
@@ -597,12 +603,20 @@ static int euler_ex_step(DATA* data, SOLVER_INFO* solverInfo)
   SIMULATION_DATA *sDataOld = (SIMULATION_DATA*)data->localData[1];
   modelica_real* stateDer = sDataOld->realVars + data->modelData.nStates;
 
+   /* adjust next time step by rounding. No rounding, 
+   * when time event is activated, since then currectStepSize is exact */
+  if (data->simulationInfo.sampleActivated){
+    solverInfo->currentTime = sDataOld->timeValue + solverInfo->currentStepSize;
+  } else {
+    solverInfo->currentTime = _omc_round(sDataOld->timeValue + solverInfo->currentStepSize, solverInfo->stepPrecision);
+    solverInfo->currentStepSize = solverInfo->currentTime - sDataOld->timeValue;
+  }
+
   for(i = 0; i < data->modelData.nStates; i++)
   {
     sData->realVars[i] = sDataOld->realVars[i] + stateDer[i] * solverInfo->currentStepSize;
   }
-  sData->timeValue = sDataOld->timeValue + solverInfo->currentStepSize;
-  solverInfo->currentTime += solverInfo->currentStepSize;
+  sData->timeValue = solverInfo->currentTime;
   return 0;
 }
 
@@ -617,6 +631,14 @@ static int rungekutta_step(DATA* data, SOLVER_INFO* solverInfo)
   modelica_real* stateDer = sData->realVars + data->modelData.nStates;
   modelica_real* stateDerOld = sDataOld->realVars + data->modelData.nStates;
 
+  /* adjust next time step by rounding. No rounding, 
+   * when time event is activated, since then currectStepSize is exact */
+  if (data->simulationInfo.sampleActivated){
+    solverInfo->currentTime = sDataOld->timeValue + solverInfo->currentStepSize;
+  } else {
+    solverInfo->currentTime = _omc_round(sDataOld->timeValue + solverInfo->currentStepSize, solverInfo->stepPrecision);
+    solverInfo->currentStepSize = solverInfo->currentTime - sDataOld->timeValue;
+  }
 
   /* We calculate k[0] before returning from this function.
    * We only want to calculate f() 4 times per call */
@@ -648,8 +670,7 @@ static int rungekutta_step(DATA* data, SOLVER_INFO* solverInfo)
     }
     sData->realVars[i] = sDataOld->realVars[i] + solverInfo->currentStepSize * sum;
   }
-  sData->timeValue = sDataOld->timeValue + solverInfo->currentStepSize;
-  solverInfo->currentTime += solverInfo->currentStepSize;
+  sData->timeValue = solverInfo->currentTime;
   return 0;
 }
 

SimulationRuntime/c/simulation/solver/solver_main.h

@@ -49,6 +49,7 @@ typedef struct SOLVER_INFO
   double laststep;
   double offset;
   int solverMethod;
+  int stepPrecision;
 
   /* set by solver if an internal root finding method is activated  */
   modelica_boolean solverRootFinding;

SimulationRuntime/c/util/simulation_options.c

@@ -50,6 +50,7 @@ const char *FLAG_NAME[FLAG_MAX+1] = {
   /* FLAG_MEASURETIMEPLOTFORMAT */ "measureTimePlotFormat",
   /* FLAG_NLS */                   "nls",
   /* FLAG_NOEMIT */                "noemit",
+  /* FLAG_NOEVENTEMIT */           "noEventEmit",
   /* FLAG_OUTPUT */                "output",
   /* FLAG_OVERRIDE */              "override",
   /* FLAG_OVERRIDE_FILE */         "overrideFile",
@@ -81,6 +82,7 @@ const char *FLAG_DESC[FLAG_MAX+1] = {
   /* FLAG_MEASURETIMEPLOTFORMAT */ "value specifies the output format of the measure time functionality",
   /* FLAG_NLS */                   "value specifies the nonlinear solver",
   /* FLAG_NOEMIT */                "do not emit any results to the result file",
+  /* FLAG_NOEVENTEMIT */           "do not emit event points to the result file",
   /* FLAG_OUTPUT */                "output the variables a, b and c at the end of the simulation to the standard output",
   /* FLAG_OVERRIDE */              "override the variables or the simulation settings in the XML setup file",
   /* FLAG_OVERRIDE_FILE */         "will override the variables or the simulation settings in the XML setup file with the values from the file",
@@ -112,6 +114,7 @@ const char *FLAG_DETAILED_DESC[FLAG_MAX+1] = {
   /* FLAG_MEASURETIMEPLOTFORMAT */ "value specifies the output format of the measure time functionality\n  svg\n  jpg\n  ps\n  gif\n  ...",
   /* FLAG_NLS */                   "value specifies the nonlinear solver",
   /* FLAG_NOEMIT */                "do not emit any results to the result file",
+  /* FLAG_NOEVENTEMIT */           "do not emit event points to the result file",
   /* FLAG_OUTPUT */                "output the variables a, b and c at the end of the simulation to the standard output\n  time = value, a = value, b = value, c = value",
   /* FLAG_OVERRIDE */              "override the variables or the simulation settings in the XML setup file\n  e.g. var1=start1,var2=start2,par3=start3,startTime=val1,stopTime=val2,stepSize=val3,\n       tolerance=val4,solver=\"see -s\",outputFormat=\"mat|plt|csv|empty\",variableFilter=\"filter\"",
   /* FLAG_OVERRIDE_FILE */         "will override the variables or the simulation settings in the XML setup file with the values from the file\n  note that: -overrideFile CANNOT be used with -override\n  use when variables for -override are too many and do not fit in command line size\n  overrideFileName contains lines of the form: var1=start1",
@@ -143,6 +146,7 @@ const int FLAG_TYPE[FLAG_MAX] = {
   /* FLAG_MEASURETIMEPLOTFORMAT */ FLAG_TYPE_OPTION,
   /* FLAG_NLS */                   FLAG_TYPE_OPTION,
   /* FLAG_NOEMIT */                FLAG_TYPE_FLAG,
+  /* FLAG_NOEVENTEMIT */           FLAG_TYPE_FLAG,
   /* FLAG_OUTPUT */                FLAG_TYPE_OPTION,
   /* FLAG_OVERRIDE */              FLAG_TYPE_OPTION,
   /* FLAG_OVERRIDE_FILE */         FLAG_TYPE_OPTION,

SimulationRuntime/c/util/simulation_options.h

@@ -57,6 +57,7 @@ enum _FLAG
   FLAG_MEASURETIMEPLOTFORMAT,
   FLAG_NLS,
   FLAG_NOEMIT,
+  FLAG_NOEVENTEMIT,
   FLAG_OUTPUT,
   FLAG_OVERRIDE,
   FLAG_OVERRIDE_FILE,