Sp-алгоритм метода решений частичных проблем собственных значений. Градиентный метод с дроблением шага

0

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

 

на тему: « Sp-алгоритм метода решений частичных проблем собственных значений. Градиентный метод с дроблением шага.

Метод Девидона-Флетчера- Пауэлла»

 

 

 

 

Содержание

 

СОБСТВЕННЫЕ ЗНАЧЕНИЯ И СОБСТВЕННЫЕ ВЕКТОРЫ МАТРИЦ. 3

  1. Основные понятия. 3
  2. Частичная проблема собственных значений. 6

2.1 Метод простой итерации. 6

2.2 Степенной метод. 7

2.3 Метод скалярных произведений. 9

  1. Градиентный метод с дроблением шага. 10

3.1 Метод градиентного спуска с дроблением шага. 11

  1. Метод сопряженных градиентов. 13

Приложение А. 15

Приложение Б. 19

Приложение В. 24

Список используемых источников. 31

 

 

СОБСТВЕННЫЕ ЗНАЧЕНИЯ И СОБСТВЕННЫЕ ВЕКТОРЫ МАТРИЦ.

1. Основные понятия

При решении многих теоретических и прикладных задач возникает необходимость определения собственных значений и собственных векторов матриц.

Собственным значением матрицы A, называется такое число , которое удовлетворяет уравнению:

 

(1)

а вектор x, соответствующий данному собственному значению и удовлетворяющий уравнению (1), называется собственным вектором матрицы A.

Перенеся неизвестные уравнения (1) в левую часть, получим:

 

(2)

Это однородная система линейных уравнений. Она имеет ненулевое решение x лишь при условии :

 

(3)

Матрица (A −λE) называется характеристической и имеет вид

 

.

Определитель называется характеристическим определителем. В развернутом виде этот определитель есть многочлен n-й степени от и называется характеристическим многочленом. Корни этого многочлена являются собственными значениями матрицы A, иначе называемые характеристическими числами многочлена

 

.

 

(4)

 

Числа,,…называются коэффициентами характеристического многочлена. Таким образом, собственные значения матрицы A определяются из уравнения (3) или (4), а соответствующие им собственные векторы x – из уравнения (2), которое можно записать в виде

 

.

 

Эта система уравнений называется характеристической. Ее решение неединственно, так как обычно одно уравнение является линейной комбинацией других. Решения получаются с точностью до множителя. Чтобы избавиться от многозначности, собственные векторы нормируют – делят компоненты вектора на какую-либо его норму – или одну какую-то его компоненту принимают равной единице. Тогда остальные компоненты определяются однозначно.

Задачи определения собственных значений и собственных векторов матриц обычно подразделяется на две:

1) задачу определения всех собственных значений и принадлежащих им собственных векторов матриц, называемую полной проблемой собственных значений;

2) задачу определения одного или нескольких собственных значений и

принадлежащих им собственных векторов, называемую частичной проблемой собственных значений.

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

Некоторые свойства собственных значений матрицы:

1) все собственные значения симметричной матрицы действительны;

2) если собственные значения матрицы действительны и различны, то соответствующие им собственные векторы ортогональны и образуют базис рассматриваемого пространства;

3) если две матрицы A и B подобны, то есть если они связаны соотношением подобия , то их собственные значения совпадают, P – некоторая матрица.

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

Самым лучшим упрощением матрицы было бы приведение ее к треугольному виду

 

Тогда определитель этой матрицы равен произведению диагональных элементов. В этом случае характеристический полином имеет вид

 

и собственные значения матрицы определяются непосредственно:

 

 

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

 

  1. Частичная проблема собственных значений

 

Для решения частичной проблемы собственных значений, когда определяют одно или несколько собственных значений и собственных векторов, обычно используют итерационные методы. Рассмотрим некоторые из них.

 

2.1 Метод простой итерации

 

Этот метод используется для вычисления наибольшего собственного значения матрицы A и соответствующего ему собственного вектора. Это следует из следующей теоремы. Если матрица A обладает простым превалирующим собственным значением , таким, что , то векторная последовательность

 

      

 

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

Для обеспечения устойчивости метода и уменьшения погрешности расчетов необходимо нормировать итерируемый вектор на каждом шаге в процессе расчетов:

где – норма вектора , – компонента вектора ;

 

– нормированный вектор, i = 1, 2,.., n . В противном случае компоненты вектора увеличиваются примерно в раз после каждой итерации, что в конечном счете может привести к переполнению разрядной сетки ЭВМ

и увеличению погрешности расчета.

Алгоритм расчета в векторном виде записывается следующим образом:

 

 

(5)

– нормированный вектор.

Процесс итераций прекращается при выполнении условия:

где – заданная погрешность вычисления.

После определения собственного вектора соответствующее ему собственное значение вычисляется по формуле Релея

, где

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

Для этого умножим систему на обратную матрицу :

,

а затем, разделив обе части уравнения на и учитывая, что , получаем:

.

Эта задача отличается от рассмотренной выше тем, что здесь вычисляется наибольшее собственное значение обратной матрицы, что достигается при наименьшем для матрицы , так как собственные значения матриц и обратны друг другу. Следовательно, рассмотренный выше итерационный процесс (5) может быть использован также для нахождения наименьшего собственного значения обратной матрицы .

 

2.2 Степенной метод.

 

Пусть матрица является матрицей простой структуры, т.е. имеет ровно n линейно независимых собственных векторов (базис) x1, x2,…,xn,, и пусть . Ставится задача приближенного вычисления наибольшего по модулю собственного числа и соответствующего ему собственного вектора х1 матрицы .

Для нахождения выберем Y(0)  — начальный вектор. Y(0)  следует выбирать так, чтобы коэффициент в разложении Y(0) =x1+x2+…xn  не был бы слишком мал. Здесь x1, x2,…,xn,— собственные векторы, так что Axi=xi, =1,…,n. Y(0) выбирается опытным путем.

Далее строится последовательность векторов Y(1), Y(2),… по формуле Y(k+1)=AY(k).

.

Пусть о вещественной матрицы известно, что она является матрицей простой структуры, т.е имеет ровно линейнонезависимых собственных векторов(базис):

 

(6)

Пусть нумерация этих векторов отвечает упорядочению соответствующих им собственных чисел по убыванию модулей: .

Возьмем произвольный ненулевой вектор и запишем его разложение по базису собственных векторов:

Далее умножим на матрицу :

в силу (1) данное уравнение можно представить в виде .

Следуя этому же принципу в результате второй итерации получим:

Таким образом на той итерации вектора с помощью матрицы получим вектор :

 

(7)

Или, с учетом представления x1, x2,…,xn,, в исходном базисе

Возьмём отношение компонент итерированного вектора:

 

(8)

Представляя вектор на основе (7) :

 

(9)

можно сделать вывод, что в силу , в фигурных скобках выражения (9) линейной комбинации векторов x1, x2,…,xn,, с ростом начнет доминировать первое слагаемое. Это означает, что вектор от итерации к итерации будет давать все большее хорошее приближение к собственному вектору x1 по направлению.

2.3 Метод скалярных произведений.

Модификация степенного метода называется методом скалярных произведений. Реализовать её можно в виде SP-алгоритма.

Алгоритм метода

Шаг 1. Дана симметричная -матрица , число для начального сравнения и произвольный начальный вектор ,где сколь угодно малое число. Положим .

Шаг 2. Вычистить скаляры :

 
 

и вектор:

 

.

Шаг 3.Вычислить (итерация нормированного вектора).

Шаг 4. Вычислить :

 
 
 
 
 

Шаг 5.Если , то положить и вернуться к шагу 3, иначе завершить работу алгоритма, считая .                        

3. Градиентный метод с дроблением шага.

 

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

где - параметр, определяющий длину и направление очередного шага;

() - антиградиент функции в текущей точке.

Вектор задает направление убывания функции в точке x, если .

Градиентные методы представляют собой методы спуска, в которых в качестве направлений убывания выбирается направление антиградиента функции.

,где   ()

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

Начиная с некоторого x0 , строится последовательность такая, что , при всех .

Критерием остановки итерационной последовательности являются обычно следующие неравенства:

     или     .

В градиентном методе с дроблением шага точка определяется по формуле

 

(10)

 

 

где величина шага находится из условия:

 

(11)

где -градиент функции в точке .

 

Алгоритм:

 

  1. Задаем начальную точку , начальную величину шага и коэффициент дробления шага . Полагаем счетчик числа итераций =0.
  2. По формуле (10) вычисляем компоненты вектора.
  3. Вычисляем величину - значение функции в точке.
  4. Если условие (11) выполнено, то переходим к следующему пункту. Иначе – переходим к пункту 6.
  5. Полагаем и переходим к пункту 2.
  6. Проверяем условие окончания поиска (см. ниже). Если условие окончания поиска выполнено, то полагаем , и завершаем итерации. Иначе – полагаем переходим к п. 2.

В качестве критерия окончания поиска можно использоваться условия (3) .

 

3.1 Метод градиентного спуска с дроблением шага.

В этом варианте градиентного метода величина шага на каждой итерации выбирается из условия выполнения неравенства:

 

где — некоторая заранее выбранная константа. Условие гарантирует (если, конечно, такие удастся найти), что получающаяся последовательность будет релаксационной. Процедуру нахождения такого обычно оформляют так.

Выбирается число и некоторый начальный шаг . Теперь для каждого полагают и делают шаг градиентного метода. Если с таким условие выполняется, то переходят к следующему . Если же условие не выполняется, то умножают на ("дробят шаг") и повторяют эту процедуру до тех пор пока не будет выполняться равенство:

 

В условиях теоремы об условной сходимости градиентного метода с постоянным шагом эта процедура для каждого за конечное число шагов приводит к нужному .

Можно показать, что в условиях теоремы (о линейной сходимости градиентного метода с постоянным щагом (Рисунок 1)) градиентный метод с дроблением шага линейно сходится. Описанный алгоритм избавляет нас от проблемы выбора  на каждом шаге, заменяя ее на проблему выбора параметров , и , к которым градиентный метод менее чувствителен. При этом, разумеется, объем вычислений возрастает (в связи с необходимостью процедуры дробления шага), впрочем, не очень сильно, поскольку в большинстве задач основные вычислительные затраты ложатся на вычисление градиента.

 

Рисунок 1- иллюстрация градиентного метода с постоянным шагом.

 

4.    Метод сопряженных градиентов

Определение сопряженности формулируется следующим образом: два вектора x и y называют -сопряженными (или сопряженными по отношению к матрице ) или –ортогональными, если скалярное произведение x и y равно нулю, то есть:

xTy = 0 (12)

Сопряженность можно считать обобщением понятия ортогональности. Действительно, когда матрица – единичная матрица, в соответствии с равенством 3.1, векторы x и y – ортогональны. Остается выяснить, каким образом вычислять сопряженные направления.

Для нахождения сопряженных направлений существует несколько способов. Один из возможных способов – использовать методы линейной алгебры, в частности, процесс ортогонализации Грамма–Шмидта. Но для этого необходимо знать матрицу , поэтому для большинства задач этот метод не годится. Так же, существуют другие, итеративные способы вычисления сопряженного направления, самый известный – формула Флетчера-Ривса:

           (13)

где:

     (14)

r(i+1) – антиградиент в точке х(i+1)

                                                

            d(i) - новое сопряженное направление

Формула 13 означает, что новое сопряженное направление получается сложением антиградиента в точке поворота и предыдущего направления движения, умноженного на коэффициент, вычисленный по формуле 14. Направления, вычисленные по формуле 13, оказываются сопряженными. То есть для квадратичных функций метод сопряженных градиентов находит минимум за n шагов (n – размерность пространства поиска). Для функций общего вида алгоритм перестает быть конечным и становится итеративным. При этом, Флетчер и Ривс предлагают возобновлять алгоритмическую процедуру через каждые n + 1 шагов.

Можно привести еще одну формулу для определения сопряженного направления, формула Полака–Райбера:

(15)

Метод Флетчера-Ривса сходится, если начальная точка достаточно близка к требуемому минимуму, тогда как метод Полака-Райбера может в редких случаях бесконечно зацикливаться . Однако последний часто сходится быстрее первого метода.

Алгоритм:

  1. Вычисляется антиградиент в произвольной точке x(0).
  2. Осуществляется спуск в вычисленном направлении пока функция уменьшается, иными словами, поиск a(i), который минимизирует
    a(i) находим методом наискорейшего спуска.
  3. Переход в точку, найденную в предыдущем пункте
  4. Вычисление антиградиента в этой точке
  5. Вычисления по формуле (14) или (15). Чтобы забыть последнее направление поиска и стартовать алгоритм заново в направлении скорейшего спуска, для формулы Флетчера–Ривса присваивается 0 через каждые n + 1 шагов. Вычисление нового сопряженного направления
  6. Переход на пункт 2.

 

 

 

 

 

 

Приложение А.

Блок-схема PM алгоритма.

Тестовый пример

 

Рассмотрим исходную систему.

При использовании РМ - алгоритма, получим следующий результат (Рис. 1):

 

 

Рис. 1. Окно программы «РМ – алгоритм»

 

 

 

 

 

 

Текст программы

 

procedure PM(A:matr;var l:real;var z:vector);

var i,k:integer;

   y,x1,x2,l1,l2:vector;

   s:real;

begin

for i:=1 to n do

   begin

   y[i]:=1;

   l2[i]:=1;

   end;

k:=1;

x2:=ch_v(1/normaV(y),y);

repeat

x1:=x2;

l1:=l2;

y:=MV(A,x1);

x2:=ch_v(1/normaV(y),y);

l2:=otn(y,x1);

until usl(l1,l2);

s:=0;

for i:=1 to n do

   s:=s+l2[i];

l:=s/n;

z:=y;

end;

function otn(y,x:vector):vector;

var i:integer;

begin

for i:=1 to n do

   if x[i]<>0 then result[i]:=y[i]/x[i];

end;

 

function usl(y,x:vector):boolean;

const eps=0.001;

var i:integer;

begin

result:=true;

for i:=1 to n do

   if abs(y[i]-x[i])>eps then result:=false;

end;

end.

 

Приложение Б.

Блок-схема градиентного метода с дроблением шага

 

 

 

 

 

Тестовый пример:

Рассмотрим систему функций.

При использовании градиентного метода с дроблением шага, получим следующий результат (Рис. 2):

 

 

 

Рис. 2. Окно программы «градиентного метода с дроблением шага»

 

 

Текст программы:

 

             unit Unit1;

 

interface

 

uses

Windows, Messages, SysUtils, Variants, Classes, Graphics, Controls, Forms,

Dialogs, StdCtrls, Grids;

 

type

TForm1 = class(TForm)

   StringGrid1: TStringGrid;

   Button1: TButton;

   Edit1: TEdit;

   Label2: TLabel;

   Label3: TLabel;

   Label6: TLabel;

   Edit3: TEdit;

   Label7: TLabel;

   Button2: TButton;

   Label9: TLabel;

   procedure Button1Click(Sender: TObject);

   procedure Button2Click(Sender: TObject);

private

   { Private declarations }

public

   { Public declarations }

end;

mas=array[1..100]of real;

function f (x:mas):real;

function pr1 (u:real):real;

function pr2 (u:real):real;

procedure drob (e:real; x0:mas;var xmin:mas; n:integer);

var

Form1: TForm1;

 

implementation

 

{$R *.dfm}

 

function f (x:mas):real;

begin

f:=sqr(x[1]-2)+sqr(x[2]-3);

end;

 

function pr1 (u:real):real;

begin

pr1:=2*(u-2);

end;

 

function pr2 (u:real):real;

begin

pr2:=2*(u-3);

end;

 

procedure drob (e:real; x0:mas;var xmin:mas; n:integer);

var i:integer;

h,x1:mas;

   a,amin,d,norma:real;

begin

a:=0.95;

d:=0.1;

h[1]:=-pr1(x0[1]);

h[2]:=-pr2(x0[2]);

for i:=1 to n do x1[i]:=x0[i]+a*h[i];

while sqrt(sqr(pr1(x1[1]))+sqr(pr2(x1[2])))>e do begin

h[1]:=-pr1(x0[1]);

h[2]:=-pr2(x0[2]);

for i:=1 to n do x1[i]:=x0[i]+a*h[i];

norma:=sqrt(sqr(pr1(x0[1]))+sqr(pr2(x0[2])));

repeat

if (f(x0)-f(x1))>=0.5*a*norma then begin

                                 a:=d*a;

                                 for i:=1 to n do x1[i]:=x0[i]+a*h[i];

                                 end;

until (f(x0)-f(x1))<0.5*a*norma;

if sqrt(sqr(pr1(x1[1]))+sqr(pr2(x1[2])))<e then xmin:=x1 else x0:=x1;

end;  

end;

 

procedure TForm1.Button1Click(Sender: TObject);

var x0,xmin:mas;

   e:real;

   i,n:integer;

begin

e:=strtofloat(Edit3.Text);

n:=strtoint(Edit1.Text);

for i:=1 to n do x0[i]:=strtofloat(Stringgrid1.Cells[0,i-1]);

drob (e,x0,xmin,n);

for i:=1 to n do Stringgrid1.Cells[1,i-1]:=floattostr(xmin[i]);

end;

 

procedure TForm1.Button2Click(Sender: TObject);

begin

StringGrid1.RowCount:=Strtoint(Edit1.Text);

end;

 

end.

 

 

 

 

                        

 

 

Приложение В.

Блок-схема метода сопряжённых градиентов

 

 

 

 

 

 

 

 

Тестовый пример:

Реализуем программно метод сопряжённых градиентов.

При использовании метода, получим следующий результат (Рис. 3):

 

 

 

Рис. 3. Окно программы «метода сопряжённых градиентов»

 


 

Тестовый пример

 

unit Unit1;

 

interface

 

uses

Windows, Messages, SysUtils, Variants, Classes, Graphics, Controls, Forms,

Dialogs, StdCtrls;

 

type

TForm1 = class(TForm)

   Button1: TButton;

   Edit1: TEdit;

   Edit2: TEdit;

   Edit3: TEdit;

   Label1: TLabel;

   Button2: TButton;

   procedure Button1Click(Sender: TObject);

   procedure Button2Click(Sender: TObject);

private

   { Private declarations }

public

   { Public declarations }

end;

 

var

Form1: TForm1;

 

implementation

 

{$R *.dfm}

 

procedure TForm1.Button1Click(Sender: TObject);

Type

Mas=Array[1..3]of real;

 

 

 

 

Function NormaGr (x:mas):real;

begin

NormaGr:=sqrt(sqr(4*(x[1]-5))+sqr(2*(x[2]-2))+sqr(2*(x[3]-3)));

end;

 

 

Function F (x:mas):real;

begin

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

end;

 

 

procedure proiz (x:mas; var gr:mas);

begin

gr[1]:=4*(x[1]-5);

gr[2]:=2*(x[2]-2);

gr[3]:=2*(x[3]-3);

end;

 

const n=3;

 

var

x,k,p,gr,w:Mas;

EPS,h,b:real;

i:integer;

 

begin

x[1]:=1;

x[2]:=3;

x[3]:=12;

h:=0.5;

eps:=0.001;

proiz(x,gr);

for i:=1 to n do

p[i]:=-gr[i];

repeat

   for i:=1 to n do

       k[i]:=x[i]-p[i]*h;

{proiz(x,gr);

h:=h+normagrk(grad)/normagrx(gr)*h;}

 

   w:=x;

if F(k)<F(x) then

     begin

         while F(k)<F(x) do

         begin

               for i:=1 to n do

               x[i]:=k[i];

             // proiz(x,gr);

               for i:=1 to n do

                       {proiz(x,gr);}

                       k[i]:=x[i]+p[i]*h;

           end;

               {h1:=proiz(x,gr)+normagr(x,gr,grad)*h0;}

     end

else

     begin

         while (F(k)>F(x)) and (h>eps) do

         h:=h/2;

             begin

                   for i:=1 to n do

                       begin

                         k[i]:=x[i]+p[i]*h;

                       end;

{ h1:=proiz(x,gr)+normagr(x,gr,grad)*h0;}

             end;

 

     b:=sqr(normagr(k)/normagr(w));

     x:=k;

     proiz(k,gr);

     for i:=1 to n do

     p[i]:=-gr[i]+b*p[i];

     end;

   until {(sqrt(sqr(gr[1])+sqr(gr[2])+sqr(gr[3])))}normagr(k)<eps;

   Edit1.Text:=floattostr(x[1]);

   Edit2.Text:=floattostr(x[2]);

   Edit3.Text:=floattostr(x[3]);

 

end;

 

procedure TForm1.Button2Click(Sender: TObject);

begin

Form1.Close;

end;

 

end.

 

 

 

 

 

 

 

Список используемых источников

 

  • Вержбицкий В.М. //Численные методы. Численные методы: Учеб. пособие для вузов. — М.: Высш. шк., 2001. — 382 с.:ил. ISBN 5-06-003982-Х
  • Бахвалов Н. С., Жидков Н. П., Кобелько Г.М. //Численные методы.-М.: Наука, 1987г.
  • Васильев Ф.П.//Численные методы решения экстремальных задач.-М.: Наука, 1990г.
  • Волков Б. А.// Численные методы.-М.: Наука, 1990г.
  • Пантелеев А.В., Летова Т.А.// Методы оптимизации в примерах и задачах: Учеб. пособие. -2-е изд.,исправл.-М.:Высш.шк., 2005. -544с.:ил. ISBN 5-06-004137-9
  • Пшеничный Б.Н., Ланилин Ю.М.//Численные методы в экстремальных задачах. –М.: Наука, 1975.

 

Скачать: Kursach.doc

Категория: Курсовые / Курсовые по математике

Уважаемый посетитель, Вы зашли на сайт как незарегистрированный пользователь.
Мы рекомендуем Вам зарегистрироваться либо войти на сайт под своим именем.