A quick question

[Yongwen Wu]

Hi there,

Recently I am looking into Prof. Hedengren's dynamic optimization material. When I tried to follow this:

http://apmonitor.com/do/index.php/Main/LevelControl

I tried to do parameter estimation using GEKKO, and apply the dynamic model as provided, it runs OK. Here is part of code:

#df(dataframe) contains all the data after step response response. 0: time, 1: u, 2: h1, 3:h3

from gekko import GEKKO

t_data = df[0].to_list()

m = GEKKO()

m.time = t_data

#parameters

u=m.MV(value=df[1].to_list()) #input

K1=m.FV(lb=0.01,ub=100)

tau1=m.FV(lb=1,ub=200)

K2=m.FV(lb=0.01,ub=100)

tau2=m.FV(lb=1,ub=200)

#variables

h1=m.CV(value=df[2].to_list())

h2=m.CV(value=df[3].to_list())

#Equations

m.Equation(tau1*h1.dt()==-h1+K1*u)

m.Equation(tau2*h2.dt()==-h2+K2*h1)

#options

m.options.IMODE = 5   # dynamic estimation

m.options.NODES = 5   # collocation nodes

m.options.EV_TYPE=2

#STATUTS

u.STATUS=0

u.FSTATUS=1

K1.STATUS=1

K1.FSTATUS=0

tau1.STATUS=1

tau1.FSTATUS=0

K2.STATUS=1

K2.FSTATUS=0

tau2.STATUS=1

tau2.FSTATUS=0

h1.STATUS=1

h1.FSTATUS=1

h2.STATUS=1

h2.FSTATUS=1

m.solve(display=True)   # display solver output

It works perfectly fine, and I tried to change the dynamic model to the original model as follows:

t_data = df[0].to_list()

g = GEKKO()

g.time = t_data

#parameters

u=g.MV(value=df[1].to_list()) #input

C1=g.FV()

C2=g.FV()

C3=g.FV()

C4=g.FV()

#variables

h1=g.CV(value=df[2].to_list())

h2=g.CV(value=df[3].to_list())

#Equations (major changes)

g.Equation(h1.dt()== C1*u-C2*g.sqrt(h1))

g.Equation(h2.dt()== C3*u-C4*g.sqrt(h2)+C2*g.sqrt(h1))

#options

g.options.IMODE = 5   # dynamic estimation

g.options.NODES = 5   # collocation nodes

g.options.EV_TYPE=2

#STATUTS

u.STATUS=0

u.FSTATUS=1

C1.STATUS=1

C1.FSTATUS=0

C2.STATUS=1

C2.FSTATUS=0

C3.STATUS=1

C3.FSTATUS=0

C4.STATUS=1

C4.FSTATUS=0

h1.STATUS=1

h1.FSTATUS=1

h2.STATUS=1

h2.FSTATUS=1

g.solve(display=True)   # display solver output

The major change is from linear to nonlinear, and I tried to run it, but looks it is stuck there forever as shown:

Anyone has any idea about this?

I suppose it supports nonlinear term in the equations as Ipopt is a nonlinear solver, but a little bit confused that "running with linear solver ma57".

Thanks for any help.

Yongwen

[Yongwen Wu]

It is running as I changed "g = GEKKO(remote=False)", seems the remote server is down or something. The issue is after I changed, it is using "APOPT" rather than "OPOPT", and basically it doesn't work properly I assume that not a good start point is given.

But I didn't figure out why IPOPT is not applicable?

and Also, I tried to use MPC to do the control for this same problem, as the provided code. Since the first tank is easily full as the pump sees a big increase in a short time, so I tried to apply DMAX setting on its MV (p variable originally), but it seems the upper bound of pump is set not the maximum delta MV in a single time step.

I only add in the original python code:

p.DMAX=0.2

Any ideal?

Thanks

[John Hedengren]

Dear Yongwen,

 

Here are a couple things to try:

 

Set a non-zero default for C1 – C4 and set a lower bound (such as 1e-5).

 

Before

  C1=g.FV()

  C2=g.FV()

  C3=g.FV()

  C4=g.FV()

 

After

  C1,C2,C3,C4 = m.Array(g.FV,4,value=0.1,lb=1e-5)

 

Also, starting with a non-zero h1 and h2 can help because they start right at a point that is almost undefined with sqrt(h1) and sqrt(h2) if the solver tries negative values. Another way to deal with this is to transform the equations to remove sqrt by taking the square of both sides. That strategy is used successfully in the reservoir simulation exercise:

 

m.Equations([Vout[i]**2 == c[i]**2 * h[i] for i in range(4)])

 

Let us know if you need any additional support. The public server is working but it has different solver options than local solve methods. You may be getting better performance with IPOPT (and linear solver MUMPS) for the local solve versus IPOPT (with linear solver ma57) for the public server. There are license restrictions for distributing the Harwell Subroutine Libraries.

 

-John Hedengren

 

From: apmo...@googlegroups.com <apmo...@googlegroups.com> On Behalf Of Yongwen Wu
Sent: Wednesday, October 6, 2021 8:32 AM
To: apmonitor <apmo...@googlegroups.com>
Subject: [APM] Re: A quick question

 

It is running as I changed "g = GEKKO(remote=False)", seems the remote server is down or something. The issue is after I changed, it is using "APOPT" rather than "OPOPT", and basically it doesn't work properly I assume that not a good start point is given.

 

But I didn't figure out why IPOPT is not applicable?

 

 

 

and Also, I tried to use MPC to do the control for this same problem, as the provided code. Since the first tank is easily full as the pump sees a big increase in a short time, so I tried to apply DMAX setting on its MV (p variable originally), but it seems the upper bound of pump is set not the maximum delta MV in a single time step.

 

I only add in the original python code:

p.DMAX=0.2

 

 

 

Any ideal?

 

Thanks

 

 

 

 

Anyone has any idea about this?

 

I suppose it supports nonlinear term in the equations as Ipopt is a nonlinear solver, but a little bit confused that "running with linear solver ma57".

 

Thanks for any help.

Yongwen

--
--
APMonitor user's group e-mail list.
- To post a message, send email to apmo...@googlegroups.com
- To unsubscribe, send email to apmonitor+...@googlegroups.com
- Visit this group at http://groups.google.com/group/apmonitor
---
You received this message because you are subscribed to the Google Groups "apmonitor" group.
To unsubscribe from this group and stop receiving emails from it, send an email to apmonitor+...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/apmonitor/d8a62d25-d9fd-4c11-a6c1-b37303e04c81n%40googlegroups.com.