From d36b059fcac25fef14dece8281cafd357be88c94 Mon Sep 17 00:00:00 2001 From: Peter Aronsson Date: Wed, 13 Sep 2006 07:09:48 +0000 Subject: [PATCH] Fixed bug with initialization (alg. vars fixed) Fixed bug with simplification of arrays resulting in wrong operators. git-svn-id: https://openmodelica.org/svn/OpenModelica/trunk@2524 f25d12d1-65f4-0310-ae8a-bbce733d8d8e --- Compiler/Exp.mo | 46 +++++++++++++++++++++----------- Compiler/Static.mo | 25 ++++++++--------- c_runtime/simulation_runtime.cpp | 45 ++++++++++++++++--------------- 3 files changed, 67 insertions(+), 49 deletions(-) diff --git a/Compiler/Exp.mo b/Compiler/Exp.mo index 1fd3688942d..6474801d71c 100644 --- a/Compiler/Exp.mo +++ b/Compiler/Exp.mo @@ -2478,9 +2478,9 @@ algorithm matchcontinue (inExp,inInteger) local Exp e_1,e,e1_1,e2_1,e1,e2,exp; - Type t,t_1; + Type t,t_1,t2; Integer indx,i_1,n; - Operator op; + Operator op,op2; Boolean b; list exps,expl_1; list> expl; @@ -2489,51 +2489,67 @@ algorithm case (UNARY(operator = UMINUS_ARR(ty = t),exp = e),indx) equation e_1 = simplifyAsub(e, indx); + t2 = typeof(e_1); + b = typeBuiltin(t2); + op2 = Util.if_(b,UMINUS(t2),UMINUS_ARR(t2)); then - UNARY(UMINUS_ARR(t),e_1); + UNARY(op2,e_1); case (UNARY(operator = UPLUS_ARR(ty = t),exp = e),indx) equation e_1 = simplifyAsub(e, indx); + t2 = typeof(e_1); + b = typeBuiltin(t2); + op2 = Util.if_(b,UPLUS(t2),UPLUS_ARR(t2)); then - UNARY(UPLUS(t),e_1); + UNARY(op2,e_1); case (BINARY(exp1 = e1,operator = SUB_ARR(ty = t),exp2 = e2),indx) + local Boolean b; Type t2; Operator op2; equation e1_1 = simplifyAsub(e1, indx); e2_1 = simplifyAsub(e2, indx); + t2 = typeof(e1_1); + b = typeBuiltin(t2); + op2 = Util.if_(b,SUB(t2),SUB_ARR(t2)); then - BINARY(e1_1,SUB_ARR(t),e2_1); + BINARY(e1_1,op2,e2_1); case (BINARY(exp1 = e1,operator = MUL_SCALAR_ARRAY(ty = t),exp2 = e2),indx) equation e2_1 = simplifyAsub(e2, indx); e1_1 = simplify(e1); - op = simplifyAsubOperator(e2_1, MUL(t), MUL_SCALAR_ARRAY(t)); + t2 = typeof(e2_1); + b = typeBuiltin(t2); + op = Util.if_(b,MUL(t2),MUL_SCALAR_ARRAY(t2)); then BINARY(e1_1,op,e2_1); case (BINARY(exp1 = e1,operator = MUL_ARRAY_SCALAR(ty = t),exp2 = e2),indx) equation e1_1 = simplifyAsub(e1, indx); e2_1 = simplify(e2); - op = simplifyAsubOperator(e2_1, MUL(t), MUL_SCALAR_ARRAY(t)); + t2 = typeof(e1_1); + b = typeBuiltin(t2); + op = Util.if_(b,MUL(t2),MUL_SCALAR_ARRAY(t2)); then BINARY(e1_1,op,e2_1); case (BINARY(exp1 = e1,operator = DIV_ARRAY_SCALAR(ty = t),exp2 = e2),indx) equation e1_1 = simplifyAsub(e1, indx); e2_1 = simplify(e2); + t2 = typeof(e1_1); + b = typeBuiltin(t2); + op = Util.if_(b,DIV(t2),DIV_ARRAY_SCALAR(t2)); then BINARY(e1_1,DIV(t),e2_1); case (BINARY(exp1 = e1,operator = ADD_ARR(ty = t),exp2 = e2),indx) + local Boolean b; Type t2; Operator op2; equation e1_1 = simplifyAsub(e1, indx); e2_1 = simplifyAsub(e2, indx); + t2 = typeof(e1_1); + b = typeBuiltin(t2); + op2 = Util.if_(b,ADD(t2),ADD_ARR(t2)); then - BINARY(e1_1,ADD_ARR(t),e2_1); - case (BINARY(exp1 = e1,operator = SUB_ARR(ty = t),exp2 = e2),indx) - equation - e1_1 = simplifyAsub(e1, indx); - e2_1 = simplifyAsub(e2, indx); - then - BINARY(e1_1,SUB_ARR(t),e2_1); + BINARY(e1_1,op2,e2_1); + case (ARRAY(ty = t,scalar = b,array = exps),indx) equation i_1 = indx - 1; @@ -3051,7 +3067,7 @@ algorithm end matchcontinue; end arrayEltType; -protected function unliftArray "function: unliftArray +public function unliftArray "function: unliftArray Converts an array type into its element type. " diff --git a/Compiler/Static.mo b/Compiler/Static.mo index d27740d4080..c3f3f41f886 100644 --- a/Compiler/Static.mo +++ b/Compiler/Static.mo @@ -476,7 +476,7 @@ algorithm then (cache,Exp.ARRAY(at,a,es_1),Types.PROP((Types.T_ARRAY(Types.DIM(SOME(l)),t),NONE),const),st); case (cache,env,Absyn.MATRIX(matrix = es),impl,st) - local list> es; + local list> es; Exp.Type mt; equation (cache,_,tps,_) = elabExpListList(cache,env, es, impl, st) "matrix expressions, e.g. {1,0;0,1} with elements of simple type." ; tps_1 = Util.listListMap(tps, Types.getPropType); @@ -2127,9 +2127,7 @@ algorithm (cache,Exp.ARRAY(tp,sc,expl),Types.PROP((Types.T_ARRAY(d1,(Types.T_ARRAY(d2,eltp),_)),_),_),_) = elabExp(cache,env, matexp, impl, NONE); dim1 = Types.arraydimInt(d1); - print("calling elab builtin_transpose2\n"); exp_2 = elabBuiltinTranspose2(expl, 1, dim1); - print("elab builtin_transpose2 done\n"); newtp = (Types.T_ARRAY(d2,(Types.T_ARRAY(d1,eltp),NONE)),NONE); prop = Types.PROP(newtp,Types.C_VAR()); then @@ -6466,8 +6464,8 @@ algorithm b1 = (ds < 20); b2 = (ds2 < 20); true = boolAnd(b1, b2); - elt_tp = Exp.arrayEltType(exptp); - e = createCrefArray2d(cr, 1, ds, ds2, elt_tp, t); + //elt_tp = Exp.arrayEltType(exptp); + e = createCrefArray2d(cr, 1, ds, ds2, exptp, t); then e; @@ -6476,8 +6474,8 @@ algorithm equation false = Types.isArray(t); (ds < 20) = true; - elt_tp = Exp.arrayEltType(exptp); - e = createCrefArray(cr, 1, ds, elt_tp, t); + //elt_tp = Exp.arrayEltType(exptp); + e = createCrefArray(cr, 1, ds, exptp, t); then e; case (e,_) then e; @@ -6537,7 +6535,7 @@ algorithm local Exp.ComponentRef cr,cr_1; Integer indx,ds,indx_1; - Exp.Type et; + Exp.Type et,elt_tp; tuple> t; list expl; Exp.Exp e_1; @@ -6554,7 +6552,8 @@ algorithm Exp.WHOLEDIM()::ss = Exp.crefLastSubs(cr); cr_1 = Exp.crefStripLastSubs(cr); cr_1 = Exp.subscriptCref(cr_1, Exp.INDEX(Exp.ICONST(indx))::ss); - e_1 = crefVectorize(Exp.CREF(cr_1,et), t); + elt_tp = Exp.unliftArray(et); + e_1 = crefVectorize(Exp.CREF(cr_1,elt_tp), t); then Exp.ARRAY(et,false,(e_1 :: expl)); case (cr,indx,ds,et,t) /* no subscript */ @@ -6563,7 +6562,8 @@ algorithm {} = Exp.crefLastSubs(cr); Exp.ARRAY(_,_,expl) = createCrefArray(cr, indx_1, ds, et, t); cr_1 = Exp.subscriptCref(cr, {Exp.INDEX(Exp.ICONST(indx))}); - e_1 = crefVectorize(Exp.CREF(cr_1,et), t); + elt_tp = Exp.unliftArray(et); + e_1 = crefVectorize(Exp.CREF(cr_1,elt_tp), t); then Exp.ARRAY(et,false,(e_1 :: expl)); case (cr,indx,ds,et,t) /* index */ @@ -6602,7 +6602,7 @@ algorithm local Exp.ComponentRef cr,cr_1; Integer indx,ds,ds2,indx_1; - Exp.Type et,tp; + Exp.Type et,tp,elt_tp; tuple> t; list>> ms; Boolean sc; @@ -6619,7 +6619,8 @@ algorithm indx_1 = indx + 1; Exp.MATRIX(_,_,ms) = createCrefArray2d(cr, indx_1, ds, ds2, et, t); cr_1 = Exp.subscriptCref(cr, {Exp.INDEX(Exp.ICONST(indx))}); - Exp.ARRAY(tp,sc,expl) = crefVectorize(Exp.CREF(cr_1,et), t); + elt_tp = Exp.unliftArray(et); + Exp.ARRAY(tp,sc,expl) = crefVectorize(Exp.CREF(cr_1,elt_tp), t); scs = Util.listFill(sc, ds2); row = Util.listThreadTuple(expl, scs); then diff --git a/c_runtime/simulation_runtime.cpp b/c_runtime/simulation_runtime.cpp index cf90897a9df..ffd08ec1c35 100644 --- a/c_runtime/simulation_runtime.cpp +++ b/c_runtime/simulation_runtime.cpp @@ -287,8 +287,8 @@ void euler (DATA * data, } } - -double *static_y; +/* This function calculates the residual value as the sum of squared residual equations. + */ void leastSquare(long *nz, double *z, double *funcValue) { @@ -301,15 +301,17 @@ void leastSquare(long *nz, double *z, double *funcValue) if (globalData->initFixed[indAct++]==0) globalData->parameters[ind] = z[indz++]; - //functionODE(x,xd,y,p,nx,ny,np,&t); - //functionDAE_output(&t,x,xd,y,p); functionODE(); functionDAE_output(); - for (ind=0,indy=0,indAct=2*globalData->nStates; indnAlgebraic; ind++) +/* for (ind=0,indy=0,indAct=2*globalData->nStates; indnAlgebraic; ind++) if (globalData->initFixed[indAct++]==1) globalData->algebraics [ind] = static_y[indy++]; - + + Comment from Bernhard: Even though algebraic variables are "fixed", they are calculated from + the states, so they should be allowed to change when states vary, + and NOT be replaced by their initial values as above. +*/ initial_residual(); for (ind=0, *funcValue=0; indnInitialResiduals; ind++) @@ -318,7 +320,6 @@ void leastSquare(long *nz, double *z, double *funcValue) if (sim_verbose) { cout << "initial residual: " << *funcValue << endl; } - } /** function reportResidualValue @@ -328,10 +329,16 @@ void leastSquare(long *nz, double *z, double *funcValue) int reportResidualValue(double funcValue) { + int i; if (funcValue > 1e-3) { std::cerr << "Error in initialization. System of initial equations are not consistent." << std::endl; std::cerr << "(Least Square function value is " << funcValue << ")" << std::endl; - return -1; + for (i=0; inInitialResiduals; i++) { + if (fabs(globalData->initialResiduals[i]) > 1e-6) { + cout << "residual[" << i << "] = " << globalData->initialResiduals[i] << endl; + } + } + return 0 /*-1*/; } return 0; } @@ -431,7 +438,7 @@ int newuoa_initialization(long& nz,double *z) int initialize(const std::string*method) { long nz; - int ind, indAct, indz, indy, n_static_y; + int ind, indAct, indz, indy; std::string init_method; if (method == NULL) { @@ -449,13 +456,14 @@ int initialize(const std::string*method) if (globalData->initFixed[ind]==0) nz++; } - for (ind=2*globalData->nStates, n_static_y=0; - ind<2*globalData->nStates+globalData->nAlgebraic; ind++){ - if (globalData->initFixed[ind]==1) - n_static_y++; - } - if (n_static_y>0) static_y = new double[n_static_y]; + if (sim_verbose) { + cout << "fixed attribute for states:" << endl; + for(int i=0;inStates; i++) { + cout << getName(&globalData->states[i]) << "(fixed=" << (globalData->initFixed[i]?"true":"false") << ")" + << endl; + } + } // No initial values to calculate. if (nz == 0) { @@ -475,13 +483,6 @@ int initialize(const std::string*method) z[indz++] = globalData->states[ind]; } } - for (ind=0,indy=0,indAct=2*globalData->nStates; indnAlgebraic; ind++) - { - if (globalData->initFixed[indAct++]==1) - { - static_y[indy++] = globalData->algebraics[ind]; - } - } if (init_method == std::string("simplex")) { return simplex_initialization(nz,z);