Методы решения систем нелинейных уравнений

СОДЕРЖАНИЕ

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 точностью :

Ответ:

нелинейное уравнение решение ньютон

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с.