Мое основное увлечение - это аэрокосмос. И "космические" задачи - выведение полезных грузов на орбиту, возврат с орбиты через атмосферу являются функциями от целого набора параметров. К примеру, управление РН на первой ступени даже в двухмерной (без учета изменения азимута и dog-leg маневра)описывается по меньшей мере через 5 переменных:
продолжительность вертикального участка набора высоты после старта;
угол атаки при начале гравитационного разворота;
продолжительность маневра заклонения;
угол атаки после выполнения гравитационного разворота;
длительности пауз между отсечкой двигателей первой ступени, отделением второй ступени и запуском двигателей второй ступени.
Еще можно добавить сюда дросселирование двигателей, ограничения по скоростному напору, требования по запасу топлива для возврата и посадки первой ступени.
Так что даже в самой простой постановке нужно рассматривать вектор из 10-15 компонентов (не забываем, есть еще 2-ья и 3-ья ступень, поля падения отработавших ступеней и обтекателей и еще много интересного). Аналитически через набор готовых формул такую задачу не решить. Нужны численные методы, среди которых известны, просты и понятны методы градиентного поиска. (я знаю про методы случайного поиска, генетические алгоритмы и swarm-методы, но мне нужно прозрачное, быстрое, грязное и надежное решение, которое не будет пытаться стать умнее меня)
Градие́нт (от лат gradiens, род. падеж gradientis — шагающий, растущий) — вектор, своим направлением указывающий направление наискорейшего возрастания некоторой величины.
![Вычисление градиента из конечных разностей Вычисление градиента из конечных разностей](https://habrastorage.org/getpro/habr/upload_files/f27/6cc/9d8/f276cc9d8c93f8688a4662055efa8334.png)
Осталось только найти величину шага, по которому мы будем продвигаться в направлении градиента вдоль каждого измерения. Постоянный шаг - простое и понятное решение, но разные переменные оказывают разное влияние на целевую функцию, и при постоянном шаге велика вероятность начать блуждание вдоль склонов "хребта" из локальных оптимумов.
Для поиска оптимального шага воспользуемся методом дихотомии. Простой, не требующий даже вычисления производных.
![Делим отрезок пополам, сопоставляем значения вокруг точки дробления, идем дальше Делим отрезок пополам, сопоставляем значения вокруг точки дробления, идем дальше](https://habrastorage.org/getpro/habr/upload_files/8ef/126/109/8ef12610900a9455eaff0a41c84889c9.png)
Объединим эти два метода, на каждом шаге сначала вычислим градиент, затем найдем оптимальные величины шага для каждой компоненты вектора аргументов целевой функции, затем перейдем в следующую точку, и если модуль разности значения целевой функции в этой точке и ранее полученного максимального значения окажется меньше достаточно малого эпсилона, то мы закончим поиск.
А чтобы не вывихнуть свой мозг сразу, начнем с простой и проверяемой задачи квадратичной регрессии. Суть проста - есть набор экспериментальных данных, описывающих зависимость некоего Y от параметра X (набор зашумлен и замусорен).
Надо найти такие коэффициенты a0, a1 и a2 для параболической функции, чтобы эта функция-аппроксимация, пропущенная через экспериментальный набор точек, имела бы наименьшее среднее отклонение между квадратичной аппроксимацией и экспериментальными данными. Такая задача была многократно решена самыми разными способами, и эти решения уже давно встроены в Excel, openOffice, MathCAD.
Поставим задачу: для фиксированной набора данных вида [X, Y] найти значения коэффициентов регрессии a0, a1, a2, чтобы сумма для всех X из нашей выборки |R(a0, a1, a2, X) - Y| была минимальна.
Начнем движение из точки [-5, -5, -5] (это не самая лучшая точка, и в идеале поиск стартовой точки должен решаться отдельно, но старт из субоптимальный условий - хорошая проверка расчетной схемы на устойчивость).
![От разности в 5540 к скромным 27 за 14 шагов От разности в 5540 к скромным 27 за 14 шагов](https://habrastorage.org/getpro/habr/upload_files/d19/033/f1b/d19033f1b5a626fbf405bf567a07883b.png)
Теперь сравним полученные результаты с квадратичной регрессией, встроенной в openOffice.
![Красная линия - встроенныая квадратичная регрессия planMaker-а, темно-малиновая - регрессия, коэффициенты которой получены градиентным поиском Красная линия - встроенныая квадратичная регрессия planMaker-а, темно-малиновая - регрессия, коэффициенты которой получены градиентным поиском](https://habrastorage.org/getpro/habr/upload_files/85b/7e7/6a1/85b7e76a175d06c659f91badffa04d8a.png)
Теперь сравним количественные характеристики.
![*-индекс значений, полученных градиентным способом, t-индекс встроенной регрессии планмейкера. Самопальный метод выигрывает на 9,46% *-индекс значений, полученных градиентным способом, t-индекс встроенной регрессии планмейкера. Самопальный метод выигрывает на 9,46%](https://habrastorage.org/getpro/habr/upload_files/1b9/91f/b85/1b991fb85ed9d84da6834dc1948309ea.png)
Выводы
Градиентный поиск с подбором длины шага методом дихотомии действительно может находить максимум функции от вектора произвольной размерности.
Работоспособность метода подтверждена на задаче поиска коэффициентов квадратичной регрессии, в качестве эталона использовалась встроенная в openOffice planmaker функция полиноминальной регрессии.
Когда-нибудь я прикручу этот градиентный движок к своему рокетсайенсу.
Исходники алгоритма живут на моём гитхабе
Комментарии (11)
thegriglat
12.02.2022 18:11+1Если я правильно понял задачу, вам может подойти программа professor (с hepforge), она осуществляет поиск глобального минимума в многопараметрической системе (если есть симуляция, естественно)
the_stucky Автор
12.02.2022 18:25Спасибо. Не знал об этом проекте, выглядит серьезно и тяжеловесно, надо вникнуть
Jury_78
12.02.2022 18:31+2"Валик-джан, я тебе один умный вещь скажу, но только ты не обижайся:"
Обратите внимание на phython там есть отличные математические пакеты, например scipy.
the_stucky Автор
12.02.2022 18:36Остается только подружить питоновские либы с написанным на ноде движком траекторного моделирования. Ну или превозмочь и переписать траекторный блок на Cython каком-нибудь, чтобы настало щастье
aleex
12.02.2022 20:19+1Квадратичную регрессию можно было не считать в OpenOffice, а взять готовые формулы из "Метода наименьших квадратов" Л.В. Коломийцева.
В качестве метода поиска глобального минимума можете взять Shuffled Complex Evolution, метод старый и не шибко быстрый, но несложный в реализации.
Кстати, вы ракету как материальную точку моделируете или же как твердое тело?
the_stucky Автор
12.02.2022 20:36В существующей модели - материальная точка с тремя степенями свободы (итого 7-мь дифуров - 6 для координат + 1 для массы), управление в тангаже/курсе задается функциями угла атаки и крена от параметров полета. Добавить моделирование твердого тела вполне реально, вопрос тут в исходных данных (моменты инерции в динамике с учетом расхода топлива + коэф. моментов тангажа/рысканья/крена + характеристики органов управления от параметров полета и углов отклонения)
masai
13.02.2022 20:19Раз уж статья называется «Градиентный поиск коэффициентов квадратической регрессии», то замечу, что искать три коэффициента градиентным спуском — это из пушки по воробьям. Для среднеквадратической ошибки (MSE) можно взять нормальное уравнение.
Впрочем, у вас оптимизируется средняя абсолютная ошибка (MAE). Поэтому, кстати, сравнение с тем, что выдает офис не совсем корректно, так как он ищет скорее всего именно MSE. Так как целевые функции разные, то и результаты разные. Для шумных данных MSE будет выдавать среднее, а MAE медиану.
the_stucky Автор
13.02.2022 21:09Для проверки метода поиска регрессия оказалась удобной в том плане, что для нее не нужно далеко тянуться за эталоном. И по результатам какую-то откровенную дичь метод гнать не стал, что показалось мне хорошим признаком. Теперь нужно проверить его на чем-нибудь поближе к теме. Решение в двухмерной (дальность-высота) постановке задачи на дальность планирования с ограничениями на программные углы атаки и посадочную скорость, например
sshikov
Это все конечно хорошо, но вы по-моему забыли про одну вещь — все эти задачи оптимизации — они как правило с ограничениями. Такой же быстрый, прозрачный, грязный и т.п. подход к ним — это штрафы за нарушение ограничений. Но они обычно приводят к тому, что гладкая функция становится овражной, и простое решение перестает работать, или работает медленно.
У меня нет для вас ответа, как нужно, это просто сомнение.
the_stucky Автор
Здесь не может быть единого ответа. Можно двигаться из области штрафа по антиградиенту функции ограничений, можно действовать через "наказание случайностью", и при попадании в неудачную точку совершать откат с изменением параметров шага. И еще много чего. А пока - только проверка концепта на задаче, для которой легко найти эталонное решение