## Level Regulation with MPC

## Main.LevelControl History

Show minor edits - Show changes to output

Changed lines 5-7 from:

to:

[[https://en.wikipedia.org/wiki/Pumped-storage_hydroelectricity|Pumped-storage hydroelectricity]] is a promising method to improve energy dispatch potential with increasing renewable sources such as wind and solar. The following application is a dual reservoir system with pumped water that either enters the upper or lower reservoir.

You are asked to develop a model predictive controller (MPC) to reach a level set point for tank 2 (e.g. 50% full, level=0.5), initially starting from empty tanks. Tune the controller for satisfactory performance with setpoint changes in the level or with disturbances to the valve position ({`\gamma_1`}). Demonstrate the control performance and discuss methods used to achieve satisfactory control during startup from empty and with valve position changes (show both of these).

You are asked to develop a model predictive controller (MPC) to reach a level set point for tank 2 (e.g. 50% full, level=0.5), initially starting from empty tanks. Tune the controller for satisfactory performance with setpoint changes in the level or with disturbances to the valve position ({`\gamma_1`}). Demonstrate the control performance and discuss methods used to achieve satisfactory control during startup from empty and with valve position changes (show both of these).

Added lines 12-15:

(:html:)

<iframe width="560" height="315" src="https://www.youtube.com/embed/YvU3Rd3JBc0" frameborder="0" allow="autoplay; encrypted-media" allowfullscreen></iframe>

(:htmlend:)

Changed line 18 from:

where ''1'' and ''2'' refer to the tanks, {`\tau_~~1`~~} and {`\tau_~~2`~~} are time constants, {`K_~~1`~~} and {`K_~~2`~~} are gains, {`h_1`} and {`h_2`} are liquid level heights, and ''p'' is the pump flow between ''0'' and ''1''.

to:

where ''1'' and ''2'' refer to the tanks, {`\tau_1=18.4`} and {`\tau_2=24.4`} are time constants, {`K_1=1.3`} and {`K_2=1.0`} are gains, {`h_1`} and {`h_2`} are liquid level heights, and ''p'' is the pump flow between ''0'' and ''1''.

Changed lines 48-49 from:

h1 = m.~~SV~~(value=0)

to:

h1 = m.Var(value=0.0)

Changed line 51 from:

h2 = m.CV(value=0)

to:

h2 = m.CV(value=0.0)

Changed lines 55-56 from:

to:

h2.TR_INIT = 1

Changed lines 63-64 from:

m.options.CV_TYPE = ~~1~~

to:

m.options.CV_TYPE = 2

Deleted line 115:

Added line 119:

h2.SP = sp[i]

Changed lines 124-125 from:

to:

Changed line 136 from:

to:

Changed line 138 from:

if (i%5==~~1~~):

to:

if (i%5==3):

Added lines 164-165:

Attach:download.png [[Attach:final_solution3_steps_python.zip|Dual Tank: Step Test for Model Identification (Python/Excel)]]

Deleted lines 172-176:

(:html:)

<iframe width="560" height="315" src="https://www.youtube.com/embed/U7uyj9BaNKg" frameborder="0" allowfullscreen></iframe>

(:htmlend:)

Added lines 10-11:

Attach:download.png [[Attach:level_simulation.zip|Dual Tank Simulation Files (Simulink and Python)]]

Deleted lines 19-20:

Added lines 157-158:

%width=550px%Attach:gravity_tank_mpc.png

Changed line 20 from:

Attach:gravity_tank_mpc.png

to:

%width=550px%Attach:gravity_tank_mpc.png

Changed line 164 from:

Attach:2tank_linear_MPC.png

to:

%width=400px%Attach:2tank_linear_MPC.png

Added lines 159-176:

!!!! Python (APM) Solution

Attach:download.png [[Attach:final_solution3_python_mpc.zip|Dual Tank: Linear MPC Solution (Python)]]

Attach:2tank_linear_MPC.png

!!!! Simulink Solution

Attach:download.png [[Attach:final_solution3_steps_simulink.zip|Dual Tank: Step Test for Model Identification (Simulink/Excel)]]

Attach:download.png [[Attach:final_solution3_steps_python.zip|Dual Tank: Step Test for Model Identification (Python/Excel)]]

(:html:)

<iframe width="560" height="315" src="https://www.youtube.com/embed/U7uyj9BaNKg" frameborder="0" allowfullscreen></iframe>

(:htmlend:)

Attach:download.png [[Attach:final_solution3_linear_mpc.zip|Dual Tank: Linear MPC Solution (MATLAB/Simulink)]]

Changed line 7 from:

%width=~~300px~~%Attach:2tank_level.png

to:

%width=200px%Attach:2tank_level.png

Added lines 7-8:

%width=300px%Attach:2tank_level.png

Changed lines 13-14 from:

{$\tau_{1} ~~* ~~\frac{h_1}{dt} = -h_1 + K_1 \, p$}

{$\tau_{2}~~ *~~ \frac{h_2}{dt} = -h_2 + K_2 \, h_1$}

{$\tau_{2}

to:

{$\tau_{1} \frac{h_1}{dt} = -h_1 + K_1 \, p$}

{$\tau_{2} \frac{h_2}{dt} = -h_2 + K_2 \, h_1$}

{$\tau_{2} \frac{h_2}{dt} = -h_2 + K_2 \, h_1$}

Added lines 1-156:

(:title Level Regulation with MPC:)

(:keywords Python, MATLAB, linear control, model predictive control, level, height, liquid, dynamic programming:)

(:description Design a model predictive controller for level control in a dual tank system. Show that the controller is able to maintain a set point that changes at frequent intervals.:)

Develop a model predictive controller (MPC) to reach a level set point for tank 2 (e.g. 50% full, level=0.5), initially starting from empty tanks. Tune the controller for satisfactory performance with setpoint changes in the level or with disturbances to the valve position ({`\gamma_1`}). Demonstrate the control performance and discuss methods used to achieve satisfactory control during startup from empty and with valve position changes (show both of these).

A pump transports water to tank 1 where it drains to tank 2. Both tanks have an opening at the bottom that allows liquid to flow out. A valve is used to divert pumped liquid into tank 1 (valve position=0), into tank 2 (valve position=1) or fractionally into both (valve position 0-1). The valve position is not available to the controller as a measurement and is considered a disturbance to the process. Only the level of tank 2 (bottom tank) is available as a measured value. The pump rate can be manipulated between 0 and 1. The pump has a maximum rate of change of 0.1 every second.

For control you can use linear or nonlinear model predictive control. For estimation, you can use moving horizon estimation (MHE) or simple bias updating. The simplest option is linear MPC with bias updating. The following is a linear second order model that approximates the level dynamics.

{$\tau_{1} * \frac{h_1}{dt} = -h_1 + K_1 \, p$}

{$\tau_{2} * \frac{h_2}{dt} = -h_2 + K_2 \, h_1$}

where ''1'' and ''2'' refer to the tanks, {`\tau_1`} and {`\tau_2`} are time constants, {`K_1`} and {`K_2`} are gains, {`h_1`} and {`h_2`} are liquid level heights, and ''p'' is the pump flow between ''0'' and ''1''.

!!!! Python (GEKKO Solution)

Attach:gravity_tank_mpc.png

(:toggle hide gekko button show="Show GEKKO (Python) Code":)

(:div id=gekko:)

(:source lang=python:)

import numpy as np

import matplotlib.pyplot as plt

from scipy.integrate import odeint

import csv

from gekko import GEKKO

# create MPC with GEKKO

m = GEKKO()

m.time = [0,1,2,4,8,12,16,20]

# empirical constants

Kp_h1 = 1.3

tau_h1 = 18.4

Kp_h2 = 1

tau_h2 = 24.4

# manipulated variable

p = m.MV(value=0,lb=1e-5,ub=1)

p.STATUS = 1

p.DCOST = 0.01

p.FSTATUS = 0

# unmeasured state

h1 = m.SV(value=0)

# controlled variable

h2 = m.CV(value=0)

h2.STATUS = 1

h2.FSTATUS = 1

h2.TAU = 20

# equations

m.Equation(tau_h1*h1.dt()==-h1 + Kp_h1*p)

m.Equation(tau_h2*h2.dt()==-h2 + Kp_h2*h1)

# options

m.options.IMODE = 6

m.options.CV_TYPE = 1

# simulated system (for measurements)

def tank(levels,t,pump,valve):

h1 = max(1.0e-10,levels[0])

h2 = max(1.0e-10,levels[1])

c1 = 0.08 # inlet valve coefficient

c2 = 0.04 # tank outlet coefficient

dhdt1 = c1 * (1.0-valve) * pump - c2 * np.sqrt(h1)

dhdt2 = c1 * valve * pump + c2 * np.sqrt(h1) - c2 * np.sqrt(h2)

if h1>=1.0 and dhdt1>0.0:

dhdt1 = 0

if h2>=1.0 and dhdt2>0.0:

dhdt2 = 0

dhdt = [dhdt1,dhdt2]

return dhdt

# Initial conditions (levels)

h0 = [0,0]

# Time points to report the solution

tf = 400

t = np.linspace(0,tf,tf+1)

# Set point

sp = np.zeros(tf+1)

sp[5:100] = 0.5

sp[100:200] = 0.8

sp[200:300] = 0.2

sp[300:] = 0.5

# Inputs that can be adjusted

pump = np.zeros(tf+1)

# Disturbance

valve = 0.0

# Record the solution

y = np.zeros((tf+1,2))

y[0,:] = h0

# Create plot

plt.figure(figsize=(10,7))

plt.ion()

plt.show()

# Simulate the tank step test

for i in range(1,tf):

#########################

# MPC ###################

#########################

# measured height

h2.MEAS = y[i,1]

print('h2: ' + str(h2.MEAS))

# set point deadband

h2.SPHI = sp[i]+0.01

h2.SPLO = sp[i]-0.01

# solve MPC

m.solve(disp=False)

# retrieve 1st pump new value

pump[i] = p.NEWVAL

print('pump: ' + str(pump[i]))

#########################

# System ################

#########################

# Specify the pump and valve

inputs = (pump[i],valve)

# Integrate the model

h = odeint(tank,h0,[0,1],inputs)

# Record the result

y[i+1,:] = h[-1,:]

# Reset the initial condition

h0 = h[-1,:]

# update plot every 5 cycles

if (i%5==1):

plt.clf()

plt.subplot(2,1,1)

plt.plot(t[0:i],sp[0:i],'k-')

plt.plot(t[0:i],y[0:i,0],'b-')

plt.plot(t[0:i],y[0:i,1],'r--')

plt.ylabel('Height (m)')

plt.legend(['Set point','Height 1','Height 2'])

plt.subplot(2,1,2)

plt.plot(t[0:i],pump[0:i],'k-')

plt.ylabel('Pump')

plt.xlabel('Time (sec)')

plt.draw()

plt.pause(0.01)

# Construct and save data file

data = np.vstack((t,pump))

data = np.hstack((np.transpose(data),y))

np.savetxt('data.txt',data,delimiter=',')

(:sourceend:)

(:divend:)

(:keywords Python, MATLAB, linear control, model predictive control, level, height, liquid, dynamic programming:)

(:description Design a model predictive controller for level control in a dual tank system. Show that the controller is able to maintain a set point that changes at frequent intervals.:)

Develop a model predictive controller (MPC) to reach a level set point for tank 2 (e.g. 50% full, level=0.5), initially starting from empty tanks. Tune the controller for satisfactory performance with setpoint changes in the level or with disturbances to the valve position ({`\gamma_1`}). Demonstrate the control performance and discuss methods used to achieve satisfactory control during startup from empty and with valve position changes (show both of these).

A pump transports water to tank 1 where it drains to tank 2. Both tanks have an opening at the bottom that allows liquid to flow out. A valve is used to divert pumped liquid into tank 1 (valve position=0), into tank 2 (valve position=1) or fractionally into both (valve position 0-1). The valve position is not available to the controller as a measurement and is considered a disturbance to the process. Only the level of tank 2 (bottom tank) is available as a measured value. The pump rate can be manipulated between 0 and 1. The pump has a maximum rate of change of 0.1 every second.

For control you can use linear or nonlinear model predictive control. For estimation, you can use moving horizon estimation (MHE) or simple bias updating. The simplest option is linear MPC with bias updating. The following is a linear second order model that approximates the level dynamics.

{$\tau_{1} * \frac{h_1}{dt} = -h_1 + K_1 \, p$}

{$\tau_{2} * \frac{h_2}{dt} = -h_2 + K_2 \, h_1$}

where ''1'' and ''2'' refer to the tanks, {`\tau_1`} and {`\tau_2`} are time constants, {`K_1`} and {`K_2`} are gains, {`h_1`} and {`h_2`} are liquid level heights, and ''p'' is the pump flow between ''0'' and ''1''.

!!!! Python (GEKKO Solution)

Attach:gravity_tank_mpc.png

(:toggle hide gekko button show="Show GEKKO (Python) Code":)

(:div id=gekko:)

(:source lang=python:)

import numpy as np

import matplotlib.pyplot as plt

from scipy.integrate import odeint

import csv

from gekko import GEKKO

# create MPC with GEKKO

m = GEKKO()

m.time = [0,1,2,4,8,12,16,20]

# empirical constants

Kp_h1 = 1.3

tau_h1 = 18.4

Kp_h2 = 1

tau_h2 = 24.4

# manipulated variable

p = m.MV(value=0,lb=1e-5,ub=1)

p.STATUS = 1

p.DCOST = 0.01

p.FSTATUS = 0

# unmeasured state

h1 = m.SV(value=0)

# controlled variable

h2 = m.CV(value=0)

h2.STATUS = 1

h2.FSTATUS = 1

h2.TAU = 20

# equations

m.Equation(tau_h1*h1.dt()==-h1 + Kp_h1*p)

m.Equation(tau_h2*h2.dt()==-h2 + Kp_h2*h1)

# options

m.options.IMODE = 6

m.options.CV_TYPE = 1

# simulated system (for measurements)

def tank(levels,t,pump,valve):

h1 = max(1.0e-10,levels[0])

h2 = max(1.0e-10,levels[1])

c1 = 0.08 # inlet valve coefficient

c2 = 0.04 # tank outlet coefficient

dhdt1 = c1 * (1.0-valve) * pump - c2 * np.sqrt(h1)

dhdt2 = c1 * valve * pump + c2 * np.sqrt(h1) - c2 * np.sqrt(h2)

if h1>=1.0 and dhdt1>0.0:

dhdt1 = 0

if h2>=1.0 and dhdt2>0.0:

dhdt2 = 0

dhdt = [dhdt1,dhdt2]

return dhdt

# Initial conditions (levels)

h0 = [0,0]

# Time points to report the solution

tf = 400

t = np.linspace(0,tf,tf+1)

# Set point

sp = np.zeros(tf+1)

sp[5:100] = 0.5

sp[100:200] = 0.8

sp[200:300] = 0.2

sp[300:] = 0.5

# Inputs that can be adjusted

pump = np.zeros(tf+1)

# Disturbance

valve = 0.0

# Record the solution

y = np.zeros((tf+1,2))

y[0,:] = h0

# Create plot

plt.figure(figsize=(10,7))

plt.ion()

plt.show()

# Simulate the tank step test

for i in range(1,tf):

#########################

# MPC ###################

#########################

# measured height

h2.MEAS = y[i,1]

print('h2: ' + str(h2.MEAS))

# set point deadband

h2.SPHI = sp[i]+0.01

h2.SPLO = sp[i]-0.01

# solve MPC

m.solve(disp=False)

# retrieve 1st pump new value

pump[i] = p.NEWVAL

print('pump: ' + str(pump[i]))

#########################

# System ################

#########################

# Specify the pump and valve

inputs = (pump[i],valve)

# Integrate the model

h = odeint(tank,h0,[0,1],inputs)

# Record the result

y[i+1,:] = h[-1,:]

# Reset the initial condition

h0 = h[-1,:]

# update plot every 5 cycles

if (i%5==1):

plt.clf()

plt.subplot(2,1,1)

plt.plot(t[0:i],sp[0:i],'k-')

plt.plot(t[0:i],y[0:i,0],'b-')

plt.plot(t[0:i],y[0:i,1],'r--')

plt.ylabel('Height (m)')

plt.legend(['Set point','Height 1','Height 2'])

plt.subplot(2,1,2)

plt.plot(t[0:i],pump[0:i],'k-')

plt.ylabel('Pump')

plt.xlabel('Time (sec)')

plt.draw()

plt.pause(0.01)

# Construct and save data file

data = np.vstack((t,pump))

data = np.hstack((np.transpose(data),y))

np.savetxt('data.txt',data,delimiter=',')

(:sourceend:)

(:divend:)