关于初始值约束的参数估计不适用于 m.connection 或 fixed_initial,对此有何建议?

问题描述

遵循Rate coefficients and initial concentration estimation using gekko - advice about model structure的建议 我修改了两个数据集的代码并运行了(我使用的是 COLDSTART=2)。

enter image description here

但是,我仍然有这些疑问,如果有一些建议就好了:

  1. 我想估计 kdoc1 和 kdoc2 以及 DOC1DOC2 的初始值 (t=0),即 DOC10DOC20,关于 { {1}} 和 DOC10[1] = DOC10[0]*0.5。我期望得到 DOC20[1] = DOC20[0]*0.5 的关系,之后 t=0DOC10 会衰减,但它没有发生,值保持不变。

  2. 我没有使用 DOC20,因为我遇到了这个错误

m.connection
  1. 并且不使用 p13 not found in results file p14 not found in results file p15 not found in results file p16 not found in results file ,因为它不可行

  2. 一个问题是我需要对变量初始化进行注释才能实现结果。

  3. 估计值等于提议的初始值。好像程序不合适,没有给出估计值。

代码如下:

fixed_initial =False

解决方法

这里有两个想法:

  1. 尝试一系列固定的初始条件 DOC10DOC20,然后对每个组合进行重新优化。这可以是 parallelized with multithreading,以便多个优化问题同时运行以加速解决方案。
import numpy as np
from gekko import GEKKO

# Optimize at mesh points
x = np.arange(20.0,30.0,2.0)
y = np.arange(30.0,50.0,4.0)
amg,bmg = np.meshgrid(x,y)

# Initialize results array
obj = np.empty_like(amg)

m = GEKKO(remote=False)
objective = float('NaN')

a,b = m.Array(m.FV,2)

# model variables,equations,objective
x1 = m.Var(1,lb=1,ub=5)
x2 = m.Var(5,ub=5)
x3 = m.Var(5,ub=5)
x4 = m.Var(1,ub=5)
m.Equation(x1*x2*x3*x4>=a)
m.Equation(x1**2+x2**2+x3**2+x4**2==b)
m.Minimize(x1*x4*(x1+x2+x3)+x3)
m.options.SOLVER = 1 # APOPT solver

# Calculate obj at all meshgrid points
for i in range(amg.shape[0]):
    for j in range(bmg.shape[1]):
        a.MEAS = amg[i,j]
        b.MEAS = bmg[i,j]

        m.solve(disp=False)

        obj[i,j] = m.options.OBJFCNVAL
        print(amg[i,j],bmg[i,obj[i,j])

# plot 3D figure of results
from mpl_toolkits.mplot3d import Axes3D
import matplotlib.pyplot as plt
from matplotlib import cm
import numpy as np

fig = plt.figure()
ax = fig.gca(projection='3d')
surf = ax.plot_surface(amg,bmg,obj,\
       rstride=1,cstride=1,cmap=cm.coolwarm,\
       vmin = 10,vmax = 25,linewidth=0,antialiased=False)
ax.set_xlabel('a')
ax.set_ylabel('b')
ax.set_zlabel('obj')
plt.show()

range of optimization values

  1. 如果优化器需要选择这些值,另一种选择是在有限元上写出正交搭配的离散化方程。这为 IMODE=3 提供了更多控制,因为所有方程和变量都是针对每个时间点单独定义的。这是带有 IMODE=4 的 Gekko 模型:
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

以及带有 IMODE=3 的等效模型:

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,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)

这将导致模型的重大重写,但如果需要更多控制,则可能是必要的。有关 Orthogonal Collocation on Finite Elements 的其他信息在动态优化课程中。