Parameter Estimation 1#
\[min \sum ((X_{1}(t_{i}) - X_{1,meas}(t_i))^{2})\]
\[s.t \;\; X_{1,dot} = X_{2} \;\;\;\;\;\;\;\; X_{1}(0) = p_{1}\]
\[\;\; X_{2,dot} = 1-2*X_{2}-X_{1} \;\;\;\;\;\;\;\; X_{2}(0) = p_{2}\]
\[-1.5 \leq p_{1}, p_{2} \leq 1.5\]
\[t_{f} = 6\]
import pyomo.environ as pyo
from pyomo.dae import ContinuousSet, DerivativeVar
measurements = {1:0.264, 2:0.594, 3:0.801, 5:0.959}
model = pyo.ConcreteModel()
model.t = ContinuousSet(initialize=measurements.keys(),bounds=(0, 6))
model.x1 = pyo.Var(model.t)
model.x2 = pyo.Var(model.t)
model.p1 = pyo.Var(bounds=(-1.5,1.5))
model.p2 = pyo.Var(bounds=(-1.5,1.5))
model.x1dot = DerivativeVar(model.x1,wrt=model.t)
model.x2dot = DerivativeVar(model.x2)
def _init_conditions(model):
yield model.x1[0] == model.p1
yield model.x2[0] == model.p2
model.init_conditions = pyo.ConstraintList(rule=_init_conditions)
# Alternate way to declare initial conditions
#def _initx1(model):
# return model.x1[0] == model.p1
#model.initx1 = pyo.Constraint(rule=_initx1)
#def _initx2(model):
# return model.x2[0] == model.p2
#model.initx2 = pyo.Constraint(rule=_initx2)
def _x1dot(model,i):
return model.x1dot[i] == model.x2[i]
model.x1dotcon = pyo.Constraint(model.t, rule=_x1dot)
def _x2dot(model,i):
return model.x2dot[i] == 1-2*model.x2[i]-model.x1[i]
model.x2dotcon = pyo.Constraint(model.t, rule=_x2dot)
def _obj(model):
return sum((model.x1[i]-measurements[i])**2 for i in measurements.keys())
model.obj = pyo.Objective(rule=_obj)
# Discretize model using Orthogonal Collocation
discretizer = pyo.TransformationFactory('dae.collocation')
discretizer.apply_to(model,nfe=8,ncp=5)
solver = pyo.SolverFactory('ipopt')
results = solver.solve(model,tee=True)
t_meas = sorted(list(measurements.keys()))
x1_meas = [pyo.value(measurements[i]) for i in sorted(measurements.keys())]
t = list(model.t)
x1 = [pyo.value(model.x1[i]) for i in model.t]
import matplotlib.pyplot as plt
plt.plot(t,x1)
plt.plot(t_meas,x1_meas,'o')
plt.xlabel('t')
plt.ylabel('x')
plt.title('Dynamic Parameter Estimation Using Collocation')
plt.show()
Ipopt 3.14.17:
******************************************************************************
This program contains Ipopt, a library for large-scale nonlinear optimization.
Ipopt is released as open source code under the Eclipse Public License (EPL).
For more information visit https://github.com/coin-or/Ipopt
******************************************************************************
This is Ipopt version 3.14.17, running with linear solver MUMPS 5.7.3.
Number of nonzeros in equality constraint Jacobian...: 769
Number of nonzeros in inequality constraint Jacobian.: 0
Number of nonzeros in Lagrangian Hessian.............: 4
Total number of variables............................: 166
variables with only lower bounds: 0
variables with lower and upper bounds: 2
variables with only upper bounds: 0
Total number of equality constraints.................: 164
Total number of inequality constraints...............: 0
inequality constraints with only lower bounds: 0
inequality constraints with lower and upper bounds: 0
inequality constraints with only upper bounds: 0
iter objective inf_pr inf_du lg(mu) ||d|| lg(rg) alpha_du alpha_pr ls
0 1.9838140e+00 1.00e+00 7.15e-02 -1.0 0.00e+00 - 0.00e+00 0.00e+00 0
1 3.7714229e-07 2.54e-15 7.19e-16 -1.0 1.00e+00 - 1.00e+00 1.00e+00h 1
2 3.6867493e-07 2.57e-15 5.91e-17 -2.5 1.45e-04 - 1.00e+00 1.00e+00h 1
3 3.5639520e-07 2.83e-15 6.61e-17 -3.8 1.72e-03 - 1.00e+00 1.00e+00h 1
4 3.5566463e-07 5.86e-15 8.20e-17 -5.7 5.27e-04 - 1.00e+00 1.00e+00h 1
5 3.5566442e-07 1.07e-14 8.02e-17 -8.6 9.04e-06 - 1.00e+00 1.00e+00h 1
Number of Iterations....: 5
(scaled) (unscaled)
Objective...............: 3.5566442170075718e-07 3.5566442170075718e-07
Dual infeasibility......: 8.0230960763927328e-17 8.0230960763927328e-17
Constraint violation....: 1.0651202142497596e-14 1.0651202142497596e-14
Variable bound violation: 0.0000000000000000e+00 0.0000000000000000e+00
Complementarity.........: 2.5137424147587802e-09 2.5137424147587802e-09
Overall NLP error.......: 2.5137424147587802e-09 2.5137424147587802e-09
Number of objective function evaluations = 6
Number of objective gradient evaluations = 6
Number of equality constraint evaluations = 6
Number of inequality constraint evaluations = 0
Number of equality constraint Jacobian evaluations = 6
Number of inequality constraint Jacobian evaluations = 0
Number of Lagrangian Hessian evaluations = 5
Total seconds in IPOPT = 0.005
EXIT: Optimal Solution Found.
