85666 (612542), страница 3

Просмтор этого файла доступен только зарегистрированным пользователям. Но у нас супер быстрая регистрация: достаточно только электронной почты!

Текст из файла (страница 3)

где значение параметра отвечает симметричному варианту Пауелла методу Бройдена, а — методу Давидона - Флечера - Пауелла.

Во втором из рассматриваемых классов квазиньютоновских методов матрица Н_к аппроксимирует матрицу . Здесь матрица Н_о должна быть симметричной, а вместо (3.29) используется формула

Где соответствует методу Бройдена-Флечера-Голвдфарба-Шенно, что является одним из наилучших (с вычислительной точки зрения), который учитывает симметричность матрицы Якоби.

Описанные выше квазиньютоновские методы сходятся лишь при достаточно хорошом начальном приближении х⁽⁰⁾. Для расширения области их сходимости можно использовать прием, который имеет название одномерного поиска.

Пусть имеем квазиньютоновское направление (или ). Используем длину шага = 1 и проверим неравенство

(3.34)

где - евклидовая норма. Если оно выполняется, то заканчиваем одномерный поиск и считаем

(3.35)

т.е. уменьшаем длину шага (устанавливая, например, ), пока не выполнится (3.34). На этом заканчиваем одномерный поиск и переходим к формуле (3.35).

Как видим, одномерный поиск (в случае успеха) обеспечивает монотонное уменьшение нормы отклонения с ростом к. Если квазиньютоновское направление сильно отличается от ньютоновского, то одномерный поиск может оказаться неудачным, и тогда необходимо возобновить матрицу В_к, (или )., приравняв ее, например, конечно-разносной аппроксимации матрицы Якоби (или ). Критерием окончания итераций для квазиньютоновских методов есть неровность

3. Другие итерационные методы решения систем нелинейных уравнений

3.1 Метод Пикара

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

Так, если в уравнениях системы можно выделить линейную l(X) и нелинейную g(X)части функций f_i(X) = l_i(X) + g_i(x), то удобней применить к ней метод Пикара.

В таком случае систему уравнений можно записать в виде

l_i(X) = - g_i(X), i=1,2,3...n;

или в векторной форме A X= - G(X);

где A- матрица коэффициентов линейных частей уравнений;

G(X) - вектор-функция нелинейных частей уравнений.

Выберем некоторый начальный вектор X⁽⁰⁾ и построим итерационный процесс в виде

A X^(k+1)=-G(X^(k)).

Для выполнения одной итерации таким методом необходимо решать систему линейных уравнений, у которой вектором свободных членов будут нелинейные части функций f_i(X). Причем поскольку матрица A остается неизменной при всех итерациях, то для решения СЛАУ можно использовать специальные алгоритмы, предусматривающие возможность преобразования только столбца свободных членов.

3.2 Метод релаксаций

Перепишем систему в виде

X=X+ F(X),

где  - некоторая константа, и построим итерационный процесс по схеме

X^(k+1) = X^(k) +  F(X^(k)).

Параметр  должен быть таким, чтобы в окрестности решения выполнялось достаточное условие сходимости

||Е+  W|| < 1,

где E- единичная матрица.

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

||X^(k)-X^(k-1)||<||X^(k-1)-X^(k-2)||.

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

3.3 Метод градиентного спуска

Пусть имеем систему уравнений (А)

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

(В)

Очевидно, что каждое решение системы (А) превращает в ноль функцию U(x); наоборот, числа , для которых функция U(x) равняется нулю, является корнем системы (А).

Предположим, что система имеет лишь изолированное решение, которое представляет собой точку строго минимума функции U(x) в n-мерном пространстве .

Пусть x - вектор системы (А) и x0 - его нулевое приближение. Через точку x0 проведем поверхность уровня функции U(x). Если точка x0 довольно близка к корню x, то при наших предположениях поверхность уровня

U(x)= U(x0)

будет похожа на эллипсоид.

Из точки x0 движемся по нормали к поверхности U(x)= U(x0) до тех пор, пока эта нормаль не затронет в некоторой точке x1 какой-то другой поверхности уровня U(x)= U(x1).

Потом, отправляясь от точки x1, снова движемся по нормали к поверхности уровня U(x)= U(x1) до тех пор, пока эта нормаль не затронет в некоторой точке x2 новой поверхности уровня U(x)= U(x2), и т.д.

Поскольку U(x0)>U(x1)>U(x2)>..., то, двигаясь таким путем, мы быстро приближаемся к точке с наименьшим значением U ("дно ямы"), что отвечает искомому решению исходной системы. Обозначим через

градиент функции U(x).

Находить нужное решение будем по формуле:

Остается определить множители . Для этого рассмотрим скалярную функцию

Функция () дает изменение уровня функции U вдоль соответствующей нормали к поверхности уровня в точке x_p. Множитель надо выбрать таким образом, чтобы () имела минимум. Беря производную по и приравнивая ее нулю, получаем уравнение

.

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

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

Раскладывая функции за степенями с точностью до линейных членов, получим:

,

где .

Отсюда

Итак,

, где

- матрица Якоби вектор- функции f.

Дальше, имеем:

.

Отсюда

,

где W'(x) - транспонированная матрица Якоби.

Поэтому окончательно

,

причем

.

3. Программная реализация итерационных методов

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

3.1 Метод простых итераций

Приведём систему к виду:

Проверим условие сходимости метода простых итераций.

Для этого построим матрицу Якоби

> f1:=0.1-x0^2+2*y0*z0;

f2:=-0.2+y0^2-3*x0*z0;

f3:=0.3-z0^2-2*x0*y0;

> f1x:=diff(f1,x0):

> f1y:=diff(f1,y0):

> f1z:=diff(f1,z0):

> f2x:=diff(f2,x0):

> f2y:=diff(f2,y0):

> f2z:=diff(f2,z0):

> f3x:=diff(f3,x0):

> f3y:=diff(f3,y0):

> f3z:=diff(f3,z0):

> A:=<,,>;

И найдём ей обратную, норму обратной матрицы сначала в общем виде:

> A1:=MatrixInverse(A);

> norma:=MatrixNorm(A1,1);

Найдём значения при которых норма обратной матрицы Якоби меньше единицы.

> x0:=1; y0:=1; z0:=1;

> norma;

Это означает, что по формулам

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

Построим итерационную последовательность

> restart;

> with(LinearAlgebra):

> x0:=0:

y0:=0:

z0:=0:

> x:=0.1-x0^2+2*y0*z0;

y:=-0.2+y0^2-3*x0*z0;

z:=0.3-z0^2-2*x0*y0;

i:=1;

> while (abs(x-x0)>0.0001)and(abs(y-y0)>0.0001)and(abs(z-z0)>0.0001) do

x0:=x:

y0:=y:

z0:=z:

x:=0.1-x0^2+2*y0*z0;

y:=-0.2+y0^2-3*x0*z0;

z:=0.3-z0^2-2*x0*y0;

i:=i+1;

end do:

Получили ответ:

Количество итераций:

Погрешность решения:

Отсюда можно получить коэффициент сжатия последовательности:

При

> P:= 0.3*q^22/(1-q)-0.0001;

> q:= fsolve(P);

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

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

procedure TForm1.Button3Click(Sender: TObject);

var i:integer;

x0,y0,z0,x,y,z,eps: real;

begin

x0:=StrToFloat(Edit1.Text);

y0:=StrToFloat(Edit2.text);

z0:=StrToFloat(Edit3.Text);

eps:=StrToFloat(Edit20.Text);

i:=1;

x:=0.1-x0*x0+2*y0*z0;

y:=-0.2+y0*y0-3*x0*z0;

z:=0.3-z0*z0-2*x0*y0;

repeat

i:=i+1;

x0:=x;

y0:=y;

z0:=z;

x:=0.1-x0*x0+2*y*z;

y:=-0.2+y0*y0-3*x0*z0;

z:=0.3-z0*z0-2*x0*y0;

until ((abs(x-x0)

Edit8.Text:=FloatToStr(x);

Edit9.Text:=FloatToStr(y);

Edit10.Text:=FloatToStr(z);

Edit11.Text:=IntToStr(i);

end;

Преобразование Эйткена на примере метода простых итереций:

> restart;

> x0:=0:

y0:=0:

z0:=0:

> f1:=0.1-x0^2+2*y0*z0;

f2:=-0.2+y0^2-3*x0*z0;

f3:=0.3-z0^2-2*x0*y0;

ff1:=0.1-f1^2+2*f2*f3;

ff2:=-0.2+f2^2-3*f1*f3;

ff3:=0.3-f3^2-2*f1*f2;

x:=(x0*ff1-f1^2)/(ff1-2*f1+x0);

y:=(y0*ff2-f2^2)/(ff2-2*f2+y0);

z:=(z0*ff3-f3^2)/(ff3-2*f3+z0);

i:=1;

while (abs(x-x0)>0.0001)do

x0:=x:

y0:=y:

z0:=z:

f1:=0.1-x0^2+2*y0*z0;

f2:=-0.2+y0^2-3*x0*z0;

f3:=0.3-z0^2-2*x0*y0;

ff1:=0.1-f1^2+2*f2*f3;

ff2:=-0.2+f2^2-3*f1*f3;

ff3:=0.3-f3^2-2*f1*f2;

x:=(x0*ff1-f1^2)/(ff1-2*f1+x0);

y:=(y0*ff2-f2^2)/(ff2-2*f2+y0);

z:=(z0*ff3-f3^2)/(ff3-2*f3+z0):

i:=i+1;

end do:

Получили ответ:

Количество итераций:

3.2 Метод градиентного спуска

Построим функцию:

> U:=(0.1-x^2+2*y*z-x)^2+(-0.2+y^2-3*x*z-y)^2+(0.3-z^2-2*x*y-z)^2;

Найдём градиент функции:

> Ux:= diff(U,x);

Uy:= diff(U,y);

Uz:= diff(U,z);

Выберем начальное приближение и построим итерационную последовательность:

> x:=0;

y:=0;

z:=0;

> N1:=2*(.1-x^2+2*y*z-x)*(-2*x-1)-6*(-.2+y^2-3*x*z-y)*z-4*(.3-z^2-2*x*y-z)*y;

> N2:=4*(.1-x^2+2*y*z-x)*z+2*(-.2+y^2-3*x*z-y)*(2*y-1)-4*(.3-z^2-2*x*y-z)*x;

> N3:=4*(.1-x^2+2*y*z-x)*y-6*(-.2+y^2-3*x*z-y)*x+2*(.3-z^2-2*x*y-z)*(-2*z-1);

> x:=x-lambda*N1;

y:=y-lambda*N2;

z:=z-lambda*N3;

i:=1;

> N1:=2*(.1-x^2+2*y*z-x)*(-2*x-1)-6*(-.2+y^2-3*x*z-y)*z-4*(.3-z^2-2*x*y-z)*y;

N2:=4*(.1-x^2+2*y*z-x)*z+2*(-.2+y^2-3*x*z-y)*(2*y-1)-4*(.3-z^2-2*x*y-z)*x;

N3:=4*(.1-x^2+2*y*z-x)*y-6*(-.2+y^2-3*x*z-y)*x+2*(.3-z^2-2*x*y-z)*(-2*z-1);

x:=x-lambda*N1;

y:=y-lambda*N2;

z:=z-lambda*N3;

> while (abs(N3)>0.0001) do

N1:=2*(.1-x^2+2*y*z-x)*(-2*x-1)-6*(-.2+y^2-3*x*z-y)*z-4*(.3-z^2-2*x*y-z)*y:

N2:=4*(.1-x^2+2*y*z-x)*z+2*(-.2+y^2-3*x*z-y)*(2*y-1)-4*(.3-z^2-2*x*y-z)*x:

N3:=4*(.1-x^2+2*y*z-x)*y-6*(-.2+y^2-3*x*z-y)*x+2*(.3-z^2-2*x*y-z)*(-2*z-1):

x:=x-lambda*N1:

y:=y-lambda*N2:

z:=z-lambda*N3:

i:=i+1:

end do:

Получили ответ:

Количество итераций и данным шагом :

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

procedure TForm1.Button1Click(Sender: TObject);

const

lambda=-0.0001;

n=3;

Характеристики

Список файлов курсовой работы