Методы нахождения условного и безусловного экстремума

Метод штрафных функций:

Исходные данные:

f(х1,х2)=(х1–4)2+(х2–4)2 +х1·х2;

x(0) = [-10,-10]T – начальная точка.

Ограничение на решение:

h(x1,x2)=2·x1 + 3·x2 – 12= 0.                           

Найдём минимум целевой функции f с заданным квадратичным штрафом:

h(x1,x2)=2·x1 + 3·x2 – 12= 0.

P(x1,x2,R) =(х1–4)2+(х2–4)2 +х1·х2 + 1/R·(2·x1 + 3·x2 – 12)2.

Совместное решение даёт

Устремляя R к нулю, получаем = [18/7;16/7]; f () = 76/7; h()=0.

Устремляя R к бесконечности, получаем = [8/3;8/3]; f () = 32/3; h()=4/3.

Т.е. при изменении R от нуля до бесконечности, решение будет изменяться от минимума задачи с учётом ограничений до минимума функции без учёта  ограничений.

Исследуем функцию при различных значениях параметра R  0 с помощью метода Ньютона:

R =100, 10, 1, 0.1  т.е. 1/R = 0.01; 0.1; 1; 10.

1.      R=100:

х(1) = x(0) – 2P(х(0) )-1·P (х(0));

х(1) =[2.6585;2.6341]T;

P (x(1),R ) = 10.6829.

2.     R=10:

х(2) = x(1) – 2P(х(1) )-1·P (х(1));

х(2) =[2.6207;2.4828]T;

P (x(2),R ) = 10.7586;

P (x(1),R )- P (x(2),R ) = 10.6829-10.7586=-0.0757.

3.     R=1:

х(3) = x(2) – 2P(х(2) ) -1·P (х(2));

х(3) =[2.5806;2.3226]T;

P (x(3),R ) =10.8387;

P (x(2),R )- P (x(3),R )=-0.0801.

4.     R=0.1:

х(4) = x(3) – 2f (х(3) ) -1·f (х(3));

х(4) =[2.5724;2.2898]T;

P (x(4),R ) =10.8551;

P (x(3),R )- P (x(4),R )=-0.0164.

Сведем все данные в таблицу.

Мы подтвердили, что при увеличении штрафного параметра вес ограничения уменьшается, что доставляет минимум задачи безусловной оптимизации. Наоборот, при уменьшении штрафного параметра до нуля, вес ограничения возрастает, что доставляет минимум задачи условной оптимизации.

Рис 9.  Графическое решение задачи

6.     Заключение

Изучение методов решения оптимизационных задач является необходимым для инженеров в настоящее время, из-за необходимости управления сложными техническими, физическими, электронными процессами, возникающими вследствие развития цивилизации. Так как в последнее время особое внимание уделяется быстрому росту и развитию в сфере разработки высокоинтеллектуальных технологий. К сожалению, до сих пор есть целый ряд оптимизационных задач, точное решение которых так и не было найдено. Например, задачи по транспортировке грузов среди большого количества городов. Так как сложность таких задач возрастает с каждым добавленным городом.

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

Из-за этого следует выбирать метод решения исходя из конкретной задачи, а также вычислительных средств, которые будут необходимы.

Большинство методов хорошо программируются (Приложение 1 –  метод Гаусса-Зайделя), но не все стоит программировать. Методы прямого поиска удобно решать при помощи компьютеров, так как необходимо произвести много однотипных вычислений. С другой стороны градиентные методы требуют информацию о производных, которую придется также вводить в программу. Затрудняется работа по программированию задачи.

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

7.     Приложение 1. Листинг программы

Файл Unit1.cpp

//---------------------------------------------------------------------------

#include <vcl.h>

#pragma hdrstop

#include "Unit1.h"

#include "math.h"

#include "Table.h"

//---------------------------------------------------------------------------

#pragma package(smart_init)

#pragma link "CSPIN"

#pragma resource "*.dfm"

TGaussZeidel *GaussZeidel;

TStringList *ListSeries=new TStringList;

TStringList *ListPoints=new TStringList;

extern TStringList *Points=new TStringList;

bool fl=false;

double fun(double x,double y);

//---------------------------------------------------------------------------

__fastcall TGaussZeidel::TGaussZeidel(TComponent* Owner)

        : TForm(Owner)

{

}

//---------------------------------------------------------------------------

void __fastcall TGaussZeidel::CountLinsChange(TObject *Sender)

{

   LevelLines->UndoZoom();

   double x, y;//координаты точек линий уровня

   double a=StrToFloat(Waa->Text),

          b=StrToFloat(Wbb->Text),

          A=StrToFloat(WA-> Text),

          B=StrToFloat(WB->Text),

          C=StrToFloat(WC->Text);//задаваемые параметры функции

   double atemp, btemp, ctemp, d;//для определения координаты y

   double Atemp, Btemp, Ctemp, D, xmin, xmax;//для определения границ

Atemp=C*C-4*A*B;//определение границ (шаг 1)

Btemp=8*A*B*a-4*B*C*b;//определение границ (шаг 2)

atemp=B;//определени координаты y (шаг 1)

if (ListSeries->Count>CountLins->Value*2)

{ delete ListSeries->Objects[ListSeries->Count-1];

   delete ListSeries->Objects[ListSeries->Count-2];

   ListSeries->Delete(ListSeries->Count-1);

   ListSeries->Delete(ListSeries->Count-1);

   return;

}

/*******************************************************/

double Xmin, Ymin, Zmin;//точное решение

Xmin=(4*A*B*a-2*B*C*b)/(4*A*B-C*C);

Ymin=b-C*Xmin/(2*B);

Zmin=fun(Xmin,Ymin);

/*******************************************************/

double z=Zmin+CountLins->Value*5;

   TFastLineSeries *Series1=new TFastLineSeries(GaussZeidel);

   ListSeries->AddObject("",Series1);

   Series1->ParentChart=LevelLines;

   Series1->SeriesColor=clBlue;

   TFastLineSeries *Series2=new TFastLineSeries(GaussZeidel);

   ListSeries->AddObject("",Series2);

   Series2->ParentChart=LevelLines;

   Series2->SeriesColor=clBlue;

   Ctemp=4*B*z-4*A*B*a*a;//определение границ (шаг 3)

   D=Btemp*Btemp-4*Atemp*Ctemp;//определение границ (шаг 4)

   xmin=(-Btemp+sqrt(D))/(2*Atemp);//определение границ (шаг 5)

   xmax=(-Btemp-sqrt(D))/(2*Atemp);//определение границ (шаг 6)

   for (x=xmin;x<=xmax;x+=0.02)

   { btemp=C*x-2*B*b;//определение координаты y (шаг 2)

     ctemp=A*x*x-2*A*a*x+A*a*a+B*b*b-z;//определение координаты y (шаг 3)

     d=btemp*btemp-4*atemp*ctemp;//определение координаты y (шаг 4)

     if (d<0) d=0;//из-за округления при x=xmin или x=xmax d м. б. <0

     y=(-btemp+sqrt(d))/(2*atemp);//определение координаты y (шаг 5)

     Series1->AddXY(x,y,"",clBlue);

     y=(-btemp-sqrt(d))/(2*atemp);//определение координаты y (шаг 6)

     Series2->AddXY(x,y,"",clBlue);

   }

   x=xmax;//в цикле x не достигает точно xmax

   btemp=C*x-2*B*b;//определение координаты y (шаг 2)

   ctemp=A*x*x-2*A*a*x+A*a*a+B*b*b-z;//определение координаты y (шаг 3)

   d=btemp*btemp-4*atemp*ctemp;//определение координаты y (шаг 4)

   if (d<0) d=0;//из-за округления при x=xmin или x=xmax d м. б. <0

   y=(-btemp+sqrt(d))/(2*atemp);//определение координаты y (шаг 5)

   Series1->AddXY(x,y,"",clBlue);

   y=(-btemp-sqrt(d))/(2*atemp);//определение координаты y (шаг 6)

   Series2->AddXY(x,y,"",clBlue);

}

//---------------------------------------------------------------------------

void __fastcall TGaussZeidel::CountLinsKeyPress(TObject *Sender, char &Key)

{

  Key=0;

}

//---------------------------------------------------------------------------

void __fastcall TGaussZeidel::FormDestroy(TObject *Sender)

{

  for (int i=0;i<ListSeries->Count;i++)

    delete ListSeries->Objects[i];

  delete ListSeries;

  for (int i=0;i<ListPoints->Count;i++)

    delete ListPoints->Objects[i];

  delete ListPoints;

  delete Points;

}

//---------------------------------------------------------------------------

void __fastcall TGaussZeidel::ParFunClick(TObject *Sender)

{

  LevelLines->UndoZoom();

  for (int i=ListSeries->Count-1;i>=0;i--)

  { delete ListSeries->Objects[i];

    ListSeries->Delete(i);

  }

  for (int i=ListPoints->Count-1;i>=0;i--)

  { delete ListPoints->Objects[i];

    ListPoints->Delete(i);

  }

   StartPoint->Clear();

   PointOpt->Clear();

   double x, y;//координаты точек линий уровня

   double a, b, A, B, C;//задаваемые параметры функции

   try

   { a=StrToFloat(Waa->Text),

     b=StrToFloat(Wbb->Text),

     A=StrToFloat(WA-> Text),

     B=StrToFloat(WB->Text),

     C=StrToFloat(WC->Text);

   }

   catch (EConvertError &)

   { Start->Enabled=false;

     Application->MessageBox("Проверте параметры функции","Ошибка!",MB_OK);

     return;

   }

   if (4*A*B-C*C==0)

   { Application->MessageBox("Функция не имеет экстремума","Предупреждение",MB_OK);

     return;

   }

   if (4*A*B-C*C<0)

   { Application->MessageBox("Функция имеет максимум","Предупреждение",MB_OK);

     return;

   }

  CountLins->Enabled=true;

  ParFun->Enabled=false;

  if (!Point->ModalResult || Point->ModalResult==2)

  { StartPoint->AddXY(StrToFloat(Wx->Text),StrToFloat(Wy->Text),"",clLime);

    Point->ModalResult=2;

  }

  if (!Step->Checked) Start->Enabled=true;

  else Stepp->Enabled=true;

  Step->Enabled=true;

   double atemp, btemp, ctemp, d;//для определения координаты y

   double Atemp, Btemp, Ctemp, D, xmin, xmax;//для определения границ

/*******************************************************/

double Xmin, Ymin, Zmin;//точное решение

Xmin=(4*A*B*a-2*B*C*b)/(4*A*B-C*C);

Ymin=b-C*Xmin/(2*B);

Zmin=fun(Xmin,Ymin);

PointOpt->AddXY(Xmin,Ymin,"",clBlue);

/*******************************************************/

Atemp=C*C-4*A*B;//определение границ (шаг 1)

Btemp=8*A*B*a-4*B*C*b;//определение границ (шаг 2)

atemp=B;//определени координаты y (шаг 1)

for (double z=Zmin+5;z<Zmin+CountLins->Value*5+1;z+=5)

{ TFastLineSeries *Series1=new TFastLineSeries(GaussZeidel);

   ListSeries->AddObject("",Series1);

   Series1->ParentChart=LevelLines;

   Series1->SeriesColor=clBlue;

   TFastLineSeries *Series2=new TFastLineSeries(GaussZeidel);

   ListSeries->AddObject("",Series2);

   Series2->ParentChart=LevelLines;

   Series2->SeriesColor=clBlue;

   Ctemp=4*B*z-4*A*B*a*a;//определение границ (шаг 3)

   D=Btemp*Btemp-4*Atemp*Ctemp;//определение границ (шаг 4)

   xmin=(-Btemp+sqrt(D))/(2*Atemp);//определение границ (шаг 5)

   xmax=(-Btemp-sqrt(D))/(2*Atemp);//определение границ (шаг 6)

   for (x=xmin;x<=xmax;x+=0.02)

   { btemp=C*x-2*B*b;//определение координаты y (шаг 2)

     ctemp=A*x*x-2*A*a*x+A*a*a+B*b*b-z;//определение координаты y (шаг 3)

     d=btemp*btemp-4*atemp*ctemp;//определение координаты y (шаг 4)

     if (d<0) d=0;//из-за округления при x=xmin или x=xmax d м. б. <0

     y=(-btemp+sqrt(d))/(2*atemp);//определение координаты y (шаг 5)

     Series1->AddXY(x,y,"",clBlue);

     y=(-btemp-sqrt(d))/(2*atemp);//определение координаты y (шаг 6)

     Series2->AddXY(x,y,"",clBlue);

   }

   x=xmax;//в цикле x не достигает точно xmax

   btemp=C*x-2*B*b;//определение координаты y (шаг 2)

   ctemp=A*x*x-2*A*a*x+A*a*a+B*b*b-z;//определение координаты y (шаг 3)

   d=btemp*btemp-4*atemp*ctemp;//определение координаты y (шаг 4)

   if (d<0) d=0;//из-за округления при x=xmin или x=xmax d м. б. <0

   y=(-btemp+sqrt(d))/(2*atemp);//определение координаты y (шаг 5)

   Series1->AddXY(x,y,"",clBlue);

   y=(-btemp-sqrt(d))/(2*atemp);//определение координаты y (шаг 6)

   Series2->AddXY(x,y,"",clBlue);

}

}

//---------------------------------------------------------------------------

void __fastcall TGaussZeidel::StartClick(TObject *Sender)

{

  Table->Enabled=false;

  LevelLines->UndoZoom();

  if (!fl)

  {

   for (int i=ListPoints->Count-1;i>=0;i--)

   { delete ListPoints->Objects[i];

     ListPoints->Delete(i);

   }

   Points->Clear();

   double e=StrToFloat(We->Text),

          h=StrToFloat(Wh->Text),

          k=StrToFloat(Wk->Text),

          x, y;

   try

   { x=StrToFloat(Wx->Text),

     y=StrToFloat(Wy->Text);

   }

   catch (EConvertError &)

   { Application->MessageBox("Проверте параметры поиска","Ошибка!",MB_OK);

     return;

   }

  StartPoint->Clear();

  if (!Point->ModalResult || Point->ModalResult==2)

  { StartPoint->AddXY(StrToFloat(Wx->Text),StrToFloat(Wy->Text),"",clLime);

    Point->ModalResult=2;

  }

  Stop->Enabled=true;

  Start->Enabled=false;

  if (!Step->Checked) Pause->Enabled=true;

  WA->Enabled=false; WB->Enabled=false; WC->Enabled=false;

  Waa->Enabled=false; Wbb->Enabled=false; We->Enabled=false;

  Wh->Enabled=false; Wk->Enabled=false; Wx->Enabled=false; Wy->Enabled=false;

  FTable->Tablica->RowCount=2;

   double f1, f=fun(x,y), htemp;

   int i=1;

   if (h<=e) h=2*e;

   while (h>e)

   {

     if (f<fun(x+h,y)) htemp=-h; else htemp=h;

     Points->Add(x); Points->Add(y);

     FTable->Tablica->RowCount++;

     FTable->Tablica->Cells[1][i]=h; i++;

     f1=fun(x+=htemp,y);

     Points->Add(x); Points->Add(y);

     FTable->Tablica->RowCount++;

     FTable->Tablica->Cells[1][i]=h; i++;

     if (f1>f) h/=k;

     while (f1<f)

     { f=f1;

       x+=htemp; f1=fun(x,y); Points->Add(x); Points->Add(y);

       FTable->Tablica->RowCount++;