Смекни!
smekni.com

Поиск максимума одной функции многих переменных методом покоординатного спуска и с помощью метода (стр. 2 из 3)

Перебирая все xi, найдем максимум функции.

Перебирая всевозможные параметры p и q, получим некоторые наборы

(в зависимости от p и q) на которых функция достигает максимума.

3. Решение задачи с использованием метода покоординатного спуска

3.1 Описание метода покоординатного спуска

Изложим этот метод на примере функции трех переменных

.

Выберем нулевое приближение

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

Затем из новой точки сделаем спуск по направлению, параллельному оси

, т. е. рассмотрим
, найдем ее минимум и обозначим его через
. Второй шаг приводит нас в точку
. Из этой точки делаем третий шаг – спуск параллельно оси
и находим минимум функции
. Приход в точку
завершает цикл спусков.

Будем повторять циклы. На каждом спуске функция не возрастает, и при этом значения функции ограничены снизу ее значением в минимуме

. Следовательно, итерации сходятся к некоторому пределу
. Будет ли здесь иметь место равенство, т. е. сойдутся ли спуски к минимуму и как быстро?

Это зависит от функции и выбора нулевого приближения.

Метод спуска по координатам несложен и легко программируется на ЭВМ. Но сходится он медленно. Поэтому его используют в качестве первой попытки при нахождении минимума.


3.2 Алгоритм решения

Будем перебирать с с шагом какой-либо малой длины.

Зададим нулевое приближение, например:

Найдем частные производные

.

Пусть при каком-то j

Тогда k-ое приближение считаем по формулам:

Шаг t будем выбирать таким образом, чтобы

и
.

В результате, закончив процесс по критерию

, где
-заданная точность мы придем к набору
, при котором функция f максимальна.

Подставим найденный набор

и соответствующее
в функцию f1=
и перебрав все с, выберем те
, при которых f1 минимальна.

Заключение

В ходе решения поставленной задачи рассмотрены случаи, когда n=4,5,6. Были найдены две основные области наборов

:

1) xi=0 или 1(количество 0 и 1 одинаково)

2) xi=0.5,

.

Причем участок 1<p<1.5 покрыт первой областью, при q>>p

–– из первой области, при q≈p
–– из второй области, а при p→∞ область делилась между ними примерно пополам.

На участке p>2 появлялись области вида:

i<k => xi=0;

i>n-k => xi=1;

k-1<i<n-k+1=> xi=0.5.

На участке 2<q<3 p

2 существует область, в которой максимум достигается при
вида:

xi=a => xn-i=1-a, 0<a<1.


Список использованных источников

1. Амосов А.А., Дубинский Ю.А., Копченова Н.В. Вычислительные методы для инженеров. М.: Высшая школа, 1994. 543с.

2. Березин И.С. и Жидков Н. П. Методы вычислений. т.1. М.: “Наука”, 1965. 633c.

3. Подбельский В.В. и Фомин С.С. Программирование на языке Си. М.: “Финансы и статистика”, 2000. 599с.


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

//вывод на экран областей максимума функции

#include "stdafx.h"

#include "KE2.h"

#include "math.h"

#include "KE2Doc.h"

#include "KE2View.h"

#ifdef _DEBUG

#define new DEBUG_NEW

#undef THIS_FILE

static char THIS_FILE[] = __FILE__;

#endif

IMPLEMENT_DYNCREATE(CKE2View, CView)

BEGIN_MESSAGE_MAP(CKE2View, CView)

//{{AFX_MSG_MAP(CKE2View)

// NOTE - the ClassWizard will add and remove mapping macros here.

// DO NOT EDIT what you see in these blocks of generated code!

//}}AFX_MSG_MAP

// Standard printing commands

ON_COMMAND(ID_FILE_PRINT, CView::OnFilePrint)

ON_COMMAND(ID_FILE_PRINT_DIRECT, CView::OnFilePrint)

ON_COMMAND(ID_FILE_PRINT_PREVIEW, CView::OnFilePrintPreview)

END_MESSAGE_MAP()

CKE2View::CKE2View()

{

}

CKE2View::~CKE2View()

{

}

BOOL CKE2View::PreCreateWindow(CREATESTRUCT& cs)

{

return CView::PreCreateWindow(cs);

}

void CKE2View::OnDraw(CDC* pDC)

{

CKE2Doc* pDoc = GetDocument();

ASSERT_VALID(pDoc);

drawPoint(pDC);

// TODO: add draw code for native data here

}

BOOL CKE2View::OnPreparePrinting(CPrintInfo* pInfo)

{

// default preparation

return DoPreparePrinting(pInfo);

}

void CKE2View::OnBeginPrinting(CDC* /*pDC*/, CPrintInfo* /*pInfo*/)

{

// TODO: add extra initialization before printing

}

void CKE2View::OnEndPrinting(CDC* /*pDC*/, CPrintInfo* /*pInfo*/)

{

// TODO: add cleanup after printing

}

#ifdef _DEBUG

void CKE2View::AssertValid() const

{

CView::AssertValid();

}

void CKE2View::Dump(CDumpContext& dc) const

{

CView::Dump(dc);

}

CKE2Doc* CKE2View::GetDocument() // non-debug version is inline

{

ASSERT(m_pDocument->IsKindOf(RUNTIME_CLASS(CKE2Doc)));

return (CKE2Doc*)m_pDocument;

}

#endif //_DEBUG

int sgn(float a)

{

int sg;

if (a>0) sg=1;

if (a==0) sg=0;

if (a<0) sg=-1;

return(sg);

}

#define n 6

void CKE2View::drawPoint(CDC *pDC)

{

double **c,*f1,*f,*x,*w,*e,max,p=2,q=2,xx,yy;

int i=0,j=0,k,m,a,b,*l,bb=0;

c=new double*[10000];

for(i=0;i<10000;i++)

{

c[i]=new double[3];

memset(c[i],0,sizeof(double)*3);

}

f=new double[10000];

e=new double[10000];

w=new double[10000];

f1=new double[10000];

x=new double[n];

l=new int[10000];

for(xx=0.5;xx<1;xx+=0.01)

for(yy=xx;yy<1);yy+=0.01)

{

p=1./(1.-xx);

q=1./(1.-yy);

memset(w,0,sizeof(double)*10000);

memset(e,0,sizeof(double)*10000);

memset(f1,0,sizeof(double)*10000);

memset(x,0,sizeof(double)*n);

x[n-1]=1;

j=0;

for(i=0;i<10000;i++)

{j=0;

f1[i]=1;c[i][0]=0;c[i][1]=1;c[i][2]=0.5;

while(fabs(f1[i])>0.00000001)

{

f1[i]=0;

for(k=0;k<n;k++)

{ f1[i]+=pow((fabs(x[k]-c[i][2])),(p-1))*sgn(x[k]-c[i][2]);}

if (f1[i]<-0.00000001)

{max=c[i][2];c[i][2]=c[i][2]-(fabs(c[i][2]-c[i][1])/2.0);c[i][1]=max;}

if (f1[i]>0.00000001)

{max=c[i][2];c[i][2]=c[i][2]+(fabs(c[i][2]-c[i][1])/2.0);c[i][1]=max;}

if (fabs(f1[i])<=0.00000001)

{c[i][0]=c[i][2];goto B;}

}

B:

c[i][0]=c[i][2];

for(a=0;a<n;a++)

{

for(b=0;b<n;b++)

w[i]+=pow((fabs(x[a]-x[b])),q);

e[i]+=pow((fabs(x[a]-c[i][0])),p);

}

f[i]=pow((e[i]/n),(1./p))/pow((w[i]/(n*n)),(1./q));

x[n-2]+=0.1;

for(k=2;k<n;k++)

{

if(x[n-k]>1.04)

{

x[n-k-1]+=0.1;

x[n-k]=x[n-k-1];

for(m=2;m<k;m++)

x[n-m]=x[n-k-1];

}

if (x[0]!=0) goto A;

}

}

A:

max=f[0];k=0;

for(m=0;m<i;m++)

{

if (fabs(max-f[m])<0.001) {k++;l[k]=m;}

if (max<f[m]) {k=0;max=f[m];l[k]=m;}

}

for(a=0;a<n-1;a++)

x[a]=0;


for(a=0;a<l[0];a++)

{

x[n-2]+=0.1;

for(k=2;k<n;k++)

if(x[n-k]>1.04)

{

x[n-k-1]+=0.1;

x[n-k]=x[n-k-1];

for(m=2;m<k;m++)

x[n-m]=x[n-k-1];

}

}

b=0;

for(k=0;k<n;k++)

{

if((x[k]==0)||(fabs(x[k]-1)<0.04)) b++;

else

{

if(fabs(x[k]-0.5)<0.04) b+=2;

else b=-n;

}

}

b-=n;

if (b<0) b=24;

if (b==0) b=58;

if(b==bb) continue;

bb=b;

c=%f&bsol;n",p,q,l[0],l[k],k+1,max,c[l[0]][0]);

COLORREF cr(RGB((b%3)*127,(b%4)*85,(b%5)*63));

CBrush r(cr);

CPen rp(PS_SOLID,0,cr);

pDC->SelectObject(&rp);

pDC->SelectObject(&r);

CPoint r1[3]={CPoint(0,360),CPoint(int(720./p),-int(720./q)+360),CPoint(int(720./p),360)};

pDC->Polygon(r1,3);

}

}

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

//Покоординатный спуск

#include<stdAfx.h>

#include<stdio.h>

#include<iostream.h>

#include<conio.h>