add further deltaX option for the solvers

Willi Braun · OpenModelica-Hudson · commit 024573ebc634 · 2016-10-27T18:28:42.000+02:00
diff --git a/SimulationRuntime/c/linearization/linearize.cpp b/SimulationRuntime/c/linearization/linearize.cpp
@@ -112,7 +112,7 @@ int functionODE_residual(DATA* data, threadData_t *threadData, double *dx, doubl
 /*  Calculate the jacobian matrix by numerical finite difference */
 int functionJacAC_num(DATA* data, threadData_t *threadData, double *matrixA, double *matrixC, double *matrixCz)
 {
-    const double delta_h = numericalDifferentiationDeltaX;
+    const double delta_h = numericalDifferentiationDeltaXlinearize;
     double delta_hh;
     double xsave;
 
@@ -203,7 +203,7 @@ int functionJacAC_num(DATA* data, threadData_t *threadData, double *matrixA, dou
 
 int functionJacBD_num(DATA* data, threadData_t *threadData, double *matrixB, double *matrixD, double *matrixDz)
 {
-    const double delta_h = numericalDifferentiationDeltaX;
+    const double delta_h = numericalDifferentiationDeltaXlinearize;
     double delta_hh;
     double usave;
     double* u;
diff --git a/SimulationRuntime/c/simulation/simulation_runtime.cpp b/SimulationRuntime/c/simulation/simulation_runtime.cpp
@@ -508,14 +508,21 @@ int startNonInteractiveSimulation(int argc, char**argv, DATA* data, threadData_t
     infoStreamPrint(LOG_STDOUT, 0, "Linearization will performed at point of time: %f", data->simulationInfo->stopTime);
   }
 
-  /* set delta x  */
-  if(omc_flag[FLAG_DELTA_X]) {
-    numericalDifferentiationDeltaX = atof(omc_flagValue[FLAG_DELTA_X]);
-    infoStreamPrint(LOG_SOLVER, 0, "Set delta x for numerical differentiation to %f", numericalDifferentiationDeltaX);
+  /* set delta x for linearization */
+  if(omc_flag[FLAG_DELTA_X_LINEARIZE]) {
+    numericalDifferentiationDeltaXlinearize = atof(omc_flagValue[FLAG_DELTA_X_LINEARIZE]);
+    infoStreamPrint(LOG_SOLVER, 0, "Set delta x for numerical differentiation of the linearization to %f", numericalDifferentiationDeltaXlinearize);
   }else{
-    numericalDifferentiationDeltaX = sqrt(DBL_EPSILON*2e1);
+    numericalDifferentiationDeltaXlinearize = sqrt(DBL_EPSILON*2e1);
   }
 
+  /* set delta x for integration methods dassl, ida */
+  if(omc_flag[FLAG_DELTA_X_SOLVER]) {
+    numericalDifferentiationDeltaXsolver = atof(omc_flagValue[FLAG_DELTA_X_SOLVER]);
+    infoStreamPrint(LOG_SOLVER, 0, "Set delta x for numerical differentiation of the integrator to %f", numericalDifferentiationDeltaXsolver);
+  }else{
+    numericalDifferentiationDeltaXsolver = sqrt(DBL_EPSILON);
+  }
 
   if(omc_flag[FLAG_S]) {
     if (omc_flagValue[FLAG_S]) {
diff --git a/SimulationRuntime/c/simulation/solver/dassl.c b/SimulationRuntime/c/simulation/solver/dassl.c
@@ -172,7 +172,6 @@ int dassl_initial(DATA* data, threadData_t *threadData, SOLVER_INFO* solverInfo,
 
   dasslData->idid = 0;
 
-  dasslData->sqrteps = sqrt(DBL_EPSILON);
   dasslData->ysave = (double*) malloc(N*sizeof(double));
   dasslData->ypsave = (double*) malloc(N*sizeof(double));
   dasslData->delta_hh = (double*) malloc(N*sizeof(double));
@@ -1137,7 +1136,7 @@ int jacA_num(DATA* data, double *t, double *y, double *yprime, double *delta, do
   TRACE_PUSH
   DASSL_DATA* dasslData = (DASSL_DATA*)(void*)((double**)rpar)[1];
 
-  double delta_h = dasslData->sqrteps;
+  double delta_h = numericalDifferentiationDeltaXsolver;
   double delta_hh,delta_hhh, deltaInv;
   double ysave;
   double ypsave;
@@ -1236,7 +1235,7 @@ int jacA_numColored(DATA* data, double *t, double *y, double *yprime, double *de
   TRACE_PUSH
   const int index = data->callback->INDEX_JAC_A;
   DASSL_DATA* dasslData = (DASSL_DATA*)(void*)((double**)rpar)[1];
-  double delta_h = dasslData->sqrteps;
+  double delta_h = numericalDifferentiationDeltaXsolver;
   double delta_hhh;
   int ires;
   double* delta_hh = dasslData->delta_hh;
diff --git a/SimulationRuntime/c/simulation/solver/dassl.h b/SimulationRuntime/c/simulation/solver/dassl.h
@@ -67,7 +67,6 @@ typedef struct DASSL_DATA{
   int *jroot;
 
   /* variables used in jacobian calculation */
-  double sqrteps;
   double *ysave;
   double *ypsave;
   double *delta_hh;
diff --git a/SimulationRuntime/c/simulation/solver/ida_solver.c b/SimulationRuntime/c/simulation/solver/ida_solver.c
@@ -187,7 +187,6 @@ ida_solver_initial(DATA* data, threadData_t *threadData, SOLVER_INFO* solverInfo
       idaData->yp);
 
   /* allocate memory for jacobians calculation */
-  idaData->sqrteps = sqrt(DBL_EPSILON);
   idaData->ysave = (double*) malloc(idaData->N*sizeof(double));
   idaData->ypsave = (double*) malloc(idaData->N*sizeof(double));
   idaData->delta_hh = (double*) malloc(idaData->N*sizeof(double));
@@ -1102,7 +1101,7 @@ int jacOwnNumColoredIDA(double tt, N_Vector yy, N_Vector yp, N_Vector rr, DlsMat
   double *ysave = idaData->ysave;
   double *ypsave = idaData->ypsave;
 
-  double delta_h = idaData->sqrteps;
+  double delta_h = numericalDifferentiationDeltaXsolver;
   double delta_hhh;
   long int i,j,l,ii;
 
@@ -1258,7 +1257,7 @@ int jacOwnNumIDA(double tt, N_Vector yy, N_Vector yp, N_Vector rr, DlsMat Jac, d
 
   double ysave, ypsave;
 
-  double delta_h = idaData->sqrteps;
+  double delta_h = numericalDifferentiationDeltaXsolver;
   double delta_hh;
   double delta_hhh;
   double deltaInv;
@@ -1384,7 +1383,7 @@ int jacobianSparseNumIDA(double tt, N_Vector yy, N_Vector yp, N_Vector rr, SlsMa
   double *ysave = idaData->ysave;
   double *ypsave = idaData->ypsave;
 
-  double delta_h = idaData->sqrteps;
+  double delta_h = numericalDifferentiationDeltaXsolver;
   double *delta_hh = idaData->delta_hh;
   double delta_hhh;
   double deltaInv;
diff --git a/SimulationRuntime/c/simulation/solver/model_help.c b/SimulationRuntime/c/simulation/solver/model_help.c
@@ -57,7 +57,8 @@ double newtonXTol = 1e-12;
 double newtonFTol = 1e-12;
 const size_t SIZERINGBUFFER = 3;
 int compiledInDAEMode = 0;
-double numericalDifferentiationDeltaX = 1e-8;
+double numericalDifferentiationDeltaXlinearize = 1e-8;
+double numericalDifferentiationDeltaXsolver = 1e-8;
 
 static double tolZC;
 
diff --git a/SimulationRuntime/c/simulation/solver/model_help.h b/SimulationRuntime/c/simulation/solver/model_help.h
@@ -81,7 +81,8 @@ extern double newtonXTol;
 extern double newtonFTol;
 extern const size_t SIZERINGBUFFER;
 extern int compiledInDAEMode;
-extern double numericalDifferentiationDeltaX;
+extern double numericalDifferentiationDeltaXlinearize;
+extern double numericalDifferentiationDeltaXsolver;
 
 void initializeDataStruc(DATA *data, threadData_t *threadData);
 
diff --git a/SimulationRuntime/c/util/simulation_options.c b/SimulationRuntime/c/util/simulation_options.c
@@ -39,7 +39,8 @@ const char *FLAG_NAME[FLAG_MAX+1] = {
   /* FLAG_CPU */                   "cpu",
   /* FLAG_CSV_OSTEP */             "csvOstep",
   /* FLAG_DAE_MODE */              "daeMode",
-  /* FLAG_DELTA_X */               "deltaX",
+  /* FLAG_DELTA_X_LINEARIZE */     "deltaXLinearize",
+  /* FLAG_DELTA_X_SOLVER */        "deltaXSolver",
   /* FLAG_EMBEDDED_SERVER */       "embeddedServer",
   /* FLAG_EMIT_PROTECTED */        "emit_protected",
   /* FLAG_F */                     "f",
@@ -119,7 +120,8 @@ const char *FLAG_DESC[FLAG_MAX+1] = {
   /* FLAG_CPU */                   "dumps the cpu-time into the result file",
   /* FLAG_CSV_OSTEP */             "value specifies csv-files for debuge values for optimizer step",
   /* FLAG_DAE_MODE */              "flag to let the integrator use daeResiduals",
-  /* FLAG_DELTA_X */               "value specifies the delta x value for numerical differentiation, used by linearization or integrator",
+  /* FLAG_DELTA_X_LINEARIZE */     "value specifies the delta x value for numerical differentiation used by linearization. The default value is 1e-5.",
+  /* FLAG_DELTA_X_SOLVER */        "value specifies the delta x value for numerical differentiation used by integrator. The default values is sqrt(DBL_EPSILON).",
   /* FLAG_EMBEDDED_SERVER */       "enables an embedded server. Valid values: none, opc-da [broken], opc-ua [experimental], or the path to a shared object.",
   /* FLAG_EMIT_PROTECTED */        "emits protected variables to the result-file",
   /* FLAG_F */                     "value specifies a new setup XML file to the generated simulation code",
@@ -207,8 +209,10 @@ const char *FLAG_DETAILED_DESC[FLAG_MAX+1] = {
   "  Value specifies csv-files for debuge values for optimizer step",
   /* FLAG_DAE_MODE */
   "  Enables daeMode simulation if the model was compiled with the omc flag --daeMode and the IDA integrator is used.",
-  /* FLAG_DELTA_X */
-  "value specifies the delta x value for numerical differentiation, used by linearization or integrator",
+  /* FLAG_DELTA_X_LINEARIZE */
+  "value specifies the delta x value for numerical differentiation used by linearization. The default value is sqrt(DBL_EPSILON*2e1).",
+  /* FLAG_DELTA_X_SOLVER */
+  "value specifies the delta x value for numerical differentiation used by integrator. The default values is sqrt(DBL_EPSILON).",
   /* FLAG_EMBEDDED_SERVER */
   "  Enables an embedded server. Valid values:\n\n"
   "  * none - default, run without embedded server\n"
@@ -407,7 +411,8 @@ const int FLAG_TYPE[FLAG_MAX] = {
   /* FLAG_CPU */                   FLAG_TYPE_FLAG,
   /* FLAG_CSV_OSTEP */             FLAG_TYPE_OPTION,
   /* FLAG_DAE_SOLVING */           FLAG_TYPE_FLAG,
-  /* FLAG_DELTA_X */               FLAG_TYPE_OPTION,
+  /* FLAG_DELTA_X_LINEARIZE */     FLAG_TYPE_OPTION,
+  /* FLAG_DELTA_X_SOLVER */        FLAG_TYPE_OPTION,
   /* FLAG_EMBEDDED_SERVER */       FLAG_TYPE_OPTION,
   /* FLAG_EMIT_PROTECTED */        FLAG_TYPE_FLAG,
   /* FLAG_F */                     FLAG_TYPE_OPTION,
diff --git a/SimulationRuntime/c/util/simulation_options.h b/SimulationRuntime/c/util/simulation_options.h
@@ -47,7 +47,8 @@ enum _FLAG
   FLAG_CPU,
   FLAG_CSV_OSTEP,
   FLAG_DAE_MODE,
-  FLAG_DELTA_X,
+  FLAG_DELTA_X_LINEARIZE,
+  FLAG_DELTA_X_SOLVER,
   FLAG_EMBEDDED_SERVER,
   FLAG_EMIT_PROTECTED,
   FLAG_F,
