Лабораторная работа №3
.docМинистерство Образования Российской Федерации
Уфимский Государственный Авиационный Технический Университет
Кафедра ТК
Отчет по лабораторной работе №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
Тесты