Линейная интерполяция в ОДУ (ГЕККО)
Я потратил некоторое время, чтобы решить эту проблему, но до сих пор не могу найти решение. По сути, у меня есть параметр, скажем, k, как кусочно-линейная функция времени, например:
Итак, мы можем записать k как массив в зависимости от времени – k(t). Дело в том, что у меня есть система ОДУ, содержащая k, скажем, простая, например: dx/dt = k(t), и я бы хотел, чтобы GEKKO линейно интерполировала k из приведенной выше кусочно-линейной функции во время интегрирования ОДУ. По сути, это проблема моделирования.
После этого я рассмотрю изменение k(t) (т. е. оптимальное управление с использованием кусочно-линейного управления) для достижения некоторой целевой функции.
Есть ли какие-нибудь предложения о том, как правильно сделать это в GEKKO?
Спасибо
Я попробовал np.interp, но он не принимал переменные GK. Я также пытался написать собственный код для линейной интерполяции, но не смог понять, как извлечь момент времени (m.time) во время интегрирования.
Ответы 1
Создайте линейную интерполяцию перед созданием m.Param().
kt = np.array([0,0.2,0.45,0.75,1.0])
kv = np.array([0.5,1.0,0.2,1.2,1.0])
k = np.interp(m.time, kt, kv)
k = m.Param(value=k)
Вот пример решения на интервале времени 0,1 от 0 до 1.
from gekko import GEKKO
import numpy as np
kt = np.array([0,0.2,0.45,0.75,1.0])
kv = np.array([0.5,1.0,0.2,1.2,1.0])
m = GEKKO()
m.time = np.linspace(0, 1, 11)
k = np.interp(m.time, kt, kv)
k = m.Param(value=k)
x = m.Var(value=0)
m.Equation(x.dt() == k)
m.options.IMODE = 4
m.solve(disp=False)
import matplotlib.pyplot as plt
plt.figure(figsize=(5,2.5))
plt.plot(kt,kv,'bo-',label='k data')
plt.plot(m.time,k.value,'rx-',label='k interp')
plt.plot(m.time,x.value,'k-',label='x')
plt.plot(); plt.legend(); plt.grid()
plt.tight_layout(); plt.savefig('p.png',dpi=300)
plt.show()
Если есть ошибка интерполяции (как показано на рисунке выше), увеличьте количество временных точек или включите временные точки kt в m.time.
# include kt time points in m.time
m.time = np.concatenate((m.time, kt))
m.time = np.sort(m.time)
m.time = np.unique(m.time)
from gekko import GEKKO
import numpy as np
kt = np.array([0,0.2,0.45,0.75,1.0])
kv = np.array([0.5,1.0,0.2,1.2,1.0])
m = GEKKO()
m.time = np.linspace(0, 1, 11)
# include kt time points in m.time
m.time = np.concatenate((m.time, kt))
m.time = np.sort(m.time)
m.time = np.unique(m.time)
k = np.interp(m.time, kt, kv)
k = m.Param(value=k)
x = m.Var(value=0)
m.Equation(x.dt() == k)
m.options.IMODE = 4
m.solve(disp=False)
import matplotlib.pyplot as plt
plt.figure(figsize=(5,2.5))
plt.plot(kt,kv,'bo-',label='k data')
plt.plot(m.time,k.value,'rx-',label='k interp')
plt.plot(m.time,x.value,'k-',label='x')
plt.plot(); plt.legend(); plt.grid()
plt.tight_layout(); plt.savefig('p.png',dpi=300)
plt.show()
Если k вычислено, установите k=m.MV() с k.STATUS=1 для линейной интерполяции временных точек внутренней коллокации . Существует глобальная опцияm.options.MV_TYPE для удержания нулевого или первого порядка (линейная интерполяция).
Это хороший вопрос, и я думаю, у меня есть решение. Не могли бы вы задать отдельный вопрос StackOverflow с кодом и вашей лучшей попыткой?
Я уже задавал отдельный вопрос. Спасибо за помощь.
Вот ссылка на дополнительный вопрос: stackoverflow.com/questions/78077443/…
Дорогой Джон, Большое спасибо за ваши подробные ответы. Я попробовал ваш пример, и все работало отлично. У меня только один уточняющий вопрос. Когда я устанавливаю k = m.MV(...) с k.STATUS = 1 для оптимизации (IMODE 6), кажется, что решатель рассчитал все 13 точек k (в соответствии с длиной k или m.time). Если я хочу, чтобы решатель вычислил только пять точек в соответствии с исходными данными, как мне изменить код? Я также не хочу уменьшать дискретизацию времени, поскольку это влияет на точность.