Solving Pendulum2 From Schittkowski Dae Test Suite?

I just tried to solve one of the DAE problems from Schittkowski DAE test suite (http://klaus-schittkowski.de/mc_dae.htm) but no success (Python 3.7). m=GEKKO(remote=False) m.time =

Solution 1:

Gekko reports the time steps that you request and does not fill in additional points. You are observing aliasing that comes with a different sampling frequency of the simulation than the natural frequency of the oscillating system. To resolve this, you can either increase the sampling frequency or else match the natural frequency so that you are always plotting the minima and maxima. Here is a solution with 5000 data points. With IMODE=7 Gekko scales linearly with increased horizon time.

from gekko import GEKKO
import numpy as np
import matplotlib.pyplot as plt
m=GEKKO(remote=False)
m.time = np.linspace(0.0,100.0,5000)
m1 = m.Param(value=1.0) 
g = m.Param(value=9.81)
s = m.Param(value=1.0)

m.options.IMODE=7
m.options.NODES=3#Initial values t=0
x = m.Var(value=0.0)
y = m.Var(value=-1.0)
v = m.Var(value=1.0)
w = m.Var(value=0.0)
l = m.Var(value=4.405)

m.Equation(x.dt() == v) 
m.Equation(y.dt() == w)
m.Equation(m1*v.dt() == -2*x*l)
m.Equation(m1*w.dt() == -m1*g - 2*y*l)

# Index-3 DAE#m.Equation(x**2+y**2==1)# Index-2 DAE
m.Equation(x*v  == -y*w)           

m.solve(disp=False)
plt.plot(m.time,x)
plt.show()

With IMODE=7 you also don't need the initialization step. If you are solving parameter optimization then you will benefit from the initialization step. One other suggestion is to use the index-3 DAE form instead of the index-2 DAE form:

# Index-3 DAE
m.Equation(x**2+y**2==1)

# Index-2 DAE# m.Equation(x*v  == -y*w) 

There is more information on the DAE index and advantages of using a higher index DAE solver like GEKKO versus DAE solvers that are limited to index-1 or index-2 Hessenberg form in the article Initialization strategies for optimization of dynamic systems.

There is no difference for this case study but numerical drift can be an issue if using index-0 (ODE) form.