Skip to content

Commit

Permalink
symbolical QR-solver: improved handling of know-vars
Browse files Browse the repository at this point in the history 
git-svn-id: https://openmodelica.org/svn/OpenModelica/trunk@25410 f25d12d1-65f4-0310-ae8a-bbce733d8d8e
  • Loading branch information
Vitalij Ruge committed Apr 7, 2015
1 parent bd14462 commit 456f638
Show file tree
Hide file tree
Showing 2 changed files with 143 additions and 27 deletions.
131 changes: 120 additions & 11 deletions Compiler/BackEnd/BackendEquation.mo
View file
Expand Up @@ -492,6 +492,89 @@ algorithm
(_, (_,cr_lst)) := traverseExpsOfEquation(e, Expression.traverseSubexpressionsHelper, (Expression.traversingComponentRefFinder, cr_lst));
end traversingEquationCrefFinder;


public function equationUnknownCrefsNotConst "author: Frenkel TUD 2012-05
From the equation and a variable array return all
variables in the equation an not in the variable array or not const"
input list<BackendDAE.Equation> inEquationLst;
input BackendDAE.Variables inVars;
input BackendDAE.Variables inKnVars;
output list<DAE.ComponentRef> cr_lst;
protected
HashTable.HashTable ht;
algorithm
ht := HashTable.emptyHashTable();
(_, (_, (_, _, ht))) := traverseExpsOfEquationList(inEquationLst, Expression.traverseSubexpressionsHelper, (checkEquationsUnknownCrefsExpNotConst, (inVars, inKnVars, ht)));
cr_lst := BaseHashTable.hashTableKeyList(ht);
end equationUnknownCrefsNotConst;

protected function checkEquationsUnknownCrefsExpNotConst
input DAE.Exp inExp;
input tuple<BackendDAE.Variables, BackendDAE.Variables, HashTable.HashTable> inTuple;
output DAE.Exp outExp;
output tuple<BackendDAE.Variables, BackendDAE.Variables, HashTable.HashTable> outTuple;
algorithm
(outExp,outTuple) := matchcontinue (inExp,inTuple)
local
DAE.Exp e, e1;
BackendDAE.Variables vars, knvars;
HashTable.HashTable ht;
DAE.ComponentRef cr;
list<DAE.Exp> expl;
list<DAE.Var> varLst;
list<BackendDAE.Var> var_lst;

// special case for records
case (e as DAE.CREF(componentRef = cr, ty= DAE.T_COMPLEX(varLst=varLst, complexClassType=ClassInf.RECORD(_))), _)
equation
expl = List.map1(varLst, Expression.generateCrefsExpFromExpVar, cr);
(_, outTuple) = Expression.traverseExpList(expl, checkEquationsUnknownCrefsExp, inTuple);
then (e, outTuple);

// special case for arrays
case (e as DAE.CREF(ty = DAE.T_ARRAY()), _)
equation
(e1, true) = Expression.extendArrExp(e, false);
(_, outTuple) = Expression.traverseExpBottomUp(e1, checkEquationsUnknownCrefsExp, inTuple);
then (e, outTuple);

// case for function pointers
case (DAE.CREF(ty=DAE.T_FUNCTION_REFERENCE_FUNC()), _)
then (inExp, inTuple);

// already there
case (DAE.CREF(componentRef = cr), (_, _, ht))
equation
_ = BaseHashTable.get(cr, ht);
then (inExp, inTuple);

// known
case (DAE.CREF(componentRef = cr), (vars, _, _))
algorithm
(var_lst, _) := BackendVariable.getVar(cr, vars);
for var in var_lst loop
true := BackendVariable.isVarConst(var) or BackendVariable.isParam(var);
end for;
then (inExp, inTuple);

case (DAE.CREF(componentRef = cr), (_, knvars, _))
algorithm
(var_lst, _) := BackendVariable.getVar(cr, knvars);
for var in var_lst loop
false := BackendVariable.isInput(var);
end for;
then (inExp, inTuple);

// add it
case (DAE.CREF(componentRef = cr), (vars, knvars, ht))
equation
ht = BaseHashTable.add((cr, 0), ht);
then (inExp, (vars, knvars, ht));

else (inExp,inTuple);
end matchcontinue;
end checkEquationsUnknownCrefsExpNotConst;

public function equationUnknownCrefs "author: Frenkel TUD 2012-05
From the equation and a variable array return all
variables in the equation an not in the variable array."
Expand Down Expand Up @@ -2220,27 +2303,51 @@ public function makeTmpEqnForExp
input Integer offset;
input BackendDAE.EquationArray ieqns;
input BackendDAE.Variables ivars;
output DAE.Exp oExp;
output BackendDAE.EquationArray oeqns;
output BackendDAE.Variables ovars;
input BackendDAE.Shared ishared;

output DAE.Exp oExp = iExp;
output BackendDAE.EquationArray oeqns = ieqns;
output BackendDAE.Variables ovars = ivars;
output BackendDAE.Shared oshared = ishared;

protected
DAE.ComponentRef cr;
BackendDAE.Var tmpvar;
String name_ = "$TMP$" + intString(offset) + "$" + name;
String name_ = "OMC__HELPER__VAR" + intString(offset) + "$" + name;
DAE.Exp x, y;
BackendDAE.Equation eqn;
list<DAE.ComponentRef> cr_lst;
BackendDAE.Variables knowVars = BackendVariable.daeKnVars(oshared);
Boolean b;

algorithm
if not (Expression.isCref(iExp) or Expression.isConst(iExp)) then

cr := ComponentReference.makeCrefIdent(name_, DAE.T_REAL_DEFAULT , {});
tmpvar := BackendVariable.makeVar(cr);
ovars := BackendVariable.addVar(tmpvar, ivars);
x := Expression.crefExp(cr);
(y, _) := ExpressionSimplify.simplify(iExp);
oeqns := BackendEquation.addEquation(BackendDAE.EQUATION(x, y, DAE.emptyElementSource, BackendDAE.EQ_ATTR_DEFAULT_UNKNOWN), ieqns);
oExp := x;
else
oExp := iExp;
oeqns := ieqns;
ovars := ivars;

eqn := BackendDAE.EQUATION(x, y, DAE.emptyElementSource, BackendDAE.EQ_ATTR_DEFAULT_UNKNOWN);
//BackendDump.printEquation(eqn);
cr_lst := equationUnknownCrefsNotConst({eqn}, ivars, knowVars);

b := match(cr_lst)
local DAE.ComponentRef cr1;
case({}) then true;
case({cr1}) then ComponentReference.crefEqualNoStringCompare(cr,cr1);
else false;
end match;

if b then
tmpvar := BackendVariable.setBindExp(tmpvar, SOME(y));
tmpvar := BackendVariable.setVarKind(tmpvar, BackendDAE.CONST());
oshared := BackendVariable.addKnVarDAE(tmpvar, oshared);
else
oeqns := BackendEquation.addEquation(eqn, oeqns);
ovars := BackendVariable.addVar(tmpvar, ovars);
end if;
end if;

end makeTmpEqnForExp;
Expand All @@ -2255,13 +2362,15 @@ public function normalizationVec
input Integer offset;
input BackendDAE.EquationArray ieqns;
input BackendDAE.Variables ivars;
input BackendDAE.Shared ishared;
output array<DAE.Exp> nvec;
output BackendDAE.EquationArray oeqns;
output BackendDAE.Variables ovars;
output BackendDAE.Shared oshared;
protected
DAE.Exp len = Expression.lenVec(vec);
algorithm
(len,oeqns,ovars) := makeTmpEqnForExp(len, name, offset, ieqns, ivars);
(len,oeqns,ovars,oshared) := makeTmpEqnForExp(len, name, offset, ieqns, ivars,ishared);
if Expression.isZero(len) then
fail();
end if;
Expand Down
39 changes: 23 additions & 16 deletions Compiler/BackEnd/ResolveLoops.mo
View file
Expand Up @@ -2216,7 +2216,7 @@ algorithm
beqs = listReverse(beqs);
n = listLength(beqs);
names = List.map(var_lst,BackendVariable.varCref);
(eqns,vars, n) = solveLinearSystem4(beqs, jac, names, var_lst, n, eqns, vars, offset);
(eqns,vars, n, shared) = solveLinearSystem4(beqs, jac, names, var_lst, n, eqns, vars, offset, shared);
//eqns = List.fold(eqn_indxs,BackendEquation.equationRemove,eqns);
then
(BackendDAE.EQSYSTEM(vars,eqns,NONE(),NONE(),matching,stateSets,partitionKind),shared,n);
Expand All @@ -2235,16 +2235,19 @@ protected function solveLinearSystem4
input BackendDAE.EquationArray ieqns;
input BackendDAE.Variables ivars;
input Integer offset;
input BackendDAE.Shared ishared;
output BackendDAE.EquationArray oeqns = ieqns;
output BackendDAE.Variables ovars = ivars;
output Integer offset_ = offset + 1;
output BackendDAE.Shared oshared = ishared;
protected
array<DAE.Exp> R;
array<DAE.Exp> Qb = arrayCreate(n,DAE.RCONST(0.0));
array<DAE.Exp> b = arrayCreate(n,DAE.RCONST(0.0));
array<DAE.Exp> A = arrayCreate(n*n,DAE.RCONST(0.0));
array<DAE.Exp> ax = arrayCreate(n,DAE.RCONST(0.0));
DAE.Exp a, x;
array<DAE.Exp> solvedX = arrayCreate(n,DAE.RCONST(0.0));
DAE.Exp a, x, eqn_exp, eqn_scalar;
Integer m, ii, jj, mm;
list<DAE.Exp> x_lst = List.map(cr_x, Expression.crefExp);
DAE.ComponentRef cr;
Expand All @@ -2262,14 +2265,14 @@ algorithm
for i in 1:mm loop
(jj,ii,BackendDAE.RESIDUAL_EQUATION(exp = a)) :: jac_ := jac_; // jac(1) = a11, jac(2)=a12,.., jac(n+1) = an1
m := ii + (jj-1)*n;
(a, oeqns, ovars) := BackendEquation.makeTmpEqnForExp(a, "QR$A$" + intString(m), offset, oeqns, ovars);
(a, oeqns, ovars, oshared) := BackendEquation.makeTmpEqnForExp(a, "QR$A$" + intString(m), offset, oeqns, ovars, oshared);
arrayUpdate(A,m,a);
end for;

// b
m := 1;
for b_ in b_lst loop
(a, oeqns, ovars) := BackendEquation.makeTmpEqnForExp(b_, "QR$b$" + intString(m), offset, oeqns, ovars);
(a, oeqns, ovars, oshared) := BackendEquation.makeTmpEqnForExp(b_, "QR$b$" + intString(m), offset, oeqns, ovars, oshared);
arrayUpdate(b, m, a);
m := m + 1;
end for;
Expand All @@ -2289,10 +2292,10 @@ algorithm
//qrDecomposition3(ax, n, false, "x");

// A*x = b -> R*x = Q'b
(R, Qb, oeqns, ovars) := qrDecomposition(A, n, b, oeqns, ovars, offset);
(R, Qb, oeqns, ovars, oshared) := qrDecomposition(A, n, b, oeqns, ovars, offset, oshared);

// R*x = Q'*b
for i in 1:n loop
for i in n:-1:1 loop
m := (i-1)*n;
a := Expression.makeSum1(list(Expression.expMul(arrayGet(R, m + j), arrayGet(ax, j)) for j in i:n));
eqn := BackendDAE.EQUATION(a, arrayGet(Qb,i), DAE.emptyElementSource, BackendDAE.EQ_ATTR_DEFAULT_UNKNOWN);
Expand All @@ -2311,11 +2314,13 @@ protected function qrDecomposition
input array<DAE.Exp> ib;
input BackendDAE.EquationArray ieqns;
input BackendDAE.Variables ivars;
input Integer offset;
input BackendDAE.Shared ishared;
output array<DAE.Exp> R = arrayCreate(n*n,DAE.RCONST(0.0));
output array<DAE.Exp> b = arrayCreate(n,DAE.RCONST(0.0));
output BackendDAE.EquationArray oeqns = ieqns;
output BackendDAE.Variables ovars = ivars;
input Integer offset;
output BackendDAE.Shared oshared;
protected
array<DAE.Exp> Q = arrayCreate(n*n,DAE.RCONST(0.0));
array<DAE.Exp> v = arrayCreate(n,DAE.RCONST(0.0));
Expand All @@ -2334,27 +2339,27 @@ protected
algorithm
//Gram–Schmidt process
v := qrDecomposition1(A,n,kk);
(u, oeqns, ovars) := BackendEquation.normalizationVec(v,"QR$NOM$" + intString(kk), offset, oeqns, ovars);
(u, oeqns, ovars, oshared) := BackendEquation.normalizationVec(v,"QR$NOM$" + intString(kk), offset, oeqns, ovars, ishared);

for j in 1:n loop
(a,_) := ExpressionSimplify.simplify(arrayGet(u,j));
(a, oeqns, ovars) := BackendEquation.makeTmpEqnForExp(a, "QR$Q$" + intString(kk + (j-1)*n), offset, oeqns, ovars);
(a, oeqns, ovars,oshared) := BackendEquation.makeTmpEqnForExp(a, "QR$Q$" + intString(kk + (j-1)*n), offset, oeqns, ovars,oshared);
arrayUpdate(Q, kk + (j-1)*n, a);
end for;

for k in 1:m loop
v := qrDecomposition1(A,n,k+1);
for j in 1:k loop
u := qrDecomposition1(Q,n,j);
(v, oeqns, ovars) := gramSchmidtProcessHelper(v,u,"QR$W$" + intString(kk) + "$" + intString(kk), offset, oeqns, ovars);
(v, oeqns, ovars,oshared) := gramSchmidtProcessHelper(v,u,"QR$W$" + intString(kk) + "$" + intString(kk), offset, oeqns, ovars,oshared);
kk := kk +1;
end for;
(u, oeqns, ovars) := BackendEquation.normalizationVec(v,"QR$NOM$" + intString(k+1), offset, oeqns, ovars);
(u, oeqns, ovars, oshared) := BackendEquation.normalizationVec(v,"QR$NOM$" + intString(k+1), offset, oeqns, ovars, oshared);
//qrDecomposition3(u, n, false, "u");
for j in 1:n loop
nn := k+1 + (j-1)*n;
(a,_) := ExpressionSimplify.simplify(arrayGet(u,j));
(a, oeqns, ovars) := BackendEquation.makeTmpEqnForExp(a, "QR$Q$" + intString(nn), offset, oeqns, ovars);
(a, oeqns, ovars, oshared) := BackendEquation.makeTmpEqnForExp(a, "QR$Q$" + intString(nn), offset, oeqns, ovars, oshared);
arrayUpdate(Q, nn, a);
end for;

Expand All @@ -2370,7 +2375,7 @@ algorithm
for j in i:n loop
y := qrDecomposition1(A,n,j);
a := Expression.makeScalarProduct(x,y);
(a, oeqns, ovars) := BackendEquation.makeTmpEqnForExp(a, "QR$R$" + intString(m + j), offset, oeqns, ovars);
(a, oeqns, ovars,oshared) := BackendEquation.makeTmpEqnForExp(a, "QR$R$" + intString(m + j), offset, oeqns, ovars,oshared);
arrayUpdate(R, m+j, a);
end for;
end for;
Expand All @@ -2382,7 +2387,7 @@ algorithm
x := qrDecomposition1(Q,n,i);
//qrDecomposition3(x, n, false, "x" + intString(i));
a := Expression.makeScalarProduct(x,ib);
(a, oeqns, ovars) := BackendEquation.makeTmpEqnForExp(a, "QR$Qb$" + intString(i), offset, oeqns, ovars);
(a, oeqns, ovars,oshared) := BackendEquation.makeTmpEqnForExp(a, "QR$Qb$" + intString(i), offset, oeqns, ovars, oshared);
arrayUpdate(b, i, a);
end for;
//qrDecomposition3(b, n, false, "Qb");
Expand Down Expand Up @@ -2446,18 +2451,20 @@ protected function gramSchmidtProcessHelper
input Integer offset;
input BackendDAE.EquationArray ieqns;
input BackendDAE.Variables ivars;
input BackendDAE.Shared ishared;
output array<DAE.Exp> v;
output BackendDAE.EquationArray oeqns;
output BackendDAE.Variables ovars;
output BackendDAE.Shared oshared;
protected
DAE.Exp h = Expression.makeScalarProduct(w,u);
Integer n = arrayLength(w);
algorithm
(h,oeqns,ovars) := BackendEquation.makeTmpEqnForExp(h, name + "_h", offset, ieqns, ivars);
(h,oeqns,ovars,oshared) := BackendEquation.makeTmpEqnForExp(h, name + "_h", offset, ieqns, ivars, ishared);
v := Array.map1(u, Expression.expMul, h);
v := Expression.subVec(w,v);
for i in 1:n loop
(h,oeqns,ovars) := BackendEquation.makeTmpEqnForExp(arrayGet(v,i), name + "_" + intString(i), offset, oeqns, ovars);
(h,oeqns,ovars,oshared) := BackendEquation.makeTmpEqnForExp(arrayGet(v,i), name + "_" + intString(i), offset, oeqns, ovars, oshared);
arrayUpdate(v,i,h);
end for;

Expand Down

0 comments on commit 456f638

Please sign in to comment.