Пусть уравнение (2.1) можно заменить эквивалентным ему уравнением

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

x = (x). (2.4)

Например, уравнение - 0.5 = 0 можно заменить эквивалентным ему уравнением x = 0.5sinx.

Выберем каким-либо образом начальное приближение x₀. Вычислим значение функции (x) при x = x₀ и найдем уточненное значение x₁ = (x₀). Подставим теперь x₁в уравнение (2.4) и получим новое приближение x₂ = (x₁) и т. д. Продолжая этот процесс неограниченно, получим последовательность приближений к корню:

x_n+1 = (x_n). (2.5)

Формула (2.5) является расчетной формулой метода простых итераций.

Если последовательность {x_n} сходится при n, т. е. существует

x^* = x_n , (2.6)

и функция (x) непрерывна, то, переходя к пределу в (2.5) и учитывая (2.6), получим:

x^* = x_n = (x _{n -}₁) = (x_{n -}₁) = (x^*).

Таким образом, x^* = (x^*), следовательно, x^* - корень уравнения (2.4).

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

Теорема 2.2. Если в интервале, содержащем корень x^*уравнения (2.4), а также его последовательные приближения x₀, x₁, …, x_n, …, вычисляемые по формуле (2.5), выполнено условие:

|(x)| q < 1, (2.7)

то x^* = x_n.

т. е. итерационный процесс сходится и справедлива следующая оценка погрешности:

|x_n - x^*| qⁿ|x₀ - x^*| (2.8)

Оценка (2.8) является априорной. Она показывает, что метод простой итерации сходится со скоростью геометрической прогрессии с знаменателем q. Чем меньше q, тем выше скорость сходимости.

Как следует из теоремы 2.2, условие (2.7) является достаточным для сходимости метода простых итераций. Его выполнение гарантирует сходимость процесса (2.5), но невыполнение условия (2.7), вообще говоря, не означает, что итерационный процесс будет расходиться.

На рис. 2.3 - 2.6 показаны четыре случая взаимного расположения линий y = x и y = (x) и соответствующие итерационные процессы.

Рис. 2.3 и 2.4 соответствуют случаю |(x)| < 1, и итерационный процесс сходится. При этом, если (x) > 0 (рис. 2.3), сходимость носит односторонний характер, а если (x)< 0 (рис. 2.4), сходимость носит двусторонний, колебательный характер. Рис. 2.5 и 2.6 соответствуют случаю |(x)| > 1 - итерационный процесс расходится. При этом может быть односторонняя (рис. 2.5) и двусторонняя (рис 2.6) расходимость.

Рис. 2.3 Рис. 2.4 Рис. 2.5

Рис. 2.6

Погрешность метода. Если известна величина q в условии (2.7), то применима следующая апостериорная оценка погрешности:

|x_n - x^*| |x_n - x_{n -}₁|, n > 1. (2.9)

Критерий окончания. Из оценки (2.9) вытекает следующий критерий окончания итерационного процесса. Вычисления следует продолжать до выполнения неравенства

|x_n - x_{n -}₁| < .

Если это условие выполнено, то можно считать, что x_n является приближением к x^* с точностью .

Если q 0.5, то можно пользоваться более простым критерием окончания:

|x_n - x_{n -}₁| < . (2.10)

Пример 2.2.

Используем метод простой итерации для решения уравнения f(x) = sin x - x² = 0с точностью = 0.001.

Преобразуем уравнение к виду (2.4):

x = , т. е. (x)= .

Нетрудно убедиться, что корень уравнения находится на отрезке [/6, /3]. Например, вычислив значения f(x)на концах отрезка, получим: f(/6)> 0, а f(/3)< 0, т. е. функция на концах отрезка имеет разные знаки, что в соответствии с теоремой 2.1 указывает на то, что внутри отрезка есть корень. Расположение корня наглядно иллюстрирует рис.2.7.

Рис. 2.7

Подсчитаем, первую и вторую производные функции (x):

(x) = , "(x) = .

Так как "(x) > 0 на отрезке [/6, /3], то производная (x) монотонно возрастает на этом отрезке и принимает максимальное значение на правом конце отрезка, т. е. в точке /3. Поэтому, справедлива оценка:

|(x)| |(/3)| 0.312.

Таким образом, условие (2.7) выполнено, q < 0.5, и можно воспользоваться критерием окончания вычислений в виде (2.10). В табл. 2.2 приведены приближения, полученные по расчетной формуле (2.5). В качестве начального приближения выбрано значение x₀ = 1.

Таблица 2.2


n	x_n
	0.8415 0.8861 0.8742 0.8774 0.8765

Критерий окончания выполняется при n = 5, |x₅ - x₄| < 0.001. Сходимость двусторонняя, качественный характер такой сходимости представлен на рис. 2.4. Приближенное значение корня с требуемой точностью x^* 0.8765.

2.5 Метод Ньютона (метод касательных)

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

|x_n - x_{n -}₁| < . (2.18)

Пример 2.3.

Применим метод Ньютона для вычисления . где a > 0, p - натуральное число. Вычисление эквивалентно решению уравнения x^p = a. Таким образом, нужно найти корень уравнения f(x) = 0, f(x) = x^p - a, f (x) = px^{p -}¹. Итерационная формула метода (2.13) примет вид:

x_{n +1}= x_n- = x_n+ . (2.19)

Используя формулу (2.19), найдем с точностью = 10^-3.

x_{n +1}= x_n+ .

Простой корень уравнения x³ - 7 = 0 расположен на отрезке [1, 2]. Действительно, на концах отрезка [1, 2] функция f(x) = x³ - 7 принимает разные знаки, f (1) < 0, f (2) > 0. Кроме того, при x = 2 выполнено достаточное условие сходимости (2.16): f (2)f" (2)0.

Поэтому в качестве начального приближения можно взять x₀ = 2.

Результаты приведены в табл. 2.3.

Таблица 2.3


n	x_n
	0.8415 0.8861 0.8742 0.8774 0.8765

2.6 Метод секущих (метод хорд)

В этом и следующем разделе рассмотрим модификации метода Ньютона.

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

f (x_n) ,

то вместо формулы (2.13) получим

x_{n +1}= x_n-. . (2.20)

Это означает, что касательные заменены секущими. Метод секущих является двухшаговым методом, для вычисления приближения x_{n +1}необходимо вычислить два предыдущих приближения x_n и x_{n - 1}, и, в частности, на первой итерации надо знать два начальных значения x₀и x₁.

Формула (2.20) является расчетной формулой метода секущих. На рис. 2.9 приведена геометрическая иллюстрация метода секущих.

Рис. 2.9

Очередное приближение x_{n +1}получается как точка пересечения с осью OX секущей, соединяющей точки графика функции f(x) с координатами (x_{n -1}, f(x_{n - 1})) и (x_n, f(x_n)).

Сходимость метода. Сходимость метода секущих устанавливает следующая теорема.

Теорема 2.4 Пусть x^*- простой корень уравнения f(x) = 0, и в некоторой окрестности этого корня функция f дважды непрерывно дифференцируема, причем f"(x) 0. Тогда найдется такая малая -окрестность корня x^*, что при произвольном выборе начальных приближений x₀ и x₁из этой окрестности итерационная последовательность, определенная по формуле (2.20) сходится и справедлива оценка:

|x_{n +}₁ - x^*| C |x_n - x^*| ^p, n 0, p = 1.618. (2.21)

Сравнение оценок (2.15) и (2.21) показывает, что p < 2, и метод секущих сходится медленнее, чем метод Ньютона. Но в методе Ньютона на каждой итерации надо вычислять и функцию, и производную, а в методе секущих - только функцию. Поэтому при одинаковом объеме вычислений в методе секущих можно сделать примерно вдвое больше итераций и получить более высокую точность.

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

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

|x_n - x_{n -}₁| < . (2.22)

Пример 2.4.

Применим метод секущих для вычисления положительного корня уравнения 4(1 - x²) - e^x = 0 с точностью = 10^-3.

Корень этого уравнения находится на отрезке [0, 1], так как f (0) = 3 > 0, а f (1) = -e < 0. Подсчитаем вторую производную функции: f "(x) = -8 - e^x. Условие f(x)f " (x)0 выполняется для точки b = 1. В качестве начального приближения возьмем x₀ = b = 1. В качестве второго начального значения возьмем x₁= 0.5. Проведем вычисления по расчетной формуле (2.20). Результаты приведены в табл. 2.4.

Таблица 2.4


n	x_n
	1.0000 0.5000 0.6660 0.7093 0.7033 0.7034

2.7 Метод ложного положения

Рассмотрим еще одну модификацию метода Ньютона.

Пусть известно, что простой корень x^*уравнения f(x) = 0 находится на отрезке [a, b] и на одном из концов отрезка выполняется условие f(x)f"(x)0. Возьмем эту точку в качестве начального приближения. Пусть для определенности это будет b. Положим x₀ = a. Будем проводить из точки B = (b, f(b)) прямые через расположенные на графике функции точки B_nс координатами (x_n, f(x_n), n = 0, 1, … . Абсцисса точки пересечения такой прямой с осью OX есть очередное приближение x_n+₁.

Геометрическая иллюстрация метода приведена на рис. 2.10.

Рис. 2.10

Прямые на этом рисунке заменяют касательные в методе Ньютона (рис. 2.8). Эта замена основана на приближенном равенстве

f (x_n) . (2.23)

Заменим в расчетной формуле Ньютона (2.13) производную f (x_n) правой частью приближенного равенства (2.23). В результате получим расчетную формулу метода ложного положения:

x_{n +1}= x_n-.. (2.24)

Метод ложного положения обладает только линейной сходимостью. Сходимость тем выше, чем меньше отрезок [a, b].

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

|x_n - x_{n -}₁| < . (2.25)

Пример 2.5.

Применим метод ложного положения для вычисления корня уравнения x³ + 2x - 11 = 0 с точностью = 10^-3.

Корень этого уравнения находится на отрезке [1, 2], так как f (1) = -8 < 0, а f (2) = 1 > 0. Для ускорения сходимости возьмем более узкий отрезок [1.9, 2], поскольку f (1.9) < 0, а f (2) > 0. Вторая производная функции f (x) = x³ + 2x - 11 равна 6x. Условие f(x)f"(x)0 выполняется для точки b = 2. В качестве начального приближения возьмем x₀ = a = 1.9. По формуле (2.24) имеем

x₁= x₀-. = 1.9 + 1.9254.

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

Таблица 2.5


n	x_n
	1.9 1.9254 1.9263 1.9263

Тема 3. Решение систем линейных алгебраических уравнений

3.1 Постановка задачи

Требуется найти решение системы линейных уравнений:

a₁₁x₁ + a₁₂x₂ + a₁₃x₃ + … + a_1nx_n= b₁

a₂₁x₁ + a₂₂x₂ + a₂₃x₃ + … + a_2nx_n= b₂

a₃₁x₁ + a₃₂x₂ + a₃₃x₃ + … + a_3nx_n= b₃(3.1)

a_n₁x₁ + a_n₂x₂ + a_n₃x₃ + … + a_nnx_n= b_n

или в матричной форме:

Ax = b, (3.2)

где

a₁₁a₁₂a₁₃ … a_1nx₁ b₁

a₂₁a₂₂a₂₃ … a_2nx₂b₂

A = a₃₁ a₃₂a₃₃ … a_3nx =x_{3 ,}b =b₃

a_n₁ a_n₂ a_n₃ a_nn x_nb_n

По правилу Крамера система n линейных уравнений имеет единственное решение, если определитель системы отличен от нуля (det A 0) и значение каждого из неизвестных определяется следующим образом:

x_j = , j = 1, …, n, (3.3)

где det A_j - определитель матрицы, получаемой заменой j-го столбца матрицы A столбцом правых частей b.

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

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

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

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

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

Эти методы будут рассмотрены в следующих разделах.

3.2 Метод исключения Гаусса. Схема единственного деления

Основная идея метода исключений Гаусса состоит в том, что система уравнений (3.1) приводится к эквивалентной ей системе с верхней треугольной матрицей (прямой ход исключений), а затем неизвестные вычисляются последовательной подстановкой (обратный ход исключений).

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

Прямой ход состоит из n - 1 шагов. На первом шаге исключается переменная x₁ из всех уравнений, кроме первого. Для этого нужно из второго, третьего, …, n-го уравнений вычесть первое, умноженное на величину

m = , i = 2, 3, …, n. (3.4)

При этом коэффициенты при x₁обратятся в нуль во всех уравнениях, кроме первого.

Введем обозначения:

a = a_ij - ma_1j , b= b_i - mb₁. (3.5)

Легко убедиться, что для всех уравнений, начиная со второго, a= 0, i = 2, 3, …, n. Преобразованная система запишется в виде:

a₁₁x₁ + a₁₂x₂ + a₁₃x₃ + … + a_1nx_n= b₁

ax₂ + ax₃ + … + ax_n= b

a x₂ + ax₃ + … + ax_n= b (3.6)

ax₂ + ax₃ + … + ax_n= b

Все уравнения (3.6), кроме первого, образуют систему (n - 1)-го порядка. Применяя к ней ту же процедуру, мы можем исключить из третьего, четвертого, …, n-го уравнений переменную x₂. Точно так же исключаем переменную x₃из последних n - 3 уравнений.

На некотором k-ом шаге в предположении, что главный элемент k-ого шага a0, переменная x_k исключается с помощью формул:

m = ,

a = a - ma ,

b= b - mb, i, j = k + 1, k +2, …, n. (3.7)

Индекс k принимает значения 1, 2, …, n - 1.

При k = n - 1 получим треугольную систему:

a₁₁x₁ + a₁₂x₂ + a₁₃x₃ + … + a_1nx_n= b₁

ax₂ + ax₃ + …+ ax_n= b

ax₃ + …+ ax_n= b (3.8)

ax_n= b

с треугольной матрицей A_n.

Приведение системы (3.1) к треугольному виду (3.8) составляет прямой ход метода Гаусса.

При использовании метода Гаусса нет необходимости в предварительном обосновании существования и единственности решения (т. е. доказательства, что det A 0). Если на k-ом шаге все элементы a (i = k, k +1, …, n) окажутся равными нулю, то система (3.1) не имеет единственного решения.

Обратный ход состоит в вычислении переменных. Из последнего уравнения (3.8) определяем x_n_... Подставляя его в предпоследнее уравнение, находим x_n-₁, и т. д. Общие формулы имеют вид:

x_n= ,

x_k= (b- a x_k+₁ - a x_k+₂ - … - a x_n), k = n - 1, n - 2, …, 1 (3.9)

Трудоемкость метода. Для реализации метода исключения Гаусса требуется примерно 2/3n³ операций для прямого хода и n² операций для обратного хода. Таким образом, общее количество операций составляет примерно 2/3n³ + n².

Пример 3.1.

Применим метод исключения Гаусса по схеме единственного деления для решения системы уравнений:

2.0x₁ + 1.0x₂- 0.1x₃ + 1.0x₄= 2.7

0.4x₁ + 0.5x₂ + 4.0x₃- 8.5x₄= 21.9

0.3x₁- 1.0x₂ + 1.0x₃ + 5.2x₄= -3.9 (3.10)

1.0x₁ + 0.2x₂ + 2.5x₃- 1.0x₄= 9.9

Будем делать округление чисел до четырех знаков после десятичной точки.

Прямой ход. 1-ый шаг. Вычислим множители:

m = = = 0.2; m = = = 0.15; m = = = 0.5.

Вычитая из второго, третьего и четвертого уравнений системы (3.10) первое уравнение, умноженное соответственно на m, m, m, получим новую систему:

2.0x₁ + 1.0x₂- 0.1x₃ + 1.0x₄= 2.7

0.3x₂ + 4.02x₃- 8.70x₄= 21.36

-1.15x₂ + 1.015x₃ + 5.05x₄= -4.305 (3. 11)

-0.30x₂ + 2.55x₃-1.50x₄= 8.55

2-ой шаг. Вычислим множители:

m = = = - 3.83333; m = = = -1.0.

Вычитая из третьего и четвертого уравнений системы (3.11) второе уравнение, умноженное соответственно на m и m, приходим к системе:

2.0x₁ + 1.0x₂- 0.1x₃ + 1.0x₄= 2.7

0.3x₂ + 4.02x₃- 8.70x₄= 21.36

16. 425x₃- 28.300x₄= 77.575 (3.12)

6.570x₃-10.200x₄= 29.910

3-ий шаг. Вычислим множитель:

m = = = 0.4.

Вычитая из четвертого уравнения системы (3.12) третье, умноженное на m, приведем систему к треугольному виду:

2.0x₁ + 1.0x₂- 0.1x₃ + 1.0x₄= 2.7

0.3x₂ + 4.02x₃- 8.70x₄= 21.36

16. 425x₃- 28.300x₄= 77.575 (3.13)

1.12x₄= -1.12

Обратный ход. Из последнего уравнения системы (3.13) находим x₄= 1.000. Подставляя значение x₄в третье уравнение, получим x₃= 2.000. Подставляя найденные значения x₄и x₃ во второе уравнение, найдем x₂ = 3.000. Наконец, из первого уравнения, подставив в него найденные значения x₄, x₃и x₂, вычислим x₁= -1.000.

Итак система (3.10) имеет следующее решение:

x₁= 1.000, x₂ = 2.000, x₃= 3.000, x₄= - 1.000.

3.3 Метод исключения Гаусса с выбором главного элемента по столбцу

Хотя метод Гаусса является точным методом, ошибки округления могут привести к существенным погрешностям результата. Кроме того исключение по формулам (3.7) нельзя проводить, если элемент главной диагонали a равен нулю. Если элемент a мал, то велики ошибки округления при делении на этот элемент. Для уменьшения ошибок округления применяют метод исключения Гаусса с выбором главного элемента по столбцу. Прямой ход так же, как и для схемы единственного деления, состоит из n - 1 шагов. На первом шаге прежде, чем исключать переменную x₁, уравнения переставляются так, чтобы в левом верхнем углу был наибольший по модулю коэффициент a_i₁, i = 1, 2, …, n. В дальнейшем, на k-м шаге, прежде, чем исключать переменную x_k, уравнения переставляются так, чтобы в левом верхнем углу был наибольший по модулю коэффициент a_ik, i = k, k + 1, …, n. После этой перестановки исключение переменной x_kпроизводят, как в схеме единственного деления.

Трудоемкость метода. Дополнительные действия по выбору главных элементов требуют примерно n² операций, что практически не влияет на общую трудоемкость метода.

Пример 3.2.

Применим метод исключения Гаусса с выбором главного элемента по столбцу для решения системы уравнений (3.10) из примера 3.1. Прямой ход. 1-ый шаг. Так как коэффициент a₁₁ = 2.0 наибольший из коэффициентов первого столбца, перестановки строк не требуется и 1-ый шаг полностью совпадает с 1-ым шагом примера 3.1. Из второго, третьего и четвертого уравнений исключается переменная x₁ и система приводится к виду (3.11).

2-ой шаг. Наибольший по модулю коэффициент при x₂ в системе (3.11) a = -1.15. Поэтому переставим уравнения следующим образом:

2.0x₁ + 1.0x₂- 0.1x₃ + 1.0x₄= 2.7

-1.15x₂ + 1.015x₃ + 5.05x₄= -4.305 (3.14)

0.3x₂ + 4.02x₃- 8.70x₄= 21.36

-0.30x₂ + 2.55x₃-1.50x₄= 8.55

Вычислим множители:

m = = = -0.26087 m = = = 0.26087.

Вычитая из третьего и четвертого уравнений системы (3.14) второе уравнение, умноженное соответственно на m и m, приходим к системе:

2.0x₁ + 1.0x₂- 0.1x₃ + 1.0x₄= 2.7

-1.15x₂ + 1.015x₃ + 5.05x₄= -4.305 (3.15)

4.28478x₃- 7.38261x₄= 20.23696

2.28522x₃-2.81739x₄= 9.67305

3-ий шаг. Вычислим множитель:

m = = = 0.53333.

Вычитая из четвертого уравнения системы (3.15) третье, умноженное на m, приведем систему к треугольному виду:

2.0x₁ + 1.0x₂- 0.1x₃ + 1.0x₄= 2.7

-1.15x₂ + 1.015x₃ + 5.05x₄= -4.305 (3.16)

4.28478x₃- 7.38261x₄= 20.23696

1.11998x₄= -1.11998

Обратный ход. Обратный ход полностью совпадает с обратным ходом примера 3.1. Решение системы имеет вид:

x₁= 1.000, x₂ = 2.000, x₃= 3.000, x₄= - 1.000.

3.4 Вычисление определителя методом исключения Гаусса

Из курса линейной алгебры известно, что определитель треугольной матрицы равен произведению диагональных элементов. В результате метода исключений Гаусса система линейных уравнений (3.2) с квадратной матрицей A приводится к эквивалентной ей системе (3.8) с треугольной матрицей A_n. Поэтому

det A = (-1)^s det A_n,

где s - число перестановок строк, (s = 0, если использовался метод Гаусса по схеме единственного деления).Таким образом,

det A = (-1)^s a₁₁ aa …a (3.17)

Итак, для вычисления определителя det A необходимо выполнить процедуру прямого хода в методе Гаусса для системы уравнений Ax = 0, затем найти произведение главных элементов, стоящих на диагонали треугольной матрицы и умножить это произведение на (-1)^s, где s - число перестановок строк.

Пример 3.3.

Вычислим определитель det A =

2.0 1.0 0.1 1.0

0.4 0.5 4.0 8.5

0.3 1.0 1.0 5.2

1.0 0.2 2.5 1.0

Данный определитель совпадает с определителем системы, рассмотренной в примере 3.1. Он равен произведению диагональных элементов треугольной матрицы (3.13):

det A = 2.0 0.30 16.425 1.12 = 11.0376.

Если же обратиться к примеру 3.2, то, учитывая, что была одна перестановка строк, т.е. s = 1, получим:

det A = (-1) 2.0 (-1.15) 4.28478 1.11998 = 11.0375.

3.5 Вычисление обратной матрицы методом исключения Гаусса

Обратной матрицей к матрице A называется матрица A^-1, для которой выполнено соотношение:

A A^-¹ = E, (3.18)

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

1 0 0 … 0

0 1 0 … 0

E = 0 0 1 … 0 . (3.19)

0 0 0 … 1

Квадратная матрица A называется невырожденной, если det A 0. Всякая невырожденная матрица имеет обратную матрицу.

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

Пусть A - квадратная невырожденная матрица порядка n:

a₁₁ a₁₂ a₁₃ … a_1n

a₂₁ a₂₂a₂₃ … a_2n

A = a₃₁ a₃₂a₃₃ … a_3n

a_n₁ a_n₂ a_n₃ … a_nn

и A^-1- ееобратная матрица:

x₁₁ x₁₂ x₁₃ … x_1n

x₂₁ x₂₂x₂₃ … x_2n

A^-1= x₃₁ x₃₂x₃₃ … x_3n

x_n₁ x_n₂ x_n₃ … x_nn

k = n - 1, n - 2, …, 1 (3.9)

Пример 3.1.

Применим метод исключения Гаусса по схеме единственного деления для решения системы уравнений:

2.0x₁ + 1.0x₂- 0.1x₃ + 1.0x₄= 2.7

0.4x₁ + 0.5x₂ + 4.0x₃- 8.5x₄= 21.9

0.3x₁- 1.0x₂ + 1.0x₃ + 5.2x₄= -3.9 (3.10)

1.0x₁ + 0.2x₂ + 2.5x₃- 1.0x₄= 9.9

Будем делать округление чисел до четырех знаков после десятичной точки.

Прямой ход. 1-ый шаг. Вычислим множители:

m = = = 0.2; m = = = 0.15; m = = = 0.5.

2.0x₁ + 1.0x₂- 0.1x₃ + 1.0x₄= 2.7

0.3x₂ + 4.02x₃- 8.70x₄= 21.36

-1.15x₂ + 1.015x₃ + 5.05x₄= -4.305 (3. 11)

-0.30x₂ + 2.55x₃-1.50x₄= 8.55

2-ой шаг. Вычислим множители:

m = = = - 3.83333; m = = = -1.0.

2.0x₁ + 1.0x₂- 0.1x₃ + 1.0x₄= 2.7

0.3x₂ + 4.02x₃- 8.70x₄= 21.36

16. 425x₃- 28.300x₄= 77.575 (3.12)

6.570x₃-10.200x₄= 29.910

3-ий шаг. Вычислим множитель:

m = = = 0.4.

Вычитая из четвертого уравнения системы (3.12) третье, умноженное на m, приведем систему к треугольному виду:

2.0x₁ + 1.0x₂- 0.1x₃ + 1.0x₄= 2.7

0.3x₂ + 4.02x₃- 8.70x₄= 21.36

16. 425x₃- 28.300x₄= 77.575 (3.13)

1.12x₄= -1.12

Итак система (3.10) имеет следующее решение:

x₁= 1.000, x₂ = 2.000, x₃= 3.000, x₄= - 1.000.

3.3 Метод исключения Гаусса с выбором главного элемента по столбцу

Пример 3.2.

2.0x₁ + 1.0x₂- 0.1x₃ + 1.0x₄= 2.7

-1.15x₂ + 1.015x₃ + 5.05x₄= -4.305 (3.14)

0.3x₂ + 4.02x₃- 8.70x₄= 21.36

-0.30x₂ + 2.55x₃-1.50x₄= 8.55

Вычислим множители:

m = = = -0.26087 m = = = 0.26087.

2.0x₁ + 1.0x₂- 0.1x₃ + 1.0x₄= 2.7

-1.15x₂ + 1.015x₃ + 5.05x₄= -4.305 (3.15)

4.28478x₃- 7.38261x₄= 20.23696

2.28522x₃-2.81739x₄= 9.67305

3-ий шаг. Вычислим множитель:

m = = = 0.53333.

Вычитая из четвертого уравнения системы (3.15) третье, умноженное на m, приведем систему к треугольному виду:

2.0x₁ + 1.0x₂- 0.1x₃ + 1.0x₄= 2.7

-1.15x₂ + 1.015x₃ + 5.05x₄= -4.305 (3.16)

4.28478x₃- 7.38261x₄= 20.23696

1.11998x₄= -1.11998

Обратный ход. Обратный ход полностью совпадает с обратным ходом примера 3.1. Решение системы имеет вид:

x₁= 1.000, x₂ = 2.000, x₃= 3.000, x₄= - 1.000.

3.4 Вычисление определителя методом исключения Гаусса

det A = (-1)^s det A_n,

det A = (-1)^s a₁₁ aa …a (3.17)

Пример 3.3.

Вычислим определитель det A =

2.0 1.0 0.1 1.0

0.4 0.5 4.0 8.5

0.3 1.0 1.0 5.2

1.0 0.2 2.5 1.0

det A = 2.0 0.30 16.425 1.12 = 11.0376.

Если же обратиться к примеру 3.2, то, учитывая, что была одна перестановка строк, т.е. s = 1, получим:

det A = (-1) 2.0 (-1.15) 4.28478 1.11998 = 11.0375.

3.5 Вычисление обратной матрицы методом исключения Гаусса

Обратной матрицей к матрице A называется матрица A^-1, для которой выполнено соотношение:

A A^-¹ = E, (3.18)

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

1 0 0 … 0

0 1 0 … 0

E = 0 0 1 … 0 . (3.19)

0 0 0 … 1

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

Пусть A - квадратная невырожденная матрица порядка n:

a₁₁ a₁₂ a₁₃ … a_1n

a₂₁ a₂₂a₂₃ … a_2n

A = a₃₁ a₃₂a₃₃ … a_3n

a_n₁ a_n₂ a_n₃ … a_nn

и A^-1- ееобратная матрица:

x₁₁ x₁₂ x₁₃ … x_1n

x₂₁ x₂₂x₂₃ … x_2n

A^-1= x₃₁ x₃₂x₃₃ … x_3n

x_n₁ x_n₂ x_n₃ … x_nn

Используя соотношения (3.18), (3. 19) и правило умножения матриц, получим систему из n² уравнений с n² переменными x_ij, i, j = 1, 2, …, n. Чтобы получить первый столбец матрицы E, нужно почленно умножить каждую строку матрицы A на первый столбец матрицы A^-1и приравнять полученное произведение соответствующему элементу первого столбца матрицы E. В результате получим систему уравнений:

a₁₁x₁₁ + a₁₂x₂₁ + a₁₃x₃₁ + … + a_1nx_n₁= 1

a₂₁x₁₁ + a₂₂x₂₁ + a₂₃x₃₁ + … + a_2nx_n₁= 0

a₃₁x₁₁ + a₃₂x₂₁ + a₃₃x₃₁ + … + a_3nx_n₁= 0(3.20)

a_n₁x₁₁ + a_n₂x₂₁ + a_n₃x₃₁ + … + a_nnx_n₁= 0

Аналогично, чтобы получить второй столбец матрицы E, нужно почленно умножить каждую строку матрицы A на второй столбец матрицы A^-1и приравнять полученное произведение соответствующему элементу второго столбца матрицы E. В результате получим систему уравнений:

a₁₁x₁₂ + a₁₂x₂₂ + a₁₃x₃₂ + … + a_1nx_n₂= 0

a₂₁x₁₂ + a₂₂x₂₂ + a₂₃x₃₂ + … + a_2nx_n₂= 1

a₃₁x₁₂ + a₃₂x₂₂ + a₃₃x₃₂ + … + a_3nx_n₂= 0(3.21)

a_n₁x₁₂ + a_n₂x₂₂ + a_n₃x₃₂ + … + a_nnx_n₂= 0

и т. д.

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

Пример 3.4.

Вычислим обратную матрицу A^-1для матрицы

A = 1.8 -3.8 0.7 -3.7

0.7 2.1 -2.6 -2.8

7.3 8.1 1.7 -4.9

1.9 -4.3 -4.3 -4.7

По формулам (3.7) за три шага прямого хода преобразуем матрицу Aв треугольную матрицу

1.8 -3.8 0.7 -3.7

0 3.57778 -2.87222 -1.36111

0 0 17.73577 19.04992

0 0 0 5.40155

Далее, применим процедуру обратного хода четыре раза для столбцов свободных членов, преобразованных по формулам (3.7) из столбцов единичной матрицы:

1 0 0 0

0 1 0 0

0 , 0 , 1 , 0

0 0 0 1

Каждый раз будем получать столбцы матрицы A^-1. Опустив промежуточные вычисления, приведем окончательный вид матрицы A^-1:

-0.21121 -0.46003 0.16248 0.26956

-0.03533 0.16873 0.01573 -0.08920

0.23030 0.04607 -0.00944 -0.19885 .

-0.29316 -0.38837 0.06128 0.18513

3.6 Метод простой итерации Якоби

Метод Гаусса обладает довольно сложной вычислительной схемой. Кроме того, при вычислениях накапливается ошибка округления, что может привести к недостаточно точному результату. Рассмотрим метод простой итерации Якоби, свободный от этих недостатков, хотя требующий приведения исходной системы уравнений к специальному виду.

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

Ax = b (3.22)

с квадратной невырожденной матрицей A привести к виду

x = Bx + c, (3.23)

где B - квадратная невырожденная матрица с элементами b_ij, i, j = 1, 2, …, n, x - вектор-столбец неизвестных x_i, c - вектор-столбец с элементами c_i, i = 1, 2, …, n.

Существуют различные способы приведения системы (3.22) к виду (3.23). Рассмотрим самый простой. Представим систему (3.22) в развернутом виде:

a₁₁x₁ + a₁₂x₂ + a₁₃x₃ + … + a_1nx_n= b₁

a₂₁x₁ + a₂₂x₂ + a₂₃x₃ + … + a_2nx_n= b₂

a₃₁x₁ + a₃₂x₂ + a₃₃x₃ + … + a_3nx_n= b₃(3.24)

a_n₁x₁ + a_n₂x₂ + a_n₃x₃ + … + a_nnx_n= b_n

Из первого уравнения системы (3.24) выразим неизвестную x₁:

x₁ = a(b₁ - a₁₂x₂- a₁₃x₃ - … - a_1nx_n),

из второго уравнения - неизвестную x₂:

x₂ = a(b₂ - a₂₁x₁- a₂₃x₃ - … - a_2nx_n),

и т. д. В результате получим систему:

x₁= b₁₂x₂ + b₁₃x₃ + … + b_1,n-1x_n-₁ + b_1nx_n + c₁

x₂= b₂₁x₁+ b₂₃x₃ + … + b_2,n-1x_n-₁+ b_2nx_n+ c₂

x₃= b₃₁x₁ + b₃₂x₂+ … + b_3,n-1x_n-₁+ b_3nx_n+ c₃(3.25)

x_n= b_n₁x₁ + b_n₂x₂ + b_n₃x₃+ b_n,n-₁x_n-₁+ c_n

Матричная запись системы (3.25) имеет вид (3.23). На главной диагонали матрицы B находятся нулевые элементы, а остальные элементы вычисляются по формулам:

b_ij = , c_i= , i, j = 1,2, …n, i j. (3.26)

Очевидно, что диагональные элементы матрицы A должны быть отличны от нуля.

Выберем произвольно начальное приближение Обычно в качестве первого приближения берут x= c_i или x= 0. Подставим начальное приближение в правую часть (3.25). Вычисляя левые части, получим значения x, x, …, x. Продолжая этот процесс дальше, получим последовательность приближений, причем (k + 1)-е приближение строится следующим образом:

x = b₁₂x + b₁₃ x + … + b_1,n-1 x + b_1n x + c₁

x = b₂₁x₁+ b₂₃ x + … + b_2,n-1 x + b_2n x + c₂

x= b₃₁x + b₃₂ x + … + b_3,n-1 x + b_3n x + c₃ … (3.27)

x= b_n₁x + b_n₂x + b_n₃ x + b_n,n-₁ x + c._n

Система (3.27) представляет собой расчетные формулы метода простой итерации Якоби.

Сходимость метода простой итерации. Известно следующее достаточное условие сходимости метода простой итерации Якоби.

Если элементы матрицы A удовлетворяют условию:

|a_ii| > , i = 1, 2, …, n. (3.28)

то итерационная последовательность x^k сходится к точному решению x^*.

Условие (3.28) называют условием преобладания диагональных элементов матрицы A, так как оно означает, что модуль диагонального элемента i-ой строки больше суммы модулей остальных элементов этой строки, i = 1, 2, …, n.

Необходимо помнить, что условие сходимости (3.28) является лишь достаточным. Его выполнение гарантирует сходимость метода простых итераций, но его невыполнение, вообще говоря, не означает, что метод расходится.

Справедлива следующая апостериорная оценка погрешности:

max|x - x| max|x- x|, i = 1, 2, …, n, (3.29)

где = max |b_ij| i, j = 1, 2, …, n.

Правую часть оценки (3.29) легко вычислить после нахождения очередного приближения.

Критерий окончания. Если требуется найти решение с точностью , то в силу (3.29) итерационный процесс следует закончить как только на (k+1)-ом шаге выполнится неравенство:

max|x- x| < , i = 1, 2, …, n. (3.30)

Поэтому в качестве критерия окончания итерационного процесса можно использовать неравенство

max|x- x| < ₁, i = 1, 2, …, n. (3.31)

где _{1 = .}

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

max|x- x| < , i = 1, 2, …, n. (3.32)

В других случаях использование критерия (3.32) неправомерно и может привести к преждевременному окончанию итерационного процесса.

Пример 3.5.

Применим метод простой итерации Якоби для решения системы уравнений

20.9x₁ + 1.2 x₂ + 2.1x₃ + 0.9x₄= 21.70

1.2x₁ + 21.2 x₂ + 1.5x₃ + 2.5x₄= 27.46

2.1x₁ + 1.5 x₂ + 19.8x₃ + 1.3x₄= 28.76 (3.33)

0.9x₁ + 2.5 x₂ + 1.3x₃ + 32.1x₄= 49.72

Заметим, что метод простой итерации сходится, т. к. выполняется условие преобладания диагональных элементов (3.28):

|20.9| > |1.2 + 2.1 + 0.9|,

|21.2| > |1.2| + |1.5| + |2.5|,

|19.8| > |2.1| + |1.5| + |1.3|,

|32.1| > |0.9| + |2.5| + |1.3|.

Пусть требуемая точность = 10^-3. Вычисления будем проводить с четырьмя знаками после десятичной точки.

Приведем систему к виду (3.25):

x₁= - 0.0574 x₂- 0.1005x₃- 0.0431x₄+ 1.0383

x₂= -0.0566x₁ - 0.0708x₃ - 0.1179x₄+ 1.2953

x₃= -0.1061x₁ - 0.0758 x₂- 0.0657x₄+ 1.4525 (3.34)

x₄= -0.0280x₁ - 0.0779 x₂ - 0.0405x₃ + 1.5489

Величина = max |b_ij|, i, j = 1, 2, 3,4 равна 0.1179, т. е. выполняется условие , и можно пользоваться критерием окончания итерационного процесса (3.32).

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

x = 1.0383, x = 1.2953, x = 1.4525, x = 1.5489. (3.35)

Вычисления будем вести до тех пор, пока все величины |x- x|, i = 1, 2, 3, 4, а следовательно, и max|x- x| не станут меньше = 10^-3.

Последовательно вычисляем:

при k = 1

x = - 0.0574x - 0.1005x - 0.0431x + 1.0383 = 0.7512

x = -0.0566x - 0.0708x - 0.1179x + 1.2953 = 0.9511

x = -0.1061x - 0.0758 x - 0.0657x + 1.4525 = 1.1423

x = -0.0280x - 0.0779x - 0.0405x + 1.5489 = 1.3601

при k = 2

x= 0.8106, x= 1.0118, x= 1.2117, x= 1.4077.

при k = 3

x= 0.7978, x= 0.9977, x= 1.1975, x= 1.3983.

при k = 4

x= 0.8004, x= 1.0005, x= 1.2005, x = 1.4003.

Вычисляем модули разностей значений xпри k = 3 и k = 4:

| x- x| = 0.026, | x- x| = 0.028, | x- x| = 0.0030, | x- x| = 0.0020.

Так как все они больше заданной точности = 10^-3, продолжаем итерации.

При k = 5

x= 0.7999, x= 0.9999, x= 1.1999, x = 1.3999.

Вычисляем модули разностей значений xпри k = 4 и k = 5:

| x- x| = 0.0005, | x- x| = 0.0006, | x - x| = 0.0006, | x- x| = 0.0004.

Все они меньше заданной точности = 10^-3, поэтому итерации заканчиваем. Приближенным решением системы являются следующие значения:

x₁0.7999, x₂0.9999, x₃1.1999, x₄1.3999.

Для сравнения приведем точные значения переменных:

x₁= 0.8, x₂= 1.0, x₃= 1.2, x₄= 1.4.

3.7 Метод Зейделя

Модификацией метода простых итераций Якоби можно считать метод Зейделя.

В методе Якоби на (k+1)-ой итерации значения x, i = 1, 2, …, n. вычисляются подстановкой в правую часть (3.27) вычисленных на предыдущей итерации значений x. В методе Зейделя при вычислении xиспользуются значения x, x, x, уже найденные на (k+1)-ой итерации, а не x, x, …, x, как в методе Якоби, т.е. (k + 1)-е приближение строится следующим образом:

x = b₁₂x + b₁₃ x + … + b_1,n-1 x + b_1n x + c₁

x = b₂₁x + b₂₃ x + … + b_2,n-1 x + b_2n x + c₂

x= b₃₁x + b₃₂ x + … + b_3,n-1 x + b_3n x + c₃(3.36)

x= b_n₁ x + b_n₂x x + b_n₃ x x+ … + b_n,n-₁ x + c._n

Формулы (3.36) являются расчетными формулами метода Зейделя.

Введем нижнюю и верхнюю треугольные матрицы:

0 0 0 … 0 0 b₁₂b₁₃ … b_1n

b₂₁0 0 … 000 b₂₃ … b_2n

B₁= b₃₁ b₃₂0 … 0 и B₂= 000 … b_{3n .}

b_n₁ b_n₂ b_n₃ …00 0 0 … 0

Матричная запись расчетных формул (3.36) имеет вид:

x^k+¹= B₁x^k+¹+ B₂x^k+ c.(3.37)

Так как B = B₁+ B₂, точное решение x^* исходной системы удовлетворяет равенству:

x^*= B₁x^*+ B₂x^*+ c.(3.38)

Сходимость метода Зейделя. Достаточным условием сходимости метода Зейделя является выполнение неравенства:

= max |b_ij|,< 1, i, j = 1, 2, …, n. (3.39)

Неравенство (3.39) означает, что для сходимости метода Зейделя достаточно, чтобы максимальный по модулю элемент матрицы B был меньше единицы.

Если выполнено условие (3.39), то справедлива следующая апостериорная оценка погрешности:

max|x - x| max|x- x| i = 1, 2, …, n, (3.40)

где - максимальный элемент матрицы B, ₂- максимальный элемент матрицы B₂.

Правую часть оценки (3.40) легко вычислить после нахождения очередного приближения.

Критерий окончания. Если требуется найти решение с точностью , то в силу (3.37) итерационный процесс следует закончить как только на (k+1)-ом шаге выполнится неравенство:

max|x- x| < , i = 1, 2, …, n. (3.41)

Поэтому в качестве критерия окончания итерационного процесса можно использовать неравенство

max|x- x| < ₁, i = 1, 2, …, n. (3.42)

где _{1 = .}

Если выполняется условие

5rik.ru

Материалы для учебы и работы

Пусть уравнение (2.1) можно заменить эквивалентным ему уравнением