Постановка задачи оптимизации

Лекция 1 Постановка задачи оптимизации. Классификация План 1. Постановка задач оптимизации. 2. Классификация задач оптимизации. 3. Классификация методов оптимизации. 1. Постановка задач оптимизации Достаточно часто мы сталкиваемся с ситуацией, когда из некоторой совокупности возможных вариантов своего поведения или принятия решения в какой-либо области деятельности необходимо выбрать один. Выбор осуществляется путем сравнения различных вариантов при помощи некоторой количественной их оценки. В этом случае говорят о необходимости решения задачи оптимизации. Оптимизация (от латинского слова «optimus» – наилучший) – поиск наилучшего варианта, при наличии множества альтернативных. Всевозможные устройства, процессы и ситуации, применительно к которым предстоит решать задачу оптимизации, называют объектом оптимизации. Методы оптимизации занимаются построением оптимальных решений для математических моделей, при этом именно вид модели определяет метод или методы, используемые для построения оптимального решения. Для решения любой задачи необходимо формализовать объект оптимизации и представить его в виде математической модели. Математическая модель – модель, которая определена с помощью математических формализмов. Математическая модель не является точной, а является идеализацией. Для применения теории оптимизации к решению конкретных задачи нужно выполнить определѐнную последовательность действий, которая называется постановкой задачи оптимизации. Она включает этапы. Этап I. Установление границ подлежащей оптимизации системы. Система – некая изолированная часть внешнего мира. Границы системы задают пределы, которые отделяют еѐ от внешнего мира. При этом предполагается, что взаимосвязи с внешним миром зафиксированы. Первоначальный выбор границ системы может оказаться слишком жѐстким. Для полу7 чения адекватного решения нужно включить в систему дополнительные подсистемы, однако это ведѐт к увеличению размерности задачи. Следует стремиться к представлению системы в виде изолированных подсистем, которые можно рассматривать независимо от других. Этап II. Выбор количественного критерия, позволяющего выявить наилучший вариант, называемого характеристическим критерием. Критерии могут быть, в зависимости от конкретной задачи, экономического или технологического характера (минимальная стоимость, максимальный крутящий момент). Независимо от того, какой критерий принят в качестве характеристического, он должен принимать максимальное (или минимальное) значение для наилучшего варианта. Критериев может быть много, тогда задача становится многокритериальной. Существуют методы решения многокритериальных задач, но можно привести многокритериальную задачу к однокритериальной. Для этого один из критериев выбирается в качестве первичного, а остальные становятся вторичными. Первичный критерий используется как характеристический, а вторичные формируют ограничения задачи. Этап III: определение внутрисистемных переменных, через которые выражается характеристический критерий. Выбор переменных осуществляется с учѐтом следующих рекомендаций. Нужно разделить переменные, которые меняются в широком диапазоне и переменные, которые фиксированы или меняются слабо. Первые определяются, как независимые переменные, а вторые – как параметры задачи. Параметры задачи разделяют на фиксированные и те, которые испытывают флуктуации под воздействием вешней среды. Нужно выбрать только те переменные, которые оказывают наибольшее влияние на характеристический критерий. Этап IV. Построение модели, которая описывает взаимосвязь внутрисистемных переменных. Модель системы описывает взаимосвязь между переменными и отражает степень влияния этих переменных на характеристический критерий. Модель включает в себя основные уравнения материальных и энергетических балансов; уравнения, описывающие физические процессы в системе; неравенства, определяющие области допустимых значений переменных. Изменяемые при оптимизации величины, входящие в математическую модель объекта оптимизации, называют параметрами оптимизации, а соотношения, устанавливающие пределы возможного изменения этих пара- 8 метров, – ограничениями. Эти ограничения могут быть заданы в форме равенств или неравенств. Таким образом, задача в виде, пригодном для решения методом оптимизации состоит в минимизации (максимизации) вещественнозначной функции f (x) n-мерного аргумента x, компоненты которого удовлетворяют системе ограничений в виде уравнений H k ( x)  0, k  1, 2,..., m или неравенств g j ( x)  0, j  m  1,...,s. Математически будем записывать следующим образом: f ( x)  min (max), Q Q  H k ( x)  0, k  1,..., m, Q:  g j ( x)  0, j  m  1,..., s. При математической формулировке задачи условной оптимизации целевую функцию выбирают с таким знаком, чтобы решение задачи соответствовало поиску минимума этой функции. Поэтому формулировку общей задачи математического программирования обычно записывают так:  (1) f ( x )  min , Q где Q  R N — множество возможных альтернатив, рассматриваемых при поиске решения задачи. Любую точку x  Q называют допустимым решением задачи математического программирования, а само множество — множеством допустимых решений или, короче, допустимым множеством. Точку x*, в которой функция f (x) достигает своего наименьшего значения, называют оптимальным решением задачи. Если множество параметров оптимизации является подмножеством конечномерного линейного пространства, то говорят о конечномерной задаче оптимизации в отличие от бесконечномерных задач, которые рассматривают в вариационном исчислении и оптимальном управлении. При этом критерием оптимальности может быть требование достижения наибольшего или наименьшего значения одной или несколькими действительными (скалярными) функциями параметров — оптимизации, выражающими количественно меру достижения цели оптимизации рассматриваемого объекта. Каждую из таких функций принято называть целевой. Задачу (1) в дальнейшем будем называть задачей минимизации целевой функции на множестве Q. Но целевая функция может и не достигать на Q 9 наименьшего значения. Тогда говорят о точной нижней грани функции f (x) на этом множестве и вместо (1) используют запись (2) f ( x)  inf, x  Q. Отличие (1) от (2) в том, что в первом случае предполагают существование точки x * , в которой целевая функция достигает своего наименьшего значения на множестве Q, а во втором случае такая точка может и не существовать. Поэтому решение общей задачи математического программирования состоит в том, чтобы в первом случае найти точные (или с некоторой заданной точностью) значения координат точки x*  {x j } j 1,..,N и значение целевой функции f ( x*)  min f ( x) . Во втором случае построить таQ кую последовательность точек xn nN , которой бы соответствовала последовательность f ( x n ) , сходящаяся к значению inf f ( x) , и вычислить Q это значение с заданной точностью. Отметим, что в большинстве прикладных задач имеет место первый случай, поэтому использование записи вида (2) будем оговаривать особо. 2. Классификация задач оптимизации Если целевая функция единственная, то задачу конечномерной оптимизации называют задачей математического программирования, а в противном случае — задачей многокритериальной (векторной) оптимизации. При отсутствии ограничений множество Q совпадает с областью определения целевой функции. Если же рассматриваемые альтернативы должны удовлетворять некоторым ограничениям, то множество допустимых решений сужается. Рассмотрим задачу оптимизации: найти минимум функции f (x) nмерного аргумента x, компоненты которого удовлетворяют системе ограничений в виде уравнений H k ( x)  0, k  1, 2,..., m или неравенств g j ( x)  0, j  m  1,...,s. Такая задача называется задачей условной оптимизации. Если задача не содержит ограничения и рассматривается на всем пространстве, то это задача безусловной оптимизации. Задачи оптимизации классифицируются в соответствии с видом функций f ( x), H k ( x), g j ( x) , и размерностью вектора х. 10 Задачи без ограничений с N  1 называются задачами одномерной оптимизации, с N  2 – многомерной оптимизации. Если в задаче функции H k ( x), g j ( x) линейны, то это задача с линейными ограничениями. При этом целевая функция f (x) может быть как линейной, так и нелинейной. Задача условной оптимизации, в которой все функции линейны, называется задачей линейного программирования. Задачи с нелинейной целевой функцией называются задачами нелинейного программирования. При этом если f (x) квадратичная функция, то задачей квадратичного программирования. Если f (x) отношение линейных функций, а ограничения – линейные, то рассматривается задача дробнолинейного программирования. Если целевая функция представима в виде n f ( x)   f j ( x j ) , функции системы ограничений имеют вид j 1 n g i ( x)   g ik ( x k )   i , i  1, 2,...,m , соответствующая задача называется k 1 задачей сепарабельного программирования. k Функция называется позиномом, если f ( x)   c j p j ( x), где c j  R  , j 1 n  p j ( x)   xi ji ,  ji  R . Задача оптимизации, для которой f (x) , i 1 H k ( x), g j ( x) – позиномы, называется задачей геометрического про- граммирования. Пример. Задача о размещении заказа между предприятиями с различными производственными и технологическими характеристиками, но взаимозаменяемыми в смысле выполнения заказа. Требуется составить план размещения заказа так, чтобы заказ был выполнен с минимальными технологическими затратами при имеющихся производственных возможностях. Рассмотрим математическую модель задачи. Имеется m однотипных предприятий, на которых нужно разместить заказ n видов продукции в объемах x *j ,   - затратная матрица, кажj  1, 2,..., n . Известна C  cij mn дый элемент которой cij – есть общие затраты i-го предприятия на выполнения j-го вида продукции. Мощность каждого предприятия ограничена временными ресурсами Ti . Обозначим X  xij план заказов, здесь  mn xij  количество времени, затраченное i-м предприятием на изготовление 11 j-го продукта. Математическая модель данной задачи имеет вид задачи линейного программирования: Целевая функция - суммарные затраты на выполнение заказа. f  n  cij xij  min . j 1 Система ограничений: по временным затратам на производство - n  xij  Ti , i  1,..., m , j 1 выпуск продукции соответствует плану m  aij xij i 1 естественные ограничения - xij  0, i  1,..., m,  x *j , j  1,..., n , j  1,..., n . 3. Классификация методов оптимизации В соответствии с классификацией задач оптимизации классифицируются и методы оптимизации. В некоторых случаях ограничения в задаче оптимизации позволяют через один из параметров выразить остальные и исключить их из целевой функции. В результате данных действий, задача будет сведена к поиску наибольшего или наименьшего значения скалярной действительной функции f (x), выражающей критерий оптимальности. Выбирая тот или иной знак перед этой функцией, всегда можно ограничиться лишь поиском ее наименьшего значения в области определения D(f), заданной с учетом ограничений на параметр оптимизации х. Для краткости будем говорить об одномерной минимизации, имея в виду нахождение наименьшего значения функции f(x) на множестве D(f) и точек, в которых это значение достигается. Изучение методов одномерной минимизации важно не только для решения задачи min {z  f ( x)} , имеющей самостоятельное значение. Эти x[ a,b] методы являются также существенной составной частью методов многомерной минимизации, при помощи которых находят наименьшее значение действительных функций многих переменных, так как целый ряд методов нелинейного программирования включает в себя в качестве составной части решение задач одномерной оптимизации. Методы одномерной оптимизации разделяются на подклассы по следующим принципам: 12  использование в процессе поиска экстремума информации о самой функции, так как в ряде задач целевая функция задана таким образом, что точных значений производных найти нельзя (только оценить).  использование в процессе поиска экстремума информации о самой функции или ее производных.  по виду целевой функции (методы решения одно- и многоэкстремальных задач). Методы оптимизации подразделяются на аналитические, и численные (приближенные). Все численные методы решения задач безусловной оптимизации со-   стоят в том, что мы строим последовательность точек x (n) зом, чтобы последовательность функций таким обра- f ( x (n) ) была убывающей  nN , удовлетворяющая этому f ( x k 1 )  f ( x k ). Последовательность x (n) условию, называется релаксационной последовательностью. Методы решения делятся на методы с использованием информации о производных функции и без использования таковой. Численное решение задач безусловной минимизации функций многих переменных, как правило, значительно сложнее, чем решение задач минимизации функций одного переменного. В самом деле, с ростом числа переменных возрастают объемы вычислений и усложняются конструкции вычислительных алгоритмов, а также более сложным становится анализ поведения целевой функции. Методы численного решения задач многомерной безусловной минимизации многочисленны и разнообразны. Условно их можно разделить на три больших класса в зависимости от информации, используемой при реализации метода. 1. Методы нулевого порядка, или прямого поиска, стратегия которых основана на использовании информации только о свойствах целевой функции. 2. Методы первого порядка, в которых при построении итерационной процедуры наряду с информацией о целевой функции используется информация о значениях первых производных этой функции. 3. Методы второго порядка, в которых наряду с информацией о значениях целевой функции и ее производных первого порядка используется информация о вторых производных функции. 13 Вопросы и задания для самопроверки 1. Что называют точкой строгого локального экстремума функции одного переменного? 2. Что такое математическая модель объекта оптимизации? Каковы этапы построения математической модели 3. Сформулируйте математическую постановку задачи оптимизации. 4. Какова классификация задач оптимизации по размерности управляемой переменной, условиям на функции, методам решения. 5. Дайте определение оптимального решения задачи оптимизации. 6. Какая последовательность называется релаксационной? 1 7. Будет ли последовательность x n  2  релаксационной для функции n x f ( x)  ( x  3)e . 8. В прямой круговой конус вписан прямой круговой цилиндр так, что основания конуса и цилиндра лежат в одной плоскости. Найти наибольшую возможную часть объема конуса, занятую таким цилиндром. К какому типу относится модель данной задачи. 9. Постройте математическую модель задачи. В заданный треугольник ABC с высотой Н и основанием b вписать параллелограмм наибольшей площади, стороны которого параллельны двум сторонам треугольника. Классифицируйте задачу. 10. К какому типу относится задача оптимизации: найти минимум функции f ( x, y )  1  e 14  x2  y2 , при ограничении x  3 y  1 2 Лекция 2 Методы одномерной оптимизации без использования информации о производной функции План 1. Релаксационная последовательность. Оптимальный пассивный поиск. 2. Методы последовательного поиска (методы интервалов) Метод дихотомии. Метод деления пополам. Метод золотого сечения. Метод Фибоначчи. 3. Методы аппроксимации 1. Релаксационная последовательность. Оптимальный пассивный поиск Рассмотрим задачу одномерной оптимизации min {z  f ( x)} . Введем x[ a,b] некоторые определения. Монотонность функции. Функция f (x) является монотонной на интервале, если для любых x1 и x2 из этого интервала, таких, что x1  x 2 выполняется неравенство f ( x1 )  f ( x 2 ) , если функция монотонно возрастающая или f ( x1 )  f ( x 2 ) , если функция монотонно убывающая. Унимодальность. Функция f (x) является унимодальной на отрезке, если она монотонна по обе стороны от единственной на отрезке точки х0, то есть функция f (x) в полуинтервале [а, х0) убывает, а в полуинтервале (х0,b] возрастает. Примеры графиков унимодальных функций приведены на рис. 1. Точка х0 может быть внутренней точкой отрезка [а, b] или совпадать с одним из его концов. Унимодальная функция не обязательно непрерывна на отрезке [а, b]. Определение глобального минимума. Функция f (x) , определѐнная на множестве D достигает глобального минимума в точке x *  D , если f ( x * )  f ( x) для всех x  D . 15 Определение локального минимума. Функция f (x) , определѐнная на множестве D имеет локальный минимум в точке x *  D , если существует такая ε-окрестность точки x*, что для всех x из этой ε-окрестности f ( x * )  f ( x) . f(x) a x0 f(x) b a x x0 b x Рис. 1 Если функция f (x) не унимодальна, то наименьший из локальных минимумов будет глобальным минимумом (аналогично – наибольший из локальных максимумов будет глобальным максимумом). Для нахождения аналитического оптимума необходимо найти стационарные или критические точки функции, использовать какое-либо из достаточных условий экстремума (изменение знака первой производной или определенного знака четной производной) и сравнить значения функции в точках локальных экстремумов и граничных значениях. Однако в прикладных задачах нередки ситуации, когда трудно вычислить производные функции (например, если функция не задана в аналитическом виде). Более того, не исключено, что значения функции известны или могут быть вычислены только в отдельных точках. В таких ситуациях использование необходимого и достаточного условий локального минимума невозможно и следует применять другие методы решения задачи оптимизации. Методы минимизации функции одного переменного, в которых используют значения функции в точках рассматриваемого промежутка и не используют значения ее производных, называют методами прямого поиска. Можно выделить две группы методов прямого поиска, соответствующие двум принципиально различным ситуациям: 1) все N точек хк., к = 1,…, N, в которых будут вычислены значения функции, выбирают заранее (до вычисления функции в этих точках); 2) точки xk выбирают последовательно (для выбора последующей точки используют значения функции, вычисленные в предыдущих точках). 16 В первом случае поиск наименьшего значения называют пассивным, а во втором — последовательным. Так как в прикладных задачах вычисление каждого значения функции может быть достаточно трудоемким, то целесообразно выбрать такую стратегию поиска, чтобы значение целевой функции f (x*) заданной точностью было найдено наиболее экономным путем. Будем считать, что стратегия поиска определена, если:  определен алгоритм выбора точек x k ., к = 1,…, N;  определено условие прекращения поиска, т.е. условие, при выполнении которого значение f (x*) считают найденным с заданной точностью. Оптимальный пассивный поиск состоит в выборе точек, равномерно расположенных на отрезке [a, b] , координаты которых (b  a) k xk  , k  1, 2,..., N . N 1 2(b  a) При этом длина l N интервала, содержащего минимум l N  N 1 дает оценку скорости сходимости пассивного поиска с ростом числа N точек, так как скорость сходимости любого метода прямого поиска можно характеризовать скоростью уменьшения интервала неопределенности с возрастанием N. Вычисляем значения целевой функции в каждой точке f k  f ( x k ) , найдем наименьшее значение f j и, тогда точка оптимума   x*  x j 1 , x j 1 и можно считать, что x*  x j   с точностью ε. Возникает вопрос о сходимости метода, далее говоря о сходимости, будем считать последовательность приближенных решений x k сходится к оптимальному x * , если для любого сколь угодно малого положительного ε, найдется номер n , такой что для любого k  n выполнено x k  x *   , p при этом, если существует M  0 такая что x k 1  x *  M x k  x * , метод сходится со скоростью p. 2. Методы последовательного поиска Предположение об унимодальности функции позволяет разработать адаптивные (или последовательные) алгоритмы поиска оптимального значения). При этом используется понятие интервала неопределенности, как интервала содержащего точку минимума, хотя ее точное положение неиз17 вестно. Важной характеристикой данных методов является количество вычислений значений оптимизируемой функции для получения конечной величины интервала неопределенности, не превышающей заданной величины. Методы ориентированы на нахождение точки оптимума внутри заданного интервала и основаны на свойстве унимодальности функции и не используют информацию о производной функции. В алгоритмах этих методов вычисляются значения функции в промежуточных точках  k и  k (k   k ) интервала неопределенности: если f (k )  f ( k ) , то в качестве границ нового интервала рассматривается интервал [k , bk ] (рис. 2) если f (k )  f ( k ) это будет интервал [a k ,  k ] (рис. 3) если f (k )  f ( k ) , то оставим интервал [k ,  k ] . f ( k ) f ( k ) f ( k ) f ( k ) ak k k Рис. 2 bk ak k k bk Рис. 3 Такой подход позволяет построить последовательность вложенных интервалов, удовлетворяющей теореме о вложенных отрезках. Теорема1. Пусть даны монотонно возрастающая переменная xn , и монотонно убывающая переменная yn , причем всегда xn  yn . Если их разность xn  yn стремится к нулю, то обе переменные имеют общий конечный предел: lim xn  lim yn  c . n  n  То есть все вложенные интервалы имеют одну общую точку, где достигается минимум функции. Что обеспечивает сходимость методов ин- 1 Доказательство теоремы см. Фихтенгольц Г. М. Основы математического анализа. В 2ч.: Учеб. для вузов. Ч.1, с.96. 18 тервалов. Все алгоритмы различаются только способом определения промежуточных точек  k и  k . Метод дихотомии Идея метода состоит в вычислении на каждой очередной итерации двух значений целевой функции в точках, отстоящих на величину  в обе стороны от середины интервала неопределенности. Величина  в этом методе называется константой различимости, такова, что, с одной стороны, величина 2 была близка к желаемому конечному значению интервала неопределенности, с другой, значения оптимизируемой функции на краях интервала 2 были различимы. Введем обозначения:  – конечная длина интервала, определяющая точность, [a1 , b1 ] – начальный интервал неопределенности, k – номер итерации. Пусть на к-м шаге x* a k , bk  . На этом отрезке ak , bk  выберем a  bk a  bk  , k  k   (рис. 4). Вычислим значедве точки  k  k 2 2 ния в этих точках. f (k ), f ( k ) и выполним процедуру исключения отрезка, то есть если f (k )  f ( k ) , то a k 1  k , bk 1  bk , иначе a k 1  a k и bk 1   k . Таким образом, построим новый отрезок ak 1, bk 1   [ak , bk ] . Длина интервала неопределенности после k-ой итерации 1 1   bk 1  a k 1  (b1  a1 )  2 1   k k  2 2  Отсюда можно вычислить число итераций для достижения необходимой ln (b1  a1 ) /   точности  , потребуется n  итераций. На каждой итераln 2 ции минимизируемая функция вычисляется дважды. Метод деления пополам. На каждой итерации исключается половина интервала. Пусть на к-м шаге x* a k , bk  , длиной lk (рис. 5). Найдем середину интервала a  bk x k , вычислим f (x ) . Найдем промежуточные точки 2 x1  ak  lk / 4 , x2  bk  lk / 4 и значения f ( x1 ), f ( x 2 ) . Далее проведем процедуру исключения интервала: если f ( x1)  f ( x ) то исключается ин19 тервал ( x , bk ) , при этом bk 1  x , ak 1  ak , x  x1 ; если f ( x2 )  f ( x ) , то исключается интервал (a, x ) , при этом ak 1  x , bk 1  bk , x  x2 . Средняя точка последовательности получаемых интервалов всегда совпадает с одной из пробных точек x1 , x 2 , или x , найденных на предыдущих итерациях. a1 λk x* μk bk a Рис. 4 x1 x* x2 b Рис. 5 На каждой итерации требуется не более 2-х вычислений значений функции. После N вычислений длина интервала равна l  1 / 21 / N длинны исходного интервала. Таким образом, для достижения необходимой точноln 2 сти нужно провести N   итераций. ln  Метод золотого сечения. Идея метода состоит в использовании на каждой итерации для сокращения интервала неопределенности одной из внутренних точек предыдущей итерации. Чтобы выполнить процедуру исключения отрезка на этом шаге, отрезок ak , bk  необходимо двумя внутренними точками  k и  k разделить на три части, при этом должны быть выполнены следующие условия:  Пробные точки на каждой итерации находятся на одинаковых расстояниях от концов интервала неопределенности k  a k  bk   k  Для новой итерации точки k 1,  k 1 выбираются так, чтобы  k 1 совпало с  k либо  k 1 совпало с  k .  Сжатие интервала неопределенности осуществляется на каждой итерации с одним и тем же коэффициентом сжатия , удовлетворяющее уравl нению  2    1  0 . При этом отношение длин отрезков k сохраняl k 1 20 l  k 1    const . Число τ называют l k 1 l k  2 сечения, которое можно вычислить как ется от шага к шагу, то есть отношением золотого lk 5 1  0,61803... 2  При выполнении этих условий  k и  k вычисляют по формулам (1) k  a k  (1   )bk  a k ,  k  a k   bk  a k , b  a k 1 здесь   k 1 , определяется из условий выше. bk  a k Рассмотрим к-й шаг последовательного поиска. Пусть x* a k , bk  , вычисляем  k и  k по формулам (1) и соответствующие значения функции. Если f (k )  f ( k ) , то положить a k 1  k , bk 1  bk , k 1   k ,  k 1  a k 1   (bk 1  a k 1 ) , если же f (k )  f ( k ) . Положить a k 1  a k , bk 1   k ,  k 1  k , k 1  a k 1  (1   )(bk 1  a k 1 ) . Если исходный интервал имеет единичную длину, длина интервала по-  сле N вычислений равна  N 1 , иначе  n1  ln ности  потребуется n   b1  a1 . Для достижения точ- b1  a1  5 1 ln 2 итераций. Метод Фибоначчи. Метод аналогичен методу золотого сечения. Отличие состоит в том, что коэффициент сжатия интервала неопределенности меняется от итерации к итерации согласно последовательности Фибоначчи. Последовательность чисел определяется следующим образом: F0  F1  1, Fk 1  Fk  Fk 1 , k  1, 2, ... Пробные точки  k и  k вычисляют по формулам F F (2)  k  a k  n  k 1 bk  a k ,  k  a k  n  k bk  a k . Fn  k 1 Fn  k 1 21 При этом число итераций выбирается до начала вычислений и обуb  a1 словлено требуемой точностью   1 . Fn На k-м шаге имеется интервал (a k , bk ) , вычисляем  k и  k по формулам (2). Если f (k )  f ( k ) , то a k 1  k , bk 1  bk , k 1   k , F  k 1  a k 1  nk 1 (bk 1  a k 1 ) , если f (k )  f ( k ) , то a k 1  a k , Fnk F bk 1   k ,  k 1  k ,  k 1  a k 1  nk 2 (bk 1  a k 1 ) , Fnk В качестве оценки скорости сходимости методов прямого поиска можно использовать скорость убывания длины интервала неопределенности в зависимости от числа N вычисленных значений минимизируемой функции в различных методах. Так для метода золотого сечения и метода Фибоначчи отношение длин интервалов неопределенности при больших N равен l n (зол.сеч)  2 lim   1,17 . Таким образом, скорость сходимости метода n l n (Фиб ) 5 Фибоначчи при больших значениях n всего примерно на 17% выше, чем скорость сходимости метода золотого сечения. Сравнивая при больших значениях п методы золотого сечения и дихоl (зол.сеч) томии, получаем lim n  0 . Таким образом, метод золотого сеn l n (дих) чения качественно „лучше" метода дихотомии, но скорость сходимости метода дихотомии при больших значениях п выше, чем скорость сходимости метода оптимального пассивного поиска. 3. Методы аппроксимации В методах прямого поиска мы не имели никакой информации о минимизируемой функции за исключением ее значений в выбранных нами точках и предположения, что она непрерывна и является унимодальной функцией на рассматриваемом отрезке. Если функцию в некоторой окрестности точки ее минимума можно достаточно точно заменить (аппроксимировать) многочленом, то для ее минимизации целесообразно использовать так называемые методы полиномиальной аппроксимации. Их общая особенность состоит в вычислении коэффициентов многочлена по известным значениям функции в отдельных точках и последующем нахождении минимума этого многочлена с использованием необходимых и достаточных 22 условий экстремума Основная идея метода: возможность аппроксимации гладкой функции полиномом достаточно высокого порядка и использование этого полинома для оценивания точки оптимума. Качество этой оценки может быть повышено двумя способами:  увеличением степени полинома;  уменьшением интервала аппроксимации. Второй способ предпочтительнее, так как построение полинома порядка более 3 – достаточно сложная задача, а сужение интервала для унимодальной функции – достаточно простая. Использование квадратичной аппроксимации для нахождения оптимума. Чтобы функция имела минимум внутри отрезка она должна быть, по крайней мере квадратичной. Для построения квадратичной функции достаточно трех точек: M 1 ( x1 , y1 ), M 2 ( x 2 , y 2 ), M 3 ( x3 , y3 ), . Можно задать аппроксимацию функции полиномом вида: P2 ( x)  a 0  a1 ( x  x1 )  a 2 ( x  x1 )( x  x 2 ) , необходимо выбрать a 0 , a1 и a 2 так, чтобы P2 ( x1 )  y1 , P2 ( x 2 )  y 2 , P2 ( x3 )  y3 . Отсюда следует, что  y  y1 y 2  y1     3  x  x x  x 1 2 1  3 Найдѐм стационарную точку x* полинома P2 ( x) : a 0  y1 , a1  y 2  y1 1 , a2  x3  x 2 x 2  x1 (3) x 2  x1 a (4)  1 2 2a 2 Так как функция унимодальна на рассматриваемом интервале и полином P2 ( x) тоже унимодальная функция, то x* является приемлемой оценкой истинного оптимума. Метод Пауэлла основан на последовательном применении процедуры оценивания с использованием квадратичной аппроксимации. Метод квадратичной аппроксимации удобно применять после локализации точки минимума методом золотого сечения или методом Фибоначчи. Это объясняется тем, что для дважды дифференцируемой функции многочлен второго порядка достаточно хорошо аппроксимирует функцию в окрестности точки минимума. x*  23 Вопросы и задания для самопроверки 1. Сформулируйте теоремы Ферма и Лагранжа, напишите формулу конечных приращений. 2. Что называют точкой строгого локального экстремума функции одного переменного? Сформулируйте необходимые и достаточные условия экстремума такой функции. В чем различие между локальным экстремумом и наибольшим (наименьшим) значением этой функции на отрезке? 3. В каких точках отрезка линейная функция достигает своих наибольшего и наименьшего значений? Как найти точку экстремума квадратного трехчлена в интервале? 4. Что называют сходимостью метода вычислений и порядком его сходимости? Запишите условия, при выполнении которых скорость сходимости метода является линейной, сверхлинейной, квадратичной, кубической. 5. В чем отличие метода пассивного поиска и метода интервалов. 6. Проведите сравнительный анализ методов дихотомии, деления пополам, золотого сечения, Фибоначчи, Пауэлла, Какова скорость сходимости методов? 7. Имеет ли функция f ( x)  xe  x экстремум в интервале (0,3)? Если имеет, то в какой точке? Имеет ли она минимум в том же интервале, минимум на отрезке [0, 3], и если да, то в какой точке? 8. Проверьте, являются ли унимодальными следующие функции: а) f ( x)  x 2  2 x  1 на отрезках [0, 2], [1,5, 2]; б) f ( x)  x  1  1 на отрезках [-3,3], [-3,1], [1,3], [0,2]. 1  4 x0  x sin 9. Имеет ли функция f ( x)   минимум в точке х = 0, выx 0, x0 полняется ли в этой точке необходимое, достаточное условия экстремума? 24 Лекция 3 Методы одномерной оптимизации с использованием информации о производной План Метод средней точки Метод Ньютона (метод касательной) Метод секущей Метод кубической аппроксимации Сравнение методов оптимизации В методах прямого поиска при вычислении значений минимизируемой функции f(x) неизбежно возникают погрешности, к которым чувствительны алгоритмы прямого поиска, основанные на сравнении значений функции в отдельных точках. При этом погрешность вычисления оптимального 2 f решения x * порядка , где  f – погрешность вычисления целевой f " ( x*) функции. Поэтому можно использовать для вычисления оптимального решения задачи значения производных целевой функции. Теорема Ферма. Пусть f (x) дифференцируема на интервале и достигает своего наименьшего (наибольшего) значения на этом интервале в некоторой внутренней точке x * этого интервала, тогда f ' ( x*)  0. Таким образом, если унимодальная функция f(x) непрерывно дифференцируема на отрезке минимизации, то точку х* наименьшего значения функции можно вычислять как корень уравнения f ' ( x)  0 с помощью тех или иных методов численного решения нелинейных уравнений. В этом случае на точность решения задачи решающее влияние оказывает погрешность вычисления производной функции. Рассмотрим некоторые методы одномерной минимизации, основанные на использовании производной минимизируемой функции. Метод средней точки. Будем искать минимум функции f (x) непрерывно дифференцируемой и строго унимодальной на отрезке [a, b] . В этом случае единственной точкой x* [a, b] минимума будет стационарная точка, в которой f ' ( x*)  0 . 25 Отметим, что непрерывно дифференцируемая унимодальная на отрезке функция может иметь на нем более одной стационарной точки. На отрезке определяются две точки a k , bk , в которых производные имеют разные знаки, f ' (a k ) f ' (bk )  0 . Искомый оптимум находится между ними. Делим a  bk интервал пополам: x k  k , если f ' ( x k )  0 ( 0) , то из двух интерва2 лов оставляем тот, на концах которого производная имеет разные знаки. Метод средней точки напоминает метод дихотомии, но сходится к искомому значению х* быстрее. Поскольку для метода средней точки после вычисления n значений производной минимизируемой на отрезке [0, 1] 1 функции f(x) для длины интервала неопределенности получаем l n  . 2n Таким образом, для одинакового уменьшения значения l n в методе средней точки нужно вычислить вдвое меньше значений производной функции по сравнению с числом значений самой функции в методе дихотомии. Метод Ньютона (метод касательной) Если строго унимодальная на отрезке [a, b] функция f(x) дважды непрерывно дифференцируема на этом отрезке, то точку x* [a, b] минимума этой функции можно найти путем решения уравнения f ' ( x)  0 методом Ньютона, иногда называемым методом касательных. Выбираем х 0 начальное приближение называемое обычно начальной точкой. Линеаризуем функцию f ' ( x) в окрестности начальной точки, приближенно заменив дугу графика этой функции касательной в точке x0 , f ' ( x0 ) . Уравнение касательной имеет вид f ' ( x)  f ' ( x0 )  f " ( x0 )  x  x0  . Выберем в качестве следующего приближения к х* точку x1 пересечения касательной f ' ( x0 ) с осью абсцисс (рис. 1). Получаем первый элемент x1  x 0  итеf " ( x0 ) рационной последовательности {x k } . На (k+1)-м шаге по найденной на предыдущем шаге точке x k можно найти точку f ' ( xk ) (1) x k 1  x k  f " ( xk ) В общем случае сходимость метода Ньютона существенно зависит от выбора начальной точки х 0 . Для надежной работы этого метода необхо26 димо, чтобы вторая производная f " ( x) в некоторой окрестности искомой точки х* сохраняла знак, а начальная точка х 0 выбиралась из такой окрестности. В противном случае второе слагаемое в правой части (1) может стать неограниченным. Рис. 1 Поскольку для дважды непрерывно дифференцируемой функции в точке минимума f " ( x*)  0 , то должно быть и f " ( x0 )  0 . Поэтому говорят, что метод Ньютона обладает локальной сходимостью в том смысле, что надо выбрать хорошее начальное приближение, попадающее в такую окрестность точки х*, где f " ( x)  0 . Однако проверка выполнения этого условия не всегда возможна. Достаточным условием монотонной сходимости метода Ньютона будут постоянство в интервале между точками х 0 и х* знака производной f'"(x) и совпадение его со знаком f ' ( x) . Теорема. Если f (x) трижды дифференцируема в окрестности точки х*, тогда релаксационная последовательность (1) сходится к х* и имеет место оценка x k 1  x *  2 max f ' ' ' ( x) U ( x *) min f " ( x) xk  x * . 2 U ( x*) Доказательство. Рассмотрим формулу Тейлора для f ' ( x*) в окрестности точки хk, учитывая, что f ' ( x*)  0 f ' ' ' ( ) x *  x k 2 ,   x k  x * 0  f ' ( x*)  f ' ( x k )  f " ( x k )x *  x k   2 Откуда получаем f ' ' ' ( )  x *  x k 2 f " ( x k )( x *  x k )   f ' ( x k )  (2) 2 27 Рассмотрим отношение разностей приближенного и точного решения на kм и (k+1)-м шаге, подставив (1) и (2) f ' ( xk ) x *  xk  x *  x k 1 f " ( xk ) f ' ( xk )  1 , x *  xk x *  xk f " ( x k )( x *  x k ) f ' ' ' ( ) Откуда получаем, что если  0 , то f ' ( xk ) x *  x k 1 2 1 1, f ' ' ' ( ) x *  xk 2 2 f ' ( xk )  ( x *  xk ) f ' ( xk ) следовательно, последовательность x k  монотонно убывает и ограничена и lim x k  x *  0 . Далее оценим скорость сходимости метода, из (2) полуk  чаем f ' ( xk ) f ' ' ' ( xk )  x *  x k 2 ,  ( x *  x k )  f " ( xk ) 2 f " ( xk ) тогда x k 1  x *  x k  x *  f ' ( xk ) f ' ' ' ( ) x k  x *2 .  f " ( xk ) 2 f " ( xk ) Из последнего равенства следует оценка теоремы. Метод секущих Метод секущих похож на метод Ньютона, но строится не касательная к графику производной, а секущая. Геометрическая интерпретация этого метода (рис. 2) состоит в том, что в качестве очередного приближения x k 1 выбирают точку пересечения с осью абсцисс секущей к графику функции f ' ( x) , то есть f ' ( x k )x k  x k 1  (3) x k 1  x k   f ' ( x k )  f ' ( x k 1 ) Таким образом, мы избегаем вычислений производной второго порядка в f ' ( x k )  f ' ( x k 1 ) формуле (1), заменяя ее конечной разностью f " ( x k )  . x k  x k 1 Метод секущих можно считать модификацией метода Ньютона. Выбор начальной точки х0 проводят следующим образом. Если на отрезке [a, b] функция f(x) имеет знакопостоянную третью производную 28 f ' ' ' ( x), то в качестве х 0 выбирают тот конец отрезка [a, b] , на котором совпадают знаки производных первого и третьего порядaf ' (b)  bf ' (a) ков: f ' ( x)  f ' ' ' ( x)  0 , тогда x1  . Точка x1 – точка пересеf ' (b)  f ' (a) чения с осью абсцисс хорды, стягивающей дугу графика функции f ' ( x) на отрезке [a, b] (рис 2). Таким образом, первый шаг метода секущих выполняют согласно методу хорд, выбран стартовый отрезок [ x0, x1 ] для после- дующих точек а последующие шаги — в соответствии (3). Рис. 2 Этот метод имеет сверхлинейную  скорость x *  x k  C x *  x k 1 , где С – const, а   сходимости, причем 1 5  1,618... — отношение 2 золотого сечения. Модификации метода Ньютона обладают только локальной сходимостью, т.е. сходятся, если начальная точка выбрана в некоторой окрестности точки минимума функции. Если же начальная точка расположена слишком далеко от точки минимума, то подобный метод может расходиться или «зацикливаться». В отличие от метода средней точки метод секущих использует информацию не только о знаке производной, но и о значениях в пробных точках. Метод кубической аппроксимации В методах полиномиальной аппроксимации при построении многочлена, аппроксимирующего минимизируемую функцию в окрестности искомой точки х* минимума, помимо значений функции в отдельных точках могут быть использованы и значения ее производных. 29 Пусть для непрерывно дифференцируемой функции f (x) , строго выпуклой на отрезке [ x1 , x 2 ] известны значения y1  f ( x1 ), y 2  f ( x 2 ) функции и значения ее производной y11  f ' ( x1 ), y12  f ' ( x 2 ) . Если y11  y12  0 , то рассматриваемая функция строго унимодальна на этом отрезке. Рассмотрим метод поиска точки х* при условии y11  y12  0 , называемый методом кубической аппроксимации, поскольку в этом случае на отрезке [х1, х2] можно построить единственный многочлен третьей степени, располагая значениями функции и производных на концах этого отрезка. Этот многочлен, называемый кубическим интерполяционным многочленом Эрмита, можно преобразовать к виду P3 ( x)  a0  a1 ( x  x1 )  a 2 ( x  x1 )(x  x 2 )  a3 ( x  x1 ) 2 ( x  x 2 ), y  y1 y  y1 y11 , a2  2 , a 0  y1 , a1  2  x 2  x1 ( x 2  x1 ) 2 x 2  x1 a2  y12  y11 2 (4) (5) y 2  y1 (6) ( x 2  x1 ) 2 ( x 2  x1 ) 3 Из необходимого условия P3 ( x)  0 экстремума этого многочлена с учетом коэффициентов (5), (6) получаем квадратное уравнение y  y  3a1 3a 3 ( x  x1 ) 2  2 11 12 ( x  x1 )  y11  0 , x 2  x1 решение которого представим в виде (7) x0  x1   ( x 2  x1 ) ,   z  y11 y  y1 , z  y11  y12  3 2 ,   z 2  y11  y12 . 2  y11  y12 x 2  x1 Покажем, что   (0, 1) и, следовательно x0  ( x1 , x 2 ) . где   Так как f ' ( x1 )  y11  0, и f ' ( x 2 )  y12  0 , то   z 2  y11  y12  z , отсюда   z  y11  0 и 2  y11  y12  0 . Учитывая эти неравенства   0 .   z  y11 Рассмотрим   , добавим в числитель   z  0 , тогда 2  y11  y12    z  y11    z 2  y11   1 , то есть 0    1 . 2  y11  y12 2  y11  y12 Если f ' ( x0 )  0 , то x0  x * — искомая точка минимума функции f(x) на отрезке [х1, х2] при этом f ( x*)  f ( x1 ) и f ( x*)  f ( x 2 ) . Если же f ' ( x0 )  0 30 то оставляем меньший отрезок [х0, х2] или [х1, х0], такой чтобы f ' ( x0 ) y1 j  0, j  1,2 (то есть на концах отрезка производные имели разные знаки) и продолжать описанным выше способом поиск точки минимума на новом отрезке по формулам (5) – (6). После каждого приближения правильность вычислений подтверждается уменьшением минимального значения многочлена по сравнению с его минимальным значением на предыдущем шаге. Вычисления можно прекратить, когда длина интервала неопределенности, в котором гарантированно находится искомая точка х*, станет меньше заданной наибольшей допустимой величины  . Сравнение методов одномерной оптимизации Для быстрого получения предварительных результатов (начальной точки для применения других методов), а также, если требуется надѐжная работа алгоритма при неизвестной заранее целевой функции, лучше использовать методы исключения интервалов. Если требуется точное решение, необходимо воспользоваться градиентными методами (особенно кубической аппроксимацией). С другой стороны, если требуются высокая точность, но функция не задана аналитически, лучше пользоваться методами точечного оценивания, так как при использовании градиентных методов накапливается погрешность при конечно-разностной аппроксимации производных. Если сравнить методы с точки зрения поставленной задачи и вида функции, то при минимуме информации о функции следует использовать метод исключения интервалов. Если функция квадратичная или близка к таковой, то следует использовать метод Пауэлла Если функция дважды дифференцируемая, непрерывная и задана аналитически, то следует использовать градиентные методы. Методы точечного оценивания при прочих равных условиях (интервалы, гладкая функция) быстрее методов исключения интервалов. Вопросы и задания для самопроверки 1. Сформулируйте теоремы Ферма и Лагранжа, напишите формулу конечных приращений. 2. Что называют точкой строгого локального экстремума функции одного переменного? Сформулируйте необходимые и достаточные условия экстремума такой функции. В чем различие между локальным экстре31 3. 4. 5. 6. 7. 8. 9. мумом и наибольшим (наименьшим) значением этой функции на отрезке? Как проверить, является ли функция одного действительного переменного выпуклой (строго выпуклой) вниз (вверх) функцией? Сколько экстремумов может иметь выпуклая (строго выпуклая) функция одного переменного на отрезке? Напишите формулу Тейлора с остаточным членом в форме Лагранжа и Пеано для функции действительной переменной. Сравните методы средней точки и метод деления пополам. Каков порядок сходимости для каждого? Сравните методы Ньютона и секущих. Каков порядок сходимости каждого? В чем состоит идея метода кубической аппроксимации. Какие многочлены называются интерполяционным многочленом Лагранжа, Эрмита? В чем их отличие? Проверьте, является ли функция f ( x)  100(1  x 2 ) 2  (1  x) 2 унимодальной. Чему равно значение первого приближения метода Ньютона, если начальное приближение x0  0 . 10. Найдите два приближения минимизации функции f ( x)  ( x  1) sin x 2 на отрезке [2, 3] методом секущих, найти минимум функции методом кубической аппроксимации, сравните результат. 32 Лекция 4 Задачи многомерной оптимизации План 1. Постановка задач многомерной оптимизации 2. Теоремы об оценках сходимости методов спуска 1. Постановка задач многомерной оптимизации При математической формулировке задачи условной оптимизации целевую функцию выбирают с таким знаком, чтобы решение задачи соответствовало поиску минимума этой функции. Поэтому формулировку общей задачи математического программирования обычно записывают так:  (1) f ( x )  min , Q где Q  R N — множество возможных альтернатив, рассматриваемых  при поиске решения задачи. Любую точку x  Q называют допустимым решением задачи математического программирования, а само множество — множеством допустимых решений или, короче, допустимым множе ством. Точку x*, в которой функция f (x ) достигает своего наименьшего значения, называют оптимальным решением задачи. При отсутствии ограничений множество Q совпадает с областью определения целевой функции. Если же рассматриваемые альтернативы должны удовлетворять некоторым ограничениям, то множество допустимых решений сужается. Задачу (1) в дальнейшем будем называть задачей минимизации целевой функции на множестве Q. Но целевая функция может и не достигать на Q наименьшего значения. Тогда говорят о точной нижней грани функции  f (x ) на этом множестве и вместо (11) используют запись  f ( x )  inf, x  Q. (2) Отличие (1) от (2) в том, что в первом случае предполагают существование точки x * , в которой целевая функция достигает своего наименьшего значения на множестве Q, а во втором случае такая точка может и не существовать. Поэтому решение общей задачи математического программирования состоит в том, чтобы в первом случае найти точные (или с некоторой заданной точностью) значения координат точки x*  {x j } j 1,..,N и зна33 чение целевой функции f ( x*)  min f ( x) . Во втором случае построить таQ кую последовательность точек xn nN , которой бы соответствовала последовательность f ( x n ) , сходящаяся к значению inf f ( x) , и вычислить Q это значение с заданной точностью. Отметим, что в большинстве прикладных задач имеет место первый случай, поэтому использование записи вида (2) будем оговаривать особо. Сформулируем задачу многомерной безусловной оптимизации: найти минимум функции f (x) , где x  R , при отсутствии ограничений на x , при этом f (x) – это скалярная целевая функция, непрерывно дифференцируемая. При решении этого класса задача необходимо учитывать следующие факторы:  характер целевой функции решаемой задачи (одноэкстремальная или многоэкстремальная);  возможность получения в процессе оптимизации информации о производных целевой функции (возможность присутствует или отсутствует);  наличие различных подходов к организации итеративной процедуры поиска оптимума (методы, основанные на итеративном движении переменных в направлении, определяемом тем или иным способом). Все методы решения задач безусловной оптимизации состоят в том, N   что мы строим последовательность точек x (n) таким образом, чтобы по- следовательность функций f ( x (n) ) была убывающей (т. е. спускаемся  вдоль функции). На k-м шаге ( k  0 ) определяем вектор S k , в направлении которого функция f ( x ) уменьшается. В этом направлении делаем шаг ве личиной  k и получаем новую точку x k 1  x k  k S k , в которой  nN , удовлетворяющая этому f ( x k 1 )  f ( x k ). . Последовательность x (n) условию, называется релаксационной последовательностью, а соответствующие методы – методами спуска. Методы решения делятся на методы с использованием информации о производных функции и без использования таковой. Различные методы спуска отличаются выбором направления и величины шага. Как правило, для нахождения  k используется процедура одномерного поиска. 34 Численное решение задач безусловной минимизации функций многих переменных, как правило, значительно сложнее, чем решение задач минимизации функций одного переменного. В самом деле, с ростом числа переменных возрастают объемы вычислений и усложняются конструкции вычислительных алгоритмов, а также более сложным становится анализ поведения целевой функции. Методы численного решения задач многомерной безусловной минимизации многочисленны и разнообразны. Условно их можно разделить на три больших класса в зависимости от информации, используемой при реализации метода. 1. Методы нулевого порядка, или прямого поиска, стратегия которых основана на использовании информации только о свойствах целевой функции. 2. Методы первого порядка, в которых при построении итерационной процедуры наряду с информацией о целевой функции используется информация о значениях первых производных этой функции. 3. Методы второго порядка, в которых наряду с информацией о значениях целевой функции и ее производных первого порядка используется информация о вторых производных функции. 2. Теоремы об оценках сходимости методов спуска При анализе сходимости релаксационной последовательности удобно рассматривать невозрастающую последовательность { k }kN , где  k  f k  f * . Здесь f k  f ( x (k ) ) , f * минимум (инфимум) функции. Обозначим f k  grad f ( x (k ) ) - градиент функции f(x) в точке x k . Для оценки сходимости этой последовательности используют следующие утверждения. Лемма 1.2 Если для элементов последовательности { k }kN выполнены условия  k 1   k   k  k21 ,  k  0,  k  0 , то имеет место оценка n  0 n . 1  0   k k 1 2 Доказательсво леммы см. Аттетков А. В. Методы оптимизации : учеб. для вузов, с.165 35 Теорема 13. Пусть минимизируемая функция f (x) выпукла и диффе-   ренцируема на множестве D f , последовательность x (n) является релакf k 1  f k сационной. Если  k  , то имеет место оценка 2 2 ( k 1) f k 1 x  x* n  f0  f * . n 1   f 0  f *   k k 1 Лемма.2.4 Если для элементов последовательности { k }kN выполнены условия  k 1   k   k  k 1 ,  k 1  0,  k  0 , то имеет место оценка   n   0 exp   n     k  .  k 1 Теорема 2.5 Пусть минимизируемая функция f (x) выпукла и диффе-   ренцируема на множестве D f ,, множество X 0  x  D f : f ( x)  f ( x (0) )   ограничено, последовательность x (n) k  является релаксационной. Если f k 1  f k , то справедлива оценка diam X 0 f k 1  n  f n  f *   f 0  f * exp    k  .  k 1  Практически можно задать некоторое число   0 , связанное с выбранной точностью вычислений, и проводить итерации до тех пор, пока не f  fk будет нарушено неравенство k 1   , а затем либо повысить точ2 f k 1 ность вычислений, либо перейти к другому методу спуска. Поиск точки х* обычно прекращают при выполнении на к-й итерации хотя бы одного из неравенств: x (k )  x (k 1)   1 , f k  f k 1   2 , где  1 ,  2 — заданные достаточно малые положительные числа, называемые 3 Доказательство теоремы см. Аттетков А. В. Методы оптимизации : учеб. для вузов, с.166 Доказательсво леммы см. Аттетков А. В. Методы оптимизации : учеб. для вузов, с.167 5 Доказательсво теоремы см. Аттетков А. В. Методы оптимизации : учеб. для вузов, с.168 4 36 параметрами точности поиска. В случае минимизации дифференцируемой функции необходимым условием того, что х* — точка ее минимума, является равенство f ( x*)  0. При этом условием прекращения поиска на k-й итерации может быть выполнение неравенства f ( x (k 1) )   3 .   Далее строим релаксационную последовательность x ( k )   k 1, 2,... , влетворяющую соотношению x ( k 1)  x ( k )  hk f x ( k ) , где hk  удо-  f ( x (k ) , ) а   0 . Процедура построения таких точек называется градиентной процедурой поиска min f ( x) . Теорема 3. Если f (x) дифференцируемая ограниченная снизу функция, градиент которой удовлетворяет условию Липшица: 2 f ( x)  f ( y)  M x  y , а шаг 0    hn  ,  2  0 , то для граM  2 2 диентной процедуры выполняются условия: f ( x (n) ) убывает и lim x n1  x n  0 , (1) n lim f ( x ( n) )  0 , (2) n если f (x) - выпуклая, то lim f ( x n )  f ( x * )  min f ( x) . (3) n Доказательство. Используем интегральное представление приращения 1 функции, то есть f ( x  y )  f ( x)   f ( x  ty) ydt . Обозначим f n  f ( x (n) ) , f n  f ( x (n) ) , тогда получаем f n1  f n  f ( x ( n)  hn f n )  f ( x ( n) 1 )  hn  f ( x ( n)  thn f n )f n dt. Добавим и вычтем f n в первый множитель подынтегрального выражения, тогда 1   f n1  f n  hn  f n  f n  f ( x ( n)  thn f n ) f n dt  1 2 1    hn  f n dt  hn  f n  f ( x ( n)  thn f n f n dt  37  hn f n 2 1  hn  f n  f ( x ( n)  thn f n )  f n dt. Используя условие Липшица для разности градиентов в последнем слагаемом, получаем t M 1  2 2 2 2 f n1  f n  hn f n  hn M f n  t dt    hn2 f n .  2 hn      Из условия x ( k 1)  x ( k )  hk f x ( k ) следует x ( k 1)  x ( k )  hk f x ( k ) , тогда x ( n1)  x ( k ) Из 2    hn2 f x ( n) hn  неравенства f n1  f n   2 x ( n1)  x ( n) 2 2 . 2 M  2 2 следует M 1    2 , 2 hn тогда 0 Следовательно, градиентная процедура обеспечивает убывание целевой функции f ( x ( n1) )  f ( x ( n) ) . Из f n1  f n   2 x ( n1)  x ( n) соотношения 2 следует что 2 f n  f n1   2 x ( n1)  x ( n) . Так как f (x) ограничена снизу и последова- тельность f n убывает, то  lim f n   , последовательность  f n nN фунn даментальна и f n  f n1  0 , тогда x ( n1)  x ( n) 2  f n  f n1 2  0, при n   . Следовательно lim x n1  x n  0 n M 1  2 Так как f n1  f n    hn2 f n ,  2 hn   1 M 2 2 то f n  f n 1    hn f n , 2  hn Из f n 2 условия  f n  f n1  2   12 0    hn , hn  1 M   1 hn 2 , получаем, что , следовательно lim f ( x ( n) )  0 n Если функция выпуклая, то градиентный процесс сходится к минимуf (x) - выпуклая функция, следовательно му. Так как 38 f ( x)  f ( x ( n) )  f ( x ( n) )( x  x ( n) ) , x  X . Это неравенство справедливо и в некоторой фиксированной точке x *  lim x ( n) : n ( n) f ( x*)  f n  f n  ( x *  x (n) ) , отсюда f n  f ( x*)  f n  ( x то есть  x*) . Рассмотрим f n  ( x ( n)  x*)  f n  ( x *  x (n)  x ( n1)  x ( n1) )      x *  x ,     f n  x *  x ( n1)  f n  x ( n)  x ( n1)  f n  x ( n)  x ( n1)    1 ( n) x  x ( n1) hn так как x ( n1) x ( n1) ( n)  , следовательно f x   hn f x Таким образом, получаем:  ( n)  f n  f ( x*)  f n  x ( n)  x ( n1)  ( n)   x ( n1)  x ( n)   hn  1 ( n) x  x ( n1)  x *  x ( n1)  hn  1  x ( n)  x ( n1)   f n  x *  x ( n1) hn   1  x ( n1)  x ( n)   f n  x *  x ( n1) 1         0  так как f n  0 , x *  x ( n1) - ограничены и lim x n1  x n  0 n Следовательно lim f n  f ( x*) . В силу того, что последовательность n f n убывает и сходится к f (x*) , то lim f n  f ( x*)  min f ( x (n) ) . Теорема n доказана. Вопросы и задания для самопроверки 1. Что такое градиент функции многих переменных, матрица Гессе? Запишите приращения дифференцируемой и дважды дифференцируемой функции, используя эти понятия. 2. . Что такое производная функции многих переменных по направлению вектора и как она связана с градиентом функции? Имеет ли дифференцируемая функция многих переменных производные по всем направлениям? Верно ли обратное? 39 3. В чем различие между точкой экстремума и критической или стационарной точками скалярной функции многих переменных? Что называют строгим (нестрогим) локальным экстремумом такой функции? 4. Сформулируйте необходимые условия экстремума скалярной функции многих переменных: а) с использованием частных производных функции; б) с использованием градиента функции. 5. Сформулируйте достаточные условия экстремума функции многих переменных с использованием: а) понятия знакоопределенности второго дифференциала функции; б) главных миноров матрицы Гессе; в) собственных чисел матрицы Гессе. 6. Какая последовательность называется релаксационной? 1  1 7. Будет ли последовательность x ( k )  1  , 2  n  релаксационной 2   n для функции f ( x1 , x 2 )  ( x1  1) 2  ( x 2  2) 4 . 8. В чем состоит идея методов спуска. 9. Покажите, что функция Розенброка f ( x, y)  100( y  x 2 ) 2  (1  x) 2 является унимодальной. 10. Найти минимум функции f ( x, y)   x 2  y 2  xy  x  2 y аналитически. 11. Найти минимум функции f ( x, y)  x  2 y 2  e x 40 2  y2 аналитически. Лекция 5 Методы прямого поиска План Общая характеристика методов нулевого порядка Метод покоординатного пуска Метод Хука – Дживса Метод вращающихся направлений Метод поиска по симплексу Метод сопряженных направлений 1. Общая характеристика методов нулевого порядка В методах прямого поиска минимума целевой функции (или методах нулевого порядка) используют информацию только о значениях этой функции. Многие из этих методов не имеют строгого теоретического обоснования и построены на основе эвристических соображений. Поэтому вопросы сходимости методов прямого поиска еще мало изучены, а оценки скорости сходимости обычно отсутствуют. Вместе с тем эти методы идейно связаны с методами первого и второго порядков, что в ряде случаев позволяет оценивать эффективность алгоритмов прямого поиска применительно к минимизации некоторых классов функций. Распространенным способом оценки эффективности методов прямого поиска являются вычислительные эксперименты и сравнительный анализ методов по результатам таких экспериментов. Однако следует учитывать, что этот анализ не во всех случаях может приводить к однозначным выводам о преимуществах одного метода перед другим. Во-первых, это связано с тем, что сравнению обычно подвергаются не только методы, но и программные реализации соответствующих алгоритмов. Хороший метод можно „загубить" плохим программированием, неудачным выбором параметров алгоритма. Вовторых, методы могут вести себя по-разному на различных этапах процесса минимизации. Удовлетворительного способа преодоления указанных трудностей не существует. Единственное, что можно сделать в подобной ситуации, — привести данные о результатах вычислений в развернутой форме, позволяющей сравнивать методы по различным критериям. Кроме того, не следует забывать, что поиск решения всегда остается искусством, которому можно научиться лишь путем проб и ошибок, применяя различные методы при решении конкретных задач. К методам нулевого порядка относятся методы, не использующие производные для выбора направления спуска: метод Гаусса, метод вращаю41 щихся направлений (Розенброка); метод деформируемого многогранника (поиска по симплексу); метод Хука – Дживса, метод Пауэлла. 2. Метод Гаусса – Зейделя В качестве направлений поиска используются координатные векторы    e1 , e2 , ..., en , где e j  {0,0,...,1,...0}. Таким образом, в поиске по направлению меняется только переменная x j , остальные переменные фиксируются.   Задаем начальную точку x 0 . Рассматриваем направление поиска S1  e1 . Находим параметр 1 из условий одномерной минимизации     min f ( x 0  e1 ). Обозначим промежуточную точку y (1)  x 0  1e1 . Пе  рейдем ко второму направлению e 2 . Находим параметр  2 из условия од  номерной оптимизации min f ( y (1)  e2 ) , обозначим y ( 2)  y (1)  2 e2 ,  проходим y ( j 1)  y ( j ) по всем направлениям координатных осей, определяя     j 1e j 1 , где  j находим из условия min f ( y ( j )  e j 1 ),  j  0,..., N  1. Следующая точка x (1)  y ( n) . Рассматриваем ее как стартовую, далее повторяем процесс поиска по направлениям координатных осей (pис.. 1). Процесс поиска останавливаем при выполнении условия x ( k 1)  x ( k )   или f ( x ( k 1) )  f ( x ( k ) )   , где  – заданная точность. В качестве оптимального решения выбираем x *  x ( k 1) , f *  f ( x*). Риc. 1 42 3. Метод Хука – Дживса Эффективность прямого поиска точки минимума ограниченной снизу целевой функции можно повысить, если на каждом k-м шаге поиска соответствующим образом выбирать направление спуска. Для этого на каждом k-м шаге выделяют предварительный этап исследующего поиска. Целью этого этапа является выбор направления спуска путем исследования пове дения целевой функции f (x ) в окрестности точки x ( k 1) , найденной на предыдущем шаге. В результате выполнения этапа исследующего поиска находится точка ~ x (k ) )  f ( x (k 1) ). Направление спусx ( k ) , для которой f ( ~ ка, завершающего k-й шаг поиска, определяется вектором ~ x (k )  x (k 1) . Такая стратегия поиска, предложенная в 1961 г., получила название метода Хука – Дживса, который состоит из двух этапов: исследующего поиска и поиска по образцу. Исследующий поиск заключается в поиске базовой точки. Перебираем точки вдоль координатных направлений, аналогично методу Гаусса Зейделя. Если значение целевой функции в пробной точке не превышает значения в исходной, то шаг поиска рассматривается как успешный. В противном случае необходимо вернуться в предыдущую точку и сделать шаг в противоположном направлении. После перебора всех N координат исследующий поиск заканчивается. Полученная точка называется базовой. Поиск по образцу: заключается в реализации единственного шага из полученной базовой точки вдоль прямой, соединяющей еѐ с предыдущей базовой точкой. На к-м шаге переходим к этапу спуска в направлении векx (k 1) )  f ( ~ x (k ) )  f ( x (k 1) ) . Следующая тора ~ x (k 1)  ~ x (k ) , при этом f ( ~ точка находится по формуле (1) x (k )  ~ x (k )  k ( ~ x (k 1)  ~ x (k ) ), подбираем так называемый ускоряющий множитель k из одномерной минимизации функции: f ( ~ x (k )   ( ~ x ( k 1)  ~ x ( k ) ))  min , при этом  f ( x (k ) )  f ( ~ x (k 1) ) , либо k задаем постоянным, обычно полагаем k = 2. На рис. 2 иллюстрируются этапы исследующего поиска и спуска для первых двух шагов поиска точки x * — минимума целевой функции двух пе- ременных при  k  2 из начальной точки x (0) . 43 Рис. 2 Рис. 3 На рис. 3 иллюстрируется первый шаг поиска оптимальной точки минимума целевой функции двух переменных с применением на этапе исследующего поиска модифицированного метода покоординатного спуска, а на этапе спуска — процедуры подбора ускоряющего множителя  k  0 в формуле (1), исходя из условия минимума целевой функции в направлении вектора ~ x ( 2)  ~ x (1) . 4. Метод вращающихся направлений Существует еще один альтернативный подход к решению задач минимизации овражных функций, основанный на итеративной перестройке системы ортогональных направлений таким образом, чтобы они определяли движение вдоль дна оврага. Такие методы получили название методов вращающихся направлений или методы Розенброка.. Метод, реализующий эту стратегию поиска, также предусматривает проведение исследующего поиска на каждом к-м шаге. Целью исследующего поиска является выбор текущего направления спуска с учетом информации о поведении целевой функции в окрестности точки x ( k 1) , найденной на предыдущем шаге. Отличие этого метода от метода Хука – Дживса состоит в способе выбора направлений исследующего поиска. Если в методе Хука – Дживса они фиксированы и коллинеарны направлениям векторов стандартного базиса, то в рассматриваемом методе выбор этих направлений проводят в процессе минимизации целевой функции путем построения на каждом к-м шаге поиска нового ортонормированного базиса методом Грамма – Шмидта. 44 Итогом выполнения этого этапа является нахождение точки ~ x ( k ) , для которой f ( ~ x (k ) )  f ( x (k 1) ). Тогда вектор ~ x (k )  x (k 1) определит направление спуска на к-м шаге.    Рассмотрим s1 , s 2 ,..., s n – линейно независимые, взаимно ортогональ0,   ные векторы, т. е. ( s j , s k )   1, j  k, jk ,  После решения задач min f ( y ( j )  s )   j определяется новая точка  n  x ( k 1)  x ( k )    j s j . Выбор новых направлений s1, s2 ,..., sn осуществляj ется на основе следующих соотношений (метод Грамма – Шмидта): s j , если  j  0, a j , при j j  1, n  j 1 aj  bj     i s i , если  j  0; a j   (a i , s i ) s i , при j  2; i  j i 1   bj  sj  . bj На рис. 4 иллюстрируются этапы одного шага поиска точки минимума целевой функции двух переменных из начальной точки. Рис. 4 45 Экспериментальное сравнение алгоритмов Хука – Дживса и Розенброка по числу вычислений целевой функции в процессе оптимизации говорит в пользу алгоритмов вращающихся направлений, это выигрыш растет при увеличении размерности. Но необходимо учитывать более существенные вычислительные затраты на пересчет системы ортогональных направлений. Для минимизации овражных функций с извилистыми оврагами, у которых крутизна склонов намного больше, чем крутизна дна оврага, возможно применение комбинированных алгоритмов. На первом этапе – простые алгоритмы (покоординатного спуска), на втором более сложные. Условием перехода на второй этап может служить соотношение f ( xk )  f ( x1 ) , где x1 - начальная точка поиска, 0    1 подбирается в зависимости от характера и размерности вектора оптимизационных переменных. Например, для n  2 ,  [0,6; 0,8] , при n  4 ,  [0,65; 0,75] , при n  6 ,  [0,1; 0,2] 5. Метод поиска по симплексу Это наиболее известный метод среди методов, не использующих стратегию движения по направлениям, названный методом поиска по симплексу или Нелдера – Мида. Метод основан на том, что экспериментальным образцом, содержащим наименьшее количество точек, является симплекс. Регулярный (правильный) симплекс в N-мерном пространстве – это многогранник, образованный N+1 равноотстоящими точками – вершинами симплекса. В пространстве R 2 правильный симплекс – равносторонний треугольник, а в R 3 – правильный тетраэдр. Если определена базовая вершина симплекса, то координаты всех остальных вершин симплекса {x (i ) , i  1,..,N } могут быть рассчитаны по формулам  x ( 0)  d1 , i  j , (i )  j (2) xj  ( 0)  x j  d 2 , i  j , a( N  1  1) a( N  1  N  1) ; d2  где d1  , а – длина ребра симN 2 N 2 плекса. Важное свойство симплекса — это то, что новый симплекс можно построить на любой грани исходного путѐм отражения какой-либо вершины 46 относительно центра тяжести всех остальных вершин симплекса: 1 ~ x (k )  2 xc  x (k ) , где xc   x (i ) — центр тяжести (рис. 5). N ik Рис. 5 Вычисляя значения целевой функции в вершинах симплекса, получаем информацию о характере изменения этой функции в области расположения симплекса. Поиск точки минимума целевой функции с помощью правильных симплексов производится следующим образом:  На каждой итерации происходит вычисление целевой функции во всех точках симплекса и их упорядочивание по возрастанию значений.  Затем осуществляется последовательная попытка построить новые симплексы с лучшими значениями целевых функций, путем отражения точек с худшими значениями.  Если последовательная попытка отражения двух худших вершин оканчивается неудачей, то производится сжатие симплекса к точке с наименьшим значением и осуществляется новая итерация.  Поиск завершается, как только разность между значениями функции в точках симплекса становится достаточно малой. При реализации алгоритма используются следующие процедуры:  Процедура отражения (рис. 6, а): ~ x (k )  2 x  x (k ) . (3) c  Процедура сжатия вовнутрь (0    1) (рис. 6, б): ~ x (k )  x   ( x  x (k ) . c c (4) 47 x(1) ~ x ( 2) ~ x ( 2) x (0) x(2) x(2) а) б) Рис 6  Процедура отражения со сжатием (0    1) (рис 7, а) ~ x (k )  x   ( x  x (k ) ). (5)  Процедура отражения с растяжением (   1) (рис.7, б) ~ x (k )  x   ( x  x (k ) ). (6) c c c c ~ x ( 2) ~ x ( 2) x(2) x(0) а) x ( 2) б) Рис. 7. На k-м шаге вычислить значения целевой функции в точках симплекса f ( x (i ) ), i  1,..., N , упорядочить точки по возрастанию значений целевой функции. Отобразить точу с максимальным значением функции относи- тельно центра тяжести xc  1 x (i ) , по формулам (3) – (6), k  1,..., 4 .  N ik x (i )  x ( 0) , i  1,...,N . Сжать симплекс в 2 раза: x  2 С теоретической точки зрения методы поиска при помощи нерегулярного симплекса изучены недостаточно полно, однако практические вычисления указывают на их работоспособность. В задачах безусловной минимизации рекомендуют выбирать α = 1/2, β = 2, или α = 1/4, β = 5/2. Вопрос выбора оптимальных параметров для каждой конкретной целевой функции (i ) 48 может быть решен после дополнительного исследования свойств метода с помощью вычислительного эксперимента. Графическая иллюстрация процесса поиска по симплексу приведена на рис. 8. При работе с деформируемым (нерегулярным) симплексом возможно его вырождение: когда вершины попадают в пространство, размерность которого меньше размерности задачи (например, при N=2 вершины симплекса попадают на одну прямую). Формально это ~ x (k )  x (k ) , при этом симплекс уже не выйдет из этого пространства и не может продолжать поиск минимума во всем пространстве. Поэтому необходимо осуществить процедуру обновления симплекса в соответствии с формулами (2). Рис. 8 Достоинства метода: простота; малое количество заранее установленных параметров; простая стратегия поиска, вычисление только значений функции, небольшой объѐм требуемой памяти. Недостатки метода: метод работает эффективно при N  6 , алгоритм основан на циклическом движении по координатам. Это может привести к вырождению алгоритма в бесконечную последовательность исследующих поисков без поиска по образцу. 49 6. Методы сопряженных направлений (Пауэлла) Разработан целый ряд методов безусловной оптимизации, использующих понятие сопряженности векторов, которые определяют направление поиска на смежных итерациях. Важность рассмотрения этого вида методов заключается в том, что для функций, имеющий квадратичный вид, они обеспечивают сходимость к точке минимума за число шагов, не более размерности задачи. Данный метод подходит также и для других функций после разложения в ряд Тейлора в окрестности точки оптимума. Основная идея: если квадратичную функцию n переменных привести к виду суммы полных квадратов, то еѐ оптимум может быть найден в результате n одномерных поисков по преобразованным координатным направлениям. Процедура преобразования квадратичной функции 1 q( x)  a  bx  (Qx, x) к виду суммы полных квадратов эквивалентна 2 нахождению такой матрицы преобразования U, которая приводит матрицу квадратичной формы (Qx, x) к диагональному виду. Квадратичная форма (Qx, x) путѐм преобразования x  U y приводится к виду: ( Dy, y) , где D — диагональная матрица. Вместо координат вектора x в стандартной координатной системе, используются его координаты в новой системе, задаваемой векторами yj. Поскольку yj совпадают с главными осями квадратичной формы, то матрица D — диагональная. Таким образом, с помощью преобразования переменных квадратичной функции строится новая система координат, совпадающая с главными осями квадратичной функции, следовательно, одномерный поиск точки оптимума в преобразованных координатах эквивалентен поиску вдоль каждой из осей квадратичной функции (рис. 16). Таким образом, для нахождения оптимума достаточно провести n одномерных поисков вдоль векторов yj. Ясно, что такой способ нахождения точки минимума нетрудно обобщить на n-мерный случай, хотя для реализации этого способа сначала нужно будет решить громоздкую задачу на собственные значения матрицы Q порядка n. Естествен вопрос: нельзя ли избежать нахождения собственных векторов матрицы Q, но все же сократить количество итераций при поиске точки минимума? Важным понятием, используемым методами этого вида, является сопряженность векторов. Векторы s1, s2 ,..., sk называют Н-сопряженными, если они линейно независимы и ( Hs j , si )  0, при j  i , где Н – симметрическая матрица порядка n  n . 50 Рассмотрим этот вопрос на примере двумерной задачи для квадратич1 ной формы f ( x)  (Qx, x). Выберем начальную точку x0 и произвольное 2 направление спуска, определяемое вектором p1 , не обязательно сонаправленным антиградиенту w1  Qx 0 функции f (x) этой точке, но составляющее с ним острый угол, т.е. (Qx 0 , p1 )  0 . Проведя из точки x 0 в этом направлении исчерпывающий спуск, получим точку x (1) на линии уровня f ( x)  f ( x (1) ) (рис. 17). Рис. 16 Рис. 17 Затем выберем точку ~ x 0 , не лежащую на прямой, проходящей через x 0 , p1 )  0 , и проведем исчерпывающий спуск x 0 и x (1) , но такую, что (Q~ из точки ~ x 0 в том же направлении p . В результате получим точку каса1 ~ (1) ния ~ x (1) на линии уровня f ( x)  f ( x ) . Оказывается, что для получения искомой точки минимума достаточно провести исчерпывающий спуск из точки x (1) (или из точки ~ x (1) ) в направлении p 2 , коллинеарном вектору x (1)  ~ x (1) . Векторы p1 и р2, называют сопряженными относительно матрицы Q (или Q-ортогональными), а направления, определяемые такими векторами, сопряженными направлениями. Понятие сопряженности векторов является обобщением понятия их ортогональности. Вопросы и задания для самопроверки 1. Сформулируйте идею методов прямого поиска нулевого порядка. 51 2. Каким образом выбирают направления и параметр шага в методе Гаусса – Зейделя? 3. В чем заключается этап исследующего поиска в методе Хука – Дживса? 4. Как выбирается ускоряющий множитель в этапе поиска по образцу в методе Хука – Дживса? 5. Симплексный метод Нелдера-Мида, идея метода, этапы алгоритма поиска по деформируемому многограннику? 6. Алгоритм Розенброка, как осуществляется поворот системы координат? Для чего используется ортогонализация Грама-Шмидта? 7. В чем отличия метода Розенброка и метода Хука – Дживса? Каковы условия применения метода Розенброка? 8. Дайте определение регулярного симплекса. 9. Как строится новый симплекс на основе базового? 10. В чем отличие процедур отражения и сжатия вовнутрь? 11. Какие векторы называются Н – сопряженными? 12. Сформулируйте основную идею метода Пауэлла. 13. Проанализируйте на примере поиска минимума квадратичной функции с положительно определенной матрицей Q, имеет ли значение тот порядок, в котором выбираются направления исчерпывающего спуска в методе сопряженных направлений. 14. Можно ли осуществлять поиск минимума квадратичной функции в одном и том же направлении более одного раза? 52 Лекция 6 Методы первого порядка многомерной оптимизации План 1. Метод наискорейшего спуска (Коши). 2. Методы сопряженных градиентов. Метод Флетчера – Ривса. Метод Поллака – Райвера. Методы 1-го порядка используют информацию о производной функции. Если ограниченная снизу целевая функция f ( x), x  R n , является дифференцируемой на множестве R n , то алгоритм поиска точки х* ее минимума можно построить, используя информацию, по крайней мере, о градиенте этой функции. Такие методы называются градиентными. Требуется найти точку x*  ( x1* ,...,xn* ) которая определяет минимум функции f ( x)  min , где x  X  En .    (k ) (k ) (k ) Строим последовательность точек x (k ) k 1,2,...  x1 , x2 ,...,xn где hk   f ( x (k ) )   x (k 1)  x (k )  hk f x (k ) ,  , k (1) ,   0 . Процедура построения таких точек называется градиентной процедурой поиска min f ( x) . Градиентные методы безусловной оптимизации используют только первые производные целевой функции и являются методами линейной аппроксимации на каждом шаге, т. е. целевая функция на каждом шаге заменяется касательной гиперплоскостью к ее графику в текущей точке. Во всех этих методах предполагается, что f (x), f существуют и непрерывны. Все эти методы основаны на итерационной процедуре, определяемой формулой x (k 1)  x (k )   k S (k ) , где  k — величина шага, S (k ) — вектор в направлении x (k 1)  x (k ) . Градиентные методы различаются только способом определения  k , и S (k ) обычно определяется путѐм решения задачи оптимизации f (x) в 53 направлении S (k ) . Направление S (k ) зависит от того, как аппроксимируется функция f (x) . 1. Метод наискорейшего спуска (Коши) Впервые такой метод рассмотрел и применил еще О. Коши в XVIII в. Идея его проста: градиент целевой функции f (x) в любой точке есть вектор в направлении наибольшего возрастания значения функции. Следовательно, антиградиент будет направлен в сторону наибольшего убывания функции и является направлением наискорейшего спуска. Антиградиент (и градиент) ортогонален поверхности уровня f (x) в точке х. Пусть в точке x требуется определить направление наискорейшего спуска (то есть направление наибольшего локального уменьшения f (x) . Разложим f (x) в ряд Тейлора в окрестности точки x и отбросим члены второго порядка по x и выше f ( x)  f ( x (k ) )  (f ( x (k ) ), x)  ... Локальное уменьшение f (x) определяется вторым слагаемым, т. е. наибольшее уменьшение f (x) будет тогда, когда (f ( x (k ) ), x) будет иметь наибольшую отрицательную величину. Этого можно добиться выбором S (k )  f ( x (k ) ) , вектор спуска коллинеарен антиградиенту. Этот случай соответствует наискорейшему локальному спуску x (k 1)  x (k )   k f ( x (k ) ) .   На рис.1 показана релаксационная последовательность x (k ) kN , построенная по методу наискорейшего спуска x (k 1)  x (k )  k S (k ) , где S (k )  f ( x ( k ) ) f ( x ( k ) ) (2) . Параметр  k найдем решением задачи одно- мерной оптимизации min f ( x (k )   S (k ) ) .  Теорема 1. Если f (x) дифференцируемая ограниченная снизу функция, градиент которой удовлетворяет условию Липшица: 2 f ( x)  f ( y)  M x  y , а шаг 0    hn  ,  2  0 , и для люM  2 2 бой точки х выполнено неравенство 54 f ( x) 2   0  f ( x)  f ( x*), где    0  0 , то для релаксационной последовательности x (k ) kN , построенной по формуле (2), справедливы оценки f ( x (k ) )  f *  q k f ( x 0 )  f * , x ( k )  x *  Cq k / 2 , где C  0 и 0  q  1. Доказательство. Обозначим f n  f ( x (n) ). В теореме 4.3 доказано, что и f n  f n1   2 hn2 f ( x ( n) ) , следовательно, учитывая условие теоремы f ( x) 2   0  f ( x)  f ( x*), получаем f n  f n1   0 h02  f n  f ( x*), отсю- да f n1  f n   0 h02  f n  f ( x*) , тогда   f n1  f ( x*)  1   0 h02  f n  f ( x*), так как 0  1   0 h02  1 , то обозначив q  1   0 h02 , получаем f n1  f ( x*)  q f n  f ( x*)  q 2  f n1  f ( x*)  ...  q n1  f 0  f ( x*) , Откуда следует первая оценка теоремы. Далее так как f n  f n1  f ( x*) , то f n  f n1  f n  f ( x*) . Из (1) следует x ( n 1)  x ( n) 2    hn2 f x ( n) 2   2 q n  f 0  f ( x*).   2  f n  f n 1    2 q n  f n  f ( x*)    Откуда следует вторая оценка теоремы, при C 2   2 f ( x (0)  f ( x*) . Недостатки: одним из недостатков этого метода является то, что он сходится к любой стационарной точке, в том числе и седловой, которая не может быть решением. Также отмечена очень медленная сходимость наискорейшего спуска в общем случае. Дело в том, что спуск является "наискорейшим" в локальном смысле, если гиперпространство поиска сильно вытянуто ("овраг"), то антиградиент направлен почти ортогонально дну "оврага", т. е. наилучшему направлению достижения минимума. Метод обладает большой надѐжностью, но медленную сходимость вблизи точки минимума устранить нельзя. Поэтому метод самостоятельно обычно не используется, а используется как предварительная процедура для более сложных методов. Достоинство: на каждой k-ой итерации f  x k 1   f  x k  выполняется свойство убывания функции. 55 Рис.1 2. Методы сопряженных градиентов В методе сопряженных направлений происходит отклонение направления наискорейшего спуска путем добавления к нему направления, используемого на предыдущем шаге. В методе сопряженного градиента строится последовательность точек x ( n) n0,1,...,   x ( k 1)  x ( k )  k S ( k ) (3) направления поиска S (k ) , являются линейными комбинациями градиента текущего направления наискорейшего спуска, и предыдущих направлений поиска, т. е. S ( k 1)  f ( x (k ) k )    i S (i ) . i 1 Направления S ( k 1) и S ( k 1) будут H сопряженными, то есть S ( k 1) , HS ( k )  0 . Модификации метода определяются выбором параметров  i и оператора Н.   Метод Флетчера – Ривса Величины  j выбираются так, чтобы новое направление S(k) было сопряжено со всеми предыдущими направлениями, относительно оператора 56 H, удовлетворяющего условию f ( x)  Hx . При этом критерием окончания поиска является выполнение условия: (f ( x ( k ) ), S ( k ) )  0 . Доказано, что 1   2  ...   k 1  0,  k  f ( x ( k ) ) f ( x 2 ( k 1) 2 , ) и это очень ценный результат, позволяющий строить быстрый и эффективный алгоритм оптимизации. Построение релаксационной последовательности методом Флетчера – Ривса состоит в следующем. Задаем начальную точку x1 , первое направление S (1)  f ( x (1) ) , каждая последующая точка строится по формуле x (k 1)  x (k )  k S (k ) S ( k 1)  f ( x ( k 1) )  , осуществляем f ( x ( k 1) ) f ( x ( k ) ) пересчет направления 2 2 S (k ) . Метод обладает квадратичной сходимостью. Теорема 2. Если f (x) дифференцируема, ограничена снизу, сильно выпукла и градиент удовлетворят условию Липшица f ( x)  f ( y)  M x  y , то lim x ( n)  x * и имеет место оценка n 2 x ( k 1)  x *  C x ( k )  x * . Теорема 3. Если f (x) квадратичная функция, x  R n , H – положительно определенная квадратичная форма, то метод Флетчера Ривса сходится к точке минимума не более чем за n шагов. Преимуществом алгоритма Флетчера – Ривса является то, что он не требует обращения матрицы и экономит память ЭВМ, так как ему не нужны матрицы, используемые в Ньютоновских методах, но в то же время почти столь же эффективен как квази-Ньютоновские алгоритмы. Так как направления поиска взаимно сопряжены, то квадратичная функция будет минимизирована не более, чем за n шагов. Алгоритм Флетчера – Ривса чувствителен к точности одномерного поиска, поэтому при его использовании необходимо устранять любые ошибки округления, которые могут возникнуть. Кроме того, алгоритм может 57 отказать в ситуациях, где Гессиан становится плохо обусловленным. Гарантии сходимости всегда и везде у алгоритма нет, хотя практика показывает, что почти всегда алгоритм приводит к приближенному оптимуму. В общем случае используется рестарт, который позволяет получать результат. Метод Поллака – Райвера Это — алгоритм, аналогичный методу Флетчера – Ривза, с модификацией направлений. Строим последовательность точек x ( j 1)  x ( j )   j S ( j ) , где направления S ( j ) строим с учетом предыдущих направлений: S ( j 1)  f ( x ( j 1) )   j S ( j ) ,  f ( x ( j ) )  f ( x ( j 1) ), f ( x ( j ) )  . j f ( x ( j ) ) 2 При выборе значения параметра ускорения на каждой к-й итерации в рекуррентном соотношении приходится минимизировать целевую функцию f ( x (k )  k S (k ) ) по параметру λ. Это приводит к неизбежным погрешностям на каждой итерации, которые могут вызвать нарушение сходимости алгоритма. Чтобы ослабить влияние погрешностей, используют процедуру „обновления" алгоритма, которая состоит в том, что периодически через заданное число итераций полагают  k  0 . Соответствующий номер итерации называют моментом обновления алгоритма, или рестартом. Использование в алгоритме рестартов позволяет избежать накопления вычислительных погрешностей и уменьшить вероятность построения после каждых п итераций линейно зависимых направлений спуска, но приводит к росту общего числа итераций, необходимых для достижения заданной точности поиска. Вопросы и задания для самопроверки 1. Как выбирается направление и параметр шага в методе наискорейшего спуска? Какие условия должны выполнятся для сходимости метода? 2. В чем отличие выбора направлений метода сопряженных градиентов и метода сопряженных направлений? Как выбирается параметр шага в методе Флетчера – Ривса? 3. Как выбирается параметр шага в методе Поллака – Райвера? 58 Лекция 7 Методы 2-го порядка (Ньютоновские методы) План 1. Методы, используюущие производные второго порядка 2. Модификации метода Ньютона 3. Метод Ньютона-Рафсона 4. Метод Марквардта 1. Методы, использующие производные второго порядка Если целевая функция f (x) является дважды дифференцируемой в R n , то эффективность процесса поиска точки x* ее минимума можно повысить, используя информацию не только о градиенте этой функции, но и о ее матрице Гессе Н(х). Направление поиска, соответствующее наискорейшему спуску (методы, рассмотренные выше), связано с линейной аппроксимацией целевой функции. Методы, использующие вторые производные, возникли из квадратичной аппроксимации целевой функции: f (x) , которую можно получить при разложении функции в ряд Тейлора 2-го порядка, 1  ( x)  f ( x k )  f ( x k )( x  x k )  ( x  x k )T H ( x k )( x  x k ), 2 k где H ( x ) — матрица Гессе (вторых производных). Минимум f(x) (если он существует) достигается там же, где и минимум квадратичной формы φ(x). Если матрица Гессе Н(х(k)) целевой функции, вычисленная в точке х(k) , является положительно определенной, то точка х(k) минимума функции φ(x) единственна и может быть найдена из условия, что ее градиент равен нулевому вектору:   f ( x (k ) )  H ( x (k ) )( x  x (k ) )  0 , следовательно, x (k )  x (k 1)  H 1 ( x (k 1) )f ( x (k 1) ) . Алгоритм оптимизации, в котором направление поиска определяется из этого соотношения, называется методом Ньютона, а направление pk   H 1 ( x (k 1) )f ( x (k 1) ) – ньютоновским направлением, составляет с вектором градиента тупой угол (рис. 1). Сходимость метода Ньютона существенно зависит от выбора начального приближения. Если целевая функция сильно выпуклая и для любых 59 точек x, ~ x  D f относительно матрицы Гессе Н(х) целевой функции выполнено неравенство H ( x)  H ( ~ x)  L x  ~ x , L > 0, а начальное приближение выбрано удачно (точка x 0 расположена достаточно близко к точке х* минимума), то алгоритм метода Ньютона обладает квадратичной скоро2 стью сходимости, т.е. справедлива оценка x (k )  x *  C x (k 1)  x * , где C  0. Рис. 1 Но в случае, когда целевая функция не является сильно выпуклой или же начальное приближение находится далеко от искомой точки х*, метод Ньютона может расходиться. В задачах поиска минимума произвольной квадратичной функции с положительной матрицей вторых производных метод Ньютона дает решение за одну итерацию независимо от выбора начальной точки. Квадратичная скорость сходимости, а также возможность контроля за соблюдением достаточного условия минимума целевой функции f(x) на каждой к-й итерации при помощи матрицы Гессе этой функции способствуют высокой эффективности рассматриваемого алгоритма. Однако при его практической реализации возникают две проблемы. Первая из них — это сохранение положительной определенности матрицы Гессе H ( x (k 1) ) целевой функции на каждой к-й итерации, так как иначе вектор p k может не соответствовать направлению спуска, вдоль которого эта функция убывает. Более того, матрица H ( x (k 1) ) может быть вырожденной и не иметь обратной матрицы. В связи с этим более серьезной проблемой является 60 необходимость вычисления и обращения на каждой итерации матрицы порядка п, что в случае большой размерности пространства является достаточно трудоемкой операцией. На основе устранения этих проблем рассматривают модификации метода Ньютона. 2. Модификации метода Ньютона Метод Ньютона считается эталонным, с ним сравнивают все разрабатываемые оптимизационные процедуры. Для применения метода Ньютона матрица Гессе должна быть положительно определенной и хорошо обусловленной (определитель ее должен быть существенно больше нуля, т. е. отношение наибольшего и наименьшего собственных чисел должно быть близко к единице). Но матрица Гессе может быть вырожденной и не иметь обратной матрицы. Эту проблему можно решить, если направление спуска задавать векто-  ром pk   k I  H ( x (k 1) ) 1f ( x(k 1) ) , где I — единичная матрица поряд- ка n, а βk — параметр, выбираемый так, чтобы в точке x ( k 1) матрица H k   k I  H ( x (k 1) ) была положительно определена. В связи с этим более серьезной проблемой является необходимость вычисления и обращения на каждой итерации матрицы порядка n, что в случае большой размерности пространства Rn является достаточно трудоемкой операцией. На практике обычно не вычисляют матрицу, обратную к положительно определенной матрице Hk, а вектор рк находят из решения системы линейных алгебраических уравнений (СЛАУ) H k pk  f ( x (k 1) ) . Эту СЛАУ можно решить различными численными методами, например прямыми и итерационными. Эти приемы называются модификацией метода Ньютона, они используют ньютоновские направления по мере возможности и уклоняющиеся от них только тогда, когда это необходимо. Общий принцип модификаций метода Ньютона состоит в следующем: на каждой итерации сначала строится некоторая "связанная" с H ( x k ) , положительно определенная матрица Hk, а затем направление спуска Sk вы- числяется по формуле H k S k  f ( x k ) . Так как Hk положительно определена, то направление Sk обязательно будет направлением спуска. Процедуру построения Hk организуют так, чтобы она совпадала с матрицей Гессе, если она является положительно определенной. Эти процедуры строятся на основе некоторых матричных разложений. 61 Опишем алгоритм варианта метода Ньютона — поиска точки х* — минимума дважды дифференцируемой целевой функции f (x) , в котором направление спуска определяется путем решения СЛАУ. Зададим начальную точку x (0) , начальную величину шага 0  1 , коэф.дробления шага 0    1, точность ε. На k-ом шаге: вычислим f k  f ( x k ), H k  H ( x k ) , находим направление спуска pk из решения СЛАУ H k pk  f ( x k ) , сли H k положительно определенная, если нет, то Подбираем βk, так чтобы H k   k I  H ( x (k ) ) была положительно определенной. Вычисляем x k 1  x k  k pk и f ( x k 1 ) . Условия останова поиска: x k 1  x k   x или f ( x k 1 )  f ( x k )   f . При произвольном выборе начальной точки x 0 метод Ньютона обладает квадратичной скоростью сходимости, и имеет место оценка 2 L ( k 1) x (k )  x *  x  x* , kN . 1 Здесь 1  оценка наименьшего собственного числа матрицы Гессе в окрестности точки x* и H ( x)  H ( y)  L x  y . 3. Метод Ньютона-Рафсона Данный метод состоит в многократном использовании ньютоновского направления при оптимизации функций, не являющихся квадратичными. Основная итерационная формула многомерной оптимизации x k 1  x k  k 1S k 1 используется в этом методе при выборе направления оптимизации из соотношения S k 1   H 1 ( x k )f ( x k 1 ) k 1  1. Реальная длина шага скрыта в ненормализованном Ньютоновском направлении Sk+1. Так как этот метод не требует значения целевой функции в текущей точке, то его иногда называют непрямым или аналитическим методом оптимизации. Его способность определять минимум квадратичной функции за одно вычисление выглядит на первый взгляд исключительно привлекательно. Однако это "одно вычисление" требует значительных затрат. Прежде всего, необходимо вычислить n частных производных первого порядка и n(n+1)/2 — второго. Кроме того, матрица Гессе должна быть обратимой. Это требует уже порядка n 3 вычислительных операций. 62 С теми же самыми затратами методы сопряженных направлений или методы сопряженного градиента могут сделать порядка n шагов, т. е. достичь практически того же результата. Таким образом, итерация метода Ньютона-Рафсона не дает преимуществ в случае квадратичной функции. Если же функция не квадратичная, то:  начальное направление уже, вообще говоря, не указывает действительную точку минимума, а значит, итерации должны повторяться неоднократно;  шаг единичной длины может привести в точку с худшим значением целевой функции, а поиск может выдать неправильное направление, если, например, гессиан не является положительно определенным;  гессиан может стать плохо обусловленным, что сделает невозможным его инвертирование, т. е. определение направления для следующей итерации. Сама по себе стратегия не различает, к какой именно стационарной точке (минимума, максимума, седловой) приближается поиск, а вычисления значений целевой функции, по которым можно было бы отследить, не возрастает ли функция, не делаются. Значит, все зависит от того, в зоне притяжения какой стационарной точки оказывается стартовая точка поиска. Стратегия Ньютона – Рафсона редко используется сама по себе без модификации того или иного рода. 4. Метод Марквардта Это комбинация методов Ньютона и Коши. Вдали от точки минимума направление определяется по методу Коши, а в окрестности точки минимума – по методу Ньютона. Строим релаксационную последовательность: x k 1  x k  k pk , --– k  1 ,  pk   k I  H ( x (k 1) ) 1f ( x(k 1) ) . (1) На начальном этапе  k  10 4 , при этом второй член в (1) много больше первого, поэтому поиск осуществляется по методу Коши. По мере приближения к точке оптимума  k уменьшается и стремится к нулю. Таким образом, вблизи точки оптимума первый член в (1) много больше второго и поиск осуществляется по методу Ньютона. 63 Если после первого шага f ( x (1) )  f ( x 0 ) , то следует выбрать 1   0 и реализовать следующий шаг, в противном случае  0   0 , где γ > 1, и повторить предыдущий шаг. Вычислительные эксперименты показали, что метод наиболее эффективен для функций вида суммы квадратов Вопросы и задания для самопроверки Какие направления поиска называются ньютоновскими? Как определяется матрица Гессе? В чем состоит основной принцип модификаций метода Ньютона? Чем отличается метод Ньютона от метода Ньютона – Рафсона? Как строится релаксационная последовательность в методе Марквардта? 6. Для каких функций эффективно применение методов второго порядка? 7. В каких случаях модифицированный метод Ньютона сходится, а в каких нет? 8. Что необходимо предусмотреть в реализуемом алгоритме для обеспечения сходимости метода? 1. 2. 3. 4. 5. 64 Лекция 8 Методы переменной метрики План 1. Квазиньютоновские методы: общая характеристика Методы Пирсона Метод Дэвидона – Флетчера – Пауэлла Метод Бройдена – Флетчера – Шенно Методы Пауэлла и Мак – Кормика 2. Сравнение методов безусловной оптимизации 1. Квазиньютоновские методы: общая характеристика Среди алгоритмов многомерной минимизации следует выделить группу алгоритмов, которые объединяют достоинства метода наискорейшего спуска и метода Ньютона. Такие алгоритмы принято относить к так называемым квазиньютоновским методам. Особенность этих алгоритмов состоит в том, что при их применении нет необходимости вычислять и обращать матрицу Гессе целевой функции f (x) и в то же время удается сохранить высокую скорость сходимости алгоритмов, присущую методу Ньютона и его модификациям. В этих методах обратная матрица Гессе аппроксимируется другой матрицей – метрикой. Метрика изменяется на каждой итерации, и поэтому методы также называются методами с переменной метрикой.   Элементы релаксационной последовательности x ( k ) k в алгоритмах квазиньютоновских методов минимизации непрерывно дифференцируемой в Rn целевой функции f(x) строят в соответствии с рекуррентным соотношением x k 1  x k  k pk , направление спуска на каждой k-й итерации ~ ~ ~ ~ задают в виде pk   H k f ( x k ) , где H k – метрика, а H k 1  H k  Ck , где Ск - корректирующая матрица. ~ Нужно построить последовательность метрик H k , которая при   k   имела бы предел, равный H 1 ( x*) , где H (x*)  матрица Гессе в ~ точке минимума x * функции f (x) . Матрица H k аппроксимирует матрицу Гессе Н k , которая должна удовлетворять уравнению: 65   (1) x ( k 1)  x ( k )  H 1 ( x k ) f ( x ( k 1) )  f ( x ( k ) ) . ~ ~ ~ Пусть H 1 ( x ( k 1) )  H k 1   H k  H k , производим оценку по предыдущему шагу. Подставим последнее соотношение в (1), тогда получаем ~ ~ ~ x ( k 1)  x (k )  H ( x k 1 ) f ( x ( k 1) )  f ( x ( k ) )   H k  H k f k 1  f k  . ~ Отсюда получаем уравнение, относительно H k         1 ~ ~ H k f k 1  f k   x ( k 1)  x ( k )  H k f k 1  f k ,  решением которого будет матрица ~ H k  x   x (k ) y ~ f k 1  f z ,  Hk z, f k 1  f    y, f k 1  f k  ( k 1) где y, z - произвольные векторы из R3. В зависимости от выбора , y, z различают модификации методов переменной метрики. Методы Пирсона Пирсон предложил несколько методов с аппроксимацией обратного гессиана без явного вычисления вторых производных, т. е. путем наблюдения за изменениями направления антиградиента. При этом получаются сопряженные направления. Эти алгоритмы отличаются только деталями. Приведем те из них, которые получили наиболее широкое распространение в прикладных областях. Они аппроксимируют матрицу Гессе или обратную к ней, но используют для этого только первые производные. В большинстве из методов применяются сопряжѐнные направления. ~ На каждом шаге x k 1  x k  k S k , где S k   H k f ( x k )  направление поиска. В этом алгоритме обратный гессиан аппроксимируется матрицей, вычисляемой на каждом шаге по формуле ~ ( x k 1  x k )  H k (f ( x k 1 )  f ( x k )) x k 1  x k ~ ~ H k 1  H k  ( x k 1  x k )(f ( x k )  f ( x k )) или в другой модификации ~ ~ ~ H k 1  H k  x ( k 1)  x (k )  H k  f ( x ( k 1)  f ( x ( k )  Ak ,    66      T ~ H k  f ( x ( k 1)  f ( x ( k ) где Ak  . ~ ( k 1) (k ) T ( k 1) (k ) f ( x  f ( x  H k  f ( x  f ( x ~ В качестве начальной матрицы H 0 выбирается произвольная положительно определенная, симметрическая матрица (обычно единичная). Дан~ ный алгоритм Пирсона часто приводит к ситуациям, когда матрица H k     становится плохо обусловленной, а именно  она начинает осциллировать, колеблясь между положительно определенной и не положительно определенной, при этом определитель матрицы близок к нулю. Для того чтобы избежать этой ситуации, необходимо через каждые n шагов В качестве ~ матрицы переменной метрики выбирать матрицу H 0 . Метод Дэвидона – Флетчера – Пауэлла Метод Дэвидона – Флетчера – Пауэлла (ДФП) основан на использовании ньютоновских направлений, но не требует вычисления обратного гессиана на каждом шаге. Направление поиска на k-м шаге определяется ~ ~ направлением p k   H k f ( x (k ) ) , где H k  положительно определенная симметричная матрица, которая обновляется на каждом шаге и в пределе ~ становится равной обратному гессиану, то есть lim H k  H 1 ( x*)  0 . k  ~ Матрица H k строится следующим образом: в качестве начальной выбира~ ~ ем единичную матрицу, далее на k-м шаге: H k 1  H k  Ak  Bk , матрица ~ H k  Ak  Bk , матрицы Ak , Bk строим следующим образом T ~ ~ T H k v k v kT H k u (k ) u (k ) Ak   ( k ) , Bk   , (2) ~ H v , v (u , vk ) k k k   где u (k )  x (k )  x (k 1) ,     vk  f ( x (k ) )  f ( x (k 1) ). Метод эффективен на практике, если ошибка вычислений градиента ~ f ( x (k ) ) невелика и матрица H k не становится плохо обусловленной. ~ Матрица Ak обеспечивает сходимость H k к H 1 ( x*) , матрица Bk обеспе~ чивает положительную определенность H k 1 на всех этапах. k k ~ ( k 1) H  I   Ai   Bi . i 1 i 1 67 В случае квадратичной функции n 1 n 1 i 1 i 1  Ai  H 1,  Bi   I , т. е. алго- ритм ДФП использует сопряженные направления. Таким образом, метод ДФП использует как идеи ньютоновского подхода, так и свойства сопряженных направлений, и при минимизации квадратичной функции сходится не более чем за n итераций. Если оптимизируемая функция имеет вид, близкий к квадратичной функции, то метод ДФП эффективен за счет хорошей аппроксимации H 1 (метод Ньютона). Если же целевая функция имеет общий вид, то метод ДФП эффективен за счет использования сопряженных направлений. На практике оказалось, что метод ДФП может давать отрицательные ~ шаги или окончиться в нестационарной точке. Это возможно, когда H k становится плохо обусловленной. Этого можно избежать путем увеличе~ ния числа получаемых значимых цифр или перезаданием матрицы H k в виде специальной диагональной матрицы, где h jj  x (jk )  x (jk 1) . f ( x ( k 1) ) Ошибки округления, и особенно неточность линейного поиска, могут послужить причиной потери устойчивости метода ДФП и даже привести к ошибочным шагам, когда значение целевой функции на некоторой итерации возрастает, вместо того, чтобы уменьшаться. Практическая проверка эффективности алгоритма показала, что он столь же эффективен, как и метод Флетчера – Ривса. Метод Бройдена – Флетчера – Шенно ~ Возможны и другие способы построения последовательности H k мат~ ~ риц с применением H k 1  H k  Ck . Так, например, в алгоритме квазиньютоновского метода, называемого в литературе методом Бройдена – Флетчера – Шенно, или БФШ-методом, который аналогичен методу ДФП, об~ ~ ~ новление матрицы H k строится по формуле H k 1  H k  Ak  Bk  C k , ~ где Ak , Bk находим по формулам (2), а C k  H k v k , v k r k (r k ) T , при этом ~ H k vk uk k r  ~ k k k  k k . В качестве начальной матрицы выбирается едиH v ,v u ,v ничная матрица. БШФ-метод обладает слабой чувствительностью к по-   68     грешности вычислений, не требует обновления матрицы переменной мет~ ~ рики, H k является положительно определенной и lim H k  H 1 ( x*)  0 . k  Методы Пауэлла и Мак – Кормика ~ В алгоритме Пауэлла обновление матрицы H k оценивающей матрицу    T w k  w k 1 w k  w k 1 ~ k 1 ~ k Гессе производится по формуле: H , где H  v k , wk ~ w k  u (k )  H k v k . k (k ) T ~ k 1 ~ k w (u ) В алгоритме Мак – Кормика H  H  k (k ) . v ,u ~ Все рассмотренные методы обновления матриц H k сохраняют свой~ ство положительной определенности, и последовательность H k kN схо-       дится при k   к матрице H 1 ( x*) . Рассмотрим сходимость методов безусловной оптимизации. Метод называется сходящимся, если на каждой итерации выполнено неравенство x ( k 1)  x * x (k )  x*  1 . Алгоритм обладает сходимостью порядка r, если вы- полнено соотношение x ( k 1)  x * x (k )  x* r  C , (С- конечно). Если r=1, то алго- ритм обладает линейной скоростью сходимости. Если ещѐ при этом C=0, то алгоритм обладает суперлинейной скоростью сходимости. Если r=2, то скорость квадратичная. Алгоритм Флетчера Флетчером был предложен алгоритм, в котором условие окончания процесса для квадратичной функции после n шагов было отброшено, но сохранено свойство, заключающееся в том, что для квадратичных функций ~ ~ H k  H 1 ( x*) в том смысле, что собственные значения H k стремятся к собственным значениям H 1 ( x*) . Полученное Флетчером соотношение для очередного приближения матрицы имеет вид 69    x ( k 1)  x ( k ) v k T ~  H k 1  E   u k , v k    ( k 1)  ~ k   x (k )   H E  vk x   u k , v k     T  u u T  k k ,  u , v  k k  где u k  x ( k 1)  x ( k ) , v k  f ( x ( k 1) )  f ( x k ) . Сравнение алгоритмов на тестовых функциях показывает, что в этом случае алгоритм Флетчера проигрывает по эффективности алгоритму Дэвидона - Флетчера - Пауэлла. Эффективность алгоритмов Дэвидона – Флетчера - Пауэлла и Флетчера в существенной мере определяют используемые методы одномерного поиска. Зачастую успех конкретной реализации определяет, именно, эффективность используемых методов одномерного поиска. При решение задачи одномерной оптимизации min f ( x k   p k ) по параметру , выбор параметра поиска выбирается на интервале 0,  ' с помощью кубической интерполяции, где    (k )   2 f (x  f 0    '  min 1,  , f 0 наименьшее значение минимизируемой (k )   f ( x ), p k    функции в результате одномерного поиска по направлению p k . Если найденное значение параметра  больше ’, можно выбирать пробные значения на начальных шагах ( r 1)  0,1  ( r ) .  2. Сравнение методов многомерной безусловной оптимизации Существуют два пути сравнения: теоретическое исследование сходимости и численные эксперименты. Так метод Пауэлла имеет суперлинейную скорость сходимости, метод Коши – линейную, методы Ньютона – квадратичную, методы сопряжѐнных градиентов - линейная скорость сходимости, а квазиньютоновские - квадратичную скорость. В результате численных экспериментов 6 методы распределяются по количеству вычислений значений функции, устойчивости, машинному времени. Устойчивость характеризует ширину круга задач (успешно решаемых). Лучшие методы: ДФП, Пауэлла, Бройдена-Флетчера-Шенно. Другие исследователи сравнивали градиентные методы. При этом учитывалось влияние параметров сходимости методов одномерного поиска, положительная определѐнность матрицы H для квазиньютоновских методов и точность определения компонент градиента. Выводы: 6 См. Химмельблау Д. Прикладное нелинейное программирование. 70 1. Существенное превосходство квазиньютоновских методов при решении задач с функциями общего вида; 2. на методы переменной метрики большее влияние оказывает точность вычислений на ЭВМ, чем на методы сопряжѐнных градиентов. 3. комбинация метода Коши и метода деления пополам обеспечивает наибольшую точность при больших затратах машинного времени. 4. самая эффективная с точки зрения вычисления значений функции комбинация Бройдена-Флетчера-Шенно и кубической аппроксимации. Вопросы и задания для самопроверки 1. Какие методы называются квазиньютоновскими? Как строится релаксационная последовательность в этих методах? 2. Что такое переменная метрика? В чем заключается общая идея методов переменной метрики? 3. Каким свойством должны обладать аппроксимирующие матрицу Гессе метрики? 4. Как аппроксимируют обратную матрицу Гессе в методах Пирсона? Какая матрица выбирается в качестве начальной? 5. В чем отличие методов Пирсона и метода Дэвидона – Флетчера – Пауэлла? 6. Сравните методы Дэвидона – Флетчера – Пауэлла и Бройдена – Флетчера – Шенно. К каким из функций применим тот или иной метод? 7. Чем отличаются методы Пуэлла и Мак – Кормика? Задания для лабораторной работы Лабораторная работа № 1 Методы одномерного поиска Цель работы Ознакомиться с методами одномерного поиска. Сравнить различные алгоритмы по эффективности на тестовых примерах. Порядок выполнения работы 1. Найти аналитическое решение задачи min f ( x) . x[ a, b] 2. Реализовать методы с различной точностью (согласно вариантам) Варианты 1-3 методы пассивного поиска, золотого сечения, метод Ньютона. Варианты 4-7 методы дихотомии, Фибоначчи, секущих. Варианты 8-11 методы деления пополам, Пауэлла, Ньютона. 71