Методы решения систем нелинейных уравнений
СОДЕРЖАНИЕ
1. Методы решения систем
нелинейных уравнений
. Векторная запись нелинейных
систем
. Метод Ньютона, его
реализации и модификации
.1 Метод Ньютона
.2 Модифицированный метод
Ньютона
.3 Метод Ньютона с
последовательной аппроксимацией матриц
.4 Разностный метод Ньютона
.5 Обобщение полюсного метода
Ньютона на многомерный случай
. Численный пример
. Реализация метода Ньютона в
среде MATLAB
Список используемой
литературы
1. Методы
решения систем нелинейных уравнений
Рассматриваем метод Ньютона в разных модификациях (в частности, n -полюсный метод Ньютона) решения
систем алгебраических и трансцендентных уравнений. Показывается связь между
данной задачей и за дачей безусловной минимизации функции нескольких
переменных. С единых позиций изучается сходимость основного и упрощенного
методов Ньютона и метода, получаемого из метода Ньютона применением
итерационного процесса Шульца для приближенного обращения матриц Якоби.
2. Векторная
запись нелинейных систем.
Пусть требуется решить систему уравнений
(2.1)
где- заданные, вообще говоря, нелинейные (среди них могут
быть и линейные) вещественнозначные функции п вещественных переменных Обозначив
, ,
данную
систему (2.1) можно записать одним уравнением
(2.1а)
относительно
векторной функции F векторного аргумента х. Таким образом, исходную
задачу можно рассматривать как задачу о нулях нелинейного отображения . В этой постановке она является прямым обобщением
основной задачи - построения методов нахождения нулей одномерных нелинейных
отображений. Фактически это та же задача, только в пространствах большей
размерности. В любом случае следует позаботиться о правомочности тех или иных
операций над векторными переменными и векторными функциями, а также о
сходимости получаемых таким способом итерационных процессов. Часто теоремы
сходимости для этих процессов являются тривиальными обобщениями соответствующих
результатов, полученных для методов решения скалярных уравнений. Однако не все
результаты и не все методы можно перенести со случая п = 1 на случай п ≥2.
Например, здесь уже не будут работать методы дихотомии, поскольку множество
векторов не упорядочено. В то же время, переход от n = 1 до n≥2
вносит в задачу нахождения нулей нелинейного отображения свою специфику, учет
которой приводит к новым методам и к различным модификациям уже имеющихся. В
частности, большая вариативность методов решения нелинейных систем связана с
разнообразием способов, которыми можно решать линейные алгебраические задачи,
возникающие при пошаговой линеаризации данной нелинейной вектор-функции F(x).
3. Метод
Ньютона, его реализации и модификации
.1 Метод
Ньютона
Пусть
() - некоторая последовательность невырожденных
вещественных n x n-матриц. Тогда, очевидно, последовательность задач
, k =
0,1,2,...
имеет
те же решения, что и исходное уравнение (2.1а), и для приближенного нахождения
этих решений можно формально записать итерационный процесс
, k =
0,1,2,... (3.1.1)
при
. В случае - это
действительно МПИ с линейной сходимостью последовательности () Если же различны
при разных k, то формула (3.1.1) определяет большое семейство итерационных
методов с матричными параметрами.
Рассмотрим некоторые из методов этого семейства.
Положим
, где
матрица
Якоби вектор-функции F(x). Подставив это в
(3.1.1), получаем явную формулу метода Ньютона
, (3.1.2)
обобщающего
на многомерный случай скалярный метод Ньютона (5.14). Эту формулу, требующую
обращения матриц на каждой итерации, можно переписать в неявном виде:
. (3.1.3)
Применение
(3.1.3) предполагает при каждом k = 0,1,2,... решение линейной алгебраической
системы
относительно
векторной поправки , а затем прибавление этой поправки к текущему
приближению для получения следующего:
.
К
решению таких линейных систем можно привлекать самые разные методы как прямые,
так и итерационные в зависимости от размерности n решаемой
задачи и специфики матриц Якоби (например,
можно учитывать их симметричность, разреженность и т.п.).
Сравнивая
(3.1.3) с формальным разложением F(x) в ряд Тейлора
,
видим,
что последовательность () в методе Ньютона получается в результате подмены при
каждом k=0,1,2,... нелинейного уравнения F(x) = 0 или, что то же (при
достаточной гладкости F(x)), уравнения
линейным
уравнением
т.
е. с пошаговой линеаризацией. Как следствие этого факта, можно рассчитывать,
что при достаточной гладкости F(x) и достаточно хорошем начальном приближении сходимость порождаемой методом Ньютона
последовательности () к решению будет
квадратичной и в многомерном случае. Имеется ряд теорем, устанавливающих это
при тех или иных предположениях. В частности, одна из таких теорем приводится
ниже.
Новым,
по сравнению со скалярным случаем, фактором, осложняющим применение метода
Ньютона к решению n-мерных систем, является необходимость решения n-мерных
линейных задач на каждой итерации (обращения матриц в (3.1.2) или решения СЛАУ
в (3.1.3)), вычислительные затраты на которые растут с ростом n,
вообще говоря, непропорционально быстро. Уменьшение таких затрат - одно из
направлений модификации метода Ньютона.
3.2
Модифицированный метод Ньютона
Если
матрицу Якоби F'(х) вычислить и обратить лишь один раз - в начальной точке , то от метода Ньютона (3.1.2) придем к
модифицированному методу Ньютона
(3.2.1)
Этот
метод требует значительно меньших вычислительных затрат на один итерационный
шаг, но итераций при этом может потребоваться значительно больше для достижения
заданной точности по сравнению с основным методом Ньютона (3.1.2), поскольку,
являясь частным случаем МПИ (), он
имеет лишь скорость сходимости геометрической прогрессии.
Компромиссный
вариант - это вычисление и обращение матриц Якоби не на каждом итерационном
шаге, а через несколько шагов (иногда такие методы называют рекурсивными).
Например,
простое чередование основного (3.1.2) и модифицированного (3.2.1) методов
Ньютона приводит к итерационной формуле
(3.2.2)
где
k = 0,1,2,… За здесь
принимается результат последовательного применения одного шага основного, а
затем одного шага модифицированного метода, т.е. двухступенчатого процесса
(3.2.3)
Доказано,
что такой процесс при определенных условиях порождает кубически сходящуюся
последовательность ().
3.3 Метод
Ньютона с последовательной аппроксимацией матриц
Задачу обращения матриц Якоби на каждом k-м шаге метода Ньютона (3.1.2)
можно попытаться решать не точно, а приближенно. Для этого можно применить,
например, итерационный процесс Шульца, ограничиваясь минимумом - всего одним
шагом процесса второго порядка, в котором за начальную матрицу принимается
матрица, полученная в результате предыдущего (k-1)-го шага. Таким образом приходим к методу Ньютона с
последовательной аппроксимацией обратных матриц:
(3.3.1)
где
k = 0,1,2,…, а и -
начальные вектор и матрица (). Этот
метод (будем называть его более коротко ААМН - аппроксимаиионный аналог метода
Ньютона) имеет простую схему вычислений - поочередное выполнение векторных в
первой строке и матричных во второй строке его записи (3.3.1) операций.
Скорость его сходимости почти так же высока, как и у метода Ньютона.
Последовательность () может квадратично сходиться к решению уравнения F(x)=0 (при этом матричная
последовательность () также квадратично сходится к , т.е. в нормально развивающемся итерационном процессе
(3.3.1) должна наблюдаться достаточно быстрая сходимость () к нулю).
Применение
той же последовательной аппроксимации обратных матриц к простейшему
рекурсивному методу Ньютона (3.2.2) или, что то же, к двухступенчатому процессу
(3.2.3) определяет его аппроксимацирнный аналог
(3.3.2)
как
и (3.2.2), также можно отнести к- методам третьего порядка. Доказательство
кубической сходимости этого метода требует уже более жестких ограничений на
свойства F(х) и близость к , к , чем в
предыдущем методе. Заметим, что к улучшению сходимости здесь может привести
повышение порядка аппроксимации обратных матриц, например, за счет добавления
еще одного слагаемого в формуле для подсчета :
.4 Разностный
метод Ньютона
На
базе метода Ньютона (3.1.2) можно построить близкий к нему по поведению
итерационный процесс, не требующий вычисления производных. Сделаем это, заменив
частные производные в матрице Якоби J(x) разностными отношениями, т.е.
подставив в формулу (3.1.1) вместо матрицу где
При
удачном задании последовательности малых векторов (постоянной или сходящейся к нулю) полученный таким
путем разностный (или иначе, дискретный) метод Ньютона имеет сверхлинейную,
вплоть до квадратичной, скорость сходимости и обобщает на многомерный случай
метод. При задании векторного параметра h - шага дискретизации - следует
учитывать точность машинных вычислений (macheps), точность вычисления значений
функций , средние значения получаемых приближений.
3.5 Обобщение
полюсного метода Ньютона на многомерный случай
Переложим
вывод одномерного полюсного метода Ньютона на векторную основу. Касательную к
кривой в точке () зададим
условием ортогональности текущего вектора этой
прямой и ее нормального вектора, в качестве которого можно взять вектор . Уравнение прямой, проходящей через полюс и связанную с уже известным приближением точку (),
получим из условия коллинеарности текущего вектора этой прямой и ее направляющего вектора . Таким образом, точку пересечения двух прямых,
проекцию которой на ось абсцисс считаем новым приближением , находим из совокупности условий
(3.5.1)
Первое
из этих условий означает равенство нулю скалярного произведения (n,u),
второе - пропорциональность соответствующих координат векторов v и l
или, иначе, равенство нулю составленного из них определителя. Следовательно,
искомое приближение есть первая компонента вектора, служащего решением
линейной системы
(3.5.2)
(вторая
компонента - ордината точки пересечения указанных прямых - после вычисления
значения может дать информацию об отклонении от функции в точке ее
локальной аффинной модели, каковой является проведенная в точке касательная). Ясно, что получаемое из системы (3.5.2)
значение тождественно его.
Рассмотренный
векторный подход к построению одномерного полюсного метода Ньютона служит
ключом для его распространения на двумерный случай на основе таких же
геометрических, но уже пространственных соображений.
Пусть
требуется найти приближенное решение двумерной нелинейной системы в
предположении непрерывной дифференцируемости входящих в нее функций f(x, у)
и g(x, у) в некоторой области G, содержащей
искомое решение х* =(х*; у*) и приближения к нему k = 0,1,2,....
z = f(x,y) и z = g(x,y).
(3.5.3)
Эти
плоскости задаются текущими векторами
и
нормалями
соответственно,
т.е. аналогично первому из условий (3.5.1) должно быть
иначе,
. (3.5.4)
Пересечение
двух касательных плоскостей, т.е. образ, определяемый уравнениями (3.5.4), есть
прямая в трехмерном пространстве, общая точка которой с координатной плоскостью
Оху является ньютоновским приближением к
решению х* сиcтемы. Наша цель - построить третью плоскость,
пересечение которой с упомянутой прямой (линией пересечения касательных
плоскостей) давало бы точку в пространстве R3 такую,
проекция которой на плоскость Оху могла бы оказаться ближе к х*, чем .
Чтобы
осуществить поставленную цель, зафиксируем в R3 две
несовпадающие между собой и с точки -
полюсы и . Через
указанные три точки можно провести единственную плоскость (которая здесь
играет роль прямой, проходящей через полюс и точку (хк; 0) в одномерной
ситуации). Взяв текущую точку М(х; у; z) и образовав текущий вектор этой третьей плоскости, можно задать ее условием
компланарности трех векторов- и (что служит аналогом второго из условий (3.5.1)).
Запишем
совокупность всех трех описанных средствами векторной алгебры плоскостей в
координатной форме. Имеем:
Первые
две координаты вектора (x;y;z), служащего решением полученной системы уравнений,
считаем искомым приближением ().Введя
поправки
, (3.5.5)
эту
систему превращаем в систему уравнений относительно неизвестных и z:
(3.5.6)
Для
исключения вспомогательной переменной z из линейной системы (3.5.6)
выразим ее из третьего уравнения. Обозначив
, , (3.5.7)
раскрываем
фигурирующий в (3.5.6) определитель по элементам первой строки:
Отсюда
находим выражение
(3.5.8)
подставляя
которое в первые два уравнения системы (3.5.6), приходим к двумерной линейной
системе
(3.5.9)
Фактически
эта система вместе с обозначениями (3.5.7) и определяет двумерный полюсный
метод Ньютона для нелинейной системы. Надя их нее поправки , в соответствии с равенствами (3.5.5) получаем
очередное приближение :
, .
Дельнейшее
преобразование полюсного метода Ньютона, т. е. переход от размерности 2 к
произвольной размерности, совершаем формально на основе предыдущего построения.
Пусть
задана нелинейная система, функции (образующими
вектор ) в точке , можно
описать матрично-векторным уравнением
,
(3.5.10)
где
- n-мерный вектор, каждой компонентой которого служит
вспомогательная переменная ,
входящая в уравнения гиперповерхностей .
Зададим
n полюсов (i=1,2,…,n)
так, чтобы они не принадлежали одной гиперплоскости пространства . Через все эти полюсы и точку
(), определяемую известным приближением к решению системы, проводим гиперплоскость, уравнение
которой аналогично двумерному случаю задаем условием равенство нулю определителя
(n+1) порядка:
.
(3.5.11)
Векторно-матричное
уравнение (3.5.10) и скалярное уравнение (3.5.11), в принципе, уже определяют
векторный n-полюсный метод Ньютона для построения приближенной к
решению системы. Чтобы записать соответствующую линейную систему относительно
поправок
(3.5.12)
(аналогичную
схеме (3.5.9) двумерного случая), введем следующие обозначения. Положим
, ,
и
образуем квадратную (n+1)-мерную матрицу следующей структуры:
.
Тогда
на основе (3.5.10), (3.5.11) имеем (n+1)-мерную систему уравнений
относительно неизвестных :
(3.5.13)
Как
и в двумерном случае, из второго уравнения этой системы выражаем
вспомогательную неизвестную :
(3.5.14)
где
, а есть
алгебраические дополнения к элементам первой
строки матрицы (что через соответствующие миноры этой матрицы можно представить так: , ).
Заменив
в (3.5.13) все компоненты вектора z найденным их значением
(3.5.14), приходим к следующему линейному векторно матричному уравнению
относительно вектора-поправки :
,
(3.5.15)
Где
.
(3.5.16)
Уравнение
(3.5.15) вместе со связью (3.5.12), согласно которой
,
(3.5.17)
является
неявной формой п -полюсного метода Ньютона для уравнения (2.1а).
Совокупности
формул (3.5.15)-(3.5.17) можно придать другой вид:
,
(3.5.18)
который
удобно трактовать как явный метод Ньютона со своеобразной коррекцией матриц
Якоби путем прибавления к ним формирующихся по заданному правилу матриц . Как и в одномерном случае, для ускорения сходимости
последовательности приближений полюсы целесообразно изменять в такт с изменением значений
функций, и в самом простом случае есть смысл фиксировать матрицу С, а вектор брать равным или -
4. Численный
пример
Начальное
приближение:
Вектор-функция:
Матрица
Якоби вектор-функции:
Вычисляем
корень по формуле метода Ньютона c точностью :
k
|
|
|
|
|
|
|
0
|
0 -1
|
-0.841 0
|
-1.06 0.54 0 -2
|
-0.944 -0.255 0 -0.5
|
-0.794 -1
|
1
|
-0.794 -1
|
0.295 0.63
|
-1.821 -0.221 -1.588 -2
|
-0.608 0.067 0.482 -0.553
|
-0.657 -0.794
|
0.247>
|
2
|
-0.657 -0.794
|
0.058 0.062
|
-1.48 0.12 -1.314 -1.588
|
-0.633 -0.048 0.524 -0.59
|
-0.617 -0.788
|
0.040>
|
3
|
-0.617 -0.788
|
-0.0000597
0.011
|
-1.441 0.159 -1.234 -1.588
|
-0.639 -0.064 0.497 -0.58
|
-0.616 -0.788
|
0.001=
|
4
|
-0.616 -0.788
|
0.000522 0.0004
|
-1.434 0.166 -1.232 -1.576
|
-0.639 -0.067 0.5 -0.582
|
-0.616 -0.788
|
0<
|
Ответ:
нелинейное уравнение
решение ньютон
5. Реализация
метода Ньютона в среде MATLAB
М-функция.
function nwt = newton(x,y,e,F0,F1,dF0x,dF0y,dF1x,dF1y)
for i = 1:1000000=[F0(x,y); F1(x,y)];=[dF0x(x,y) dF0y(x,y);
dF1x(x,y) dF1y(x,y)];= [x;y] - dF^(-1)*F;= sqrt((x-Z(1))^2+(y-Z(2))^2);b >
e= Z(1);= Z(2);break;('Количество итераций'); disp(i);
nwt = Z;
Функция
в качестве входных параметров принимает начальное приближение (), функцию (),
частные производные функции и точность e.
Пятая
строка находит новую точку приближения. Шестая строка вычисляет норму между текущим и следующим приближением. Строки восемь
и девять запоминают точку начального приближения.
Процесс
нужно продолжать до тех пор пока . Если , процесс завершить и получим решение .
Теперь
напишем скрипт который покажет работу нашей М-функции.
Найдем решение заданной системы нелинейных уравнений
при
начальном приближении x=0, y=-1, с точностью до 0.001:
Скрипт.
% Решить систему уравненийу методом Ньютона
% sin(x+y)-1.6x
% x^2+y^2-1
% Введём функцию F(x)(систему функций)= inline('sin(x+y)-1.6*x');
F1 = inline('x^2+y^2-1');
disp(F0); disp(F1);
% Их производны
dF0x = inline('cos(x+y)-1.6');y = inline('cos(x+y)');x =
inline('2*x+0*y');y = inline('0*x+2*y');(dF0x); disp(dF0y);(dF1x); disp(dF1y);
% Начальное приближение=0; y=-1; e=0.000000001;
% Значение функций= newton(x,y,e,F0,F1,dF0x,dF0y,dF1x,dF1y);('Решение
системы');(root);
После выполнения программы получим следующее:
>>function:(x,y) = sin(x+y)-1.6*xfunction:(x,y) =
x^2+y^2-1function:x(x,y) = cos(x+y)-1.6function:y(x,y) =
cos(x+y)function:x(x,y) = 2*x+0*yfunction:y(x,y) = 0*x+2*y
Количество итераций
Решение системы
.6163
.7875
Полученное решение совпадает с рассчитанным.
Список используемой
литературы
1. Н. Бахвалов, Н. Жидков, Г. Кобельков. Численные методы.
М., 2002, 632 с.
2. Н. Калиткин. Численные методы. М., 1972,
. А. Самарский. Введение в численные методы. М., ,
270с.
. М. Лапчик, М. Рагулина, Е. Хеннер. Численные методы.
М., 2004, 384с.
. В. Потемкин. Система MATLAB. Справочное пособие. М.,
1997, 350с.
. Е. Алексеев, О. Чеснокова. MATLAB 7. М., 2006, 464с.