Добавил:
Опубликованный материал нарушает ваши авторские права? Сообщите нам.
Вуз: Предмет: Файл:

Лабораторная работа №3

.doc
Скачиваний:
17
Добавлен:
02.05.2014
Размер:
215.04 Кб
Скачать

Министерство Образования Российской Федерации

Уфимский Государственный Авиационный Технический Университет

Кафедра ТК

Отчет по лабораторной работе №3

по предмету «Методы оптимизации»

на тему: Методы Ньютона и сопряжённых градиентов

Выполнил: студент гр.Т28-320

Проверил: Хасанов А.Ю.

-Уфа, 2004 -

Лабораторная работа №3

Тема: Методы Ньютона и сопряжённых градиентов

Методы: Метод Ньютона-Рафсона с оптимальным шагом;

Модификация II метода Ньютона

Задание: найти точку минимума функции f(x,y) = При a=15;b=-0,5;c=2,25;d=2,5;

Начальное приближение (0,0) Точность E = 0,0002

Описание алгоритма в псевдокоде:

структура Pnt

начало

вещественные X, Y - вещественные;

алгоритм Pnt operator =(Pnt A)

начало

X=A.X;

Y=A.Y;

конец

конец

алгоритм F(x,y)

начало

скаляр ans - вещественное;

ans=a*x+b*y+exp(c*x*x+d*y*y);

возврат(ans);

конец

алгоритм F(A)

начало

скаляры ans - вещественное;

ans = F(A.X,A.Y);

возврат(ans);

конец

алгоритм dFx(x, y)

начало

скаляры ans - вещественное;

ans=a+2*c*x*exp(c*x*x+d*y*y);

возврат(ans);

конец

алгоритм dFy(x, y)

начало

скаляры ans - вещественное;

ans=b+2*d*y*exp(c*x*x+d*y*y);

возврат(ans);

конец

алгоритм dFxx(x,y)

начало

скаляры ans - вещественное;

ans=2*c*exp(c*x*x+d*y*y)+4*c*c*x*x*exp(c*x*x+d*y*y);

возврат ans;

конец

алгоритм dFyy(x,y)

начало

скаляры ans - вещественное;

ans=2*d*exp(c*x*x+d*y*y)+4*d*d*y*y*exp(c*x*x+d*y*y);

возврат ans;

конец

скаляры dFxy(x,y)

начало

скаляры ans;

ans=4*c*d*x*y*exp(c*x*x+d*y*y);

возврат ans;

конец

скаляры dFx(X)

начало

скаляры ans - вещественное;

ans=dFx(X.X,X.Y);

возврат(ans);

конец

алгоритм dFy(X)

начало

скаляры ans;

ans=dFy(X.X,X.Y);

возврат(ans);

конец

алгоритм dFxx(X)

начало

скаляры ans;

ans=dFxx(X.X,X.Y);

возврат ans;

конец

скаляры dFyy(X)

начало

скаляры ans;

ans=dFyy(X.X,X.Y);

возврат ans;

конец

скаляры dFxy(X)

начало

скаляры ans;

ans=dFxy(X.X,X.Y);

возврат ans;

конец

алгоритм EndSearchNewton(X, E)

начало

dFI+=2;

если(dFy(X)<E/2)&&(dFx(X)<E/2)) возврат 1;

иначе возврат 0;

все-если

конец

алгоритм Newton(X0,E)

начало

структуры X=X0,P - точка;

скаляры A,B,C,D,F - вещественные;

пока(!EndSearchNewton(X,E))

начало

A=dFx(X);

B=dFy(X);

C=dFxx(X);

D=dFyy(X);

F=dFxy(X);

P.Y=((F*A)/C-B)/(D-F*F/C);

P.X=(-A-F*P.Y)/C;

X=X+P;

I++;ddFI+=3;dFI+=2;

вывод X;

конец

возврат X;

конец

алгоритм Fi(a, X, P)

начало

скаляры Ans - вещественное;

Ans=F(X.X+a*P.X,X.Y+a*P.Y);

FI++;

возврат Ans;

конец

алгоритм GetMinFi(E, X, P)

начало

скаляры Lmb=1.618033989,a=0,b=0 - вещественное;

пока (Fi(b+0.1,X,P)<Fi(b,X,P))b+=0.1;

a=b-0.1;

скаляры x1=b-(b-a)/Lmb;

x2=a+(b-a)/Lmb;

y1=Fi(x1,X,P);

y2=Fi(x2,X,P); - вещественные

пока((fabs(b-a)/2)>E)

начало

если (y1>y2)

начало

a=x1;

x1=x2;y1=y2;

x2=a+(b-a)/Lmb;y2=Fi(x2,X,P);

конец

иначе

начало

b=x2;

x2=x1;y2=y1;

x1=b-(b-a)/Lmb;y1=Fi(x1,X,P);

конец

все-если

конец

возврат (b+a)/2;

конец

алгоритм NewtonOpt(X0, E)

начало

структура X=X0, P - точки;

скаляры A,B,C,D,F,alfa - вещественные;

пока(!EndSearchNewton(X,E))

начало

A=dFx(X);

B=dFy(X);

C=dFxx(X);

D=dFyy(X);

F=dFxy(X);

P.Y=((F*A)/C-B)/(D-F*F/C);

P.X=(-A-F*P.Y)/C;

alfa=GetMinFi(0.0001,X,P);

X=X+P*alfa;

I++;ddFI+=3;dFI+=2;

вывод X

конец

возврат X;

конец

алгоритм Newton_mod2(X0, E, M)

начало

структуры X=X0, P - вещественное;

скаляры A,B,C,D,F,alfa - вещественное;

t=M - целое;

пока(!EndSearchNewton(X,E))

начало

A=dFx(X);

B=dFy(X);

если t>=M

начало

C=dFxx(X);

D=dFyy(X);

F=dFxy(X);

t=0;

ddFI+=3;

конец

все-если

P.Y=((F*A)/C-B)/(D-F*F/C);

P.X=(-A-F*P.Y)/C;

alfa=GetMinFi(0.0001,X,P);

X=X+P*alfa;

I++;dFI+=2;

t++;

вывод X;

конец

возврат X;

конец

алгоритм btnSrchClick(TObject *Sender)

начало

cтруктура X0 - вещественное;

X0.X=StrToFloat(edtX1->Text);

X0.Y=StrToFloat(edtX2->Text);

E=StrToFloat(edtE->Text);

a=StrToFloat(edtA->Text);

b=StrToFloat(edtB->Text);

c=StrToFloat(edtC->Text);

d=StrToFloat(edtD->Text);

M=StrToInt(edtM->Text);

структура Ans - вещественное;

если (rbNewton->Checked)

начало I=0; dFI=0;ddFI=0;FI=0;Ans=Newton(X0,E);конец все-если

если (rbNewton_opt->Checked)

начало I=0; dFI=0;ddFI=0;FI=0;Ans=NewtonOpt(X0,E);конец все-если

если (rbNewton_mod2->Checked)

начало I=0; dFI=0;ddFI=0;FI=0;Ans=Newton_mod2(X0,E,M);конец все-если

вывод Ans;

вывод ddFI;

вывод dFI;

вывод FI;

вывод I;

вывод (ddFI*39+dFI*14+FI*6);

вывод (ddFI*7+dFI*4+FI*3);

вывод (ddFI*5+dFI*2+FI);

конец

Листинг программы

Main.cpp

#include <vcl.h>

#include <math.h>

#pragma hdrstop

#include "Main.h"

#pragma package(smart_init)

#pragma resource "*.dfm"

long double F(long double x,long double y)

{

long double ans;

ans=a*x+b*y+exp(c*x*x+d*y*y);

return(ans);

}

long double F(Pnt A)

{

long double ans;

ans = F(A.X,A.Y);

return(ans);

}

long double dFx(long double x,long double y)

{

long double ans;

ans=a+2*c*x*exp(c*x*x+d*y*y);

return(ans);

}

long double dFy(long double x,long double y)

{

long double ans;

ans=b+2*d*y*exp(c*x*x+d*y*y);

return(ans);

}

long double dFxx(long double x, long double y)

{

long double ans;

ans=2*c*exp(c*x*x+d*y*y)+4*c*c*x*x*exp(c*x*x+d*y*y);

return ans;

}

long double dFyy(long double x, long double y)

{

long double ans;

ans=2*d*exp(c*x*x+d*y*y)+4*d*d*y*y*exp(c*x*x+d*y*y);

return ans;

}

long double dFxy(long double x, long double y)

{

long double ans;

ans=4*c*d*x*y*exp(c*x*x+d*y*y);

return ans;

}

long double dFx(Pnt X)

{

long double ans;

ans=dFx(X.X,X.Y);

return(ans);

}

long double dFy(Pnt X)

{

long double ans;

ans=dFy(X.X,X.Y);

return(ans);

}

long double dFxx(Pnt X)

{

long double ans;

ans=dFxx(X.X,X.Y);

return ans;

}

long double dFyy(Pnt X)

{

long double ans;

ans=dFyy(X.X,X.Y);

return ans;

}

long double dFxy(Pnt X)

{

long double ans;

ans=dFxy(X.X,X.Y);

return ans;

}

TfrmMain *frmMain;

__fastcall TfrmMain::TfrmMain(TComponent* Owner)

: TForm(Owner)

{

}

void __fastcall TfrmMain::Button1Click(TObject *Sender)

{

Close();

}

int __fastcall TfrmMain::EndSearchNewton(Pnt X, long double E)

{

dFI+=2;

if((dFy(X)<E/2)&&(dFx(X)<E/2)) return 1;

else return 0;

}

Pnt __fastcall TfrmMain::Newton(Pnt X0, long double E)

{

Pnt X=X0;

Pnt P;

long double A,B,C,D,F;

while(!EndSearchNewton(X,E))

{

A=dFx(X);

B=dFy(X);

C=dFxx(X);

D=dFyy(X);

F=dFxy(X);

P.Y=((F*A)/C-B)/(D-F*F/C);

P.X=(-A-F*P.Y)/C;

X=X+P;

I++;ddFI+=3;dFI+=2;

Series1->AddXY(X.X,X.Y,"",clWhite);

}

return X;

}

long double Fi(long double a,Pnt X, Pnt P)

{

long double Ans;

Ans=F(X.X+a*P.X,X.Y+a*P.Y);

FI++;

return Ans;

}

long double GetMinFi(long double E,Pnt X, Pnt P)

{

long double Lmb=1.618033989;

long double a=0;

long double b=0;

while (Fi(b+0.1,X,P)<Fi(b,X,P))b+=0.1;

a=b-0.1;

long double x1=b-(b-a)/Lmb;

long double x2=a+(b-a)/Lmb;

long double y1=Fi(x1,X,P);

long double y2=Fi(x2,X,P);

while((fabs(b-a)/2)>E)

{

if(y1>y2)

{

a=x1;

x1=x2;y1=y2;

x2=a+(b-a)/Lmb;y2=Fi(x2,X,P);

}

else

{

b=x2;

x2=x1;y2=y1;

x1=b-(b-a)/Lmb;y1=Fi(x1,X,P);

}

}

return (b+a)/2;

}

Pnt __fastcall TfrmMain::NewtonOpt(Pnt X0, long double E)

{

Pnt X=X0;

Pnt P;

long double A,B,C,D,F,alfa;

while(!EndSearchNewton(X,E))

{

A=dFx(X);

B=dFy(X);

C=dFxx(X);

D=dFyy(X);

F=dFxy(X);

P.Y=((F*A)/C-B)/(D-F*F/C);

P.X=(-A-F*P.Y)/C;

alfa=GetMinFi(0.0001,X,P);

X=X+P*alfa;

I++;ddFI+=3;dFI+=2;

Series2->AddXY(X.X,X.Y,"",clRed);

}

return X;

}

Pnt __fastcall TfrmMain::Newton_mod2(Pnt X0, long double E, int M)

{

Pnt X=X0;

Pnt P;

long double A,B,C,D,F,alfa;

int t=M;

while(!EndSearchNewton(X,E))

{

A=dFx(X);

B=dFy(X);

if (t>=M)

{

C=dFxx(X);

D=dFyy(X);

F=dFxy(X);

t=0;

ddFI+=3;

}

P.Y=((F*A)/C-B)/(D-F*F/C);

P.X=(-A-F*P.Y)/C;

alfa=GetMinFi(0.0001,X,P);

X=X+P*alfa;

I++;dFI+=2;

t++;

Series3->AddXY(X.X,X.Y,"",clGreen);

}

return X;

}

void __fastcall TfrmMain::btnSrchClick(TObject *Sender)

{

Pnt X0;

X0.X=StrToFloat(edtX1->Text);

X0.Y=StrToFloat(edtX2->Text);

E=StrToFloat(edtE->Text);

a=StrToFloat(edtA->Text);

b=StrToFloat(edtB->Text);

c=StrToFloat(edtC->Text);

d=StrToFloat(edtD->Text);

M=StrToInt(edtM->Text);

Pnt Ans;

if (rbNewton->Checked) {I=0; dFI=0;ddFI=0;FI=0;Ans=Newton(X0,E);}

if (rbNewton_opt->Checked) {I=0; dFI=0;ddFI=0;FI=0;Ans=NewtonOpt(X0,E);}

if (rbNewton_mod2->Checked) {I=0; dFI=0;ddFI=0;FI=0;Ans=Newton_mod2(X0,E,M);}

edt_x->Text=FloatToStr(Ans.X);

edt_y->Text=FloatToStr(Ans.Y);

edtddFI->Text=FloatToStr(ddFI);

edtdFI->Text=FloatToStr(dFI);

edtFI->Text=FloatToStr(FI);

edtI->Text=FloatToStr(I);

edtMlt->Text=FloatToStr(ddFI*39+dFI*14+FI*6);

edtSum->Text=FloatToStr(ddFI*7+dFI*4+FI*3);

edtExp->Text=FloatToStr(ddFI*5+dFI*2+FI);

}

void __fastcall TfrmMain::btnClearClick(TObject *Sender)

{

Series1->Clear();

Series2->Clear();

Series3->Clear();

}

Main.h

#ifndef MainH

#define MainH

#include <Classes.hpp>

#include <Controls.hpp>

#include <StdCtrls.hpp>

#include <Forms.hpp>

#include <ExtCtrls.hpp>

#include <Chart.hpp>

#include <Series.hpp>

#include <TeEngine.hpp>

#include <TeeProcs.hpp>

struct Pnt

{

long double X;

long double Y;

Pnt operator =(Pnt A)

{

X=A.X;

Y=A.Y;

}

};

Pnt operator +(Pnt A, Pnt B)

{

Pnt C;

C.X=A.X+B.X;

C.Y=A.Y+B.Y;

return C;

};

Pnt operator -(Pnt A, Pnt B)

{

Pnt C;

C.X=A.X-B.X;

C.Y=A.Y-B.Y;

return C;

};

Pnt operator *(Pnt A, long double c)

{

Pnt C;

C.X=A.X*c;

C.Y=A.Y*c;

return C;

};

long double a,b,c,d,E,M,I,dFI,ddFI,FI;

class TfrmMain : public TForm

{

__published: // IDE-managed Components

TEdit *edtA;

TLabel *Label1;

TEdit *edtB;

TEdit *edtC;

TEdit *edtD;

TLabel *Label2;

TLabel *Label3;

TLabel *Label4;

TLabel *Label5;

TEdit *edtX1;

TLabel *Label6;

TEdit *edtX2;

TLabel *Label7;

TEdit *edtE;

TLabel *Label8;

TLabel *Label9;

TRadioGroup *RadioGroup1;

TRadioButton *rbNewton;

TRadioButton *rbNewton_opt;

TRadioButton *rbNewton_mod2;

TButton *btnSrch;

TButton *Button1;

TEdit *edt_x;

TLabel *Label10;

TEdit *edt_y;

TLabel *Label11;

TLabel *Label12;

TEdit *edtI;

TEdit *edtddFI;

TEdit *edtdFI;

TChart *Chart1;

TLabel *Label13;

TLabel *Label15;

TEdit *edtM;

TLabel *Label16;

TLabel *Label18;

TEdit *edtFI;

TLabel *Label14;

TEdit *edtMlt;

TLabel *Label17;

TEdit *edtSum;

TEdit *edtExp;

TButton *btnClear;

TLineSeries *Series1;

TLineSeries *Series2;

TLineSeries *Series3;

TLabel *Label19;

TLabel *Label20;

TLabel *Label21;

void __fastcall Button1Click(TObject *Sender);

void __fastcall btnSrchClick(TObject *Sender);

void __fastcall btnClearClick(TObject *Sender);

private: // User declarations

public: // User declarations

__fastcall TfrmMain(TComponent* Owner);

Pnt __fastcall Newton(Pnt X0, long double E);

int __fastcall EndSearchNewton(Pnt X, long double E);

Pnt __fastcall NewtonOpt(Pnt X0, long double E);

Pnt __fastcall Newton_mod2(Pnt X0, long double E, int M);

};

extern PACKAGE TfrmMain *frmMain;

#endif

Тесты