0

I started with Pyomo and have some specific questions about the package. I'm working with DAE-Toolbox and want to use this toolbox for simulation and parameter estimation. Here my code:

DAE-Model-script(from pyomo-examples):

from pyomo.environ import *
from pyomo.dae import *

model = ConcreteModel()

time_vec = [0,1,2,3,4,5,6,7,8,9,10,11,12,13]
model.t = ContinuousSet(initialize= time_vec)

meas_time_vec = [1,2,3,5]
model.MEAS_t = Set(within=model.t,initialize =meas_time_vec)    # Measurement times, must be subset of t

meas_data_vec ={1:0.264,2:0.594,3: 0.801,5: 0.959}
model.x1_meas = Param(model.MEAS_t, initialize=meas_data_vec)

model.x1 = Var(model.t, initialize =0)
model.x2 = Var(model.t, initialize=1)

model.p1 = Var(bounds=(-1.5,1.5))
model.p2 = 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 = ConstraintList(rule=_init_conditions)

# Alternate way to declare initial conditions
#def _initx1(model):
#   return model.x1[0] == model.p1      
#model.initx1 = Constraint(rule=_initx1)

#def _initx2(model):
#   return model.x2[0] == model.p2
#model.initx2 = Constraint(rule=_initx2)

def _x1dot(model,i):
    return model.x1dot[i] == model.x2[i]
model.x1dotcon = Constraint(model.t, rule=_x1dot)

def _x2dot(model,i):
    return model.x2dot[i] == 1-2*model.x2[i]-model.x1[i]
model.x2dotcon = Constraint(model.t, rule=_x2dot)

def obj(model):
    return sum((model.x1[i]-model.x1_meas[i])**2 for i in model.MEAS_t)
model.obj = Objective(rule=obj)

The run-script:

rom pyomo.environ import *
from pyomo.dae import *
from Parameter_Estimation2 import model
import copy


#model copy
model_sim_window = copy.deepcopy(model)
model_sim_loop = copy.deepcopy(model)



#(1) Parameter Estimation
discretizer = TransformationFactory('dae.collocation')
discretizer.apply_to(model,nfe=8,ncp=5)
solver=SolverFactory('ipopt')
results = solver.solve(model,tee= True)
p_1 = value(model.p1)
p_2 = value(model.p2)
print('p_1',p_1)
print('p_2',p_2)


#(2) Simulation-window
model_sim_window.obj.deactivate()
model_sim_window.t.value = (0,20)
model_sim_window.p1.value = p_1
model_sim_window.p2.value = p_2
sim = Simulator(model_sim_window, package='casadi')
tsim, profiles = sim.simulate(numpoints=100, integrator='cvodes')


#(3)Simulation-loop
model_sim_loop.obj.deactivate()
model_sim_loop.p1.value = p_1
model_sim_loop.p2.value = p_2
sim = Simulator(model_sim_loop, package='casadi')
x0=[0,0]
result =[x0]
t_vec =[0]

for t in range(0,20):

    model_sim_loop.t.value = (t,t+1)
    tsim, profiles = sim.simulate(numpoints=10, integrator='cvodes',initcon=x0)
    result.append(profiles[-1])
    t_vec.append(tsim[-1])
    x0= profiles[-1]


print('result',result)
print('time',t_vec)

Now, the questions:

  1. Is there a way to reuse the model-instance for parameter estimation
    (comment 1) and for simulation (comment 2)? I solved this problem with a
    "ugly" deepcopy.

  2. Has pyomo.dae simulator a step method, that can be used for a step wise
    integration of the dae-system in a loop. I want to the change the inputs (measurements, control signals) between the steps without a reinit of my
    model.I know, that cvodes has such a method.

  3. How can pyomo be used for model predictive control. Are there any
    examples? What is a good starting point?

Thanks & Bye

Hendrixon

Hendrixon
  • 1

2 Answers2

0

1) There is a bug in the Simulator which prevents it from being applied to discretized models so the approach you're using with two models is the way to do it currently (note, you should use model.clone() to copy the model instead of deepcopy). This bug in the Simulator will be fixed in the next Pyomo release.

2) See this section from the Simulator documentation which I think is exactly what you're asking for.

3) See this paper on a framework for nonlinear model predictive control built on Pyomo.

Bethany Nicholson
  • 2,723
  • 1
  • 11
  • 18
0

Thanks for your response!

1) I solved this problem using an abstract model.

2) I have read the section from Simulator documentation and I made some tests. But I want to implement the pyomo.dae Simulator in a "real-time" application. My workflow:

  1. Read measurement data from plant

  2. set inputs of model

  3. do a simulation step

  4. wait until next data comes in

  5. set inputs.....

  6. do a simulation step (start with the states from the last simulation step,no init is needed)

  7. ...and so on

I know that the package pyfmi (from Jmodelica) has such a step method for the Simulator/Integration, they use cvode. Currently I'm using the IDAS-Solver because I have DAE-System. Has the pyomos IDAS interface such a step-by-step integration method?

3) very interesting!

Thanks!

Hendrixon
  • 1