Skip to content

Commit

Permalink
Implement state estimation algorithm (#9848)
Browse files Browse the repository at this point in the history
* implement state estimation prototype

* generate jacobian matrix H for state estimation

* generate html report for state estimation and runtime code for unmeasured variables

* generate html report for state estimation

* remove debug logs

* improve debug headers

* expected output

* fix expected output

* generate code to count relatedBoundaryConditions

* improve debugging messages

* update tests

* Trigger build

* copy result data only for state estimation

* improve function comments

* Trigger build

* improve comments

* add and update test
  • Loading branch information
arun3688 committed Dec 8, 2022
1 parent 24e3683 commit 43cf471
Show file tree
Hide file tree
Showing 39 changed files with 2,929 additions and 317 deletions.
5 changes: 4 additions & 1 deletion OMCompiler/Compiler/BackEnd/BackendDAE.mo
Expand Up @@ -189,9 +189,12 @@ end BackendDAEType;

uniontype DataReconciliationData
record DATA_RECON
Jacobian symbolicJacobian "SET_S w.r.t ...";
Jacobian symbolicJacobian "jacobians for set-C and set-S";
Variables setcVars "setc solved vars";
Variables datareconinputs;
Option<Variables> setBVars "setB solved vars which computes boundary conditions";
Option<Jacobian> symbolicJacobianH "For solving state estimation we need two Jacobians F for data Reconciliation and H for boundary conditions set-B and set-Sprime";
Integer relatedBoundaryConditions "count number of boundary conditions which failed the extraction algorithm";
// ... maybe more DATA for the code generation
end DATA_RECON;
end DataReconciliationData;
Expand Down
1 change: 1 addition & 0 deletions OMCompiler/Compiler/BackEnd/BackendDAEUtil.mo
Expand Up @@ -8337,6 +8337,7 @@ public function allPreOptimizationModules
(BackendDAEUtil.introduceOutputAliases, "introduceOutputAliases"),
(DataReconciliation.newExtractionAlgorithm, "dataReconciliation"),
(DataReconciliation.extractBoundaryCondition, "dataReconciliationBoundaryConditions"),
(DataReconciliation.stateEstimation, "dataReconciliationStateEstimation"),
(DynamicOptimization.createDynamicOptimization,"createDynamicOptimization"),
(BackendInline.normalInlineFunction, "normalInlineFunction"),
(EvaluateParameter.evaluateParameters, "evaluateParameters"),
Expand Down
12 changes: 12 additions & 0 deletions OMCompiler/Compiler/BackEnd/BackendVariable.mo
Expand Up @@ -800,6 +800,18 @@ algorithm
end match;
end varHasUncertainValueRefine;

public function varHasUncertainValuePropagate
"Returns true if the specified variable has the attribute uncertain and the
value of it is Uncertainty.propagate, false otherwise."
input BackendDAE.Var var;
output Boolean b;
algorithm
b := match (var)
case (BackendDAE.VAR(values=SOME(DAE.VAR_ATTR_REAL(uncertainOption=SOME(DAE.PROPAGATE()))))) then true;
else false;
end match;
end varHasUncertainValuePropagate;

public function varDistribution "author: Peter Aronsson, 2012-05
Returns Distribution record of a variable."
input BackendDAE.Var var;
Expand Down
455 changes: 449 additions & 6 deletions OMCompiler/Compiler/BackEnd/DataReconciliation.mo

Large diffs are not rendered by default.

2 changes: 1 addition & 1 deletion OMCompiler/Compiler/BackEnd/Uncertainties.mo
Expand Up @@ -492,7 +492,7 @@ algorithm
(simcodejacobian,shared)=SymbolicJacobian.getSymbolicJacobian(outDiffVars,outResidualEqns,outResidualVars,outOtherEqns,outOtherVars,shared,BackendVariable.listVar(List.map1r(extractedvars,BackendVariable.getVarAt,allVars)),"F",false);
// put the jacobian also into shared object
setcVars=BackendVariable.listVar(List.map1r(getRemovedEquationSolvedVariables(tempsetC,var),BackendVariable.getVarAt,allVars));
shared.dataReconciliationData = SOME(BackendDAE.DATA_RECON(symbolicJacobian=simcodejacobian,setcVars=outResidualVars,datareconinputs=outDiffVars));
shared.dataReconciliationData = SOME(BackendDAE.DATA_RECON(symbolicJacobian=simcodejacobian,setcVars=outResidualVars,datareconinputs=outDiffVars,setBVars=NONE(), symbolicJacobianH=NONE(), relatedBoundaryConditions=0));
//BackendDump.dumpVariables(setcVars,"SET_C_SOLVEDVARS");

// Prepare the final DAE System with Set-c equations as residual equations
Expand Down
145 changes: 76 additions & 69 deletions OMCompiler/Compiler/NSimCode/NSimCode.mo
Expand Up @@ -261,7 +261,7 @@ public
list<SimJacobian> jacobians;
UnorderedMap<ComponentRef, SimVar> simcode_map;
Option<DaeModeData> daeModeData;
SimJacobian jacA, jacB, jacC, jacD, jacF;
SimJacobian jacA, jacB, jacC, jacD, jacF, jacH;
list<SimStrongComponent.Block> inlineEquations; // ToDo: what exactly is this?

case BackendDAE.MAIN(varData = varData as BVariable.VAR_DATA_SIM(), eqData = eqData as BEquation.EQ_DATA_SIM())
Expand Down Expand Up @@ -349,10 +349,11 @@ public
(jacC, simCodeIndices) := SimJacobian.empty("C", simCodeIndices);
(jacD, simCodeIndices) := SimJacobian.empty("D", simCodeIndices);
(jacF, simCodeIndices) := SimJacobian.empty("F", simCodeIndices);
(jacH, simCodeIndices) := SimJacobian.empty("H", simCodeIndices);
//jacobians := jacA :: jacB :: jacC :: jacD :: jacF :: jacobians;
jacobians := listReverse(jacF :: jacD :: jacC :: jacB :: jacA :: jacobians);
jacobians := listReverse(jacH :: jacF :: jacD :: jacC :: jacB :: jacA :: jacobians);
// jacobian blocks only from simulation jacobians
jac_blocks := SimJacobian.getJacobiansBlocks({jacA, jacB, jacC, jacD, jacF});
jac_blocks := SimJacobian.getJacobiansBlocks({jacA, jacB, jacC, jacD, jacF, jacH});
(jac_blocks, simCodeIndices) := SimStrongComponent.Block.fixIndices(jac_blocks, {}, simCodeIndices);

generic_loop_calls := list(SimGenericCall.fromEquation(tpl) for tpl in UnorderedMap.toList(simCodeIndices.generic_call_map));
Expand Down Expand Up @@ -652,6 +653,8 @@ public
Integer numSetcVars;
Integer numDataReconVars;
Integer numRealIntputVars;
Integer numSetbVars;
Integer numRelatedBoundaryConditions;
end VAR_INFO;

function create
Expand All @@ -661,79 +664,83 @@ public
output VarInfo varInfo;
algorithm
varInfo := VAR_INFO(
numZeroCrossings = sum(StateEvent.size(se) for se in eventInfo.stateEvents),
numTimeEvents = listLength(eventInfo.timeEvents),
numRelations = sum(StateEvent.size(se) for se in eventInfo.stateEvents),
numMathEventFunctions = eventInfo.numberMathEvents,
numStateVars = listLength(vars.stateVars),
numAlgVars = listLength(vars.algVars),
numDiscreteReal = listLength(vars.discreteAlgVars),
numIntAlgVars = listLength(vars.intAlgVars),
numBoolAlgVars = listLength(vars.boolAlgVars),
numAlgAliasVars = listLength(vars.aliasVars),
numIntAliasVars = listLength(vars.intAliasVars),
numBoolAliasVars = listLength(vars.boolAliasVars),
numParams = listLength(vars.paramVars),
numIntParams = listLength(vars.intParamVars),
numBoolParams = listLength(vars.boolParamVars),
numOutVars = listLength(vars.outputVars),
numInVars = listLength(vars.inputVars),
numExternalObjects = listLength(vars.extObjVars),
numStringAlgVars = listLength(vars.stringAlgVars),
numStringParamVars = listLength(vars.stringParamVars),
numStringAliasVars = listLength(vars.stringAliasVars),
numEquations = simCodeIndices.equationIndex,
numLinearSystems = simCodeIndices.linearSystemIndex,
numNonLinearSystems = simCodeIndices.nonlinearSystemIndex,
numMixedSystems = 0,
numStateSets = 0,
numJacobians = simCodeIndices.nonlinearSystemIndex + 5, // #nonlinSystems + 5 simulation jacs (add state sets later!)
numOptimizeConstraints = 0,
numOptimizeFinalConstraints = 0,
numSensitivityParameters = 0,
numSetcVars = 0,
numDataReconVars = 0,
numRealIntputVars = 0);
numZeroCrossings = sum(StateEvent.size(se) for se in eventInfo.stateEvents),
numTimeEvents = listLength(eventInfo.timeEvents),
numRelations = sum(StateEvent.size(se) for se in eventInfo.stateEvents),
numMathEventFunctions = eventInfo.numberMathEvents,
numStateVars = listLength(vars.stateVars),
numAlgVars = listLength(vars.algVars),
numDiscreteReal = listLength(vars.discreteAlgVars),
numIntAlgVars = listLength(vars.intAlgVars),
numBoolAlgVars = listLength(vars.boolAlgVars),
numAlgAliasVars = listLength(vars.aliasVars),
numIntAliasVars = listLength(vars.intAliasVars),
numBoolAliasVars = listLength(vars.boolAliasVars),
numParams = listLength(vars.paramVars),
numIntParams = listLength(vars.intParamVars),
numBoolParams = listLength(vars.boolParamVars),
numOutVars = listLength(vars.outputVars),
numInVars = listLength(vars.inputVars),
numExternalObjects = listLength(vars.extObjVars),
numStringAlgVars = listLength(vars.stringAlgVars),
numStringParamVars = listLength(vars.stringParamVars),
numStringAliasVars = listLength(vars.stringAliasVars),
numEquations = simCodeIndices.equationIndex,
numLinearSystems = simCodeIndices.linearSystemIndex,
numNonLinearSystems = simCodeIndices.nonlinearSystemIndex,
numMixedSystems = 0,
numStateSets = 0,
numJacobians = simCodeIndices.nonlinearSystemIndex + 5, // #nonlinSystems + 5 simulation jacs (add state sets later!)
numOptimizeConstraints = 0,
numOptimizeFinalConstraints = 0,
numSensitivityParameters = 0,
numSetcVars = 0,
numDataReconVars = 0,
numRealIntputVars = 0,
numSetbVars = 0,
numRelatedBoundaryConditions = 0);
end create;

function convert
input VarInfo varInfo;
output OldSimCode.VarInfo oldVarInfo;
algorithm
oldVarInfo := OldSimCode.VARINFO(
numZeroCrossings = varInfo.numZeroCrossings,
numTimeEvents = varInfo.numTimeEvents,
numRelations = varInfo.numRelations,
numMathEventFunctions = varInfo.numMathEventFunctions,
numStateVars = varInfo.numStateVars,
numAlgVars = varInfo.numAlgVars,
numDiscreteReal = varInfo.numDiscreteReal,
numIntAlgVars = varInfo.numIntAlgVars,
numBoolAlgVars = varInfo.numBoolAlgVars,
numAlgAliasVars = varInfo.numAlgAliasVars,
numIntAliasVars = varInfo.numIntAliasVars,
numBoolAliasVars = varInfo.numBoolAliasVars,
numParams = varInfo.numParams,
numIntParams = varInfo.numIntParams,
numBoolParams = varInfo.numBoolParams,
numOutVars = varInfo.numOutVars,
numInVars = varInfo.numInVars,
numExternalObjects = varInfo.numExternalObjects,
numStringAlgVars = varInfo.numStringAlgVars,
numStringParamVars = varInfo.numStringParamVars,
numStringAliasVars = varInfo.numStringAliasVars,
numEquations = varInfo.numEquations,
numLinearSystems = varInfo.numLinearSystems,
numNonLinearSystems = varInfo.numNonLinearSystems,
numMixedSystems = varInfo.numMixedSystems,
numStateSets = varInfo.numStateSets,
numJacobians = varInfo.numJacobians,
numOptimizeConstraints = varInfo.numOptimizeConstraints,
numOptimizeFinalConstraints = varInfo.numOptimizeFinalConstraints,
numSensitivityParameters = varInfo.numSensitivityParameters,
numSetcVars = varInfo.numSetcVars,
numDataReconVars = varInfo.numDataReconVars,
numRealInputVars = varInfo.numRealIntputVars);
numZeroCrossings = varInfo.numZeroCrossings,
numTimeEvents = varInfo.numTimeEvents,
numRelations = varInfo.numRelations,
numMathEventFunctions = varInfo.numMathEventFunctions,
numStateVars = varInfo.numStateVars,
numAlgVars = varInfo.numAlgVars,
numDiscreteReal = varInfo.numDiscreteReal,
numIntAlgVars = varInfo.numIntAlgVars,
numBoolAlgVars = varInfo.numBoolAlgVars,
numAlgAliasVars = varInfo.numAlgAliasVars,
numIntAliasVars = varInfo.numIntAliasVars,
numBoolAliasVars = varInfo.numBoolAliasVars,
numParams = varInfo.numParams,
numIntParams = varInfo.numIntParams,
numBoolParams = varInfo.numBoolParams,
numOutVars = varInfo.numOutVars,
numInVars = varInfo.numInVars,
numExternalObjects = varInfo.numExternalObjects,
numStringAlgVars = varInfo.numStringAlgVars,
numStringParamVars = varInfo.numStringParamVars,
numStringAliasVars = varInfo.numStringAliasVars,
numEquations = varInfo.numEquations,
numLinearSystems = varInfo.numLinearSystems,
numNonLinearSystems = varInfo.numNonLinearSystems,
numMixedSystems = varInfo.numMixedSystems,
numStateSets = varInfo.numStateSets,
numJacobians = varInfo.numJacobians,
numOptimizeConstraints = varInfo.numOptimizeConstraints,
numOptimizeFinalConstraints = varInfo.numOptimizeFinalConstraints,
numSensitivityParameters = varInfo.numSensitivityParameters,
numSetcVars = varInfo.numSetcVars,
numDataReconVars = varInfo.numDataReconVars,
numRealInputVars = varInfo.numRealIntputVars,
numSetbVars = varInfo.numSetbVars,
numRelatedBoundaryConditions = varInfo.numRelatedBoundaryConditions);
end convert;
end VarInfo;

Expand Down
1 change: 1 addition & 0 deletions OMCompiler/Compiler/NSimCode/NSimCodeUtil.mo
Expand Up @@ -84,6 +84,7 @@ public
addListSimCodeMap(simVars.sensitivityVars, simcode_map);
addListSimCodeMap(simVars.dataReconSetcVars, simcode_map);
addListSimCodeMap(simVars.dataReconinputVars, simcode_map);
addListSimCodeMap(simVars.dataReconSetBVars, simcode_map);
end createSimCodeMap;

function addListSimCodeMap
Expand Down
13 changes: 9 additions & 4 deletions OMCompiler/Compiler/NSimCode/NSimVar.mo
Expand Up @@ -710,6 +710,7 @@ public
list<SimVar> sensitivityVars "variable used to calculate sensitivities for parameters nSensitivitityParameters + nRealParam*nStates";
list<SimVar> dataReconSetcVars;
list<SimVar> dataReconinputVars;
list<SimVar> dataReconSetBVars;
end SIMVARS;

function toString
Expand Down Expand Up @@ -751,6 +752,7 @@ public
list<SimVar> sensitivityVars = {};
list<SimVar> dataReconSetcVars = {};
list<SimVar> dataReconinputVars = {};
list<SimVar> dataReconSetBVars = {};
algorithm
_ := match varData
local
Expand Down Expand Up @@ -807,7 +809,8 @@ public
realOptimizeFinalConstraintsVars = realOptimizeFinalConstraintsVars,
sensitivityVars = sensitivityVars,
dataReconSetcVars = dataReconSetcVars,
dataReconinputVars = dataReconinputVars
dataReconinputVars = dataReconinputVars,
dataReconSetBVars = dataReconSetBVars
);
end create;

Expand Down Expand Up @@ -863,7 +866,8 @@ public
+ listLength(simVars.realOptimizeFinalConstraintsVars)
+ listLength(simVars.sensitivityVars)
+ listLength(simVars.dataReconSetcVars)
+ listLength(simVars.dataReconinputVars);
+ listLength(simVars.dataReconinputVars)
+ listLength(simVars.dataReconSetBVars);
end size;

function convert
Expand Down Expand Up @@ -899,7 +903,8 @@ public
realOptimizeFinalConstraintsVars = SimVar.convertList(simVars.realOptimizeFinalConstraintsVars),
sensitivityVars = SimVar.convertList(simVars.sensitivityVars),
dataReconSetcVars = SimVar.convertList(simVars.dataReconSetcVars),
dataReconinputVars = SimVar.convertList(simVars.dataReconinputVars));
dataReconinputVars = SimVar.convertList(simVars.dataReconinputVars),
dataReconSetBVars = SimVar.convertList(simVars.dataReconSetBVars));
end convert;

function createSimVarLists
Expand Down Expand Up @@ -1060,7 +1065,7 @@ public
end SimVars;

constant SimVars emptySimVars = SIMVARS({}, {}, {}, {}, {}, {}, {}, {}, {}, {}, {},
{}, {}, {}, {}, {}, {}, {}, {}, {}, {}, {}, {}, {}, {}, {}, {}, {}, {}, {});
{}, {}, {}, {}, {}, {}, {}, {}, {}, {}, {}, {}, {}, {}, {}, {}, {}, {}, {}, {});

type SplitType = enumeration(NONE, TYPE);
type VarType = enumeration(SIMULATION, PARAMETER, ALIAS, RESIDUAL); // ToDo: PRE, OLD, RELATIONS...
Expand Down
6 changes: 3 additions & 3 deletions OMCompiler/Compiler/SimCode/ReduceDAE.mo
Expand Up @@ -2174,14 +2174,14 @@ algorithm
local
list<SimCodeVar.SimVar> states,derVar,alg,disAlg,intAlg,boolAlg,inVar,outVar,algAlias,intAlias,boolAlias,param,
intParam,boolParam,stringAlg,stringParam,stringAlias,extObjVar,const,intConst,boolConst,stringConst,jacobianVar,
seedVar,realOptConst,realOptFinalConst,sensVar,setcVar,datareconinputvar;
seedVar,realOptConst,realOptFinalConst,sensVar,setcVar,datareconinputvar,setBVar;
SimCodeVar.SimVar simVar_1,simVar_2;
list<SimCodeVar.SimVar> param_1,param_2;
Integer i,p;
String name, name1, name2, indexStr;
case (SimCodeVar.SIMVARS(states,derVar,alg,disAlg,intAlg,boolAlg,inVar,outVar,algAlias,intAlias,boolAlias,param,
intParam,boolParam,stringAlg,stringParam,stringAlias,extObjVar,const,intConst,boolConst,stringConst,jacobianVar,
seedVar,realOptConst,realOptFinalConst,sensVar,setcVar,datareconinputvar),p,i)
seedVar,realOptConst,realOptFinalConst,sensVar,setcVar,datareconinputvar,setBVar),p,i)

equation
indexStr = intString(i);
Expand All @@ -2207,7 +2207,7 @@ algorithm

(SimCodeVar.SIMVARS(states,derVar,alg,disAlg,intAlg,boolAlg,inVar,outVar,algAlias,intAlias,boolAlias,param_2,
intParam,boolParam,stringAlg,stringParam,stringAlias,extObjVar,const,intConst,boolConst,stringConst,jacobianVar,
seedVar,realOptConst,realOptFinalConst,sensVar,setcVar,datareconinputvar),name);
seedVar,realOptConst,realOptFinalConst,sensVar,setcVar,datareconinputvar,setBVar),name);

end matchcontinue;
end createLabelVar;
Expand Down
2 changes: 2 additions & 0 deletions OMCompiler/Compiler/SimCode/SimCode.mo
Expand Up @@ -327,6 +327,8 @@ uniontype VarInfo "Number of variables of various types in a Modelica model."
Integer numSetcVars;
Integer numDataReconVars;
Integer numRealInputVars "for fmi cs to interpolate inputs";
Integer numSetbVars "for data reconciliation setB vars";
Integer numRelatedBoundaryConditions "for data reconciliation count number of boundary conditions which failed the extraction algorithm";
end VARINFO;
end VarInfo;

Expand Down

0 comments on commit 43cf471

Please sign in to comment.