- bugfix mixed equation systems

- move Modelica.Mechanics.MultiBody.Examples.Rotational3DEffects.GearConstraint.mos Modelica.Electrical.Analog.Examples.ShowVariableResistor.mos Modelica.Electrical.Machines.Examples.DCSE_Start.mos Modelica.Electrical.Analog.Examples.HeatingResistor.mos Modelica.Electrical.Analog.Examples.CharacteristicIdealDiodes.mos to working tests git-svn-id: https://openmodelica.org/svn/OpenModelica/trunk@12556 f25d12d1-65f4-0310-ae8a-bbce733d8d8e
OpenModelica · Aug 17, 2012 · 49600db · 49600db
1 parent 259a8bd
commit 49600db
Show file tree

Hide file tree

Showing 4 changed files with 42 additions and 150 deletions.
diff --git a/Compiler/BackEnd/BackendDAETransform.mo b/Compiler/BackEnd/BackendDAETransform.mo
@@ -1023,7 +1023,7 @@ protected function analyseStrongComponentBlock"function: analyseStrongComponentB
   input BackendDAE.EqSystem isyst;
   input BackendDAE.Shared ishared;  
   input array<Integer> inAss1;
-  input array<Integer> inAss2; 
+  input array<Integer> inAss2;
   input Boolean inLoop; //true if the function call itself
   output BackendDAE.StrongComponent outComp;
 algorithm
@@ -4263,7 +4263,7 @@ algorithm
 
         eqn_1 = Derive.differentiateEquationTime(eqn, v, shared);
         // print( "differentiated equation " +& intString(e) +& " " +& BackendDump.equationStr(eqn) +& "\n");
-        // print( "differentiated equation " +& intString(e) +& " " +& BackendDump.equationStr(eqn) +& "\n to \n");
+        // print( "differentiated equation " +& intString(e) +& " " +& BackendDump.equationStr(eqn_1) +& "\n to \n");
         (eqn_1,so) = traverseBackendDAEExpsEqn(eqn_1, replaceStateOrderExp,inStateOrd);
         // print(BackendDump.equationStr(eqn_1) +& "\n");
         Debug.fcall(Flags.BLT_DUMP, debugdifferentiateEqns,(eqn,eqn_1));

diff --git a/Compiler/BackEnd/SimCode.mo b/Compiler/BackEnd/SimCode.mo
@@ -5138,19 +5138,15 @@ algorithm
   (valuesRet, value_dims) :=
   match (inBackendDAEEquationLst1,inBackendDAEVarLst2,inBackendDAEEquationLst3,inBackendDAEVarLst4)
     local
-      list<DAE.Exp> rels;
       list<list<Integer>> values,values_1;
-      list<Integer> values_2;
       list<BackendDAE.Equation> cont_e,disc_e;
       list<BackendDAE.Var> cont_v,disc_v;
 
     case (cont_e,cont_v,disc_e,disc_v)
       equation
-        rels = mixedCollectRelations(cont_e, disc_e);
-        (values,value_dims) = generateMixedDiscretePossibleValues2(rels, disc_v, 0);
+        (values,value_dims) = generateMixedDiscretePossibleValues2(disc_v, 0);
         values_1 = generateMixedDiscreteCombinationValues(values);
-        values_2 = List.flatten(values_1);
-        valuesRet = values_2;
+        valuesRet = List.flatten(values_1);
       then
         (valuesRet, value_dims);
     // failure
@@ -8540,81 +8536,27 @@ algorithm
   str := stringAppendList({"if (change(",crStr,")) { needToIterate=1; }"});
 end buildDiscreteVarChangesAddEvent;
 
-protected function mixedCollectRelations "function: mixedCollectRelations
-  author: PA
-"
-  input list<BackendDAE.Equation> c_eqn;
-  input list<BackendDAE.Equation> d_eqn;
-  output list<DAE.Exp> res;
-protected
-  list<DAE.Exp> l1,l2;
-algorithm
-  l1 := mixedCollectRelations2(c_eqn);
-  l2 := mixedCollectRelations2(d_eqn);
-  res := listAppend(l1, l2);
-end mixedCollectRelations;
-
-protected function mixedCollectRelations2 "function: mixedCollectRelations2
-  author: PA
-
-  Helper function to mixed_collect_functions.
-"
-  input list<BackendDAE.Equation> inBackendDAEEquationLst;
-  output list<DAE.Exp> outExpExpLst;
-algorithm
-  outExpExpLst:=
-  matchcontinue (inBackendDAEEquationLst)
-    local
-      list<DAE.Exp> l1,l2,l3,res;
-      DAE.Exp e1,e2;
-      list<BackendDAE.Equation> es;
-      Expression.ComponentRef cr;
-    case ({}) then {};
-    case ((BackendDAE.EQUATION(exp = e1,scalar = e2) :: es))
-      equation
-        l1 = Expression.getRelations(e1);
-        l2 = Expression.getRelations(e2);
-        l3 = mixedCollectRelations2(es);
-        res = List.flatten({l1,l2,l3});
-      then
-        res;
-    case ((BackendDAE.SOLVED_EQUATION(componentRef = cr,exp = e1) :: es))
-      equation
-        l1 = Expression.getRelations(e1);
-        l2 = mixedCollectRelations2(es);
-        res = listAppend(l1, l2);
-      then
-        res;
-    case (_) then {};
-  end matchcontinue;
-end mixedCollectRelations2;
-
 protected function generateMixedDiscreteCombinationValues "function generateMixedDiscreteCombinationValues
   author: PA
   Generates all combinations of the values given as argument"
-  input list<list<Integer>> inStringLstLst;
-  output list<list<Integer>> outStringLstLst;
+  input list<list<Integer>> iValues;
+  output list<list<Integer>> oValues;
 algorithm
-  outStringLstLst := matchcontinue (inStringLstLst)
+  oValues := match (iValues)
     local
       list<Integer> value;
-      list<list<Integer>> values_1,values;
-
     case ({value}) then {value};
-
-    case (values)
-      equation
-        values_1 = generateMixedDiscreteCombinationValues1(values);
+    else
       then
-        values_1;
-  end matchcontinue;
+        generateMixedDiscreteCombinationValues1(iValues);
+  end match;
 end generateMixedDiscreteCombinationValues;
 
 protected function generateMixedDiscreteCombinationValues1
-  input list<list<Integer>> inStringLstLst;
-  output list<list<Integer>> outStringLstLst;
+  input list<list<Integer>> iValues;
+  output list<list<Integer>> oValues;
 algorithm
-  outStringLstLst := matchcontinue (inStringLstLst)
+  oValues := matchcontinue (iValues)
     local
       list<list<Integer>> values_1,values_2,values;
       list<Integer> value;
@@ -8659,13 +8601,12 @@ end generateMixedDiscreteCombinationValues2;
 
 protected function generateMixedDiscretePossibleValues2 "function: generateMixedDiscretePossibleValues2
   Helper function to generateMixedDiscretePossibleValues."
-  input list<DAE.Exp> inExpExpLst;
   input list<BackendDAE.Var> inBackendDAEVarLst;
   input Integer inInteger;
   output list<list<Integer>> outStringLstLst;
   output list<Integer> outIntegerLst;
 algorithm
-  (outStringLstLst,outIntegerLst) := matchcontinue (inExpExpLst,inBackendDAEVarLst,inInteger)
+  (outStringLstLst,outIntegerLst) := matchcontinue (inBackendDAEVarLst,inInteger)
     local
       Integer cg_id;
       list<list<Integer>> values;
@@ -8674,16 +8615,16 @@ algorithm
       BackendDAE.Var v;
       list<BackendDAE.Var> vs;
 
-    case (_,{},cg_id) then ({},{});  /* discrete vars cg var_id values value dimension */
+    case ({},cg_id) then ({},{});  /* discrete vars cg var_id values value dimension */
 
-    case (rels,(v :: vs),cg_id) /* booleans, generate true (1.0) and false (0.0) */
+    case ((v :: vs),cg_id) /* booleans, generate true (1.0) and false (0.0) */
       equation
         DAE.T_BOOL(source = _) = BackendVariable.varType(v);
-        (values,dims) = generateMixedDiscretePossibleValues2(rels, vs, cg_id);
+        (values,dims) = generateMixedDiscretePossibleValues2(vs, cg_id);
       then
         ({1,0} :: values, 2 :: dims);
 
-    case (rels,(v :: vs),_)
+    case ((v :: vs),_)
       equation
         DAE.T_INTEGER(source = _) = BackendVariable.varType(v);
         Error.addMessage(Error.INTERNAL_ERROR,

diff --git a/SimulationRuntime/c/math-support/matrix.h b/SimulationRuntime/c/math-support/matrix.h
@@ -88,45 +88,6 @@ void * _omc_hybrj_(void(*) (int*, double*, double*, double *, int*, int*, void*
   printf("}\n"); \
 } while(0)
 
-#define solve_nonlinear_system_mixed(residual, no, userdata) do { \
-   int giveUp = 0; \
-   int retries = 0; \
-   int retries2 = 0; \
-   while(!giveUp) { \
-     giveUp = 1; \
-     _omc_hybrd_(residual,&n, nls_x,nls_fvec,&xtol,&maxfev,&ml,&mu,&epsfcn, \
-          nls_diag,&mode,&factor,&nprint,&info,&nfev,nls_fjac,&ldfjac, \
-        nls_r,&lr,nls_qtf,nls_wa1,nls_wa2,nls_wa3,nls_wa4, userdata); \
-      if (info == 0) { \
-        DEBUG_INFO2(LOG_NONLIN_SYS,"improper input parameters to nonlinear eq. syst %s:%d.\n", __FILE__, __LINE__); \
-      } \
-      if ((info == 4 || info == 5) && retries < 3) { /* first try to decrease factor*/ \
-      retries++; giveUp = 0; \
-      factor = factor / 10.0; \
-      DEBUG_INFO2(LOG_NONLIN_SYS,"Solving nonlinear system: iteration not making progress, trying to decrease factor to %f",factor);  \
-      if ((info == 4 || info == 5) && retries < 5) { /* Then, try with different starting point*/  \
-        int i = 0; \
-        for (i = 0; i < n; i++) { nls_x[i]+=0.1; }; \
-        retries++; giveUp = 0; \
-        DEBUG_INFO2(LOG_NONLIN_SYS,"Solving nonlinear system: iteration not making progress, trying with different starting points (+1e-6)"); \
-        if ((info == 4 || info == 5) && retries2 < 1) { /*Then try with old values (instead of extrapolating )*/ \
-          retries = 0; retries2++; giveUp = 0; \
-          int i = 0; \
-          for (i = 0; i < n; i++) { nls_x[i] = nls_xold[i]; } \
-        } else if (info >= 2 && info <= 5) { \
-          int i = 0; \
-          modelErrorCode=ERROR_NONLINSYS; \
-          DEBUG_INFO2(LOG_NONLIN_SYS, "error solving nonlinear system nr. %d at time %f", no, data->localData[0]->timeValue); \
-          if (DEBUG_FLAG(LOG_NONLIN_SYS)) { \
-            for (i = 0; i < n; i++) { \
-               DEBUG_INFO_AL2(LOG_NONLIN_SYS," residual[%d] = %f",i,nls_fvec[i]); \
-               DEBUG_INFO_AL2(LOG_NONLIN_SYS," x[%d] = %f",i,nls_x[i]); \
-            } \
-          } \
-        } \
-      }\
-} while(0) /* (no trailing ;)*/
-
 #define solve_nonlinear_system(residual, no, userdata) do { \
    int giveUp = 0; \
    int retries = 0; \
@@ -142,8 +103,10 @@ void * _omc_hybrj_(void(*) (int*, double*, double*, double *, int*, int*, void*
      _omc_hybrd_(residual,&n, nls_x,nls_fvec,&xtol,&maxfev,&ml,&mu,&epsfcn, \
           nls_diag,&mode,&factor,&nprint,&info,&nfev,nls_fjac,&ldfjac, \
         nls_r,&lr,nls_qtf,nls_wa1,nls_wa2,nls_wa3,nls_wa4, userdata); \
-      if (info == 0) \
+      if (info == 0) {\
           printErrorEqSyst(IMPROPER_INPUT, no, data->localData[0]->timeValue); \
+          data->found_solution = -1; \
+      }\
       if (info == 1){ \
         if (DEBUG_FLAG(LOG_NONLIN_SYS)) { \
           INFO_AL("### System solved! ###"); \
@@ -197,6 +160,7 @@ void * _omc_hybrj_(void(*) (int*, double*, double*, double *, int*, int*, void*
         if (DEBUG_FLAG(LOG_NONLIN_SYS)) \
         INFO_AL(" - iteration making no progress:\tremove scaling factor at all!"); \
       } else if (info >= 2 && info <= 5) { \
+        data->found_solution = -1; \
         modelErrorCode=ERROR_NONLINSYS; \
         printErrorEqSyst(ERROR_AT_TIME, no, data->localData[0]->timeValue); \
         if (DEBUG_FLAG(LOG_DEBUG)) { \
@@ -255,30 +219,15 @@ assert(ipiv != 0); \
 _omc_dgesv_(&n,&nrhs,&A[0],&lda,ipiv,&b[0],&ldb,&info); \
  if (info < 0) { \
    DEBUG_INFO3(LOG_NONLIN_SYS,"Error solving linear system of equations (no. %d) at time %f. Argument %d illegal.\n",id,data->localData[0]->timeValue,info); \
+   data->found_solution = -1; \
  } \
  else if (info > 0) { \
    DEBUG_INFO2(LOG_NONLIN_SYS,"Error solving linear system of equations (no. %d) at time %f, system is singular.\n",id,data->localData[0]->timeValue); \
+   data->found_solution = -1; \
  } \
 free(ipiv); \
 } while (0) /* (no trailing ; ) */
 
-#define solve_linear_equation_system_mixed(A,b,size,id) do { integer n=size; \
-integer nrhs = 1; /* number of righthand sides*/\
-integer lda = n /* Leading dimension of A */; integer ldb=n; /* Leading dimension of b*/\
-integer * ipiv = (integer*) calloc(n,sizeof(integer)); /* Pivott indices */ \
-assert(ipiv != 0); \
-integer info = 0; /* output */ \
-_omc_dgesv_(&n,&nrhs,&A[0],&lda,ipiv,&b[0],&ldb,&info); \
- if (info < 0) { \
-   if (sim_verbose >= LOG_NONLIN_SYS) \
-     printf("Error solving linear system of equations (no. %d) at time %f. Argument %d illegal.\n",id,localData->timeValue,info); fflush(NULL); \
- } \
- else if (info > 0) { \
-     found_solution = -1; \
- } \
- free(ipiv);\
-} while (0) /* (no trailing ; ) */
-
 #define start_nonlinear_system(size) { double nls_x[size]; \
 double nls_xEx[size] = {0}; \
 double nls_xold[size] = {0}; \
@@ -336,33 +285,34 @@ int ldfjac = size;
 (data->localData[1]->timeValue-data->localData[2]->timeValue))
 
 #define mixed_equation_system(size) do { \
-    int found_solution = 0; \
     int cur_value_indx = 0; \
+    data->found_solution = 0; \
     do { \
         double discrete_loc[size] = {0}; \
         double discrete_loc2[size] = {0};
 
-#define mixed_equation_system_end(size) } while (!found_solution); \
+#define mixed_equation_system_end(size) } while (!data->found_solution); \
  } while(0)
 
 #define check_discrete_values(size,numValues) \
 do { \
   int i = 0; \
-  if (found_solution == -1) { \
+  if (data->found_solution == -1) { \
   /*system of equations failed */ \
-      found_solution = 0; \
+      data->found_solution = 0; \
   } else { \
-      found_solution = 1; \
+      data->found_solution = 1; \
       for (i = 0; i < size; i++) { \
           if (fabs((discrete_loc[i] - discrete_loc2[i])) > 1e-12) {\
-              found_solution=0;\
+          data->found_solution=0;\
+              break;\
           }\
       }\
   }\
-  if (!found_solution ) { \
+  if (!data->found_solution ) { \
       cur_value_indx++; \
       if (cur_value_indx >= numValues/size) { \
-          found_solution = -1; \
+          data->found_solution = -1; \
       } else {\
       /* try next set of values*/ \
           for (i = 0; i < size; i++) { \
@@ -371,7 +321,7 @@ do { \
       } \
   } \
   /* we found a solution*/ \
-  if (found_solution && DEBUG_FLAG(LOG_NONLIN_SYS)){ \
+  if (data->found_solution && DEBUG_FLAG(LOG_NONLIN_SYS)){ \
       int i = 0; \
       printf("Result of mixed system discrete variables:\n"); \
       for (i = 0; i < size; i++) { \

diff --git a/SimulationRuntime/c/simulation_data.h b/SimulationRuntime/c/simulation_data.h
@@ -7,16 +7,16 @@
 *
 * All rights reserved.
 *
-* THIS PROGRAM IS PROVIDED UNDER THE TERMS OF GPL VERSION 3 
-* AND THIS OSMC PUBLIC LICENSE (OSMC-PL). 
-* ANY USE, REPRODUCTION OR DISTRIBUTION OF THIS PROGRAM CONSTITUTES RECIPIENT'S  
+* THIS PROGRAM IS PROVIDED UNDER THE TERMS OF GPL VERSION 3
+* AND THIS OSMC PUBLIC LICENSE (OSMC-PL).
+* ANY USE, REPRODUCTION OR DISTRIBUTION OF THIS PROGRAM CONSTITUTES RECIPIENT'S
 * ACCEPTANCE OF THE OSMC PUBLIC LICENSE.
 *
 * The OpenModelica software and the Open Source Modelica
 * Consortium (OSMC) Public License (OSMC-PL) are obtained
 * from Link�ping University, either from the above address,
-* from the URLs: http://www.ida.liu.se/projects/OpenModelica or  
-* http://www.openmodelica.org, and in the OpenModelica distribution. 
+* from the URLs: http://www.ida.liu.se/projects/OpenModelica or
+* http://www.openmodelica.org, and in the OpenModelica distribution.
 * GNU version 3 is obtained from: http://www.gnu.org/copyleft/gpl.html.
 *
 * This program is distributed WITHOUT ANY WARRANTY; without
@@ -30,7 +30,7 @@
 */
 
 /*! \file simulation_data.h
- * Description: This is the C header file to provide all information 
+ * Description: This is the C header file to provide all information
  * for simulation
  */
 
@@ -141,7 +141,7 @@
   /* Alias data with various types*/
   typedef struct DATA_REAL_ALIAS
   {
-    int negate;  
+    int negate;
     int nameID;  /* Pointer to Alias */
     char aliasType; /* 0 variable, 1 parameter, 2 time */
     VAR_INFO info;
@@ -385,6 +385,7 @@
     SIMULATION_DATA **localData;
     MODEL_DATA modelData;           /* static stuff */
     SIMULATION_INFO simulationInfo;
+    int found_solution;             /* helper for mixed systems */
   }DATA;
 
 #endif