Максимизация высоты взлёта ракеты
Оптимизация параметров управления ракетой
Введение
Актуальной задачей в области прикладной математики и теории управления является решение нелинейных задач оптимального управления динамическими системами. Данный пример рассматривает частный, но показательный случай такой задачи — максимизацию конечной высоты вертикального полёта ракеты с учётом её динамики.
Целью настоящего исследования является демонстрация методологии формулировки и численного решения подобных задач с использованием современных инструментов вычислительной оптимизации.
Представленная модель учитывает ключевые физические факторы: переменную массу вследствие расхода топлива, гравитацию, зависящую от высоты, и аэродинамическое сопротивление, что делает её подходящей для демонстрации принципов оптимального управления.
Присоединим необходимые библиотеки.
neededLibs = ["JuMP", "Ipopt", "Measures"]
for lib in neededLibs
try
eval(Meta.parse("using $lib"))
catch ex
Pkg.add(lib)
eval(Meta.parse("using $lib"))
end
end
Необходимо подобрать такие оптимальные параметры, чтобы обеспечить максимально возможную высоту вертикально запускаемой ракеты.
Мы можем управлять тягой ракеты с учётом её массы, расхода топлива, гравитации и аэродинамическое сопротивления.
Рассмотрим базовое описание модели.
В нашей модели присутствуют три переменных состояния:
-
Скорость:
-
Высота:
-
Масса ракеты и оставшегося топлива:
и одна управляющая переменная:
- Тяга: .
Динамика ракеты описывается тремя уравнениями:
-
Скорость подъёма:
-
Ускорение:
-
Скорость потери массы: ,
где сопротивление является функцией высоты и скорости, гравитация — функцией высоты, а — константа.
Эти силы определяются следующим образом:
и
.
Мы используем дискретизированную по времени модель с фиксированным количеством шагов .
Таким образом, наша цель — максимизировать .
Исходные данные
В данной модели использованы нормализованные безразмерные параметры.
Инициализируем исходные данные.
h_0 = 1 # Начальная высота
v_0 = 0 # Начальная скорость
m_0 = 1.0 # Начальная масса
m_T = 0.6 # Конечная масса
g_0 = 1 # Гравитация на поверхности
h_c = 500 # Параметр высоты для сопротивления
c = 0.5 * sqrt(g_0 * h_0) # Удельный импульс (тяга к расходу топлива)
D_c = 0.5 * 620 * m_0 / g_0 # Коэффициент сопротивления
u_t_max = 3.5 * g_0 * m_0 # Максимальная тяга
T_max = 0.2 # Время полёта (в секундах)
T = 1_000 # Количество временных шагов
Δt = T_max / T # Длительность шага по времени
Создадим модель и выберем метод оптимизации. Поскольку это нелинейная задача, нам необходимо использовать нелинейный метод, такой как Ipopt.
model = Model(Ipopt.Optimizer)
set_silent(model)
Далее создадим переменные, которые описывают состояние и управление ракетой. Все они меняются со временем. Для того, чтобы метод быстрее нашёл решение, мы задаём для каждой переменной её приблизительное начальное значение.
@variable(model, x_v[1:T] >= 0, start = v_0) # Скорость
@variable(model, x_h[1:T] >= 0, start = h_0) # Высота
@variable(model, x_m[1:T] >= m_T, start = m_0) # Масса
@variable(model, 0 <= u_t[1:T] <= u_t_max, start = 0); # Тяга
Зададим начальные и конечные условия.
fix(x_v[1], v_0; force = true)
fix(x_h[1], h_0; force = true)
fix(x_m[1], m_0; force = true)
fix(u_t[T], 0.0; force = true)
Цель оптимизации — сделать так, чтобы к концу манёвра ракета достигла наибольшей возможной высоты.
@objective(model, Max, x_h[T])
Сопротивление и гравитация в модели задаются с помощью функций:
D(x_h, x_v) = D_c * x_v^2 * exp(-h_c * (x_h - h_0) / h_0)
g(x_h) = g_0 * (h_0 / x_h)^2
Ограничения оптимизации определяются в виде уравнений, описывающих движение ракеты.
ddt(x::Vector, t::Int) = (x[t] - x[t-1]) / Δt
@constraint(model, [t in 2:T], ddt(x_h, t) == x_v[t-1])
@constraint(model, [t in 2:T], ddt(x_v, t) == (u_t[t-1] - D(x_h[t-1], x_v[t-1])) / x_m[t-1] - g(x_h[t-1]))
@constraint(model, [t in 2:T], ddt(x_m, t) == -u_t[t-1] / c);
Выполним задачу оптимизации и проверим решение:
optimize!(model)
assert_is_solved_and_feasible(model)
solution_summary(model)
Определим функцию для построения графиков высоты, массы, скорости и тяги ракеты с оптимальными параметрами.
function траектория(y; kwargs...)
return Plots.plot(
(1:T) * Δt,
value.(y);
xlabel = "Время в секундах",
legend = false,
linewidth = 2,
kwargs...,
)
end
Построим график высоты полёта ракеты.
Plots.plot(траектория(x_h; ylabel = "Высота"))
Построим график изменения массы ракеты.
Plots.plot(траектория(x_m; ylabel = "Масса"))
Построим график скорости ракеты.
Plots.plot(траектория(x_v; ylabel = "Скорость"))
Построим график оптимальной тяги ракеты.
Plots.plot(траектория(u_t; ylabel = "Тяга"))
Заключение
В рамках данного исследования была успешно решена задача нелинейного оптимального управления для модели вертикального полёта ракеты. Поставленная цель — максимизация конечной высоты — достигнута путём построения дискретной аппроксимации исходной непрерывной задачи и её последующей оптимизации.
Методология, представленная в работе, является универсальной и может быть распространена на более сложные случаи, такие как задачи трёхмерной динамики, многостадийные полёты или системы с дополнительными фазовыми ограничениями. Таким образом, исследование вносит вклад в практику применения методов нелинейного программирования для решения инженерных задач оптимального управления.