Московский физико-технический институт
Опубликован: 25.10.2007 | Доступ: свободный | Студентов: 3209 / 642 | Оценка: 4.50 / 4.33 | Длительность: 24:00:00
ISBN: 978-5-9556-0065-9
Специальности: Программист, Математик
Лекция 3:

Численное решение систем линейных алгебраических уравнений

Аннотация: Рассматриваются наиболее употребительные приближенные методы решения систем линейных алгебраических уравнений. Вводятся согласованные нормы векторов и матриц. Вычисляется число обусловленности в различных нормах. Анализируется влияние ошибок округления на погрешность результата. Дается понятие о спектральных задачах. Для самосопряженной матрицы рассматривается метод вращений поиска собственных значений
Ключевые слова: системы линейных алгебраических уравнений, алгебра, ПО, матрица, вектор, компонент, определитель матрицы, прямые методы, итерационные методы, обусловленность системы, норма вектора, пространство, норма, согласованные нормы вектора и матрицы, выражение, грани, максимум, собственный вектор, неравенство, равенство, число обусловленности матрицы, погрешность, параметр, место, прямой, обратный, вычисление, метод Гаусса, LU-разложение, метод Гаусса с выбором главного элемента, задача Штурма-Лиувилля, численное решение уравнений, метод прогонки, Метод простых итераций, метод итерации, методы Якоби, Зейделя, верхней релаксации, матрица перехода, метод градиентного спуска, метод наискорейшего спуска, градиентный метод, метод сопряженных градиентов, спектральные задачи, собственное число, действительный, алгоритм, значение, отношение, метод вращений, определение, барьер, расширенная матрица, разложение вектора

К численному решению систем линейных алгебраических уравнений (СЛАУ) сводятся многие задачи математической физики. Математические модели, представляющие собой СЛАУ большой размерности, встречаются в математической экономике, биологии и т.п. Теория получения приближенных решений СЛАУ — часть вычислительной линейной алгебры. Сама вычислительная линейная алгебра, по-видимому, является наиболее обширной темой во всем курсе вычислительной математики. По прикладной линейной алгебре существует обширная литература (например, [2.1, 2.2, 2.3, 2.4, 2.5], а программы, реализующие наиболее популярные алгоритмы вычислительной линейной алгебры, являются неотъемлемой частью прикладного программного обеспечения, в частности, современных математических пакетов.

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

Рассмотрим СЛАУ вида

$ \mathbf{Au}= \mathbf{f},$ ( 2.1)

где \mathbf{A} — невырожденная ( \det \mathbf{A}\ne 0 ) квадратная матрица размером n x n

\mathbf{A}= \left( \begin{array}{cccc}
   a_{11} & a_{12} &  \ldots  & a_{1n}\\
   a_{21} & a_{22} &  \ldots  & a_{2n}\\
    \ldots  &  \ldots  &  \ldots  &  \ldots   \\
   a_{n1} & a_{n2} &  \ldots  & a_{nn}\\
 \end{array}\right),

\mathbf{u}={\{u_1, \ldots , u_n \}}^Tвектор-столбец решения, \mathbf{f}={\{f_1, \ldots , f_n \}}^Tвектор-столбец правой части.

Так как матрица системы — невырожденная, \Delta  = \det\mathbf{A}\ne 0, то решение системы (2.1) существует и единственно.

Из курса линейной алгебры [2.6] известно правило Крамера нахождения решения. Так, каждый компонент вектора неизвестных может быть вычислен как

$ u_i  = \frac{\Delta_i}{\Delta }, $

где \Delta _{i}определитель матрицы, получаемой из \mathbf{A} заменой i столбца столбцом правых частей. Однако несложные арифметические оценки позволяют понять, что использование этой формулы приводит к неоправданно большим затратам машинного времени [2.3]. Так, например, если одно слагаемое в \Delta вычисляется за 10 -6 с, то время расчета для n = 100 на существующих в момент написания книги компьютерах будет измеряться годами.

На самом деле в настоящее время с помощью компьютеров численно решаются СЛАУ намного более высокого порядка (примерно до n \approx  10^{6} ). Такие решения осуществляются при помощи прямых или итерационных численных методов. Прямые методы позволяют в предположении отсутствия ошибок округления (при проведении расчетов на идеальном, т.е. бесконечноразрядном компьютере) получить точное решение задачи за конечное число арифметических действий Итерационные методы, или методы последовательных приближений, позволяют вычислить последовательность \{{\mathbf{u}}_k\}, сходящуюся к решению задач при k \to \infty (на практике, разумеется, ограничиваются конечным k, в зависимости от требуемой точности).

Однако неточность в задании правых частей и элементов матрицы \mathbf{A} может приводить к значительным погрешностям при вычислении решения (2.1). В "Предмет вычислительной математики. Обусловленность задачи, устойчивость алгоритма, погрешности вычислений. Задача численного дифференцирования" на примере было показано, что такое явление наблюдается в случае плохо обусловленной системы. Остановимся подробнее на важном вопросе оценки погрешности решения СЛАУ.

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

2.2. Согласованные нормы векторов и матриц

В векторном n -мерном линейном нормированном пространстве введем следующие нормы вектора:

кубическая:

{\|\mathbf{u}\|}_1  = \max\limits_{1 \le i \le n}|u_i|, ( 2.2а)

октаэдрическая:

{\|\mathbf{u}\|}_2  = \sum\limits_{i = 1}^n|u_i|, ( 2.2б)

евклидова (в комплексном случае — эрмитова):

{\|\mathbf{u}\|}_3  ={\left(\sum\limits_{i = 1}^n{|u_i|}^2\right)}^{1/2}={(\mathbf{u}, \mathbf{u})}^{1/2}. ( 2.2в)

Рассмотрим квадратную матрицу \mathbf{A} и связанное с ней линейное преобразование \mathbf{v}= \mathbf{Au}, где \mathbf{v}, \mathbf{u} \in L^n ( Lnn -мерное линейное нормированное пространство). Норма матрицы определяется как действительное неотрицательное число, характеризующее это преобразование и определяющееся как

$ \|\mathbf{A}\|= \sup\limits_{\|\mathbf{u}\| \ne 0}
\frac{\|\mathbf{Au}\|}{\|\mathbf{u}\|}$ ( 2.3)

Укажем некоторые свойства нормы матрицы:

\begin{gather*}
\|\mathbf{A + B}\|= \|\mathbf{A}\|+ \|\mathbf{B}\|, \\ 
\|\mathbf{\lambda A}\|=|\lambda| \|\mathbf{A}\|, \\ 
\|\mathbf{AB}\| \le \|\mathbf{A}\| \|\mathbf{B}\|, \\ 
\|\mathbf{A}\|= 0 \text{ тогда и только тогда, когда }\mathbf{A}= 0.
\end{gather*}

Заметим, что норму матрицы (2.3) называют подчиненной норме вектора. Говорят, что норма матрицы \mathbf{A} согласована с нормой вектора \mathbf{u}, если выполнено условие

\|\mathbf{Au}\| \le \|\mathbf{A}\|\|\mathbf{u}\|.

Нетрудно видеть, что подчиненная норма согласована с соответствующей метрикой векторного пространства. В самом деле

$ \|\mathbf{A}\|= \sup\limits_{\|\mathbf{u}\| \ne 0} \frac{\|\mathbf{Au}\|}{\|\mathbf{u}\|}\ge \frac{\|\mathbf{Au}\|}{\|\mathbf{u}\|}\mbox{ откуда }\|\mathbf{Au}\| \le \|\mathbf{A}\|\cdot \|\mathbf{u}\|. $

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

\begin{gather*}
{\|\mathbf{A}\|}_1 = \max\limits_{1 \le i \le n}\sum\limits_{j = 1}^n|a_{ij}|, \\ 
{\|\mathbf{A}\|}_2 = \max\limits_{1 \le j \le n}\sum\limits_{i = 1}^n{|a_{ij}|}, \\ 
{\|\mathbf{A}\|}_3 = \sqrt{\max\limits_{1 \le i \le n} \lambda^i (\mathbf{A}^*  \cdot \mathbf{A})}.
\end{gather*}

Покажем, как получается выражение для согласованной нормы матрицы {\|\mathbf{A}\|}_1 , соответствующей норме вектора {\|\mathbf{u}\|}_1.

Вычислим норму вектора {\|\mathbf{Au}\|}_1:

\begin{multline*}
{\|\mathbf{Au}\|}_1 = \max\limits_i \left|{\sum\limits_j{a_{ij}u_j}}\right| \le \max\limits_i (\sum\limits_j{|{a_{ij}}|}|{u_j}|) \le \\ 
\le (\max\limits_i \sum\limits_j{|{a_{ij}}|}) \max\limits_j|{u_j}| = (\max\limits_i \sum\limits_j{|{a_{ij}}|}){\|\mathbf{u}\|}_1 ,
\end{multline*}

откуда

\frac{{\|\mathbf{Au}\|}_1}{{\|\mathbf{u}\|}_1}\le \max\limits_i \sum\limits_j |{a_{ij}}|.
По определению нормы матрицы как точной верхней грани отношения
\frac{\|\mathbf{Au}\|}{\|\mathbf{u}\|}, \max\limits_j \sum\limits_i{|{a_{ij}}|} = {\|\mathbf{A}\|}_1 ,
если существует вектор, на котором точная верхняя грань достигается.

Покажем, что таким вектором является, например,

{\mathbf{v}}_k = \{sign a_{k1}, \ldots , sign a_{kn}\}^T,
при этом допустим, что максимум в последнем неравенстве достигается при i = k.

Поскольку \|{\mathbf{v}}_k\|= k, то \sum\limits_j{a_{kj}v_j} = \sum\limits_j {|a_{kj}|} =  \max\limits_i \sum\limits_j{|{a_{ij}}|}.

Тогда, в соответствии с выражением для первой нормы вектора, получаем

\|{\mathbf{Av}}\|}_1 = \max\limits_i \sum\limits_j{|{a_{ij}}|}.

Таким образом, точная верхняя грань в рассмотренном неравенстве достижима и действительно {\|\mathbf{A}\|}_1  = \max \sum\limits_j{|{a_{ij}}|}.

Для третьей нормы (2.2в)

$ {\|\mathbf{A}\|}_3  = \sup\limits_\mathbf{u}\frac{{\|\mathbf{Au}\|}_3}{{\|\mathbf{u}\|}_3}= \sup\limits_\mathbf{u} \sqrt{\frac{(\mathbf{Au,Au})}{(\mathbf{u,u})}}= \\ 
 = \sup\limits_\mathbf{u}\sqrt{\frac{(\mathbf{A}^*\mathbf{Au,u})}{(\mathbf{u,u})}} $

Заметим, что матрица \mathbf{B}= \mathbf{A}^*\mathbf{A} — симметричная. Без ограничения общности предположим, что все собственные числа матрицы различны. Матрица обладает всеми действительными собственными значениями, и каждому собственному числу соответствует собственный вектор. Все собственные векторы взаимно ортогональны. Можно рассмотреть ортонормированную систему собственных векторов \omega _{1}, \dots  , \omega _{n}; \lambda _{1}, \dots  , \lambda _{n} — соответствующие им собственные значения. Любой вектор \mathbf{u} можно представить в виде своего разложения по базису из собственных векторов: \sum\limits_i{\xi_i\omega_i}. Кроме того, (\mathbf{A}^*\mathbf{A})\omega_i  = \lambda_i \omega_i. Поэтому

$ \sup\limits_{\|\mathbf{u}\| \ne 0}\sqrt{\frac{(\mathbf{A}^* \mathbf{Au,Au})}{(\mathbf{u,u})}} = \sup\limits_{\|\mathbf{u}\| \ne 0}\sqrt{\sum\limits_i
 \frac{(\lambda_i\xi_i \omega_i, \xi_i \omega_i)}{(\xi_i\omega_i, \xi_i\omega_i)}} = \\ 
= \sup\limits_{\|\mathbf{u}\| \ne 0}\sqrt{\frac{\sum \lambda_i{\xi}^2_i}{\sum {\xi}^2_i}} 
= \sqrt{\max \lambda_i (\mathbf{A}^* \mathbf{A})}, $

причем точная верхняя грань достигается при \mathbf{u}= \omega_i. Действительно,

$  \sup\limits_u \sqrt{\frac{(\mathbf{A}^* \mathbf{A}\omega_i, \omega_i)}{(\omega_i, \omega_i)}} = \sup\limits_u \sqrt{\lambda^i (\mathbf{A}^*\mathbf{A})} = \sqrt{\max\limits_i {\lambda_i (\mathbf{A}^* \mathbf{A})}},
$

т.к. \mathbf{A}^*\mathbf{A}\omega_i  = \lambda_i \omega_i, откуда (\mathbf{A}^*  \mathbf{A}\omega_i, \omega_i) = \lambda_i (\omega_i, \omega_i),

$  \frac{(\mathbf{A}^*\mathbf{A}\omega_i, \omega_i)}{(\omega_i, 
\omega_i)} = \lambda_i. $

В важном частном случае симметричной (самосопряженной) матрицы \mathbf{A} имеем \lambda_{\mathbf{A}^* \mathbf{A}}^i  = \lambda_{{\mathbf{A}}^2}^i  = {|{\lambda_{\mathbf{A}}^i}|}^2, поэтому {\|\mathbf{A}\|}_3  = \max\limits_i |{\lambda_{\mathbf{A}}^i}|.

Сергей Постнов
Сергей Постнов
Россия, Москва, МГУ им. М.В. Ломоносова, 2004