5.7. РЕФАКТОРИНГ АЛГОРИТМОВ И ЭВРОРИТМОВ

5.7. РЕФАКТОРИНГ АЛГОРИТМОВ И ЭВРОРИТМОВ

Алгоритм Нелдера — Мида является широко известным и применяется в качестве алгоритма прямого поиска локального экстремума вещественных функций от 2 до 6 вещественных переменных.

Следующий абзац содержит фрагмент текста из книги Д. Химмельблау [26], в котором содержится часть описания алгоритма Нелдера — Мида (метода деформируемого многогранника).

В методе Нелдера и Мида минимизируется функция n независимых переменных с использованием n + 1 вершин деформируемого многогранника в Eⁿ. Каждая вершина может быть идентифицирована вектором x. Вершина (точка) в Eⁿ, в которой значение f(x) максимально, проектируется через центр тяжести (центроид) оставшихся вершин. Улучшенные (более низкие) значения целевой функции находятся последовательной заменой точки с максимальным значением f(x) на более "хорошие" точки, пока не будет найден минимум f(x).

Начальный многогранник обычно выбирается в виде регулярного симплекса (но это не обязательно) с точкой в начале координат. Процедура отыскания вершины в Eⁿ, в которой f(x) имеет лучшее значение, состоит из следующих операций:

Отражение — проектирование x^(k)_h через центр тяжести в соответствии с соотношением x^(k)_n+3 = x^(k)_n+2 + ?(x^(k)_n+2 — x^(k)_h), где ? > 0 является коэффициентом отражения; x^(k)_n+2 — центр тяжести, x^(k)_h — вершина, в которой функция f(x) принимает наибольшее из n + 1 ее значений на k-м этапе.

Растяжение выполняется, если f(x^(k)_n+3) ? f(x^(k)_i), то вектор (x^(k)_n+3 — x^(k)_n+2) растягивается в соответствии с соотношением x^(k)_n+4 = x^(k)_n+2 + ?(x^(k)_n+3 — x^(k)_n+2), где ? — коэффициент растяжения. Если f(x^(k)_n+4) < f(x^(k)_i), то x^(k)_h заменяется на x^(k)_n+4 и процедура продолжается с шага 1 (k = k + 1), иначе x^(k)_n заменяется на x^(k)_n+3, процедура продолжается с шага 1 (k = k + 1).

Сжатие — если f(x^(k)_n+3) > f(x^(k)_i) для всех i ? h, то вектор (x^(k)_h — x^(k)_n+2) сжимается в соответствии с формулой x^(k)_n+5 = x^(k)_n+2 + ?(x^(k)_h — x^(k)_n+2), где 0 < ? < 1 представляет собой коэффициент сжатия. Затем x^(k)_h заменяется на x^(k)_n+5 и переходит на шаг 1 (k = k + 1).

Редукция — если f(x^(k)_n+3) > f(x^(k)_h), все векторы (x^(k)_i — x^(k)₁), i = 1…n + 1 уменьшаются в 2 раза с отсчетом от x^(k)₁ в соответствии с формулой x^(k)_i = x^(k)_i + 0,5(x^(k)_i — x^(k)₁), i = 1…n + 1. Затем возвращаемся к шагу 1 для продолжения поиска на k + 1 шаге.

Критерий окончания поиска, использованный Нелдером и Мидом, состоял в проверке условия среднего квадратичного отклонения функций f(x^(k)_i) от произвольного малого числа ?.

Химмельблау воспользовался традиционным математическим стилем для изложения описания алгоритма. Алгоритм имеет структуру вида "спагетти", что видно из схемы алгоритма, разработанной автором (рис. 5.15). Схема в виде "спагетти" — это скрытые go to, которые мы обнаруживаем в тексте программы, разработанной автором. Он же приводит графические рисунки принципа расчета точек и графический рисунок последовательности продвижения лучших точек симплекса по шагам метода при решении тестовой задачи. Затратив 11 страниц текста, приведя текст неструктурированной программы на 12 страницах, автор книги [26] так и не сумел понятно описать свой алгоритм.

Рис. 5.15. Схема алгоритма, разработанная Д. Химмельблау

Сравните способ подачи описания алгоритма, составленный автором [22], и функциональное описание алгоритма после рефакторинга. Функциональное описание алгоритма приведено ниже.

Особенностью алгоритма является продвижение к экстремуму целым облаком пробных точек, который условно назван симплексом (деформируемым многогранником). Общее число точек в симплексе равно увеличенному на единицу числу переменных поиска. Продвижение к экстремуму не по линии от одной точки к другой, а внутри какого-то множества пробных точек обеспечило нереагирование на мелкие дефекты целевой функции и прохождение относительно "широких" оврагов.

Вводимая информация:

n — число переменных минимизируемой функции;

X⁰ = (x⁰₁, x⁰₂…, x⁰_n) — точка первого начального приближения;

h — шаг регуляризации симплекса;

?_x — погрешность нахождения экстремума по параметрам.

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

Процедура регуляризации симплекса состоит в копировании во все поля значений аргументов симплекса значения вектора X⁰, и далее производятся изменения значений диагональных компонент векторов с номерами от 2 до n + 1 путем их увеличения на шаг регуляризации симплекса h.

В симплексе Q₁, Q₂, Q₃…, Q_n+1 — это рассчитанные значения минимизируемой функции при соответствующих по строке значениях аргументов минимизируемой функции.

Далее до выполнения условия окончания поиска осуществляются итерации (шаги) поиска. Условием окончания поиска этого метода является неудаленность от лучшей точки симплекса остальных точек симплекса более чем на e²_x.

Каждая итерация начинается с нахождения номеров l и k, соответственно лучшей и худшей по значению функции точек симплекса. Далее осуществляется расчет точки X_ц — положения центра масс всех точек симплекса за исключением наихудшей точки 1. Это позднее улучшение алгоритма. Д. Химмельблау не исключал наихудшей точки.

где n — количество точек в симплексе; i — номер компоненты вектора X(i = 1, 2…, n); j — номер точки в симплексе.

Считаются нулевыми значения x_kj слагаемых сумм, где k — номер наихудшей точки.

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

X_? = X_ц + ?(X_ц — X_k), X_? = X_ц + ?(X_ц — X_k), X_? = X_ц + ?(X_ц — X_k).

Целесообразно эти похожие формулы реализовать одной функцией.

Значения коэффициентов ? = 1, ? = 0,5 и ? = 2 подобраны экспериментальным путем (в оригинале описания алгоритма эти значения можно выявить лишь из текста программ). Расчет точек ?, ?, ? целесообразно осуществлять одной процедурой с параметром в виде значений коэффициентов ?, ?, ?.

После расчета точки X_ц вычисляется пробная точка ?. Далее выполняется одно из альтернативных действий.

Если Q_? ? Q₁, то выполняются действия достижения сильного успеха. Если Q₁ ? Q_? < Q_k, то имеется слабый успех, а точка ? записывается на место k-точки. Если Q_? ? Q_k, то выполняются действия отсутствия успеха. Рассчитывается X_? и Q_?. Далее, если Q_? ? Q_k, то точка ? записывается на место k-точки, иначе, если точка ? хуже точки k, выполняется процедура редукции симплекса и процедура расчета значения функции в точках симплекса.

Действия достижения сильного успеха. Рассчитывается X_? и Q_?. Наилучшая из точек ? или ? записывается на место наихудшей k-точки симплекса.

Действия отсутствия успеха. Рассчитывается X_? и Q_?. Далее выполняется действие по изменению симплекса при отсутствии успеха.

Действие по изменению симплекса при отсутствии успеха представляет собой альтернативу: если Q_? ? Q_k, то точка ? записывается на место k-точки, иначе, если точка ? хуже точки k, выполняется процедура редукции симплекса.

При выполнении процедуры редукции симплекса все точки симплекса стягиваются к лучшей точке симплекса на половину своего прежнего удаления и далее выполняется процедура расчета значений целевой функции во всех точках симплекса.