Поиск минимума функции многих переменных методом Нелдера - Мида

НИЖЕГОРОДСКИЙ ГОСУДАРСТВЕННЫЙ  ТЕХНИЧЕСКИЙ УНИВЕРСИТЕТ

 

 

 

 

                                Курсовая работа 

 

по информатике

 

 

Тема:

"Поиск минимума функции многих переменных

методом Нелдера-Мида"

 

 

 

 

 

 

 

 

                                               Выполнил: студент группы 10-АХ-1 АМФ

     Чернов Д.С.                     

                                                   Проверил:                         Сухов В.И.

 

 

 

 

 

 

Нижний Новгород

2011 г.

Содержание

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Введение

Впервые метод деформируемого многогранника был предложен  Нелдером и Мидом. Они предложили метод поиска, оказавшийся весьма эффективным и легко осуществляемым на ЭВМ. Чтобы можно было оценить стратегию Нелдера и Мида, кратко опишем симплексный поиск Спендли, Хекста и Химсворта, разработанный в связи со статистическим планированием эксперимента. Вспомним, что регулярные многогранники в eⁿ являются симплексами. Например, как видно из рисунка 1, для случая двух переменных регулярный симплекс представляет собой равносторонний треугольник (три точки); в случае трёх переменных регулярный симплекс представляет собой тетраэдр (четыре точки) и т.д.

                                         

а                           б

Рисунок 1. 
Регулярные симплексы для случая двух (а) и трёх (б) независимых переменных. Стрелка указывает направление  наискорейшего улучшения.

 

При поиске минимума целевой функции f(x) пробные векторы x могут быть выбраны в точках eⁿ, находящихся в вершинах симплекса, как было первоначально предложено Спендли, Хекстом и Химсвортом. Из аналитической геометрии известно, что координаты вершин регулярного симплекса определяются следующей матрицей d, в которой столбцы представляют собой вершины, пронумерованные от 1 до (n+1), а строчки – координаты, i принимает значения от 1 до n:

 

– матрица n x (n+1),

где 

t – расстояние между двумя  вершинами

 

 

1. Постановка задачи

Целевая функция может быть вычислена  в каждой из вершин симплекса; из вершины, где целевая функция максимальна, проводится проектирующая прямая через  центр тяжести симплекса. Затем точка a исключается и строится новый симплекс, называемый отражённым, из оставшихся прежних точек и одной новой точки b, расположенной на проектирующей прямой на надлежащем расстоянии от центра тяжести. Продолжение этой процедуры, в которой каждый раз вычёркивается вершина, где целевая функция максимальна, а также использование правил уменьшения размера симплекса и предотвращения циклического движения в окрестности экстремума позволяют осуществить поиск, не использующий производные и в котором величина шага на любом этапе k фиксирована, а направление поиска можно изменять. На рисунке 2 приведены последовательные симплексы, построенные в двумерном пространстве с «хорошей» целевой функцией.

 

 

Рисунок 2. 
Последовательность регулярных симплексов, полученных при минимизации f(x). 
----- проекция

Определённые практические трудности, встречающиеся при использовании  регулярных симплексов, а именно отсутствие ускорения поиска и трудности  при проведении поиска на искривлённых «оврагах» и «хребтах», привели к необходимости некоторых улучшений методов. Далее будет изложен метод Нелдера и Мида, в котором симплекс может изменять свою форму и таким образом уже не будет оставаться симплексом. Именно поэтому здесь использовано более подходящее название «деформируемый многогранник».

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

 

 

 

 

 

 

 

 

 

2. Методика решения задачи

  В методе Нелдера–Мида вокруг начальной точки поиска в пространстве переменных размещается начальный симплекс – конфигурация из (n+1)-й точки (в пространстве R ² они образуют вершины треугольника, а в R ³ – вершины пирамиды). Затем происходит перемещение симплекса путем отражения вершины с наибольшим значением функции относительно центра тяжести противолежащего основания симплекса. При этом используются специальные операции, связанные с растяжением симплекса в направлении убывания функции и операции сжатия при неудачных пробных перемещениях

С помощью операции растяжения и  сжатия размеры и форма деформируемого многогранника адаптируются к топографии целевой функции. В результате многогранник вытягивается вдоль длинных наклонных поверхностей, изменяет направление в изогнутых впадинах, сжимается в окрестности минимума, что определяет эффективность рассмотренного метода.

Метод Нелдера–Мида имеет тот недостаток, что для сильно овражных функций  может происходить вырождение симплекса, особенно при числе переменных N > 2.

Термин «вырождение» означает, что  все точки симплекса с некоторого шага размещаются в многообразии размерности меньшей, чем N, или же попадают в малую его окрестность, величина которой много меньше расстояния между точками симплекса.

В диалоговой системе оптимизации  выход из подпрограммы, реализующей  метод деформируемого многогранника, осуществляется при предельном сжатии многогранника.

 

3. Алгоритм метода Нелдера–Мида.

Шаг 0. Задаем векторы h¹,h²,…,h^N+1, определяющие положение вершин стандартного симплекса с центром в начале координат, и числа S₁,S ₂,… ,S_N+1, определяющие размеры начального симплекса; параметры останова; –параметры отражения, растяжения, сжатия к основанию, сжатия к лучшей вершине ( >0, >1, 0< <1, ). Задаем также начальную точку y ⁰.

Шаг 1. Формируем начальный симплекс с координатами вершин y¹,…, y ^N+1

y ^j = y ⁰+ S _j h ^j ; (j = 1,…,N+1) .

Вычисляем f ^j = f(y ^j). (При этом в y ⁰ вычисление не выполняется).

Шаг 2. Определяем номера худшей и лучшей вершины

f ^h = max{f_j : j=1,...N+1} ; f ^e = min{ f_j: j=1,...N+1}.

Шаг 3. Определяем центр тяжести основания

Шаг 4. Проверяем критерий останова

    

Если  и , то выполняем останов, выдаем оценку решения y ^e, f ^e. Если условия останова не выполнены, переходим на шаг 5.

 

Шаг 5. Выполняем отражение с коэффициентом a>0

и вычисляем f^* = f(y^*) .

Шаг 6. Если , то выполняем растяжение

при заменяем y ^h:= y^**, f ^h:= f^** и переходим на шаг 2;

при f^** > f^* заменяем y ^h:=y^*, f ^h:=f^* и переходим на шаг 2;

если f^* > f ^e, то переходим на шаг 7.

Шаг 7. Если для любого j=1,…,N+1, но , выполняется f ^e < f ^* < f ^j, то заменяем

y ^h:=y^*, f ^h:=f^* и переходим на шаг 2, иначе — на шаг 8.

Шаг 8. Если f^* < f ^h, то выполняем сжатие к основанию. Для этого вычисляем

   

заменяем y ^h:=y ^{^}, f ^h:=f ^{^} и переходим на шаг 2.

Если , то выполняем сжатие к лучшей вершине

Переходим на шаг 2.

4. Блок-схема

4.1 Блок-схема основной программы

 

4.2 Блок-схема процедуры  calcMinimum

 

4.3 Блок-схема функции

5. Программа на языке Pascal

5.1 Основная программа

 

program MNM;

uses crt, unit1;

 

function f(x:SimplexPoint):double; far;

var q:integer;   {использование "дальней"  адресации для того, чтобы}

{эту функцию можно было передавать  как параметр в процедуру}

begin            {нахождения минимума}

  F := (x[1]-3)*(x[1]-3)+(x[2]-2)*(x[2]-2);

end;

 

const a:double = 1;

      beta:double = 0.5;

      gamma:double = 2.9;

 

var  e, y0: SimplexPoint;

     i,q,n:integer;

     x0:simplex;

     l:real;

 

begin

clrscr;

writeln('************ Метод Нелдера-Мида ************');

 write('Введите размерность пространства: '); readln(n);

write('Введите коэффициент a (a>0, а=1): '); readln(a);

write('Введите коэффициент гамма  (2,8<=gamma<=3):');  readln(gamma);

write('Введите коэффициент бетта  (0,4<=beta<=0,6):'); readln(beta);

 

 for i:=1 to n+1 do

for q:=1 to n do x0[i][q] :=0;

 

{ вводим вектора определяющие  положение вершин первого симплекса  и его размеры }

     writeln('Введите  начальную точку y0:');

     for q:=1 to n do

      begin

        write('   y0[',q,']: ');

        readln(Y0[q]);

      end;

     writeln('Введите вектора вершин симплекса и их длины :');

     for i:=1 to n+1 do

      begin

           writeln('вектор x',i,':');

           for q:=1 to n do

            begin

             write('   x',i,'[',q,']: ');

             readln(x0[i][q]);

            end;

           write('   Введите длину вектора x',i,': '); readln(L);

           for q:=1 to n do

           x0[i][q] := x0[i][q]*L+Y0[q];

      end;

{ вводим вектор останова}

     write('Введите  желаемую точность: ');

     read(E[1]);

     for q:=2 to n do E[q]:=E[1];

 

{********************************************************}

  calcMinimum(f, n, x0, e, a, beta, gamma);

{********************************************************}

 

{ вывод результатов  вычисления}

  writeln('==========================================');

  writeln('Результат: ');

  for i:=1 to n+1 do

   begin

     writeln(' Точка  ', i, ':  Значение функции: ',f(x0[i]):3:6);

     write('   [');

     for q:=1 to n do

      begin

        write(x0[i][q]:3:5);

        if q<>n then write(', ');

      end;

       writeln(']');

   end;

  readkey;

end.

5.2 Unit

unit unit1;

 

interface

 

const n_max=20; {размерность пространства}

 

type

   SimplexPoint = array[1..n_max] of double;    { точка }

   Simplex = array[1..n_max+1] of SimplexPoint; { симплекс }

   tFunc=function(point:SimplexPoint):double;

 

procedure CalcMinimum( formula:tFunc; {функция}

                       n:integer; {размерность пространства}

                       var N_kub:Simplex; {начальное приближение}

                       E:SimplexPoint {необходимая точность решения};

                       a,beta,gamma:real);

 

implementation

 

uses crt;

 

procedure CalcMinimum( formula:tFunc; {функция}

                       n:integer; {размерность пространства}

                       var N_kub:Simplex; {начальное приближение}

                       E:SimplexPoint; {необходимая точность решения}

                       a,beta,gamma:real);

var q:integer;

    F, Fmax, Fmin, len, len_max :double;

    i, j, num_H, num_L, pr :integer;

    Point, Point_H, Point_new, Point_new2 :SimplexPoint;

    F_new, F_new2:double;

    k:integer;

 

begin

{ бесконечный цикл}

    k:=0;

    while TRUE do

     begin

       k:=k+1;

       Writeln('==Shag ',k);

       for i:=1 to n+1 do

        begin

         write('  x[', i, ']=[');

         for q:=1 to n do

          begin

           write(n_kub[i][q]:3:5);

           if q<>n then write(', ');

          end;

         writeln('],  f(x[', i, '])=',formula(n_kub[i]):3:6);

        end;

        if readkey=#27 then break;

 

 

     { вычисляем номер самой удаленной(худшей) точки симплекса}

     { вычисляем  номер самой лучшей точки симплекса}

          Fmax := Formula(N_kub[1]);

          Fmin := Fmax;

          num_H := 1;

          num_L := 1;

          for i:=2 to n+1 do

           begin

             F := Formula(N_kub[i]);

             if F>Fmax then begin Fmax := F;num_H := i;end;

             if F<Fmin then begin Fmin := F;num_L := i;end;