Main

Optimal Control Problems

Main.DynamicOptimizationBenchmarks History

Hide minor edits - Show changes to markup

Added lines 1-4:

(:title Optimal Control Problems:) (:keywords nonlinear control, optimal control, dynamic optimization, engineering optimization, MATLAB, Python, GEKKO, differential, algebraic, modeling language, university course:) (:description Optimal control problems solved with Dynamic Optimization in MATLAB, Excel, and Python.:)

Added lines 446-449:

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

Added lines 367-370:

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

Deleted lines 235-238:

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

Added lines 282-285:

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

Added lines 236-239:

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

Deleted lines 285-288:

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

Added lines 282-285:

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

Added lines 211-214:

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

Changed lines 383-431 from:
  1. Show source
to:

import numpy as np import matplotlib.pyplot as plt from gekko import GEKKO

m = GEKKO()

nt = 101 m.time = np.linspace(0,12,nt)

  1. Parameters

u = m.MV(value=1,ub=1,lb=0) u.STATUS = 1 u.DCOST = 0

  1. Variables

x1 = m.Var(value=1) x2 = m.Var(value=0)

p = np.zeros(nt) p[-1] = 1.0 final = m.Param(value=p)

  1. Equations

m.Equation(x1.dt()==u*(10*x2-x1)) m.Equation(x2.dt()==-u*(10*x2-x1)-(1-u)*x2)

  1. Objective Function

m.Obj(-final*(1-x1-x2))

m.options.IMODE = 6 m.solve()

print('Objective: ' + str(1-x1[-1]-x2[-1]))

plt.figure(1)

plt.subplot(2,1,1) plt.plot(m.time,x1.value,'k:',LineWidth=2,label=r'$x_1$') plt.plot(m.time,x2.value,'b-',LineWidth=2,label=r'$x_2$') plt.ylabel('Value') plt.legend(loc='best')

plt.subplot(2,1,2) plt.plot(m.time,u.value,'r-',LineWidth=2,label=r'$u$') plt.legend(loc='best') plt.xlabel('Time') plt.ylabel('Value')

plt.show()

Changed lines 304-356 from:
  1. Show source
to:

import numpy as np import matplotlib.pyplot as plt from gekko import GEKKO

m = GEKKO()

nt = 101 m.time = np.linspace(0,1,nt)

  1. Parameters

T = m.MV(value=362,ub=398,lb=298) T.STATUS = 1 T.DCOST = 0

  1. Variables

x1 = m.Var(value=1) x2 = m.Var(value=0)

p = np.zeros(nt) p[-1] = 1.0 final = m.Param(value=p)

  1. Intermediates

k1 = m.Intermediate(4000*m.exp(-2500/T)) k2 = m.Intermediate(6.2e5*m.exp(-5000/T))

  1. Equations

m.Equation(x1.dt()==-k1*x1**2) m.Equation(x2.dt()==k1*x1**2 - k2*x2)

  1. Objective Function

m.Obj(-x2*final)

m.options.IMODE = 6 m.solve()

print('Objective: ' + str(x2[-1]))

plt.figure(1)

plt.subplot(2,1,1) plt.plot(m.time,x1.value,'k:',LineWidth=2,label=r'$x_1$') plt.plot(m.time,x2.value,'b-',LineWidth=2,label=r'$x_2$') plt.ylabel('Value') plt.legend(loc='best')

plt.subplot(2,1,2) plt.plot(m.time,T.value,'r--',LineWidth=2,label=r'$T$') plt.legend(loc='best') plt.xlabel('Time') plt.ylabel('Value')

plt.show()

Changed lines 235-275 from:
  1. Show source
to:

import numpy as np import matplotlib.pyplot as plt from gekko import GEKKO

m = GEKKO()

nt = 101 m.time = np.linspace(0,1,nt)

  1. Parameters

u = m.MV(value=1,ub=5,lb=0) u.STATUS = 1

  1. Variables

x1 = m.Var(value=1) x2 = m.Var(value=0)

p = np.zeros(nt) p[-1] = 1.0 final = m.Param(value=p)

  1. Equations

m.Equation(x1.dt()==-(u+0.5*u**2)*x1) m.Equation(x2.dt()==u*x1)

  1. Objective Function

m.Obj(-x2*final)

m.options.IMODE = 6 m.solve()

print('Objective: ' + str(x2[-1]))

plt.figure(1) plt.plot(m.time,x1.value,'k:',LineWidth=2,label=r'$x_1$') plt.plot(m.time,x2.value,'b-',LineWidth=2,label=r'$x_2$') plt.plot(m.time,u.value,'r--',LineWidth=2,label=r'$u$') plt.legend(loc='best') plt.xlabel('Time') plt.ylabel('Value') plt.show()

Changed line 244 from:

$$\min_{T(t)} x_2 \left( t_f \right)$$

to:

$$\max_{T(t)} x_2 \left( t_f \right)$$

Changed line 217 from:

$$\min_{u(t)} x_2 \left( t_f \right)$$

to:

$$\max_{u(t)} x_2 \left( t_f \right)$$

Changed lines 68-70 from:

plt.plot(m.time,x1,'k:',LineWidth=2,label=r'$x_1$') plt.plot(m.time,x2,'b-',LineWidth=2,label=r'$x_2$') plt.plot(m.time,u,'r--',LineWidth=2,label=r'$u$')

to:

plt.plot(m.time,x1.value,'k:',LineWidth=2,label=r'$x_1$') plt.plot(m.time,x2.value,'b-',LineWidth=2,label=r'$x_2$') plt.plot(m.time,u.value,'r--',LineWidth=2,label=r'$u$')

Changed lines 111-113 from:

plt.plot(m.time,x1,'k:',LineWidth=2,label=r'$x_1$') plt.plot(m.time,x2,'b-',LineWidth=2,label=r'$x_2$') plt.plot(m.time,u,'r--',LineWidth=2,label=r'$u$')

to:

plt.plot(m.time,x1.value,'k:',LineWidth=2,label=r'$x_1$') plt.plot(m.time,x2.value,'b-',LineWidth=2,label=r'$x_2$') plt.plot(m.time,u.value,'r--',LineWidth=2,label=r'$u$')

Changed lines 201-204 from:

plt.plot(m.time,x1,'r--',LineWidth=2,label=r'$x_1$') plt.plot(m.time,x2,'g:',LineWidth=2,label=r'$x_2$') plt.plot(m.time,x3,'k-',LineWidth=2,label=r'$x_3$') plt.plot(m.time,x4,'b-',LineWidth=2,label=r'$x_4$')

to:

plt.plot(m.time,x1.value,'r--',LineWidth=2,label=r'$x_1$') plt.plot(m.time,x2.value,'g:',LineWidth=2,label=r'$x_2$') plt.plot(m.time,x3.value,'k-',LineWidth=2,label=r'$x_3$') plt.plot(m.time,x4.value,'b-',LineWidth=2,label=r'$x_4$')

Deleted lines 29-31:

(:toggle hide solution1 button show="Show Solutions for Problem 1a and 1b":) (:div id=solution1:)

Deleted line 123:

(:divend:)

Changed lines 28-32 from:

Solution to Benchmarks 1a and 1b

to:

Solutions to Benchmarks 1a and 1b

(:toggle hide solution1 button show="Show Solutions for Problem 1a and 1b":) (:div id=solution1:)

Added line 127:

(:divend:)

Added lines 119-122:

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

Changed line 268 from:

$$\min_{u(t)} 1 - x_1 \left( t_f \right) - x_2 \left( t_f \right)$$

to:

$$\max_{u(t)} \left(1 - x_1 \left( t_f \right) - x_2 \left( t_f \right) \right)$$

Changed lines 267-275 from:
to:

$$\min_{u(t)} 1 - x_1 \left( t_f \right) - x_2 \left( t_f \right)$$ $$\mathrm{subject \; to}$$ $$\frac{dx_1}{dt}=u \left(10 \, x_2 - x_1 \right)$$ $$\frac{dx_2}{dt}=-u \left(10 \, x_2 - x_1 \right)-\left(1-u\right) x_2$$ $$x(0) = [1 \; 0]^T$$ $$0 \le u \le 1$$ $$t_f=12$$

Changed lines 238-248 from:
to:

$$\min_{T(t)} x_2 \left( t_f \right)$$ $$\mathrm{subject \; to}$$ $$\frac{dx_1}{dt}=-k_1 \, x_1^2$$ $$\frac{dx_2}{dt}=k_1 \, x_1^2 - k_2 \, x_2$$ $$k_1 = 4000 \, \exp{\left(-\frac{2500}{T}\right)}$$ $$k_2 = 6.2e5 \, \exp{\left(-\frac{5000}{T}\right)}$$ $$x(0) = [1 \; 0]^T$$ $$298 \le T \le 398$$ $$t_f=1$$

Changed lines 211-219 from:
to:

$$\min_{u(t)} x_2 \left( t_f \right)$$ $$\mathrm{subject \; to}$$ $$\frac{dx_1}{dt}=-\left(u+0.5u^2\right) x_1$$ $$\frac{dx_2}{dt}=u \, x_1$$ $$x(0) = [1 \; 0]^T$$ $$0 \le u \le 5$$ $$t_f=1$$

February 07, 2018, at 05:43 AM by 10.37.117.169 -
Changed lines 145-203 from:
  1. Show source
to:

import numpy as np import matplotlib.pyplot as plt from gekko import GEKKO

m = GEKKO()

nt = 101 m.time = np.linspace(0,1,nt)

  1. Parameters

u = m.MV(value=9,lb=-4,ub=10) u.STATUS = 1 u.DCOST = 0

  1. Variables

t = m.Var(value=0) x1 = m.Var(value=0) x2 = m.Var(value=-1) x3 = m.Var(value=-np.sqrt(5)) x4 = m.Var(value=0)

p = np.zeros(nt) p[-1] = 1.0 final = m.Param(value=p)

  1. Equations

m.Equation(t.dt()==1) m.Equation(x1.dt()==x2) m.Equation(x2.dt()==-x3*u+16*t-8) m.Equation(x3.dt()==u) m.Equation(x4.dt()==x1**2+x2**2 +0.005*(x2+16*t-8-0.1*x3*(u**2))**2)

  1. Objective Function

m.Obj(x4*final)

m.options.IMODE = 6 m.options.NODES = 4 m.options.MV_TYPE = 1 m.options.SOLVER = 3 m.solve()

print(m.path)

print('Objective = min x4(tf): ' + str(x4[-1]))

plt.figure(1) plt.subplot(2,1,1) plt.plot(m.time,u,'r-',LineWidth=2,label=r'$u$') plt.legend(loc='best') plt.subplot(2,1,2) plt.plot(m.time,x1,'r--',LineWidth=2,label=r'$x_1$') plt.plot(m.time,x2,'g:',LineWidth=2,label=r'$x_2$') plt.plot(m.time,x3,'k-',LineWidth=2,label=r'$x_3$') plt.plot(m.time,x4,'b-',LineWidth=2,label=r'$x_4$') plt.legend(loc='best') plt.xlabel('Time') plt.ylabel('Value') plt.show()

February 07, 2018, at 04:53 AM by 10.37.117.169 -
Changed line 44 from:

nt = 501

to:

nt = 101

Changed lines 124-133 from:
to:

$$\min_{u(t)} x_4 \left( t_f \right)$$ $$\mathrm{subject \; to}$$ $$\frac{dx_1}{dt}=x_2$$ $$\frac{dx_2}{dt}=-x_3 \, u + 16 \, t - 8$$ $$\frac{dx_3}{dt}=u$$ $$\frac{dx_4}{dt}=x_1^2+x_2^2+0.0005 \left(x_2 + 16 \, t -8 -0.1x_3\,u^2\right)^2$$ $$x(0) = [0 \; -1 \; -\sqrt{5} \; 0]^T$$ $$-4 \le u \le 10$$ $$t_f=1$$

February 07, 2018, at 04:48 AM by 10.37.117.169 -
Changed line 14 from:

$$x(0) = [1 0]^T$$

to:

$$x(0) = [1 \; 0]^T$$

Changed lines 17-35 from:

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

to:

Example 1b

  • Nonlinear, unconstrained, minimize final state with terminal constraint

$$\min_{u(t)} x_2 \left( t_f \right)$$ $$\mathrm{subject \; to}$$ $$\frac{dx_1}{dt}=u$$ $$\frac{dx_2}{dt}=x_1^2 + u^2$$ $$x(0) = [1 \; 0]^T$$ $$x_1 \left( t_f \right)=1$$ $$t_f=1$$

Solution to Benchmarks 1a and 1b

(:html:) <iframe width="560" height="315" src="https://www.youtube.com/embed/mmCFF3-6sGg" frameborder="0" allowfullscreen></iframe> (:htmlend:)

(:toggle hide gekko1a button show="Show GEKKO (Python) Code for 1a":)

Changed lines 77-87 from:

Example 1b

  • Nonlinear, unconstrained, minimize final state with terminal constraint

Solution to Benchmarks 1a and 1b

(:html:) <iframe width="560" height="315" src="https://www.youtube.com/embed/mmCFF3-6sGg" frameborder="0" allowfullscreen></iframe> (:htmlend:)

to:

(:toggle hide gekko1b button show="Show GEKKO (Python) Code for 1b":) (:div id=gekko1b:) (:source lang=python:) import numpy as np import matplotlib.pyplot as plt from gekko import GEKKO

m = GEKKO()

nt = 101 m.time = np.linspace(0,1,nt)

  1. Variables

x1 = m.Var(value=1) x2 = m.Var(value=0) u = m.Var(value=-0.48)

p = np.zeros(nt) p[-1] = 1.0 final = m.Param(value=p)

  1. Equations

m.Equation(x1.dt()==u) m.Equation(x2.dt()==x1**2 + u**2) m.Equation(final*(x1-1)==0)

  1. Objective Function

m.Obj(x2*final)

m.options.IMODE = 6 m.solve()

plt.figure(1) plt.plot(m.time,x1,'k:',LineWidth=2,label=r'$x_1$') plt.plot(m.time,x2,'b-',LineWidth=2,label=r'$x_2$') plt.plot(m.time,u,'r--',LineWidth=2,label=r'$u$') plt.legend(loc='best') plt.xlabel('Time') plt.ylabel('Value') plt.show() (:sourceend:) (:divend:)

Added line 125:
Added lines 133-139:

(:toggle hide gekko2 button show="Show GEKKO (Python) Code":) (:div id=gekko2:) (:source lang=python:)

  1. Show source

(:sourceend:) (:divend:)

Added lines 152-158:

(:toggle hide gekko3 button show="Show GEKKO (Python) Code":) (:div id=gekko3:) (:source lang=python:)

  1. Show source

(:sourceend:) (:divend:)

Added lines 171-177:

(:toggle hide gekko4 button show="Show GEKKO (Python) Code":) (:div id=gekko4:) (:source lang=python:)

  1. Show source

(:sourceend:) (:divend:)

Added lines 189-195:

(:toggle hide gekko5 button show="Show GEKKO (Python) Code":) (:div id=gekko5:) (:source lang=python:)

  1. Show source

(:sourceend:) (:divend:)

February 07, 2018, at 04:44 AM by 10.37.117.169 -
Changed lines 9-15 from:
to:

$$\min_{u(t)} x_2 \left( t_f \right)$$ $$\mathrm{subject \; to}$$ $$\frac{dx_1}{dt}=u$$ $$\frac{dx_2}{dt}=x_1^2 + u^2$$ $$x(0) = [1 0]^T$$ $$t_f=1$$

February 06, 2018, at 03:24 PM by 10.37.117.169 -
Changed line 3 from:

Objective: Set up and solve five dynamic optimization benchmark problems2. Create a program to optimize and display the results. Estimated Time (each): 30 minutes

to:

Objective: Set up and solve three of the five dynamic optimization benchmark problems2. Create a program to optimize and display the results. Estimated Time (each): 30 minutes

Added lines 10-52:

(:toggle hide gekko1a button show="Show GEKKO (Python) Code":) (:div id=gekko1a:) (:source lang=python:) import numpy as np import matplotlib.pyplot as plt from gekko import GEKKO

m = GEKKO()

nt = 501 m.time = np.linspace(0,1,nt)

  1. Variables

x1 = m.Var(value=1) x2 = m.Var(value=0) u = m.Var(value=-0.75)

p = np.zeros(nt) p[-1] = 1.0 final = m.Param(value=p)

  1. Equations

m.Equation(x1.dt()==u) m.Equation(x2.dt()==x1**2 + u**2)

  1. Objective Function

m.Obj(x2*final)

m.options.IMODE = 6 m.solve()

plt.figure(1) plt.plot(m.time,x1,'k:',LineWidth=2,label=r'$x_1$') plt.plot(m.time,x2,'b-',LineWidth=2,label=r'$x_2$') plt.plot(m.time,u,'r--',LineWidth=2,label=r'$u$') plt.legend(loc='best') plt.xlabel('Time') plt.ylabel('Value') plt.show() (:sourceend:) (:divend:)

Changed line 3 from:

Objective: Set up and solve several dynamic optimization benchmark problems2. Create a program to optimize and display the results. Estimated Time (each): 30 minutes

to:

Objective: Set up and solve five dynamic optimization benchmark problems2. Create a program to optimize and display the results. Estimated Time (each): 30 minutes

Added lines 64-67:

(:html:) <iframe width="560" height="315" src="https://www.youtube.com/embed/YBlOF9ATaHI" frameborder="0" allowfullscreen></iframe> (:htmlend:)

Added lines 52-55:

(:html:) <iframe width="560" height="315" src="https://www.youtube.com/embed/yFprG0iJQUE" frameborder="0" allowfullscreen></iframe> (:htmlend:)

Added lines 40-43:

(:html:) <iframe width="560" height="315" src="https://www.youtube.com/embed/sYBE3-PVS9g" frameborder="0" allowfullscreen></iframe> (:htmlend:)

Added lines 60-61:
  1. Hedengren, J. D. and Asgharzadeh Shishavan, R., Powell, K.M., and Edgar, T.F., Nonlinear Modeling, Estimation and Predictive Control in APMonitor, Computers and Chemical Engineering, Volume 70, pg. 133–148, 2014. Article
Added lines 28-31:

(:html:) <iframe width="560" height="315" src="https://www.youtube.com/embed/JP2rH7IPYvc" frameborder="0" allowfullscreen></iframe> (:htmlend:)

Added lines 4-5:

Deleted line 0:
Changed line 12 from:

Solution

to:

Solution to Benchmarks 1a and 1b

Changed lines 19-20 from:
to:

Changed line 24 from:

Solution

to:

Solution to Benchmark 2

Added lines 27-28:

Changed line 32 from:

Solution

to:

Solution to Benchmark 3

Added lines 35-36:

Changed line 40 from:

Solution

to:

Solution to Benchmark 4

Added lines 43-44:

Changed line 48 from:

Solution

to:

Solution to Benchmark 5

Added lines 50-51:

Deleted line 1:
Changed lines 6-7 from:
  • Example 1a - Nonlinear, unconstrained, minimize final state
to:

Example 1a

  • Nonlinear, unconstrained, minimize final state
Changed lines 9-10 from:
  • Example 1b - Nonlinear, unconstrained, minimize final state with terminal constraint
to:

Example 1b

  • Nonlinear, unconstrained, minimize final state with terminal constraint
Changed lines 12-20 from:
  • Example 2 - Nonlinear, constrained, minimize final state
  • Example 3 - Tubular reactor with parallel reaction
  • Example 4 - Batch reactor with consecutive reactions A->B->C

Example 5 - Catalytic reactor with A<->B->C

to:
Changed lines 14-16 from:
to:
Added lines 19-43:

Example 2

  • Nonlinear, constrained, minimize final state

Solution

Example 3

  • Tubular reactor with parallel reaction

Solution

Example 4

  • Batch reactor with consecutive reactions A->B->C

Solution

Example 5

  • Catalytic reactor with A<->B->C

Solution

Added lines 1-30:

Exercise

Objective: Set up and solve several dynamic optimization benchmark problems2. Create a program to optimize and display the results. Estimated Time (each): 30 minutes

  • Example 1a - Nonlinear, unconstrained, minimize final state
  • Example 1b - Nonlinear, unconstrained, minimize final state with terminal constraint
  • Example 2 - Nonlinear, constrained, minimize final state
  • Example 3 - Tubular reactor with parallel reaction
  • Example 4 - Batch reactor with consecutive reactions A->B->C

Example 5 - Catalytic reactor with A<->B->C

Solution

(:html:) <iframe width="560" height="315" src="https://www.youtube.com/embed/mmCFF3-6sGg" frameborder="0" allowfullscreen></iframe> (:htmlend:)

References

  1. M. Čižniar, M. Fikar, M.A. Latifi: A MATLAB Package for Dynamic Optimisation of Processes, 7th International Scientific – Technical Conference – Process Control 2006, June 13 – 16, 2006, Kouty nad Desnou, Czech Republic. Article