## Orthogonal Collocation on Finite Elements

## Main.OrthogonalCollocation History

Hide minor edits - Show changes to output

Changed line 86 from:

%width=30px%Attach:github.png [[https://github.com/APMonitor/dynopt/tree/master/CollocationActivity|~~Both~~ Notebooks on GitHub]]

to:

%width=30px%Attach:github.png [[https://github.com/APMonitor/dynopt/tree/master/CollocationActivity|Collocation Notebooks on GitHub]]

Changed lines 9-14 from:

%width=50px%Attach:colab.png [[https://colab.research.google.com/github/APMonitor/dynopt/blob/master/CollocationActivity/Collocation_Activity.ipynb|Collocation Activity on Google Colab]]

%width=50px%Attach:colab.png[[https://colab.research.google.com/github/APMonitor/dynopt/blob/master/CollocationActivity/Collocation_Activity_Solution.ipynb|Solution~~ on Google Colab]]~~

%width=30px%Attach:github.png [[https://github.com/APMonitor/dynopt/tree/master/CollocationActivity|Both Notebooks on GitHub]]

%width=50px%Attach:colab.png

%width=30px%Attach:github.png [[https://github.com/APMonitor/dynopt/tree/master/CollocationActivity|Both Notebooks on GitHub

to:

%width=50px%Attach:colab.png [[https://colab.research.google.com/github/APMonitor/dynopt/blob/master/CollocationActivity/Collocation_Activity.ipynb|Collocation Activity on Google Colab]] and [[https://colab.research.google.com/github/APMonitor/dynopt/blob/master/CollocationActivity/Collocation_Activity_Solution.ipynb|Solution]]

Added lines 85-86:

%width=30px%Attach:github.png [[https://github.com/APMonitor/dynopt/tree/master/CollocationActivity|Both Notebooks on GitHub]]

Changed lines 9-11 from:

%width=~~30px~~%Attach:colab.png [[https://colab.research.google.com/github/APMonitor/dynopt/blob/master/CollocationActivity/Collocation_Activity.ipynb|Collocation Activity on Google Colab]]

%width=~~30px~~%Attach:colab.png [[https://colab.research.google.com/github/APMonitor/dynopt/blob/master/CollocationActivity/Collocation_Activity_Solution.ipynb|Solution on Google Colab]]

%width=

to:

%width=50px%Attach:colab.png [[https://colab.research.google.com/github/APMonitor/dynopt/blob/master/CollocationActivity/Collocation_Activity.ipynb|Collocation Activity on Google Colab]]

%width=50px%Attach:colab.png [[https://colab.research.google.com/github/APMonitor/dynopt/blob/master/CollocationActivity/Collocation_Activity_Solution.ipynb|Solution on Google Colab]]

%width=50px%Attach:colab.png [[https://colab.research.google.com/github/APMonitor/dynopt/blob/master/CollocationActivity/Collocation_Activity_Solution.ipynb|Solution on Google Colab]]

Changed lines 9-10 from:

to:

%width=30px%Attach:colab.png [[https://colab.research.google.com/github/APMonitor/dynopt/blob/master/CollocationActivity/Collocation_Activity.ipynb|Collocation Activity on Google Colab]]

%width=30px%Attach:colab.png [[https://colab.research.google.com/github/APMonitor/dynopt/blob/master/CollocationActivity/Collocation_Activity_Solution.ipynb|Solution on Google Colab]]

%width=30px%Attach:github.png [[https://github.com/APMonitor/dynopt/tree/master/CollocationActivity|Both Notebooks on GitHub]]

%width=30px%Attach:colab.png [[https://colab.research.google.com/github/APMonitor/dynopt/blob/master/CollocationActivity/Collocation_Activity_Solution.ipynb|Solution on Google Colab]]

%width=30px%Attach:github.png [[https://github.com/APMonitor/dynopt/tree/master/CollocationActivity|Both Notebooks on GitHub]]

Added lines 11-85:

(:toggle hide ex1 button show="Automation Collocation with Python Gekko":)

(:div id=ex1:)

(:source lang=python:)

from gekko import GEKKO

m = GEKKO(remote=False) # create GEKKO model

y = m.Var(5.0,name='y') # create GEKKO variable

m.Equation(y.dt()==-y) # create GEKKO equation

m.time = [0,1,2,3] # time points

m.options.IMODE = 4 # simulation mode

m.options.NODES = 4 # set NODES=4

m.options.CSV_WRITE=2 # write results_all.json

# with internal nodes

m.solve(disp=False) # solve

(:sourceend:)

(:divend:)

(:toggle hide ex2 button show="Manual Collocation with Python Gekko":)

(:div id=ex2:)

(:source lang=python:)

import numpy as np

N = np.array([[0.436,-0.281, 0.121], \

[0.614, 0.064, 0.0461], \

[0.603, 0.230, 0.167]])

time = np.array([0.0, \

0.5-np.sqrt(5)/10.0, \

0.5+np.sqrt(5)/10.0, \

1.0])

from gekko import GEKKO

m = GEKKO()

y0 = 5

y = m.Array(m.Var,3,value=0)

dy = m.Array(m.Var,3,value=0)

Ndy = np.dot(N,dy)

m.Equations([Ndy[i]==(y[i]-y0) for i in range(3)])

m.Equations([dy[i]+y[i]==0 for i in range(3)])

m.solve(disp=False)

print(y)

(:sourceend:)

(:divend:)

(:toggle hide ex3 button show="Manual Collocation with Python Scipy":)

(:div id=ex3:)

(:source lang=python:)

import numpy as np

N = np.array([[0.436,-0.281, 0.121], \

[0.614, 0.064, 0.0461], \

[0.603, 0.230, 0.167]])

time = np.array([0.0, \

0.5-np.sqrt(5)/10.0, \

0.5+np.sqrt(5)/10.0, \

1.0])

from scipy.optimize import fsolve

y0 = 5

zGuess = np.zeros(6) # y1-y3 and dydt1-dydt3

def myFunction(z):

y = z[0:3]

dy = z[3:6]

F = np.empty(6)

F[0:3] = np.dot(N,dy) - (y-y0)

F[3:7] = dy + y

return F

z = fsolve(myFunction,zGuess)

print(z)

(:sourceend:)

(:divend:)

Added lines 9-12:

(:html:)

<iframe width="560" height="315" src="https://www.youtube.com/embed/Avwn53krkgY" frameborder="0" allow="accelerometer; autoplay; clipboard-write; encrypted-media; gyroscope; picture-in-picture" allowfullscreen></iframe>

(:htmlend:)

Changed lines 8-9 from:

* [[https://colab.research.google.com/github/APMonitor/dynopt/blob/master/CollocationActivity/Collocation_Activity_Solution.ipynb|Solution on Google Colab]]

~~* ~~[[https://github.com/APMonitor/dynopt/tree/master/CollocationActivity|Notebooks on GitHub]]

to:

* [[https://colab.research.google.com/github/APMonitor/dynopt/blob/master/CollocationActivity/Collocation_Activity_Solution.ipynb|Solution on Google Colab]] with [[https://github.com/APMonitor/dynopt/tree/master/CollocationActivity|Both Notebooks on GitHub]]

Changed lines 8-9 from:

* [[https://colab.research.google.com/github/APMonitor/dynopt/blob/master/CollocationActivity/Collocation_Activity_Solution.ipynb|~~Collocation Activity Solution on Google ~~Colab]]

* [[https://github.com/APMonitor/dynopt/tree/master/CollocationActivity|~~Collocation Activity ~~Notebooks on GitHub]]

* [[https://github.com/APMonitor/dynopt/tree/master/CollocationActivity|

to:

* [[https://colab.research.google.com/github/APMonitor/dynopt/blob/master/CollocationActivity/Collocation_Activity_Solution.ipynb|Solution on Google Colab]]

* [[https://github.com/APMonitor/dynopt/tree/master/CollocationActivity|Notebooks on GitHub]]

* [[https://github.com/APMonitor/dynopt/tree/master/CollocationActivity|Notebooks on GitHub]]

Changed line 8 from:

* [[https://colab.research.google.com/github/APMonitor/dynopt/blob/master/CollocationActivity/Collocation_Activity_Solution.ipynb|Collocation Activity on Google Colab]]

to:

* [[https://colab.research.google.com/github/APMonitor/dynopt/blob/master/CollocationActivity/Collocation_Activity_Solution.ipynb|Collocation Activity Solution on Google Colab]]

Changed lines 7-9 from:

* [[https://colab.research.google.com/github/APMonitor/dynopt/blob/master/CollocationActivity/Collocation_Activity.ipynb|Collocation Activity on Google Colab]] ~~with ~~[[https://colab.research.google.com/github/APMonitor/dynopt/blob/master/CollocationActivity/Collocation_Activity_Solution.ipynb|Collocation Activity on Google Colab]] ~~or ~~[[https://github.com/APMonitor/dynopt/tree/master/CollocationActivity|~~source~~]]~~.~~

to:

* [[https://colab.research.google.com/github/APMonitor/dynopt/blob/master/CollocationActivity/Collocation_Activity.ipynb|Collocation Activity on Google Colab]]

* [[https://colab.research.google.com/github/APMonitor/dynopt/blob/master/CollocationActivity/Collocation_Activity_Solution.ipynb|Collocation Activity on Google Colab]]

* [[https://github.com/APMonitor/dynopt/tree/master/CollocationActivity|Collocation Activity Notebooks on GitHub]]

* [[https://colab.research.google.com/github/APMonitor/dynopt/blob/master/CollocationActivity/Collocation_Activity_Solution.ipynb|Collocation Activity on Google Colab]]

* [[https://github.com/APMonitor/dynopt/tree/master/CollocationActivity|Collocation Activity Notebooks on GitHub]]

Changed lines 5-9 from:

Discretization of a continuous time representation allow large-scale nonlinear programming (NLP) solvers to find solutions at specified intervals in a time horizon. There are many names and related techniques for obtaining mathematical relationships between derivatives and non-derivative values. Some of the terms that are relevant to this discussion include orthogonal collocation on finite elements, direct transcription, Gauss pseudospectral method, Gaussian quadrature, Lobatto quadrature, Radau collocation, Legendre polynomials, Chebyshev polynomials, Jacobi polynomials, Laguerre polynomials, any many more. There are many papers that discuss the details of the derivation and theory behind these methods'^1-5^'. The purpose of this section is to give a practical introduction to orthogonal collocation on finite elements with Lobatto quadrature for the numerical solution of differential algebraic equations. See the documentation on [[https://apmonitor.com/wiki/index.php/Main/OptionApmNodes|Nodes]] for additional details on displaying the internal nodes.

to:

Discretization of a continuous time representation allow large-scale nonlinear programming (NLP) solvers to find solutions at specified intervals in a time horizon.

* [[https://colab.research.google.com/github/APMonitor/dynopt/blob/master/CollocationActivity/Collocation_Activity.ipynb|Collocation Activity on Google Colab]] with [[https://colab.research.google.com/github/APMonitor/dynopt/blob/master/CollocationActivity/Collocation_Activity_Solution.ipynb|Collocation Activity on Google Colab]] or [[https://github.com/APMonitor/dynopt/tree/master/CollocationActivity|source]].

There are many names and related techniques for obtaining mathematical relationships between derivatives and non-derivative values. Some of the terms that are relevant to this discussion include orthogonal collocation on finite elements, direct transcription, Gauss pseudospectral method, Gaussian quadrature, Lobatto quadrature, Radau collocation, Legendre polynomials, Chebyshev polynomials, Jacobi polynomials, Laguerre polynomials, any many more. There are many papers that discuss the details of the derivation and theory behind these methods'^1-5^'. The purpose of this section is to give a practical introduction to orthogonal collocation on finite elements with Lobatto quadrature for the numerical solution of differential algebraic equations. See the documentation on [[https://apmonitor.com/wiki/index.php/Main/OptionApmNodes|Nodes]] for additional details on displaying the internal nodes.

* [[https://colab.research.google.com/github/APMonitor/dynopt/blob/master/CollocationActivity/Collocation_Activity.ipynb|Collocation Activity on Google Colab]] with [[https://colab.research.google.com/github/APMonitor/dynopt/blob/master/CollocationActivity/Collocation_Activity_Solution.ipynb|Collocation Activity on Google Colab]] or [[https://github.com/APMonitor/dynopt/tree/master/CollocationActivity|source]].

There are many names and related techniques for obtaining mathematical relationships between derivatives and non-derivative values. Some of the terms that are relevant to this discussion include orthogonal collocation on finite elements, direct transcription, Gauss pseudospectral method, Gaussian quadrature, Lobatto quadrature, Radau collocation, Legendre polynomials, Chebyshev polynomials, Jacobi polynomials, Laguerre polynomials, any many more. There are many papers that discuss the details of the derivation and theory behind these methods'^1-5^'. The purpose of this section is to give a practical introduction to orthogonal collocation on finite elements with Lobatto quadrature for the numerical solution of differential algebraic equations. See the documentation on [[https://apmonitor.com/wiki/index.php/Main/OptionApmNodes|Nodes]] for additional details on displaying the internal nodes.

Added lines 369-372:

----

Attach:download.png [[https://colab.research.google.com/github/APMonitor/dynopt/blob/master/OrthogonalCollocation.ipynb|Collocation Solutions on Google Colab]]

Changed line 5 from:

Discretization of a continuous time representation allow large-scale nonlinear programming (NLP) solvers to find solutions at specified intervals in a time horizon. There are many names and related techniques for obtaining mathematical relationships between derivatives and non-derivative values. Some of the terms that are relevant to this discussion include orthogonal collocation on finite elements, direct transcription, Gauss pseudospectral method, Gaussian quadrature, Lobatto quadrature, Radau collocation, Legendre polynomials, Chebyshev polynomials, Jacobi polynomials, Laguerre polynomials, any many more. There are many papers that discuss the details of the derivation and theory behind these methods'^1-5^'. The purpose of this section is to give a practical introduction to orthogonal collocation on finite elements with Lobatto quadrature for the numerical solution of differential algebraic equations.

to:

Discretization of a continuous time representation allow large-scale nonlinear programming (NLP) solvers to find solutions at specified intervals in a time horizon. There are many names and related techniques for obtaining mathematical relationships between derivatives and non-derivative values. Some of the terms that are relevant to this discussion include orthogonal collocation on finite elements, direct transcription, Gauss pseudospectral method, Gaussian quadrature, Lobatto quadrature, Radau collocation, Legendre polynomials, Chebyshev polynomials, Jacobi polynomials, Laguerre polynomials, any many more. There are many papers that discuss the details of the derivation and theory behind these methods'^1-5^'. The purpose of this section is to give a practical introduction to orthogonal collocation on finite elements with Lobatto quadrature for the numerical solution of differential algebraic equations. See the documentation on [[https://apmonitor.com/wiki/index.php/Main/OptionApmNodes|Nodes]] for additional details on displaying the internal nodes.

Changed line 142 from:

plt.plot(tc(n),yc,'o',markersize=~~20~~)

to:

plt.plot(tc(n),yc,'o',markersize=10)

Changed lines 151-163 from:

from APMonitor import *

except:

# Automatically install APMonitor

import pip

pip.main(['install','APMonitor'])

try:

from APMonitor import *

except:

print('pip installation unsuccessful')

print('Visit apmonitor.com to download apm.py')

from

to:

from APMonitor.apm import *

Deleted line 187:

Changed line 35 from:

(:toggle hide python button show="APM Python Source":)

to:

(:toggle hide python button show="Show APM Python Source":)

Deleted lines 330-339:

u = m.Param(value=4)

x = m.Var(value=0)

m.Equation(5*x.dt() == -x**2 + u)

m.time = [0,tf]

m.options.imode = 4

Changed lines 334-343 from:

to:

m = GEKKO()

u = m.Param(value=4)

x = m.Var(value=0)

m.Equation(5*x.dt() == -x**2 + u)

m.time = [0,tf]

m.options.imode = 4

m.options.time_shift = 0

u = m.Param(value=4)

x = m.Var(value=0)

m.Equation(5*x.dt() == -x**2 + u)

m.time = [0,tf]

m.options.imode = 4

m.options.time_shift = 0

Added line 345:

Changed line 319 from:

plt.plot(tc(n),yc,'o',markersize=~~20~~)

to:

plt.plot(tc(n),yc,'o',markersize=10,label='Nodes = '+str(i))

Changed lines 322-323 from:

to:

plt.legend(loc='best')

Changed line 344 from:

m.clear()#apm(s,a,'clear all') # clear prior application

to:

#m.clear()#apm(s,a,'clear all') # clear prior application

Added lines 6-9:

(:html:)

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

(:htmlend:)

Deleted lines 26-27:

(:div id=solution1:)

Changed lines 359-360 from:

to:

Deleted lines 367-368:

(:div id=solution2:)

Deleted line 375:

Added lines 27-28:

(:toggle hide solution1 button show="Show Solution 1":)

(:div id=solution1:)

(:div id=solution1:)

Deleted lines 208-211:

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

(:htmlend:)

Added lines 358-362:

(:html:)

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

(:htmlend:)

(:divend:)

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

(:htmlend:)

(:divend:)

Added lines 371-372:

(:toggle hide solution2 button show="Show Solution 2":)

(:div id=solution2:)

(:div id=solution2:)

Added line 381:

(:divend:)

Changed line 31 from:

(:toggle hide python button show="Python Source":)

to:

(:toggle hide python button show="APM Python Source":)

Added lines 210-358:

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

(:div id=gekko:)

(:source lang=python:)

from __future__ import division

import numpy as np

from scipy.optimize import fsolve

from scipy.integrate import odeint

import matplotlib.pyplot as plt

# final time

tf = 1.0

# solve with ODEINT (for comparison)

def model(x,t):

u = 4.0

return (-x**2 + u)/5.0

t = np.linspace(0,tf,20)

y0 = 0

y = odeint(model,y0,t)

plt.figure(1)

plt.plot(t,y,'r-',label='ODEINT')

# ----------------------------------------------------

# Approach #1 - Write the model equations in Python

# ----------------------------------------------------

# define collocation matrices

def colloc(n):

if (n==2):

NC = np.array([[1.0]])

if (n==3):

NC = np.array([[0.75,-0.25], \

[1.00, 0.00]])

if (n==4):

NC = np.array([[0.436,-0.281, 0.121], \

[0.614, 0.064, 0.0461], \

[0.603, 0.230, 0.167]])

if (n==5):

NC = np.array([[0.278, -0.202, 0.169, -0.071], \

[0.398, 0.069, 0.064, -0.031], \

[0.387, 0.234, 0.278, -0.071], \

[0.389, 0.222, 0.389, 0.000]])

if (n==6):

NC = np.array([[0.191, -0.147, 0.139, -0.113, 0.047],

[0.276, 0.059, 0.051, -0.050, 0.022],

[0.267, 0.193, 0.252, -0.114, 0.045],

[0.269, 0.178, 0.384, 0.032, 0.019],

[0.269, 0.181, 0.374, 0.110, 0.067]])

return NC

# define collocation points from Lobatto quadrature

def tc(n):

if (n==2):

time = np.array([0.0,1.0])

if (n==3):

time = np.array([0.0,0.5,1.0])

if (n==4):

time = np.array([0.0, \

0.5-np.sqrt(5)/10.0, \

0.5+np.sqrt(5)/10.0, \

1.0])

if (n==5):

time = np.array([0.0,0.5-np.sqrt(21)/14.0, \

0.5,0.5+np.sqrt(21)/14.0, 1])

if (n==6):

time = np.array([0.0, \

0.5-np.sqrt((7.0+2.0*np.sqrt(7.0))/21.0)/2.0, \

0.5-np.sqrt((7.0-2.0*np.sqrt(7.0))/21.0)/2.0, \

0.5+np.sqrt((7.0-2.0*np.sqrt(7.0))/21.0)/2.0, \

0.5+np.sqrt((7.0+2.0*np.sqrt(7.0))/21.0)/2.0, \

1.0])

return time*tf

# solve with SciPy fsolve

def myFunction(z,*param):

n = param[0]

m = param[1]

# rename z as x and xdot variables

x = np.empty(n-1)

xdot = np.empty(n-1)

x[0:n-1] = z[0:n-1]

xdot[0:n-1] = z[n-1:m]

# initial condition (x0)

x0 = 0.0

# input parameter (u)

u = 4.0

# final time

tn = tf

# function evaluation residuals

F = np.empty(m)

# nonlinear differential equations at each node

# 5 dx/dt = -x^2 + u

F[0:n-1] = 5.0 * xdot[0:n-1] + x[0:n-1]**2 - u

# collocation equations

# tn * NC * xdot = x - x0

NC = colloc(n)

F[n-1:m] = tn * np.dot(NC,xdot) - x + x0 * np.ones(n-1)

return F

sol_py = np.empty(5) # store 5 results

for i in range(2,7):

n = i

m = (i-1)*2

zGuess = np.ones(m)

z = fsolve(myFunction,zGuess,args=(n,m))

# add to plot

yc = np.insert(z[0:n-1],0,0)

plt.plot(tc(n),yc,'o',markersize=20)

# store just the last x[n] value

sol_py[i-2] = z[n-2]

# ----------------------------------------------------

# Approach #2 - Write model in APMonitor and let

# modeling language create the collocation equations

# ----------------------------------------------------

# load GEKKO

from gekko import GEKKO

m = GEKKO()

u = m.Param(value=4)

x = m.Var(value=0)

m.Equation(5*x.dt() == -x**2 + u)

m.time = [0,tf]

m.options.imode = 4

sol_apm = np.empty(5) # store 5 results

i = 0

for nodes in range(2,7):

m.clear()#apm(s,a,'clear all') # clear prior application

m.options.nodes = nodes

m.solve() # solve problem

sol_apm[i] = x.value[-1] # store solution (last point)

i += 1

# print the solutions

print(sol_py)

print(sol_apm)

# show plot

plt.ylabel('x(t)')

plt.xlabel('time')

plt.show()

(:sourceend:)

(:divend:)

Changed lines 31-35 from:

(:~~html:)~~

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

(:htmlend:)

(:toggle hide python button show="Buy or Build TCLab":)

<iframe width="560" height

(:htmlend:)

(:toggle hide python button show="Buy or Build TCLab

to:

(:toggle hide python button show="Python Source":)

Added lines 206-209:

(:html:)

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

(:htmlend:)

Added lines 29-36:

Attach:download.png [[Attach:collocation_exercise.zip|Collocation in MATLAB and Python (APM Solution)]]

(:html:)

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

(:htmlend:)

(:toggle hide python button show="Buy or Build TCLab":)

(:div id=python:)

(:html:)

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

(:htmlend:)

(:toggle hide python button show="Buy or Build TCLab":)

(:div id=python:)

Changed lines 209-214 from:

(:html:)

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

(:htmlend

to:

(:divend:)

Changed line 19 from:

'''Objective:''' Solve a differential equation with orthogonal collocation on finite elements. Create a MATLAB or Python script to simulate and display the results. ''Estimated Time: ~~1 hour~~''

to:

'''Objective:''' Solve a differential equation with orthogonal collocation on finite elements. Create a MATLAB or Python script to simulate and display the results. ''Estimated Time: 2-3 hours''

Added lines 28-200:

(:source lang=python:)

import numpy as np

from scipy.optimize import fsolve

from scipy.integrate import odeint

import matplotlib.pyplot as plt

# final time

tf = 1.0

# solve with ODEINT (for comparison)

def model(x,t):

u = 4.0

return (-x**2 + u)/5.0

t = np.linspace(0,tf,20)

y0 = 0

y = odeint(model,y0,t)

plt.figure(1)

plt.plot(t,y,'r-',label='ODEINT')

# ----------------------------------------------------

# Approach #1 - Write the model equations in Python

# ----------------------------------------------------

# define collocation matrices

def colloc(n):

if (n==2):

NC = np.array([[1.0]])

if (n==3):

NC = np.array([[0.75,-0.25], \

[1.00, 0.00]])

if (n==4):

NC = np.array([[0.436,-0.281, 0.121], \

[0.614, 0.064, 0.0461], \

[0.603, 0.230, 0.167]])

if (n==5):

NC = np.array([[0.278, -0.202, 0.169, -0.071], \

[0.398, 0.069, 0.064, -0.031], \

[0.387, 0.234, 0.278, -0.071], \

[0.389, 0.222, 0.389, 0.000]])

if (n==6):

NC = np.array([[0.191, -0.147, 0.139, -0.113, 0.047],

[0.276, 0.059, 0.051, -0.050, 0.022],

[0.267, 0.193, 0.252, -0.114, 0.045],

[0.269, 0.178, 0.384, 0.032, 0.019],

[0.269, 0.181, 0.374, 0.110, 0.067]])

return NC

# define collocation points from Lobatto quadrature

def tc(n):

if (n==2):

time = np.array([0.0,1.0])

if (n==3):

time = np.array([0.0,0.5,1.0])

if (n==4):

time = np.array([0.0, \

0.5-np.sqrt(5)/10.0, \

0.5+np.sqrt(5)/10.0, \

1.0])

if (n==5):

time = np.array([0.0,0.5-np.sqrt(21)/14.0, \

0.5,0.5+np.sqrt(21)/14.0, 1])

if (n==6):

time = np.array([0.0, \

0.5-np.sqrt((7.0+2.0*np.sqrt(7.0))/21.0)/2.0, \

0.5-np.sqrt((7.0-2.0*np.sqrt(7.0))/21.0)/2.0, \

0.5+np.sqrt((7.0-2.0*np.sqrt(7.0))/21.0)/2.0, \

0.5+np.sqrt((7.0+2.0*np.sqrt(7.0))/21.0)/2.0, \

1.0])

return time*tf

# solve with SciPy fsolve

def myFunction(z,*param):

n = param[0]

m = param[1]

# rename z as x and xdot variables

x = np.empty(n-1)

xdot = np.empty(n-1)

x[0:n-1] = z[0:n-1]

xdot[0:n-1] = z[n-1:m]

# initial condition (x0)

x0 = 0.0

# input parameter (u)

u = 4.0

# final time

tn = tf

# function evaluation residuals

F = np.empty(m)

# nonlinear differential equations at each node

# 5 dx/dt = -x^2 + u

F[0:n-1] = 5.0 * xdot[0:n-1] + x[0:n-1]**2 - u

# collocation equations

# tn * NC * xdot = x - x0

NC = colloc(n)

F[n-1:m] = tn * np.dot(NC,xdot) - x + x0 * np.ones(n-1)

return F

sol_py = np.empty(5) # store 5 results

for i in range(2,7):

n = i

m = (i-1)*2

zGuess = np.ones(m)

z = fsolve(myFunction,zGuess,args=(n,m))

# add to plot

yc = np.insert(z[0:n-1],0,0)

plt.plot(tc(n),yc,'o',markersize=20)

# store just the last x[n] value

sol_py[i-2] = z[n-2]

# ----------------------------------------------------

# Approach #2 - Write model in APMonitor and let

# modeling language create the collocation equations

# ----------------------------------------------------

# load APM package

try:

from APMonitor import *

except:

# Automatically install APMonitor

import pip

pip.main(['install','APMonitor'])

try:

from APMonitor import *

except:

print('pip installation unsuccessful')

print('Visit apmonitor.com to download apm.py')

from apm import *

s = 'https://byu.apmonitor.com'

a = 'node_test'

fid = open('collocation.apm','w')

# write model

fid.write('Parameters \n')

fid.write(' u = 4 \n')

fid.write('Variables \n')

fid.write(' x = 0 \n')

fid.write('Equations \n')

fid.write(' 5 * $x = -x^2 + u \n')

# close file

fid.close()

fid = open('collocation.csv','w')

# write data file

fid.write('time \n')

fid.write('0 \n')

fid.write(str(tf)+'\n')

# close file

fid.close()

sol_apm = np.empty(5) # store 5 results

i = 0

for nodes in range(2,7):

apm(s,a,'clear all') # clear prior application

apm_load(s,a,'collocation.apm') # load model

csv_load(s,a,'collocation.csv') # load data

apm_option(s,a,'nlc.imode',4) # imode = 4, simulation

apm_option(s,a,'nlc.nodes',nodes) # nodes (2-6)

apm(s,a,'solve') # solve problem

y = apm_sol(s,a) # retrieve solution

sol_apm[i] = y['x'][-1] # store solution (last point)

i += 1

# print the solutions

print(sol_py)

print(sol_apm)

# show plot

plt.ylabel('x(t)')

plt.xlabel('time')

plt.show()

(:sourceend:)

Added lines 6-7:

Attach:download.png [[Attach:orthogonal_collocation_python.zip|Orthogonal Collocation in Python]]

Changed lines 37-39 from:

Compare orthogonal collocation on finite elements with 3 nodes with a numerical integrator (e.g. ODE15s in MATLAB or ODEINT in Python). Calculate the error at each of the solution points~~.~~

to:

Compare orthogonal collocation on finite elements with 3 nodes with a numerical integrator (e.g. ODE15s in MATLAB or ODEINT in Python). Calculate the error at each of the solution points for the equation (same as for Exercise 1):

5 dx/dt = -x'^2^' + u

5 dx/dt = -x'^2^' + u

Added lines 43-46:

(:html:)

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

(:htmlend:)

Changed line 41 from:

Attach:download.png [[Attach:collocation_compare.zip|Collocation Solution (APM, ODE15s (MATLAB), ODEINT (Python)~~ Compared~~)]]

to:

Attach:download.png [[Attach:collocation_compare.zip|Collocation Solution (APM, ODE15s (MATLAB), ODEINT (Python))]]

Added lines 13-14:

----

Added lines 33-34:

----

Added lines 42-43:

----

Changed lines 13-14 from:

!!!! Exercise

to:

!!!! Exercise 1

Changed lines 23-26 from:

!!!! Solution

Attach:download.png [[Attach:collocation_exercise.zip|Collocation~~Exercise ~~in MATLAB and Python]]

Attach:download.png [[Attach:collocation_exercise.zip|Collocation

to:

!!!! Solution 1

Attach:download.png [[Attach:collocation_exercise.zip|Collocation in MATLAB and Python (APM Solution)]]

Attach:download.png [[Attach:collocation_exercise.zip|Collocation in MATLAB and Python (APM Solution)]]

Added lines 30-37:

!!!! Exercise 2

Compare orthogonal collocation on finite elements with 3 nodes with a numerical integrator (e.g. ODE15s in MATLAB or ODEINT in Python). Calculate the error at each of the solution points.

!!!! Solution 2

Attach:download.png [[Attach:collocation_compare.zip|Collocation Solution (APM, ODE15s (MATLAB), ODEINT (Python) Compared)]]

Added lines 42-44:

# Kameswaran, S., Biegler, L.T., Convergence rates for direct transcription of optimal control problems using collocation at Radau points, Computational Optimization and Applications, 41 (1), pp. 81-126, 2008, DOI: 10.1007/s10589-007-9098-9.

Changed lines 5-19 from:

Discretization of a continuous time representation allow large-scale nonlinear programming (NLP) solvers to find solutions at specified intervals in a time horizon. There are many names and related techniques for obtaining mathematical relationships between derivatives and non-derivative values. Some of the terms that are relevant to this discussion include~~:~~

* Orthogonal collocation on finite elements

~~* Direct~~ transcription

~~*~~ Gauss pseudospectral method

~~*~~ Gaussian quadrature

~~*~~ Lobatto quadrature

~~*~~ Radau collocation

~~*~~ Legendre polynomials

~~*~~ Chebyshev polynomials

~~*~~ Jacobi polynomials

~~*~~ Laguerre polynomials

~~* Any~~ many more.~~..~~

There are many papers that discuss the details of the derivation and theory behind these methods'^1-5^'. The purpose of this section is to give a practical introduction to orthogonal collocation on finite elements with Lobatto quadrature for the numerical solution of differential algebraic equations.

* Orthogonal

There are many papers that discuss the details of the derivation and theory behind these methods'^1-5^'. The purpose of this section is to give a practical introduction to orthogonal collocation on finite elements with Lobatto quadrature for the numerical solution of differential algebraic equations.

to:

Discretization of a continuous time representation allow large-scale nonlinear programming (NLP) solvers to find solutions at specified intervals in a time horizon. There are many names and related techniques for obtaining mathematical relationships between derivatives and non-derivative values. Some of the terms that are relevant to this discussion include orthogonal collocation on finite elements, direct transcription, Gauss pseudospectral method, Gaussian quadrature, Lobatto quadrature, Radau collocation, Legendre polynomials, Chebyshev polynomials, Jacobi polynomials, Laguerre polynomials, any many more. There are many papers that discuss the details of the derivation and theory behind these methods'^1-5^'. The purpose of this section is to give a practical introduction to orthogonal collocation on finite elements with Lobatto quadrature for the numerical solution of differential algebraic equations.

Changed line 42 from:

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

to:

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

Changed line 39 from:

Attach:download.png [[Attach:collocation_exercise.zip|Collocation Exercise in MATLAB]]

to:

Attach:download.png [[Attach:collocation_exercise.zip|Collocation Exercise in MATLAB and Python]]

Added lines 25-42:

(:htmlend:)

!!!! Exercise

'''Objective:''' Solve a differential equation with orthogonal collocation on finite elements. Create a MATLAB or Python script to simulate and display the results. ''Estimated Time: 1 hour''

Solve the following differential equation from time 0 to 1 with orthogonal collocation on finite elements with 4 nodes for discretization in time.

5 dx/dt = -x'^2^' + u

Specify the initial condition for ''x'' as 0 and the value of the input, ''u'', as 4. Compare the solution result with 2-6 time points (nodes). Report the solution at the final time for each and comment on how the solution changes with an increase in the number of nodes.

!!!! Solution

Attach:download.png [[Attach:collocation_exercise.zip|Collocation Exercise in MATLAB]]

(:html:)

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

!!!! Exercise

'''Objective:''' Solve a differential equation with orthogonal collocation on finite elements. Create a MATLAB or Python script to simulate and display the results. ''Estimated Time: 1 hour''

Solve the following differential equation from time 0 to 1 with orthogonal collocation on finite elements with 4 nodes for discretization in time.

5 dx/dt = -x'^2^' + u

Specify the initial condition for ''x'' as 0 and the value of the input, ''u'', as 4. Compare the solution result with 2-6 time points (nodes). Report the solution at the final time for each and comment on how the solution changes with an increase in the number of nodes.

!!!! Solution

Attach:download.png [[Attach:collocation_exercise.zip|Collocation Exercise in MATLAB]]

(:html:)

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

Changed lines 19-20 from:

There are many papers that discuss the details of the derivation and theory behind these methods'^1^'. The purpose of this section is to give a practical introduction to orthogonal collocation on finite elements with Lobatto quadrature for the numerical solution of differential algebraic equations.

to:

There are many papers that discuss the details of the derivation and theory behind these methods'^1-5^'. The purpose of this section is to give a practical introduction to orthogonal collocation on finite elements with Lobatto quadrature for the numerical solution of differential algebraic equations.

Changed lines 29-38 from:

# ~~Nonlinear Modeling~~, ~~Estimation~~ and ~~Predictive Control in APMonitor~~, ~~Hedengren, J~~.~~ D~~. and ~~Asgharzadeh Shishavan~~, ~~R~~., ~~Powell, K.M.~~, and ~~Edgar, T.F., Computers and Chemical~~ Engineering, Volume 70, pg. 133–148, 2014. Available at: [[https://dx.doi.org/10.1016/j.compchemeng.2014.04.013|ScienceDirect]] or [[https://www.researchgate.net/publication/262191158_Nonlinear_Modeling_Estimation_and_Predictive_Control_in_APMonitor|ResearchGate]]

to:

# Hedengren, J. D. and Asgharzadeh Shishavan, R., Powell, K.M., and Edgar, T.F., Nonlinear Modeling, Estimation and Predictive Control in APMonitor, Computers & Chemical Engineering, Volume 70, pg. 133–148, 2014. Available at: [[https://dx.doi.org/10.1016/j.compchemeng.2014.04.013|ScienceDirect]] or [[https://www.researchgate.net/publication/262191158_Nonlinear_Modeling_Estimation_and_Predictive_Control_in_APMonitor|ResearchGate]]

# Cizniar, M., Salhi, D., Fikar, M., and Latifi, M.A., A Matlab Package for Orthogonal Collocations on Finite Elements in Dynamic Optimization, 15th Int. Conference Process Control, June 7-10, 2015. [[https://www.kirp.chtf.stuba.sk/~fikar/research/dynopt/fik05oc.pdf | Article]]

# Renfro, J., Morshedi, A., Asbjornsen, O., Simultaneous Optimization and Solution of Systems Described by Differential/Algebraic Equations, Computers & Chemical Engineering, Volume 11, Number 5, pp. 503-517, 1987.

# Finlayson, B., [[https://faculty.washington.edu/finlayso/ebook/pde/OC/OCFE.htm|Orthogonal Collocation on Finite Elements for PDE Discretization and Solution]]

# Biegler, L.T., Solution of Dynamic Optimization Problems by Successive Quadratic Programming and Orthogonal Collocation, Computers & Chemical Engineering, Volume 8, Number 3/4, pp. 243-248, 1984.

# Cizniar, M., Salhi, D., Fikar, M., and Latifi, M.A., A Matlab Package for Orthogonal Collocations on Finite Elements in Dynamic Optimization, 15th Int. Conference Process Control, June 7-10, 2015. [[https://www.kirp.chtf.stuba.sk/~fikar/research/dynopt/fik05oc.pdf | Article]]

# Renfro, J., Morshedi, A., Asbjornsen, O., Simultaneous Optimization and Solution of Systems Described by Differential/Algebraic Equations, Computers & Chemical Engineering, Volume 11, Number 5, pp. 503-517, 1987.

# Finlayson, B., [[https://faculty.washington.edu/finlayso/ebook/pde/OC/OCFE.htm|Orthogonal Collocation on Finite Elements for PDE Discretization and Solution]]

# Biegler, L.T., Solution of Dynamic Optimization Problems by Successive Quadratic Programming and Orthogonal Collocation, Computers & Chemical Engineering, Volume 8, Number 3/4, pp. 243-248, 1984.

Changed line 29 from:

# Nonlinear Modeling, Estimation and Predictive Control in APMonitor, Hedengren, J. D. and Asgharzadeh Shishavan, R., Powell, K.M., and Edgar, T.F., Computers and Chemical Engineering, Volume 70, pg. 133–148, 2014. Available at: https://dx.doi.org/10.1016/j.compchemeng.2014.04.~~013~~

to:

# Nonlinear Modeling, Estimation and Predictive Control in APMonitor, Hedengren, J. D. and Asgharzadeh Shishavan, R., Powell, K.M., and Edgar, T.F., Computers and Chemical Engineering, Volume 70, pg. 133–148, 2014. Available at: [[https://dx.doi.org/10.1016/j.compchemeng.2014.04.013|ScienceDirect]] or [[https://www.researchgate.net/publication/262191158_Nonlinear_Modeling_Estimation_and_Predictive_Control_in_APMonitor|ResearchGate]]

Changed line 24 from:

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

to:

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

Changed line 24 from:

<iframe width="560" height="315" src="https://www.youtube.com/embed/~~2FeOaGUQwKA~~?rel=0" frameborder="0" allowfullscreen></iframe>

to:

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

Changed line 19 from:

There are many papers that discuss the details of the derivation and theory behind these methods'^1^'. The purpose of this section is to give a practical introduction to orthogonal collocation on finite elements for the numerical solution of differential algebraic equations.

to:

There are many papers that discuss the details of the derivation and theory behind these methods'^1^'. The purpose of this section is to give a practical introduction to orthogonal collocation on finite elements with Lobatto quadrature for the numerical solution of differential algebraic equations.

Changed lines 19-20 from:

There are many ~~excellent ~~papers that discuss the details of the derivation and theory behind these methods. The purpose of this section is to give a practical introduction to orthogonal collocation on finite elements for the numerical solution of differential algebraic equations.

to:

There are many papers that discuss the details of the derivation and theory behind these methods'^1^'. The purpose of this section is to give a practical introduction to orthogonal collocation on finite elements for the numerical solution of differential algebraic equations.

Added lines 26-29:

!!!! References

# Nonlinear Modeling, Estimation and Predictive Control in APMonitor, Hedengren, J. D. and Asgharzadeh Shishavan, R., Powell, K.M., and Edgar, T.F., Computers and Chemical Engineering, Volume 70, pg. 133–148, 2014. Available at: https://dx.doi.org/10.1016/j.compchemeng.2014.04.013

Added lines 1-25:

(:title Orthogonal Collocation on Finite Elements:)

(:keywords direct transcription, orthogonal collocation, modeling language, differential, algebraic, tutorial:)

(:description Collocation methods are used to convert a differential algebraic model into a nonlinear programming problem that can be solved with NLP solvers:)

Discretization of a continuous time representation allow large-scale nonlinear programming (NLP) solvers to find solutions at specified intervals in a time horizon. There are many names and related techniques for obtaining mathematical relationships between derivatives and non-derivative values. Some of the terms that are relevant to this discussion include:

* Orthogonal collocation on finite elements

* Direct transcription

* Gauss pseudospectral method

* Gaussian quadrature

* Lobatto quadrature

* Radau collocation

* Legendre polynomials

* Chebyshev polynomials

* Jacobi polynomials

* Laguerre polynomials

* Any many more...

There are many excellent papers that discuss the details of the derivation and theory behind these methods. The purpose of this section is to give a practical introduction to orthogonal collocation on finite elements for the numerical solution of differential algebraic equations.

Attach:download.png [[Attach:orthogonal_collocation.zip|Orthogonal Collocation in MATLAB]]

(:html:)

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

(:htmlend:)

(:keywords direct transcription, orthogonal collocation, modeling language, differential, algebraic, tutorial:)

(:description Collocation methods are used to convert a differential algebraic model into a nonlinear programming problem that can be solved with NLP solvers:)

Discretization of a continuous time representation allow large-scale nonlinear programming (NLP) solvers to find solutions at specified intervals in a time horizon. There are many names and related techniques for obtaining mathematical relationships between derivatives and non-derivative values. Some of the terms that are relevant to this discussion include:

* Orthogonal collocation on finite elements

* Direct transcription

* Gauss pseudospectral method

* Gaussian quadrature

* Lobatto quadrature

* Radau collocation

* Legendre polynomials

* Chebyshev polynomials

* Jacobi polynomials

* Laguerre polynomials

* Any many more...

There are many excellent papers that discuss the details of the derivation and theory behind these methods. The purpose of this section is to give a practical introduction to orthogonal collocation on finite elements for the numerical solution of differential algebraic equations.

Attach:download.png [[Attach:orthogonal_collocation.zip|Orthogonal Collocation in MATLAB]]

(:html:)

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

(:htmlend:)