Как оптимизировать траекторию для прохождения путевых точек в GEKKO? ⇐ Python

[Гость]

Моя цель — рассчитать оптимальную траекторию дрона с неподвижным крылом, чтобы он достиг ряда заданных координат x (например, [150, 0]).

Я попробовал несколько подходов для обеспечения соблюдения ограничений путевых точек. При оптимизации только для одной путевой точки мне удалось добиться успеха с обеими:

m.fix_final(x, val = 150.0) Это также фиксирует производную x до 0:

m.Minimize(final*10**5*(x - 150)**2) Это сработало, но не казалось применимым к случаю с несколькими путевыми точками, поскольку время достижения путевых точек затрудняет указание того, какой узел должен минимизировать расстояние до путевой точки, и:

waypoint = m.Var(value=0) m.Equation(путевая точка == m.if3(x - 150, 0, 1)) m.fix_final(путевая точка, значение=1) Это оставляет производную свободной, но получается неуклюжей (я также не совсем понимаю, как оптимизатор с этим справляется, и он кажется очень чувствительным к начальным условиям).

Все подходы не увенчались успехом при расширении до двух путевых точек, в моей текущей попытке меня вдохновил этот подход.

импортировать numpy как np импортировать matplotlib.pyplot как plt из импорта GEKKO GEKKO из mpl_toolkits.mplot3d импортировать Axes3D импортировать matplotlib.pyplot как plt путевые точки = [150, 0] # создаем модель GEKKO m = ГЕККО (удаленный = Ложь) num_time_steps = 31 # масштабируем время 0-1 с помощью tf m.time = np.linspace(0,1,num_time_steps) # параметры m.options.NODES = 4 # узла коллокации m.options.SOLVER = 1 # решатель APOPT # 3 # решатель (IPOPT) m.options.IMODE = 6 m.options.MAX_ITER = 1000 m.options.MV_TYPE = 0 m.options.DIAGLEVEL = 1 # последний раз tf = [m.FV(value=3.0, lb = 0.1) for i in range(len(waypoints))] # lb вовремя, чтобы предотвратить путешествие во времени # Управляемые переменные CL = [m.MV(value=1, lb = 0, ub = 10) for i in range(len(waypoints))] # коэффициент подъемной силы phi = [m.MV(value=0, lb = np.deg2rad(-45), ub = np.deg2rad(45)) for i in range(len(waypoints))] # угол крена # установите STATUS в 1, чтобы позволить оптимизатору изменить для меня в диапазоне (len(путевые точки)): tf.СТАТУС = 1 CL.СТАТУС = 1 phi_dot.STATUS = 1 масса = м.Const(значение=9,5) г = m.Const (значение = 9,81) CD0 = m.Const(значение=0,01) ро = m.Const (значение = 1,225) S = m.Const(значение=0,65) f_max = m.Const (значение = 40) # Переменные состояния V = [m.Var(value=100, lb = 0, ub = 200, fix_initial=False) for i in range(len(waypoints))] # скорость psi = [m.Var(value=0, fix_initial=False) for i in range(len(waypoints))] # угол курса gamma = [m.Var(value=0, lb = np.deg2rad(-10), ub = np.deg2rad(10), fix_initial=False) for i in range(len(waypoints))] # угол траектории полета x = [m.Var(value=0, fix_initial=False) for i in range(len(waypoints))] # позиция x y = [m.Var(value=0, fix_initial=False) for i in range(len(waypoints))] # позиция y z = [m.Var(value=0, fix_initial=False) for i in range(len(waypoints))] # позиция z # указываем начальные условия m.fix_initial(x[0], val=0.0) m.fix_initial(y[0], val=0.0) m.fix_initial(z[0], значение=150.0) m.fix_initial(V[0], значение=60.0) m.fix_initial(psi[0], val=0,0) m.fix_initial(гамма[0], значение=0,0) # Аэродинамическая модель k = m.Intermediate(4 * f_max ** 2 * CD0) CD = [m.Intermediate(CD0 + CL ** 2 / k) для i в диапазоне (len(путевые точки))] D = [m.Intermediate(0,5 * rho * V ** 2 * S * CD) для i в диапазоне (len(путевые точки))] L = [m.Intermediate(0,5 * rho * V ** 2 * S * CL) для i в диапазоне (len(путевые точки))] для меня в диапазоне (len(путевые точки)): # дифференциальные уравнения, масштабированные по tf m.Equation(V.dt()==tf*(- D[i] / масса - g*m.sin(gamma[i]))) m.Equation(gamma[i].dt()== tf[i] * (L[i] * m.cos(phi[i]) - масса * g * m.cos(gamma[i])) / (масса * V[i])) m.Equation(psi[i].dt()==tf[i]*(L[i] * m.sin(phi[i]) + m.cos(psi[i])) / (масса * V [i] * m.cos(гамма[i]))) m.Equation(x[i].dt()==tf[i]*(V[i] * m.cos(gamma[i]) * m.cos(psi[i]))) m.Equation(y[i].dt()==tf[i]*(V[i] * m.cos(gamma[i]) * m.sin(psi[i]))) m.Equation(z[i].dt()==tf[i]*(V[i] * m.sin(gamma[i]))) для i в диапазоне (len(путевые точки) - 1): m.Connection(phi[i+1], phi[i], pos1 = 1, pos2 = 'конец', node1 = 1, node2 = 'конец') m.Connection(x[i+1], x[i], pos1 = 1, pos2 = 'конец', node1 = 1, node2 = 'конец') m.Connection(y[i+1], y[i], pos1 = 1, pos2 = 'конец', node1 = 1, node2 = 'конец') m.Connection(z[i+1], z[i], pos1 = 1, pos2 = 'конец', node1 = 1, node2 = 'конец') m.Connection(psi[i+1], psi[i], pos1 = 1, pos2 = 'конец', node1 = 1, node2 = 'конец') m.Connection(gamma[i+1], gamma[i], pos1 = 1, pos2 = 'конец', node1 = 1, node2 = 'конец') m.Connection(V[i+1], V[i], pos1 = 1, pos2 = 'конец', node1 = 1, node2 = 'конец') m.Connection(phi[i+1],'вычислено', pos1=1, node1=1) m.Connection(x[i+1],'вычислено', pos1=1, node1=1) m.Connection(y[i+1],'вычислено', pos1=1, node1=1) m.Connection(z[i+1],'вычислено', pos1=1, node1=1) m.Connection(psi[i+1],'вычислено', pos1=1, node1=1) m.Connection(гамма[i+1],'вычислено', pos1=1, node1=1) m.Connection(V[i+1],'вычислено', pos1=1, node1=1) f = np.zeros(num_time_steps); ж[-1]=1; окончательный = m.Param (f) # минимизировать конечное время при достижении путевых точек для меня в диапазоне (len(путевые точки)): m.Minimize(final*10**5*(x[i] - waypoints[i])**2) m.Minimize((m.sum(tf))**2) м.решить() Этот подход решает, что он не может достичь обеих путевых точек, и просто пытается расположить обе конечные точки фазы как можно ближе к средней точке между двумя целями.





Существует ли рекомендуемый способ решения проблемы с несколькими маршрутными точками такого типа? (/есть ли где-то в моем динамике глупая ошибка, запрещающая поворот?)

Заранее благодарим за любую помощь.