Системы линейных уравнений

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

Система линейных уравнений

Рассмотрим систему mm линейных уравнений с nn неизвестными x1,,xnx_1, \dots, x_n над полем KK.

{a1,1x1+a1,2x2++a1,nxn=b1a2,1x1+a2,2x2++a2,nxn=b2am,1x1+am,2x2++am,nxn=bm    Ax=b\begin{cases} a_{1,1} x_1 + a_{1,2} x_2 + \dots + a_{1,n} x_n = b_1\\ a_{2,1} x_1 + a_{2,2} x_2 + \dots + a_{2,n} x_n = b_2\\ \vdots \\ a_{m,1} x_1 + a_{m,2} x_2 + \dots + a_{m,n} x_n = b_m \end{cases} \iff A \cdot x = b

В матричной форме эта система записывается как Ax=bA \cdot x = b, где AKm×nA \in K^{m \times n} — матрица коэффициентов, xKnx \in K^n — вектор неизвестных, bKmb \in K^m — вектор правых частей системы.

Матрица AA задаёт линейное отображение A:KnKmA : K^n \to K^m.

Матрицей системы называется матрица коэффициентов при неизвестных.

A=(a1,1a1,2a1,na2,1a2,2a2,nam,1am,2am,n)A = \begin{pmatrix} a_{1,1} & a_{1,2} & \cdots & a_{1,n} \\ a_{2,1} & a_{2,2} & \cdots & a_{2,n} \\ \vdots & \vdots & \ddots & \vdots \\ a_{m,1} & a_{m,2} & \cdots & a_{m,n} \end{pmatrix}

Расширенной матрицей системы называется матрица, полученная путём приписывания к матрице коэффициентов столбца правых частей.

A=(a1,1a1,2a1,nb1a2,1a2,2a2,nb2am,1am,2am,nbm)A' = \left( \begin{array}{cccc|c} a_{1,1} & a_{1,2} & \cdots & a_{1,n} & b_1 \\ a_{2,1} & a_{2,2} & \cdots & a_{2,n} & b_2 \\ \vdots & \vdots & \ddots & \vdots & \vdots \\ a_{m,1} & a_{m,2} & \cdots & a_{m,n} & b_m \end{array} \right)

Пример. Рассмотрим систему линейных уравнений с двумя неизвестными A=(1234)A = \begin{pmatrix}1 & 2\\3 & 4\end{pmatrix}, x=(x1x2)x = \begin{pmatrix}x_1\\x_2\end{pmatrix}, b=(511)b = \begin{pmatrix}5\\11\end{pmatrix} над полем RR. Тогда система Ax=bA \cdot x = b имеет вид

{x1+2x2=53x1+4x2=11\begin{cases} x_1 + 2x_2 = 5\\ 3x_1 + 4x_2 = 11 \end{cases}

Для решения этой системы можно использовать матричный подход. Мы можем записать её в матричной форме как Ax=bA \cdot x = b, где

A=(1234)x=(x1x2)b=(511)A = \begin{pmatrix} 1 & 2 \\ 3 & 4 \end{pmatrix} \quad x = \begin{pmatrix} x_1 \\ x_2 \end{pmatrix} \quad b = \begin{pmatrix} 5 \\ 11 \end{pmatrix}

Пример. Рассмотрим систему линейных уравнений

{2x1x2+3x3=74x1+5x2x3=2\begin{cases} 2x_1 - x_2 + 3x_3 = 7\\ 4x_1 + 5x_2 - x_3 = -2 \end{cases}

Эта система имеет вид матричного уравнения с коэффициентами из матрицы AA и столбцом правых частей bb

A=(213451)b=(72)A = \begin{pmatrix} 2 & -1 & 3 \\ 4 & 5 & -1 \end{pmatrix} \quad b = \begin{pmatrix} 7 \\ -2 \end{pmatrix}

Тогда расширенная матрица этой системы будет выглядеть так

A=(21374512)A' = \left( \begin{array}{ccc|c} 2 & -1 & 3 & 7 \\ 4 & 5 & -1 & -2 \end{array} \right)

Совместная, несовместная и эквивалентные системы

Рассмотрим систему линейных уравнений Ax=bA \cdot x = b, где AKm×nA \in K^{m \times n} — матрица коэффициентов, xKnx \in K^n — вектор неизвестных, а bKmb \in K^m — вектор правых частей.

Система называется совместной, если существует хотя бы один вектор xKnx \in K^n, при котором выполняется равенство Ax=bA \cdot x = b.

Система называется несовместной, если не существует такого вектора xKnx \in K^n, при котором выполняется равенство Ax=bA \cdot x = b.

Две системы Ax=bA \cdot x = b и Ax=bA' \cdot x = b' с одинаковыми неизвестными xKnx \in K^n называются эквивалентными, если их множества решений совпадают.

Доказательство в общем виде

Элементарные преобразования строк, то есть перестановка строк, умножение строки на ненулевой скаляр и прибавление к строке другой строки, умноженной на скаляр, не меняют множество решений системы Ax=bA \cdot x = b.

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

Будем работать с расширенной матрицей A=(Ab)A' = (A \mid b) системы Ax=bA x = b, где AKm×nA \in K^{m \times n} и bKmb \in K^m.

Напомним три типа элементарных преобразований строк.

  1. Перестановка строк RiRjR_i \leftrightarrow R_j — меняем местами i-ю и j-ю строки.

  2. Умножение строки на ненулевое число RiλRiR_i \to \lambda \cdot R_i, где λK\lambda \in K и λ0\lambda \ne 0.

  3. Прибавление кратной другой строки RiRi+λRjR_i \to R_i + \lambda \cdot R_j означает, что к i-й строке прибавляем j-ю строку, умноженную на число λK\lambda \in K.

Докажем в общем виде, что каждое из этих преобразований не меняет множество решений системы.

Пусть xKnx \in K^n — вектор неизвестных.

1. Перестановка строк. Меняем местами i-ю и j-ю строки матрицы AA' . Это лишь меняет порядок уравнений, но не их содержание.

Вектор xx удовлетворяет всем уравнениям до перестановки тогда и только тогда, когда он удовлетворяет тем же уравнениям после перестановки, поэтому множество решений не меняется.

2. Умножение строки на ненулевое число. Пусть i-я строка соответствует уравнению

ai,1x1++ai,nxn=bia_{i,1} x_1 + \dots + a_{i,n} x_n = b_i

После преобразования RiλRiR_i \to \lambda \cdot R_i уравнение принимает вид

λai,1x1++λai,nxn=λbi\lambda \cdot a_{i,1} x_1 + \dots + \lambda \cdot a_{i,n} x_n = \lambda \cdot b_i

Для λ0\lambda \ne 0 эти два уравнения равносильны достаточно домножить или разделить обе части на λ\lambda. Значит множество векторов xx, удовлетворяющих системе, не меняется.

3. Прибавление кратной другой строки. Пусть i-я и j-я строки задают уравнения

ai,1x1++ai,nxn=biaj,1x1++aj,nxn=bj\begin{aligned} a_{i,1} x_1 + \dots + a_{i,n} x_n &= b_i \\ a_{j,1} x_1 + \dots + a_{j,n} x_n &= b_j \end{aligned}

После преобразования RiRi+λRjR_i \to R_i + \lambda \cdot R_j i-я строка соответствует уравнению

(ai,1+λaj,1)x1++(ai,n+λaj,n)xn=bi+λbj(a_{i,1} + \lambda \cdot a_{j,1}) x_1 + \dots + (a_{i,n} + \lambda \cdot a_{j,n}) x_n = b_i + \lambda \cdot b_j

Если вектор xx удовлетворяет исходным двум уравнениям, то левая часть нового уравнения равна bi+λbjb_i + \lambda \cdot b_j, следовательно новое уравнение тоже выполнено.

Обратно, если xx удовлетворяет новому уравнению и уравнению для строки j, то, вычитая λ\lambda-кратное j-е уравнение из нового, получаем исходное i-е уравнение, значит, оно тоже выполняется.

Следовательно, вектор xx является решением системы до преобразования тогда и только тогда, когда он является решением системы после преобразования.

Мы показали, что каждое элементарное преобразование строк задаёт биекцию между множествами решений исходной и преобразованной систем, поэтому композиция любых таких преобразований тоже не меняет множество решений.

Теперь рассмотрим конкретный пример и проследим все шаги вычислений вручную.

Возьмём систему из четырёх уравнений с тремя неизвестными x1,x2,x3Kx_1, x_2, x_3 \in K

{x1+2x2x3=12x1+4x2x3=3x12x2+2x3=03x1+6x2x3=5\begin{cases} x_1 + 2x_2 - x_3 = 1 \\ 2x_1 + 4x_2 - x_3 = 3 \\ -x_1 - 2x_2 + 2x_3 = 0 \\ 3x_1 + 6x_2 - x_3 = 5 \end{cases}

Ей соответствует расширенная матрица A=(Ab)A' = (A \mid b)

A=(1211241312203615)A' = \left( \begin{array}{ccc|c} 1 & 2 & -1 & 1 \\ 2 & 4 & -1 & 3 \\ -1 & -2 & 2 & 0 \\ 3 & 6 & -1 & 5 \end{array} \right)

Сначала применим преобразования R2R22R1R_2 \to R_2 - 2R_1 и R4R43R1R_4 \to R_4 - 3R_1

(1211241312203615)R2R22R1 , R4R43R1(1211001112200022)\left( \begin{array}{ccc|c} 1 & 2 & -1 & 1 \\ 2 & 4 & -1 & 3 \\ -1 & -2 & 2 & 0 \\ 3 & 6 & -1 & 5 \end{array} \right) \xrightarrow[]{R_2 \to R_2 - 2R_1\ ,\ R_4 \to R_4 - 3R_1} \left( \begin{array}{ccc|c} 1 & 2 & -1 & 1 \\ 0 & 0 & 1 & 1 \\ -1 & -2 & 2 & 0 \\ 0 & 0 & 2 & 2 \end{array} \right)

Это соответствует замене уравнений

  • второе уравнение становится равносильно (2x1+4x2x3)2(x1+2x2x3)=321(2x_1 + 4x_2 - x_3) - 2(x_1 + 2x_2 - x_3) = 3 - 2 \cdot 1, то есть x3=1x_3 = 1

  • четвёртое уравнение становится равносильно (3x1+6x2x3)3(x1+2x2x3)=531(3x_1 + 6x_2 - x_3) - 3(x_1 + 2x_2 - x_3) = 5 - 3 \cdot 1, то есть 2x3=22x_3 = 2, снова x3=1x_3 = 1.

Теперь заменим третью строку на R3R3+R1R_3 \to R_3 + R_1

(1211001112200022)R3R3+R1(1211001100110022)\left( \begin{array}{ccc|c} 1 & 2 & -1 & 1 \\ 0 & 0 & 1 & 1 \\ -1 & -2 & 2 & 0 \\ 0 & 0 & 2 & 2 \end{array} \right) \xrightarrow[]{R_3 \to R_3 + R_1} \left( \begin{array}{ccc|c} 1 & 2 & -1 & 1 \\ 0 & 0 & 1 & 1 \\ 0 & 0 & 1 & 1 \\ 0 & 0 & 2 & 2 \end{array} \right)

Получившееся третье уравнение по-прежнему выражает условие x3=1x_3 = 1.

Наконец, выполним преобразования R3R3R2R_3 \to R_3 - R_2 и R4R42R2R_4 \to R_4 - 2R_2

(1211001100110022)R3R3R2 , R4R42R2(1211001100000000)\left( \begin{array}{ccc|c} 1 & 2 & -1 & 1 \\ 0 & 0 & 1 & 1 \\ 0 & 0 & 1 & 1 \\ 0 & 0 & 2 & 2 \end{array} \right) \xrightarrow[]{R_3 \to R_3 - R_2\ ,\ R_4 \to R_4 - 2R_2} \left( \begin{array}{ccc|c} 1 & 2 & -1 & 1 \\ 0 & 0 & 1 & 1 \\ 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 \end{array} \right)

Эквивалентная система имеет вид

{x1+2x2x3=1x3=10=00=0\begin{cases} x_1 + 2x_2 - x_3 = 1 \\ x_3 = 1 \\ 0 = 0 \\ 0 = 0 \end{cases}

Найдём множество решений явно.

Из второго уравнения x3=1x_3 = 1. Подставляя в первое, получаем x1+2x21=1x_1 + 2x_2 - 1 = 1, то есть x1+2x2=2x_1 + 2x_2 = 2.

Положим x2=tx_2 = t, где tKt \in K — произвольный параметр, тогда

x2=tx3=1x1=22tx_2 = t \quad x_3 = 1 \quad x_1 = 2 - 2t

Множество решений можно записать как

{(22tt1)tK}\left\{ \begin{pmatrix} 2 - 2t \\ t \\ 1 \end{pmatrix} \Bigg| t \in K \right\}

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

(22t)+2t1=1(2 - 2t) + 2t - 1 = 1
2(22t)+4t1=32(2 - 2t) + 4t - 1 = 3
(22t)2t+21=0-(2 - 2t) - 2t + 2 \cdot 1 = 0
3(22t)+6t1=53(2 - 2t) + 6t - 1 = 5

Все равенства выполняются при любом tKt \in K, поэтому множество решений исходной системы совпадает с множеством решений всех систем, полученных из неё элементарными преобразованиями строк, а значит, такие системы действительно эквивалентны.

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

Метод Гаусса

Идея метода

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

В матричной записи это выглядит так. Мы рассматриваем расширенную матрицу A=(Ab)A' = (A \mid b) и с помощью элементарных преобразований строк стараемся получить матрицу в треугольном виде, где хорошо видно, какие компоненты вектора xx можно выразить через какие. Эта процедура по сути и есть метод Гаусса.

Чтобы формализовать структуру, введём понятие треугольной матрицы. Интуитивно это такой вид, когда ненулевые строки «сдвигаются» вправо как ступени лестницы, а под выбранными опорными элементами стоят нули.

Матрица RKm×nR \in K^{m \times n} называется треугольной, если выполняются следующие условия

  1. все нулевые строки расположены ниже всех ненулевых

  2. в каждой ненулевой строке первый слева ненулевой элемент стоит правее, чем первый слева ненулевой элемент предыдущей ненулевой строки

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

Пусть RKm×nR \in K^{m \times n} — треугольная матрица.

Ненулевой элемент матрицы RR, стоящий в ведущей позиции некоторой строки, называется ведущим элементом или пивотом этой строки.

Столбец, в котором стоит ведущий элемент, называется ведущим столбцом, а соответствующая неизвестная системы Ax=bA \cdot x = b называется ведущей или главной переменной.

Неизвестные, которым не соответствует ни один ведущий столбец, называются свободными переменными именно они остаются параметрами в общем решении.

Теперь можно сформулировать идею метода Гаусса более точно. Мы хотим при помощи элементарных преобразований строк превратить расширенную матрицу A=(Ab)A' = (A \mid b) в треугольную матрицу, выделить по ней ведущие переменные, а затем выразить все ведущие переменные через свободные.

Алгоритм метода Гаусса

Опишем метод Гаусса в терминах шагов над расширенной матрицей A=(Ab)A' = (A \mid b)

  1. Выбираем первый по счету ещё не обработанный столбец, соответствующий неизвестным. В этом столбце ищем ненулевой элемент в строках, начиная с текущей и ниже. Если такой элемент есть, переставляем строки так, чтобы он оказался в текущей строке. Этот элемент становится будущим ведущим элементом этой строки.

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

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

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

Выбор главного элемента и численная устойчивость

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

Частичный выбор главного элемента состоит в том, что на каждом шаге мы переставляем строки так, чтобы в текущем столбце в качестве ведущего оказался элемент наибольшего модуля. Это не меняет множество решений, но существенно улучшает устойчивость алгоритма.

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

Чтобы увидеть метод Гаусса в более алгоритмическом виде, запишем его в форме псевдокода. Ниже приведён вариант, который строит частное решение или фиксирует, что система несовместна.

Описанные шаги удобно оформить в виде алгоритма. Ниже приведён псевдокод метода Гаусса, который реализует прямой и обратный ход.

Псевдокод алгоритма

function gaussian_solve(matrix A[m, n], vector b[m]) -> result matrix M[m, n + 1] for i from 0 to m - 1 for j from 0 to n - 1 M[i, j] = A[i, j] M[i, n] = b[i] pivot_row = 0 for pivot_col from 0 to n - 1 pivot_row_candidate = -1 for i from pivot_row to m - 1 if M[i, pivot_col] != 0 pivot_row_candidate = i break if pivot_row_candidate == -1 continue swap_rows(M, pivot_row, pivot_row_candidate) pivot_value = M[pivot_row, pivot_col] for j from pivot_col to n M[pivot_row, j] = M[pivot_row, j] / pivot_value for i from pivot_row + 1 to m - 1 factor = M[i, pivot_col] for j from pivot_col to n M[i, j] = M[i, j] - factor * M[pivot_row, j] pivot_row = pivot_row + 1 if pivot_row == m break for i from 0 to m - 1 all_zero = true for j from 0 to n - 1 if M[i, j] != 0 all_zero = false break if all_zero == true and M[i, n] != 0 return ("несовместна", нет_решения) vector x[n] for j from 0 to n - 1 x[j] = 0 for i from m - 1 down to 0 leading_col = -1 for j from 0 to n - 1 if M[i, j] != 0 leading_col = j break if leading_col == -1 continue sum_known = 0 for j from leading_col + 1 to n - 1 sum_known = sum_known + M[i, j] * x[j] x[leading_col] = M[i, n] - sum_known return ("совместна", x) function swap_rows(matrix M[m, k], integer i, integer j) if i == j return for col from 0 to k - 1 temp = M[i, col] M[i, col] = M[j, col] M[j, col] = temp

Метод Гаусса даёт нам универсальный алгоритм решения систем линейных уравнений, по расширенной матрице (Ab)(A \mid b) мы можем выяснить, совместна ли система, и, если да, найти хотя бы одно её решение. Однако сам по себе алгоритм не объясняет, как устроено множество всех решений и почему в одних случаях решение единственно, в других их бесконечно много, а иногда не существует вовсе.

Чтобы понять общую структуру решений, удобно смотреть на систему Ax=bA x = b через призму линейного оператора A:KnKmA : K^n \to K^m и особенно выделить важный частный случай b=0b = 0 — однородные системы.

Однородные и неоднородные системы

Однородная и неоднородная система

Пусть задан линейный оператор A:KnKmA : K^n \to K^m, соответствующий матрице AKm×nA \in K^{m \times n}, и вектор bKmb \in K^m . Рассмотрим систему Ax=bA \cdot x = b с неизвестным вектором xKnx \in K^n.

Система называется однородной, если вектор правых частей равен нулю b=0b = 0, то есть если она имеет вид Ax=0A \cdot x = 0.

Если b0b \ne 0 и система записывается как Ax=bA \cdot x = b, то такая система называется неоднородной.

Линейный оператор AA задаёт отображение A:KnKmA : K^n \to K^m. Множество решений однородной системы Ax=0A \cdot x = 0 совпадает с ядром оператора AA

kerA={xKnAx=0}\ker A = \{x \in K^n \mid A \cdot x = 0\}

Так как AA линейно, то для любых x1,x2kerAx_1, x_2 \in \ker A и любых чисел α,βK\alpha, \beta \in K имеем A(αx1+βx2)=αAx1+βAx2=0A(\alpha \cdot x_1 + \beta \cdot x_2) = \alpha \cdot A x_1 + \beta \cdot A x_2 = 0. Следовательно, kerA\ker A является линейным подпространством в KnK^n.

Пусть теперь система Ax=bA \cdot x = b неоднородна и совместна, то есть существует некоторое x0Knx_0 \in K^n такое, что Ax0=bA \cdot x_0 = b.

Тогда вектор xKnx \in K^n является решением неоднородной системы тогда и только тогда, когда разность xx0x - x_0 лежит в ядре AA

Ax=b    A(xx0)=0A \cdot x = b \iff A \cdot (x - x_0) = 0

Поэтому множество всех решений имеет вид

{xKnAx=b}=x0+kerA={x0+yykerA}\{x \in K^n \mid A \cdot x = b\} = x_0 + \ker A = \{x_0 + y \mid y \in \ker A\}

То есть решения неоднородной системы получаются как сдвиг подпространства kerA\ker A на вектор x0x_0. Геометрически это множество можно рассматривать как прямую, плоскость или более общую «аффинную» подпространственную фигуру, которая в общем случае не проходит через нуль.

Пример. Рассмотрим систему двух уравнений с двумя неизвестными

{x1+x2=12x1+2x2=2\begin{cases} x_1 + x_2 = 1\\ 2x_1 + 2x_2 = 2 \end{cases}

Матрица коэффициентов и столбец правых частей равны

A=(1122)b=(12)A = \begin{pmatrix} 1 & 1\\ 2 & 2 \end{pmatrix} \quad b = \begin{pmatrix} 1\\ 2 \end{pmatrix}

Сначала найдём решения однородной системы Ax=0A \cdot x = 0. Уравнения имеют вид

{x1+x2=02x1+2x2=0\begin{cases} x_1 + x_2 = 0\\ 2x_1 + 2x_2 = 0 \end{cases}

Второе уравнение является удвоенным первым, поэтому достаточно одного условия x1+x2=0x_1 + x_2 = 0. Получаем x2=x1x_2 = -x_1, и любое решение однородной системы имеет вид

x=(tt)tKx = \begin{pmatrix} t\\ -t \end{pmatrix} \quad t \in K

то есть kerA={(t,t)TtK}\ker A = \{(t, -t)^{\T} \mid t \in K\}.

Теперь вернёмся к неоднородной системе Ax=bA \cdot x = b. Из первого уравнения имеем x1+x2=1x_1 + x_2 = 1. Возьмём, например, x1=0x_1 = 0, тогда x2=1x_2 = 1 и

x0=(01)x_0 = \begin{pmatrix} 0\\ 1 \end{pmatrix}

Проверим, что это решение неоднородной системы

Ax0=(1122)(01)=(12)=bA \cdot x_0 = \begin{pmatrix} 1 & 1\\ 2 & 2 \end{pmatrix} \begin{pmatrix} 0\\ 1 \end{pmatrix} = \begin{pmatrix} 1\\ 2 \end{pmatrix} = b

Общее решение неоднородной системы имеет вид x=x0+yx = x_0 + y, где ykerAy \in \ker A. Подставляя найденные выражения, получаем

x=(01)+(tt)=(t1t)tKx = \begin{pmatrix} 0\\ 1 \end{pmatrix} + \begin{pmatrix} t\\ -t \end{pmatrix} = \begin{pmatrix} t\\ 1 - t \end{pmatrix} \quad t \in K

Легко проверить обратную подстановку при любом tKt \in K выполнено x1+x2=1x_1 + x_2 = 1 и 2x1+2x2=22x_1 + 2x_2 = 2, значит множество решений действительно совпадает с x0+kerAx_0 + \ker A.

Обобщим полученный результат. Для системы Ax=bA \cdot x = b множество всех решений либо пусто, либо имеет вид x0+kerAx_0 + \ker A для некоторого решения x0x_0. Отсюда сразу следуют три возможных случая:

  • Если вектор bb не лежит в образе оператора AA, другими словами не существует x0x_0 такого что Ax0=bA \cdot x_0 = b, то система несовместна — решений нет.

  • Если система совместна и kerA={0}\ker A = \{0\}, то множество решений имеет вид x0+{0}={x0}x_0 + \{0\} = \{x_0\}, значит решение единственно.

  • Если система совместна и kerA\ker A содержит нетривиальные векторы, то множество решений x0+kerAx_0 + \ker A является аффинным подпространством, тогда решений бесконечно много.

Мы описали, как устроено множество решений системы Ax=bA x = b. В однородном случае это линейное подпространство kerA\ker A, а в неоднородном — его сдвиг вида x0+kerAx_0 + \ker A. Это даёт хорошую геометрическую картинку и качественное понимание того, почему решение может быть единственным, бесконечным множеством или не существовать вовсе.

Однако при практических вычислениях нас интересует не только структура множества решений, но и эффективные алгоритмы их нахождения. В квадратном случае AKn×nA \in K^{n \times n} особенно удобно «запоминать» прямой ход метода Гаусса в виде разложения матрицы AA на произведение треугольных матриц. Такое представление называется LU-разложением и позволяет быстро решать системы с разными правыми частями, вычислять определитель и обратную матрицу.

LU-разложение

Определение и идея

LU-разложение квадратной матрицы

Пусть AKn×nA \in K^{n \times n} — квадратная матрица. Представление

A=LUA = L \cdot U

называется LU-разложением матрицы AA, если LL — нижнетреугольная матрица с единицами на диагонали, а UU — верхнетреугольная матрица.

Интуитивно LU-разложение фиксирует те же преобразования, которые мы делаем в методе Гаусса. Матрица UU играет роль матрицы после прямого хода, а матрица LL накапливает коэффициенты, которыми мы вычитали строки.

Существование и единственность

Существование и единственность LU-разложения без перестановок

Пусть поле скаляров KK, nNn \in N и AKn×nA \in K^{n \times n} — квадратная матрица. Для каждого k=1,,nk = 1, \dots, n обозначим через AkKk×kA_k \in K^{k \times k} её верхний левый блок размера k×kk \times k.

Предположим, что все главные угловые миноры матрицы AA ненулевые

Ak0для всех k=1,,n\lvert A_k \rvert \ne 0 \quad \text{для всех } k = 1, \dots, n

Тогда существует единственная пара матриц L,UKn×nL, U \in K^{n \times n} такая, что

A=LUA = L U

где LL — нижнетреугольная матрица с единицами на диагонали, а UU — верхнетреугольная матрица.

Существование. Применим к матрице AA прямой ход метода Гаусса без перестановок строк. На шаге kk (1kn11 \le k \le n - 1) в качестве ведущего используется элемент ak,kka^{k}_{k,k} в позиции (k,k)(k, k). Ненулевость главного минора Ak|A_k| означает, что подматрица AkA_k невырождена, а значит ak,kk0a^{k}_{k,k} \ne 0 и шаг можно выполнить без перестановки строк.

На шаге kk для каждой строки i>ki > k вычитаем кратную строку kk с множителем

i,k=ai,kkak,kk\ell_{i,k} = \frac{a^{k}_{i,k}}{a^{k}_{k,k}}

Все эти множители i,k\ell_{i,k} записываем как элементы матрицы LL под диагональю, а результат прямого хода обозначаем через UU. Тогда произведение LUL U воспроизводит последовательность элементарных преобразований, возвращающую исходную матрицу: A=LUA = L U. По построению на диагонали LL стоят единицы, а UU верхнетреугольна.

Единственность. Пусть существуют два разложения

A=L1U1=L2U2A = L_1 U_1 = L_2 U_2

где L1,L2L_1, L_2 — нижнетреугольные матрицы с единичной диагональю, а U1,U2U_1, U_2 — верхнетреугольные матрицы. Тогда

L21L1=U2U11L_2^{-1} L_1 = U_2 U_1^{-1}

Левая часть — произведение нижнетреугольных матриц с единичной диагональю, значит это тоже нижнетреугольная матрица с единицами на диагонали. Правая часть — произведение верхнетреугольных матриц, то есть верхнетреугольная матрица.

Следовательно, матрица L21L1L_2^{-1} L_1 одновременно нижнетреугольна и верхнетреугольна, а значит является диагональной. Так как её диагональные элементы равны единице, получаем L21L1=InL_2^{-1} L_1 = I_n, то есть L1=L2L_1 = L_2. Из равенства A=L1U1=L1U2A = L_1 U_1 = L_1 U_2 следует U1=U2U_1 = U_2. Тем самым единственность LU-разложения доказана.

Теорема о существовании и единственности LU-разложения отвечает на вопрос «когда» и «в каком виде» матрица AA представима как произведение A=LUA = L \cdot U с нижнетреугольными и верхнетреугольными множителями. При этом в доказательстве мы опирались на метод Гаусса, где прямой ход даёт матрицу UU, а коэффициенты, которыми вычитались строки, образуют матрицу LL.

Полезно, однако, уметь строить матрицы LL и UU по явным формулам для их элементов. Такое описание удобно как для теории, так и для реализации алгоритма в виде кода. Из равенства A=LUA = L \cdot U можно поэлементно вывести рекуррентные формулы, позволяющие последовательно вычислять строки UU и столбцы LL — именно их мы сейчас выпишем.

Рекуррентные формулы построения

Пусть A=ai,jA = a_{i,j} и мы ищем матрицы LL и UU вида

L=(1000l2,1100l3,1l3,210ln,1ln,2ln,31)U=(u1,1u1,2u1,3u1,n0u2,2u2,3u2,n00u3,3u3,n000un,n)L = \begin{pmatrix} 1 & 0 & 0 & \cdots & 0 \\ l_{2,1} & 1 & 0 & \cdots & 0 \\ l_{3,1} & l_{3,2} & 1 & \cdots & 0 \\ \vdots & \vdots & \vdots & \ddots & \vdots \\ l_{n,1} & l_{n,2} & l_{n,3} & \cdots & 1 \end{pmatrix} \quad U = \begin{pmatrix} u_{1,1} & u_{1,2} & u_{1,3} & \cdots & u_{1,n} \\ 0 & u_{2,2} & u_{2,3} & \cdots & u_{2,n} \\ 0 & 0 & u_{3,3} & \cdots & u_{3,n} \\ \vdots & \vdots & \vdots & \ddots & \vdots \\ 0 & 0 & 0 & \cdots & u_{n,n} \end{pmatrix}

Из равенства A=LUA = L \cdot U поэлементно получаются формулы для элементов UU и LL

ui,j=ai,jk=1i1li,kuk,j(ij)u_{i,j} = a_{i,j} - \sum\limits_{k=1}^{i-1} l_{i,k} u_{k,j} \quad (i \le j)
li,j=ai,jk=1j1li,kuk,juj,j(i>j)l_{i,j} = \frac{a_{i,j} - \sum\limits_{k=1}^{j-1} l_{i,k} u_{k,j}}{u_{j,j}} \quad (i > j)

Сначала по этим формулам вычисляют элементы первой строки UU и первого столбца LL, затем переходят ко второй строке и второму столбцу и так далее. Условие теоремы обеспечивает uj,j0u_{j,j} \ne 0, поэтому деления на ноль не возникает.

Пример. Пусть A=(2143)A = \begin{pmatrix} 2 & 1\\ 4 & 3 \end{pmatrix} над полем RR. Ищем LL и UU такого вида

L=(10l2,11)U=(u1,1u1,20u2,2)L = \begin{pmatrix} 1 & 0 \\ l_{2,1} & 1 \end{pmatrix} \quad U = \begin{pmatrix} u_{1,1} & u_{1,2} \\ 0 & u_{2,2} \end{pmatrix}

Из равенства A=LUA = L \cdot U поэлементно получаем

u1,1=a1,1=2u1,2=a1,2=1u_{1,1} = a_{1,1} = 2 \qquad u_{1,2} = a_{1,2} = 1
l2,1=a2,1u1,1=42=2l_{2,1} = \frac{a_{2,1}}{u_{1,1}} = \frac{4}{2} = 2
u2,2=a2,2l2,1u1,2=321=1u_{2,2} = a_{2,2} - l_{2,1} u_{1,2} = 3 - 2 \cdot 1 = 1

Таким образом, L=(1021)L = \begin{pmatrix} 1 & 0\\ 2 & 1 \end{pmatrix} и U=(2101)U = \begin{pmatrix} 2 & 1\\ 0 & 1 \end{pmatrix} . Проверим, что действительно A=LUA = L \cdot U

LU=(1021)(2101)=(12+0011+0122+1021+11)=(2143)=AL \cdot U = \begin{pmatrix} 1 & 0\\ 2 & 1 \end{pmatrix} \begin{pmatrix} 2 & 1\\ 0 & 1 \end{pmatrix} = \begin{pmatrix} 1 \cdot 2 + 0 \cdot 0 & 1 \cdot 1 + 0 \cdot 1 \\ 2 \cdot 2 + 1 \cdot 0 & 2 \cdot 1 + 1 \cdot 1 \end{pmatrix} = \begin{pmatrix} 2 & 1 \\ 4 & 3 \end{pmatrix} = A

Решение системы Ax = b через LU-разложение

Пусть матрица AKn×nA \in K^{n \times n} допускает LU-разложение A=LUA = L \cdot U, где LL — нижнетреугольная матрица с единицами на диагонали, а UU — верхнетреугольная матрица с ненулевыми диагональными элементами uk,k0u_{k,k} \ne 0 для всех k=1,,nk = 1, \dots, n. Рассмотрим систему Ax=bA \cdot x = b с неизвестным xKnx \in K^n и правой частью bKnb \in K^n. Подставим разложение

LUx=bL \cdot U \cdot x = b

Обозначим Ux=yU \cdot x = y, где yKny \in K^n. Тогда сначала решаем систему Ly=bL \cdot y = b с нижнетреугольной матрицей, а затем систему Ux=yU \cdot x = y с верхнетреугольной матрицей.

Для нижнетреугольной системы Ly=bL \cdot y = b решения находятся прямой подстановкой сверху вниз. Компоненты yky_k вычисляются по формулам

yk=bkj=1k1lk,jyjk=1,2,,ny_k = b_k - \sum\limits_{j=1}^{k-1} l_{k,j} y_j \quad k = 1, 2, \dots, n

Затем для системы Ux=yU \cdot x = y выполняем обратный ход снизу вверх

xk=ykj   ⁣=   ⁣k+1nuk,j xjuk,kk=n,n1,,1x_k = \dfrac{y_k - \sum\limits_{j \;\! = \;\! k+1}^{n} u_{k,j} \ x_j}{u_{k,k}} \quad k = n, n-1, \dots, 1

Поскольку uk,k0u_{k,k} \ne 0, каждое деление корректно, и так решение системы Ax=bA \cdot x = b разбивается на два простых шага с треугольными матрицами. Когда нужно решить много систем с одной и той же матрицей AA и разными правыми частями bb, LU-разложение позволяет один раз найти LL и UU, а затем быстро пересчитывать решения для новых векторов bb.

Связь с методом Гаусса и определителя

LU-метод по сути совпадает с методом Гаусса без перестановок строк, но записан в матричной форме. Прямой ход метода Гаусса даёт верхнетреугольную матрицу UU и модифицированный столбец правых частей, а коэффициенты, которыми мы вычитали строки, образуют матрицу LL. Тем самым LU-разложение просто «запоминает» сделанные элементарные преобразования.

Если разложение A=LUA = L \cdot U построено, вычисление определителя упрощается. Поскольку диагональные элементы LL равны единице, имеем

A=LU=1i=1nui,i=i=1nui,i\lvert A \rvert = \lvert L \rvert \cdot \lvert U \rvert = 1 \cdot \prod\limits_{i=1}^{n} u_{i,i} = \prod\limits_{i=1}^{n} u_{i,i}

То есть определитель матрицы AA равен произведению диагональных элементов матрицы UU. В случае разложения с перестановками строк PA=LUP \cdot A = L \cdot U знак определителя меняется в зависимости от матрицы перестановки PP.

LU-разложение также можно использовать для вычисления обратной матрицы. Разложив A=LUA = L \cdot U, последовательно решаем системы Axi=eiA \cdot x_i = e_i для стандартных базисных векторов eie_i, i=1,,ni = 1, \dots, n. Для каждого ii сначала решаем Lyi=eiL \cdot y_i = e_i, затем Uxi=yiU \cdot x_i = y_i. Столбцы xix_i образуют матрицу A1A^{-1}.

Таким образом, LU-разложение служит универсальным инструментом один раз «дорого» разложив матрицу AA, мы можем эффективно решать системы Ax=bA \cdot x = b, находить детерминант и обратную матрицу. В следующих разделах мы увидим, как идеи метода Гаусса и LU-разложения обобщаются на переопределённые и недоопределённые системы.

В случае, когда матрица AA допускает разложение A=LUA = L \cdot U, многие важные численные характеристики AA выражаются через треугольные множители. В частности, определитель можно записать как произведение диагональных элементов матрицы UU и по его значению отличать невырождённые матрицы от вырожденных.

Поэтому имеет смысл отдельно ввести понятие определителя, разобрать его основные свойства и геометрический смысл. Это даст компактные критерии обратимости матриц и позволит аккуратно формулировать утверждения вида A=LU\lvert A \rvert = \lvert L \rvert \cdot \lvert U \rvert и эквивалентность условий A0\lvert A \rvert \ne 0 и существования обратной матрицы A1A^{-1}.

Определитель

Определитель квадратной матрицы

Пусть AKn×nA \in K^{n \times n} — квадратная матрица. Числовым отображением, сопоставляющим каждой такой матрице число A\lvert A \rvert, называется определитель, если для любой матрицы AA число A\lvert A \rvert удовлетворяет трём свойствам

  1. Линейность по строке. Определитель линейно зависит от элементов каждой строки при фиксированных остальных строках.

  2. Обращение в ноль на вырожденных матрицах. Если строки матрицы линейно зависимы, то A=0\lvert A \rvert = 0.

  3. Нормировка. Для единичной матрицы InI_n выполняется In=1\lvert I_n \rvert = 1.

Пример. Для матрицы A=(abcd)A = \begin{pmatrix} a & b\\ c & d \end{pmatrix} по известной формуле имеем A=adbc\lvert A \rvert = ad - bc.

Возьмём конкретную матрицу A=(1234)A = \begin{pmatrix} 1 & 2\\ 3 & 4 \end{pmatrix}.

A=1423=46=2\lvert A \rvert = 1 \cdot 4 - 2 \cdot 3 = 4 - 6 = -2

Для матриц порядка два существует простая явная формула abcd=adbc\left\lvert \begin{smallmatrix} a & b\\ c & d \end{smallmatrix} \right\rvert = ad - bc. Для матрицы порядка три определитель часто находят либо по формуле Саррюса, выписывая дополнительные диагонали, либо по правилу разложения по строке или столбцу разлагают A\lvert A \rvert по выбранной строке через миноры и алгебраические дополнения.

Для матриц более высокого порядка n4n \ge 4 прямые формулы становятся громоздкими, поэтому на практике определитель обычно вычисляют, приводя матрицу к треугольному виду методом Гаусса или используя LU-разложение, в треугольном случае A\lvert A \rvert равен произведению диагональных элементов с учётом знака перестановок строк.

Пример. Определитель через LU-разложение. Рассмотрим матрицу порядка три A=(211433234)A = \begin{pmatrix} 2 & 1 & 1 \\ 4 & 3 & 3 \\ 2 & 3 & 4 \end{pmatrix} над полем RR. По рекуррентным формулам из раздела про LU-разложение найдём матрицы LL и UU такие, что A=LUA = L \cdot U.

Сначала вычисляем элементы первой строки UU и первого столбца LL

u1,1=a1,1=2u1,2=a1,2=1u1,3=a1,3=1u_{1,1} = a_{1,1} = 2 \quad u_{1,2} = a_{1,2} = 1 \quad u_{1,3} = a_{1,3} = 1
l2,1=a2,1u1,1=42=2l3,1=a3,1u1,1=22=1l_{2,1} = \frac{a_{2,1}}{u_{1,1}} = \frac{4}{2} = 2 \quad l_{3,1} = \frac{a_{3,1}}{u_{1,1}} = \frac{2}{2} = 1

Переходим ко второй строке и второму столбцу. По формулам для элементов UU и LL

u2,2=a2,2l2,1u1,2=321=1u_{2,2} = a_{2,2} - l_{2,1} u_{1,2} = 3 - 2 \cdot 1 = 1
u2,3=a2,3l2,1u1,3=321=1u_{2,3} = a_{2,3} - l_{2,1} u_{1,3} = 3 - 2 \cdot 1 = 1
l3,2=a3,2l3,1u1,2u2,2=3111=2l_{3,2} = \frac{a_{3,2} - l_{3,1} u_{1,2}}{u_{2,2}} = \frac{3 - 1 \cdot 1}{1} = 2

Наконец, найдём элемент u3,3u_{3,3}

u3,3=a3,3l3,1u1,3l3,2u2,3=41121=1u_{3,3} = a_{3,3} - l_{3,1} u_{1,3} - l_{3,2} u_{2,3} = 4 - 1 \cdot 1 - 2 \cdot 1 = 1

Получаем разложение

L=(100210121)U=(211011001)A=LUL = \begin{pmatrix} 1 & 0 & 0 \\ 2 & 1 & 0 \\ 1 & 2 & 1 \end{pmatrix} \quad U = \begin{pmatrix} 2 & 1 & 1 \\ 0 & 1 & 1 \\ 0 & 0 & 1 \end{pmatrix} \quad A = L \cdot U

Определитель тогда легко находится по диагонали матрицы UU

A=LU=1(211)=2\lvert A \rvert = \lvert L \rvert \cdot \lvert U \rvert = 1 \cdot (2 \cdot 1 \cdot 1) = 2

Здесь L=1\lvert L \rvert = 1, поскольку LL — нижнетреугольная матрица с единицами на диагонали, а U\lvert U \rvert равен произведению диагональных элементов UU. Для матриц большего порядка тот же алгоритм Гаусса и LU-разложения позволяет вычислять определитель столь же систематично и удобно, в отличие от прямого разложения по строке, где число слагаемых быстро растёт.

Геометрический смысл определителя

В пространствах R2\RR^2 и R3\RR^3 у определителя есть наглядная геометрическая интерпретация. Пусть столбцы матрицы AA — координаты векторов v1,,vnv_1, \dots, v_n.

  • В R2\RR^2 число A\lvert A \rvert равно площади параллелограмма, построенного на векторах v1,v2v_1, v_2.

  • В R3\RR^3 A\lvert A \rvert равно объёму параллелепипеда, построенного на v1,v2,v3v_1, v_2, v_3.

  • В общем случае A\lvert A \rvert показывает, во сколько раз линейное отображение AA изменяет объём и меняет ли оно ориентацию знак определителя отвечает за ориентацию, а модуль — за коэффициент изменения объёма.

Если A=0\lvert A \rvert = 0 , то объём параллелепипеда равен нулю векторы v1,,vnv_1, \dots, v_n линейно зависимы и лежат в некотором подпространстве меньшей размерности. Это согласуется с тем, что такая матрица вырождена.

Линейность по строке

Свойство линейности по строке удобно формулировать так. Пусть матрица AA отличается от матрицы A1A_1 только одной строкой и в этой строке стоит вектор αp+βq\alpha p + \beta q вместо вектора pp или qq.

Тогда определитель можно разложить

A=αAp+βAq\lvert A \rvert = \alpha \cdot \lvert A_p \rvert + \beta \cdot \lvert A_q \rvert

где ApA_p и AqA_q получены из AA заменой рассматриваемой строки на pp и qq соответственно.

  • Общий множитель в одной строке можно вынести за знак определителя

    (αa1,1αa1,n)=α(a1,1a1,n)\left\lvert \begin{pmatrix} \alpha \cdot a_{1,1} & \dots & \alpha \cdot a_{1,n}\\ \ast & \dots & \ast\\ \vdots & & \vdots \end{pmatrix} \right\rvert = \alpha \cdot \left\lvert \begin{pmatrix} a_{1,1} & \dots & a_{1,n}\\ \ast & \dots & \ast\\ \vdots & & \vdots \end{pmatrix} \right\rvert
  • Если строка матрицы является суммой нескольких строк, то определитель равен сумме определителей матриц, в которых эта строка по очереди заменяется каждым слагаемым.

Пример. Рассмотрим матрицу C=(2431)C = \begin{pmatrix} 2 & 4\\ 3 & 1 \end{pmatrix} и заметим, что первая строка представима как 2(1,2)2 (1, 2).

C=(2431)=2(1231)\lvert C \rvert = \left\lvert \begin{pmatrix} 2 & 4\\ 3 & 1 \end{pmatrix} \right\rvert = 2 \left\lvert \begin{pmatrix} 1 & 2\\ 3 & 1 \end{pmatrix} \right\rvert
(1231)=1123=5\left\lvert \begin{pmatrix} 1 & 2\\ 3 & 1 \end{pmatrix} \right\rvert = 1 \cdot 1 - 2 \cdot 3 = -5
C=2(5)=10\lvert C \rvert = 2 \cdot (-5) = -10

Сложение строки с другой строкой

Прибавление к строке другой строки

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

Пусть в матрице AA заменили ii-ю строку на ai+λaja_i + \lambda \cdot a_j, где iji \ne j, и получили матрицу AA'.

По линейности по строке

A=A+λAj\lvert A' \rvert = \lvert A \rvert + \lambda \cdot \lvert A_j \rvert

где AjA_j — матрица, в которой ii-я строка заменена на aja_j . В матрице AjA_j строки aja_j и aja_j совпадают, значит они линейно зависимы и Aj=0\lvert A_j \rvert = 0. Следовательно, A=A\lvert A' \rvert = \lvert A \rvert.

Пример. Пусть A=(1234)A = \begin{pmatrix} 1 & 2\\ 3 & 4 \end{pmatrix}, и из второй строки вычитаем первую умноженную на 33

A=(1234)R2R23R1(1202)A' = \begin{pmatrix} 1 & 2\\ 3 & 4 \end{pmatrix} \xrightarrow[]{R_2 \to R_2 - 3 R_1} \begin{pmatrix} 1 & 2\\ 0 & -2 \end{pmatrix}

Вычислим оба определителя

A=1423=2\lvert A \rvert = 1 \cdot 4 - 2 \cdot 3 = -2
A=1(2)20=2\lvert A' \rvert = 1 \cdot (-2) - 2 \cdot 0 = -2

В соответствии с теоремой значение определителя не изменилось.

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

Мы увидели, что определитель отражает сразу несколько важных свойств матрицы: по нему можно судить о вырожденности и удобно вычислять его с помощью метода Гаусса или LU-разложения. Естественно возникает вопрос, как это связано с возможностью «обратить» линейное преобразование, то есть восстановить входной вектор по его образу.

Оказывается, определитель как раз и даёт критерий обратимости: квадратная матрица обратима тогда и только тогда, когда её определитель отличен от нуля. В следующем разделе мы сформулируем это более строго, введём понятие обратной матрицы и покажем, как на практике находить A1A^{-1} с помощью метода Гаусса.

Обратная матрица

Обратная матрица

Пусть AKn×nA \in K^{n \times n} — квадратная матрица над полем KK. Матрица A1Kn×nA^{-1} \in K^{n \times n} называется обратной к AA, если выполняется

AA1=A1A=In=1A A^{-1} = A^{-1} A = I_n = 1

где InI_n — единичная матрица порядка nn. Матрица AA называется невырождённой, если обратная матрица существует.

Определитель даёт удобный критерий невырожденности, а именно квадратная матрица обратима тогда и только тогда, когда её определитель отличен от нуля.

Критерий обратимости через определитель

Пусть AKn×nA \in K^{n \times n}. Тогда эквивалентны следующие утверждения

  1. Матрица AA невырождена, то есть существует A1A^{-1} такая, что AA1=A1A=In=1A A^{-1} = A^{-1} A = I_n = 1

  2. Определитель матрицы AA ненулевой A0\lvert A \rvert \ne 0

Предположим, что матрица AA обратима. Тогда отображение A:KnKnA : K^n \to K^n является биекцией, для любого bKnb \in K^n уравнение Ax=bA x = b имеет единственное решение x=A1bx = A^{-1} b.

Если бы строки матрицы AA были линейно зависимы, то существовал бы ненулевой вектор x0x \ne 0 и Ax=0A x = 0, то есть однородная система имела бы нетривиальное решение. Это противоречит биективности. Значит строки линейно независимы, матрица невырождена, и по свойству определителя получаем A0\lvert A \rvert \ne 0.

Пусть теперь A0\lvert A \rvert \ne 0. По определению определителя это означает, что строки матрицы AA линейно независимы. Аналогично, линейно независимы и столбцы матрицы AA определитель не меняется при транспонировании матрицы.

Рассмотрим линейный оператор A:KnKnA : K^n \to K^n, действующий на столбцы-веторы. Его столбцы a1,,ana_1, \dots, a_n — это образы стандартного базиса e1,,ene_1, \dots, e_n: Aei=aiA e_i = a_i.

Так как столбцы a1,,ana_1, \dots, a_n линейно независимы и их ровно nn, они образуют базис в пространстве KnK^n. Значит, любой вектор yKny \in K^n можно единственным образом представить в виде линейной комбинации этих столбцов

y=x1a1++xnany = x_1 a_1 + \dots + x_n a_n

Для вектора коэффициентов x=(x1,,xn)Tx = (x_1, \dots, x_n)^{\T} это равенство равносильно уравнению Ax=yA x = y. Таким образом, для любого yy существует единственный вектор xx такой, что Ax=yA x = y.

Итак, оператор AA биективен, он он взаимно однозначно отображает KnK^n на KnK^n. Для биективного линейного оператора существует обратный оператор A1:KnKnA^{-1} : K^n \to K^n, удовлетворяющий A1A=AA1=In=1A^{-1} A = A A^{-1} = I_n = 1. Матрица этого оператора в стандартном базисе и есть обратная матрица A1A^{-1}.

Таким образом, условие A0\lvert A \rvert \ne 0 полностью эквивалентно существованию обратной матрицы и даёт удобный числовой критерий обратимости.

Пример. Пусть A=(1234)A = \begin{pmatrix} 1 & 2\\ 3 & 4 \end{pmatrix} над полем RR.

Сначала найдём определитель

A=1423=46=2\lvert A \rvert = 1 \cdot 4 - 2 \cdot 3 = 4 - 6 = -2

Так как A0\lvert A \rvert \ne 0, матрица обратима. Для матрицы размера 2×22 \times 2 обратную можно найти по явной формуле

A1=1A(dbca)для A=(abcd)A^{-1} = \frac{1}{\lvert A \rvert} \begin{pmatrix} d & -b\\ -c & a \end{pmatrix} \quad \text{для } A = \begin{pmatrix} a & b\\ c & d \end{pmatrix}

В нашем случае a=1, b=2, c=3, d=4a = 1,\ b = 2,\ c = 3,\ d = 4, поэтому

A1=12(4231)=12(4231)A^{-1} = \frac{1}{-2} \begin{pmatrix} 4 & -2\\ -3 & 1 \end{pmatrix} = -\frac{1}{2} \begin{pmatrix} 4 & -2\\ -3 & 1 \end{pmatrix}

Проверим одно из равенств AA1=I2A A^{-1} = I_2

AA1=(1234)(12(4231))=12(462+212126+4)A A^{-1} = \begin{pmatrix} 1 & 2\\ 3 & 4 \end{pmatrix} \cdot \left( -\frac{1}{2} \begin{pmatrix} 4 & -2\\ -3 & 1 \end{pmatrix} \right) = -\frac{1}{2} \begin{pmatrix} 4 - 6 & -2 + 2\\ 12 - 12 & -6 + 4 \end{pmatrix}
AA1=12(2002)=(1001)=I2A A^{-1} = -\frac{1}{2} \begin{pmatrix} -2 & 0\\ 0 & -2 \end{pmatrix} = \begin{pmatrix} 1 & 0\\ 0 & 1 \end{pmatrix} = I_2

Значит найденная матрица действительно является обратной к AA.

Теорема даёт критерий существования обратной матрицы, но на практике нужно уметь её вычислять. Один из стандартных способов основан на тех же элементарных преобразованиях, что и метод Гаусса.

Метод Жордана–Гаусса для вычисления обратной матрицы

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

TmT2T1A=In=1T_m \dots T_2 T_1 A = I_n = 1

то

A1=TmT2T1A^{-1} = T_m \dots T_2 T_1

На практике удобно рассматривать расширенную матрицу

D=(AIn)D = (A \mid I_n)

и выполнять элементарные преобразования строк одновременно над обеими половинами, пока левая часть не превратится в единичную матрицу. Тогда правая часть становится A1A^{-1}

(AIn)(InA1)(A \mid I_n) \sim (I_n \mid A^{-1})

Пример. Найдём обратную матрицу для A=(1234)A = \begin{pmatrix} 1 & 2\\ 3 & 4 \end{pmatrix} с помощью метода Гаусса.

Рассмотрим расширенную матрицу

(AI2)=(12103401)(A \mid I_2) = \begin{pmatrix} 1 & 2 & 1 & 0\\ 3 & 4 & 0 & 1 \end{pmatrix}

Вычтем из второй строки первую, умноженную на 33

(12103401)R2R23R1(12100231)\begin{pmatrix} 1 & 2 & 1 & 0\\ 3 & 4 & 0 & 1 \end{pmatrix} \xrightarrow[]{R_2 \to R_2 - 3 R_1} \begin{pmatrix} 1 & 2 & 1 & 0\\ 0 & -2 & -3 & 1 \end{pmatrix}

Теперь умножим вторую строку на 12-\tfrac{1}{2}, чтобы получить единицу на диагонали

(12100231)R212R2(1210013212)\begin{pmatrix} 1 & 2 & 1 & 0\\ 0 & -2 & -3 & 1 \end{pmatrix} \xrightarrow[]{R_2 \to -\tfrac{1}{2} R_2} \begin{pmatrix} 1 & 2 & 1 & 0\\ 0 & 1 & \tfrac{3}{2} & -\tfrac{1}{2} \end{pmatrix}

Вычтем из первой строки вторую, умноженную на 22, чтобы занулить элемент над единицей

(1210013212)R1R12R2(1021013212)\begin{pmatrix} 1 & 2 & 1 & 0\\ 0 & 1 & \tfrac{3}{2} & -\tfrac{1}{2} \end{pmatrix} \xrightarrow[]{R_1 \to R_1 - 2 R_2} \begin{pmatrix} 1 & 0 & -2 & 1\\ 0 & 1 & \tfrac{3}{2} & -\tfrac{1}{2} \end{pmatrix}

Левая часть превратилась в I2I_2, правая часть и есть матрица A1A^{-1}

A1=(213212)A^{-1} = \begin{pmatrix} -2 & 1\\ \tfrac{3}{2} & -\tfrac{1}{2} \end{pmatrix}

Квадратные невырождённые системы

В предыдущем разделе мы показали, что для квадратной матрицы AKn×nA \in K^{n \times n} условия A0\lvert A \rvert \ne 0 и существование обратной матрицы A1A^{-1} эквивалентны. Теперь используем это, чтобы описать поведение системы Ax=bA x = b в невырождённом случае.

Система с квадратной невырождённой матрицей

Пусть AKn×nA \in K^{n \times n} и A0\lvert A \rvert \ne 0. Тогда для любого bKnb \in K^n система Ax=bA \cdot x = b имеет единственное решение.

Из эквивалентности A0\lvert A \rvert \ne 0 и обратимости следует, что существует матрица A1A^{-1} и линейный оператор A:KnKnA : K^n \to K^n биективен.

Умножим систему Ax=bA \cdot x = b слева на A1A^{-1}

A1Ax=A1bA^{-1} A \cdot x = A^{-1} b
Inx=A1bI_n \cdot x = A^{-1} b
x=A1bx = A^{-1} b

Правая часть здесь определена единственным образом, поэтому решение существует и единственно.

Пример. Возьмём матрицу A=(1234)A = \begin{pmatrix} 1 & 2\\ 3 & 4 \end{pmatrix} и вектор b=(511)b = \begin{pmatrix} 5\\ 11 \end{pmatrix}, как в самом начале статьи.

Мы уже нашли A1A^{-1}, поэтому решение можно записать как

x=A1b=12(4231)(511)x = A^{-1} b = -\frac{1}{2} \begin{pmatrix} 4 & -2\\ -3 & 1 \end{pmatrix} \begin{pmatrix} 5\\ 11 \end{pmatrix}
x=12(45+(2)1135+111)=12(202215+11)=12(24)x = -\frac{1}{2} \begin{pmatrix} 4 \cdot 5 + (-2) \cdot 11\\ -3 \cdot 5 + 1 \cdot 11 \end{pmatrix} = -\frac{1}{2} \begin{pmatrix} 20 - 22\\ -15 + 11 \end{pmatrix} = -\frac{1}{2} \begin{pmatrix} -2\\ -4 \end{pmatrix}
x=(12)x = \begin{pmatrix} 1\\ 2 \end{pmatrix}

Подстановка в исходную систему показывает, что это действительно решение.

Мы уже знаем, что для невырождённой матрицы AKn×nA \in K^{n \times n} существует обратная матрица A1A^{-1} и система Ax=bA x = b имеет единственное решение x=A1bx = A^{-1} b. Однако в явном виде элементы A1A^{-1} можно выразить через определители. Из этих формул получается классическое правило Крамера, позволяющее выписывать компоненты решения непосредственно через определители.

Правило Крамера

Пусть AKn×nA \in K^{n \times n} невырождена (A0\lvert A \rvert \neq 0) и система Ax=bA \cdot x = b совместна. Тогда каждую компоненту решения можно вычислить по формуле Крамера

xi=ΔiAΔi=Aix_i = \frac{\Delta_i}{\lvert A \rvert} \qquad \Delta_i = \lvert A_i \rvert

где AiA_i — матрица, полученная заменой ii-го столбца матрицы AA столбцом правых частей bb.

Пример. Для системы

{x1+2x2=53x1+4x2=11\begin{cases} x_1 + 2x_2 = 5\\ 3x_1 + 4x_2 = 11 \end{cases}

матрица коэффициентов и столбец правых частей равны

A=(1234)b=(511)A = \begin{pmatrix} 1 & 2\\ 3 & 4 \end{pmatrix} \quad b = \begin{pmatrix} 5\\ 11 \end{pmatrix}

Найдём определитель

A=1423=2\lvert A \rvert = 1 \cdot 4 - 2 \cdot 3 = -2

Теперь матрицы A1A_1 и A2A_2

A1=(52114)A2=(15311)A_1 = \begin{pmatrix} 5 & 2\\ 11 & 4 \end{pmatrix} \qquad A_2 = \begin{pmatrix} 1 & 5\\ 3 & 11 \end{pmatrix}
Δ1=A1=54211=2022=2\Delta_1 = \lvert A_1 \rvert = 5 \cdot 4 - 2 \cdot 11 = 20 - 22 = -2
Δ2=A2=11153=1115=4\Delta_2 = \lvert A_2 \rvert = 1 \cdot 11 - 5 \cdot 3 = 11 - 15 = -4

По правилу Крамера

x1=Δ1A=22=1x_1 = \frac{\Delta_1}{\lvert A \rvert} = \frac{-2}{-2} = 1
x2=Δ2A=42=2x_2 = \frac{\Delta_2}{\lvert A \rvert} = \frac{-4}{-2} = 2

Получили то же решение, что и через обратную матрицу.

В невырождённом квадратном случае поведение системы полностью понятно — для любого b есть единственное решение. Чтобы описать все остальные случаи, нужна более общая характеристика — ранг матрицы.

Ранг и теорема Кронекера–Капелли

Ранг матрицы

Рангом матрицы AA называется максимальное число линейно независимых строк или столбцов этой матрицы. Обозначение rankA\rank A.

Пример. Рассмотрим матрицу

A=(123246101)A = \begin{pmatrix} 1 & 2 & 3\\ 2 & 4 & 6\\ 1 & 0 & 1 \end{pmatrix}

Вторая строка равна первой, умноженной на 22, поэтому они линейно зависимы. Строки (1,2,3)(1, 2, 3) и (1,0,1)(1, 0, 1) линейно независимы. Значит rankA=2\operatorname{rank} A = 2.

То же можно увидеть методом Гаусса

(123246101)R2R22R1(123000101)R3R3R1(123000022)\begin{pmatrix} 1 & 2 & 3\\ 2 & 4 & 6\\ 1 & 0 & 1 \end{pmatrix} \xrightarrow[]{R_2 \to R_2 - 2 R_1} \begin{pmatrix} 1 & 2 & 3\\ 0 & 0 & 0\\ 1 & 0 & 1 \end{pmatrix} \xrightarrow[]{R_3 \to R_3 - R_1} \begin{pmatrix} 1 & 2 & 3\\ 0 & 0 & 0\\ 0 & -2 & -2 \end{pmatrix}

В треугольном виде видно две ненулевые строки, значит ранг равен 22.

Теорема Кронекера–Капелли

Система Ax=bA \cdot x = b совместна тогда и только тогда, когда

rankA=rankA\rank A = \rank A'

где AA' — расширенная матрица системы.

Критерий по приведённой матрице

Пусть расширенная матрица AA' приведена элементарными преобразованиями строк к треугольному виду.

Тогда система несовместна тогда и только тогда, когда в приведённой матрице появляется строка вида (0,0,,01)(0, 0, \dotsc, 0 \mid 1), то есть все коэффициенты при известных равны нулю, а свободный член нет.

Пример. Рассмотрим системы

{x1+x2=12x1+2x2=2{x1+x2=12x1+2x2=3\begin{cases} x_1 + x_2 = 1\\ 2x_1 + 2x_2 = 2 \end{cases} \qquad \begin{cases} x_1 + x_2 = 1\\ 2x_1 + 2x_2 = 3 \end{cases}

Для первой системы расширенная матрица и её треугольный вид таковы

A1=(111222)R2R22R1(111000)A'_1 = \begin{pmatrix} 1 & 1 & 1\\ 2 & 2 & 2 \end{pmatrix} \xrightarrow[]{R_2 \to R_2 - 2 R_1} \begin{pmatrix} 1 & 1 & 1\\ 0 & 0 & 0 \end{pmatrix}

Строки вида (0,01)(0, 0 \mid 1) нет, значит система совместна.

Для второй системы

A2=(111223)R2R22R1(111001)A'_2 = \begin{pmatrix} 1 & 1 & 1\\ 2 & 2 & 3 \end{pmatrix} \xrightarrow[]{R_2 \to R_2 - 2 R_1} \begin{pmatrix} 1 & 1 & 1\\ 0 & 0 & 1 \end{pmatrix}

Строка (0,01)(0, 0 \mid 1) соответствует уравнению 0=10 = 1, невозможному для любых x1,x2x_1, x_2, поэтому система несовместна.

Теорема Кронекера–Капелли даёт критерий совместности системы любой формы. Теперь посмотрим, как она проявляется в двух важных конфигурациях — когда уравнений больше неизвестных и наоборот.

Переопределённые и недоопределённые системы

Переопределённая и недоопределённая системы

Пусть над некоторым полем KK, где m,nNm, n \in N и AKm×nA \in K^{m \times n} — матрица линейного оператора A:KnKmA : K^n \to K^m. Рассмотрим систему линейных уравнений Ax=bA \cdot x = b с неизвестным xKnx \in K^n и правой частью bKmb \in K^m.

  1. Система Ax=bA \cdot x = b называется переопределённой, если число уравнений больше числа неизвестных, то есть m>nm > n. В таком случае для произвольного bb система, вообще говоря, может не иметь точного решения, и естественно рассматривать задачи поиска приближённого решения.

  2. Система Ax=bA \cdot x = b называется недоопределённой, если число уравнений меньше числа неизвестных, то есть m<nm < n. При совместности такой системы множество решений, как правило, бесконечно, и возникает задача выбрать наиболее простой вектор решения.

В случае m=nm = n систему называют квадратной, и она не относится ни к переопределённым, ни к недоопределённым.

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

Псевдообратная матрица Мура—Пенроуза

Псевдообратная матрица Мура—Пенроуза

Пусть AKm×nA \in K^{m \times n}. Матрица A+Kn×mA^+ \in K^{n \times m} называется псевдообратной в смысле Мура—Пенроуза, если выполняются четыре условия:

AA+A=AA \cdot A^+ \cdot A = A
A+AA+=A+A^+ \cdot A \cdot A^+ = A^+
(AA+)T=AA+(A \cdot A^+)^{\T} = A \cdot A^+
(A+A)T=A+A(A^+ \cdot A)^{\T} = A^+ \cdot A

Здесь T^{\T} обозначает транспонирование (в комплексном случае обычно берут сопряжённое транспонирование). Если такая матрица существует, она единственна.

В частных случаях псевдообратная записывается через обычную обратную матрицу.

  • Если AKm×nA \in K^{m \times n} имеет полный столбцовый ранг (rankA=nm\rank A = n \le m), то

    A+=(ATA)1ATA^+ = (A^{\T} \cdot A)^{-1} \cdot A^{\T}
  • Если AKm×nA \in K^{m \times n} имеет полный строковый ранг (rankA=mn\rank A = m \le n), то

    A+=AT(AAT)1A^+ = A^{\T} \cdot (A \cdot A^{\T})^{-1}

Метод наименьших квадратов

Задача наименьших квадратов

Пусть ARm×nA \in \RR^{m \times n} и bRmb \in \RR^m, обычно m>nm > n (переопределённая система). Задача наименьших квадратов состоит в поиске вектора xRnx' \in \RR^n, минимизирующего норму остатка

x=arg minx   ⁣   ⁣RnAxb2x' = \argmin\limits_{x \;\! \in \;\! \RR^n} \lVert A \cdot x - b \rVert_2

то есть мы подбираем xx так, чтобы AxA \cdot x как можно лучше приближало вектор наблюдений bb в евклидовой норме.

Условие оптимальности можно записать в виде нормальных уравнений.

Нормальные уравнения

Пусть ARm×nA \in \RR^{m \times n}, bRmb \in \RR^m и матрица AA имеет полный столбцовый ранг rankA=n\rank A = n. Тогда задача наименьших квадратов имеет единственное решение xx', которое удовлетворяет системе

ATAx=ATbA^{\T} \cdot A \cdot x' = A^{\T} \cdot b

Рассмотрим функцию f(x)=Axb22f(x) = \lVert A \cdot x - b \rVert_2^2. Раскроем квадрат:

f(x)=(Axb)T(Axb)=xTATAx2xTATb+bTbf(x) = (A \cdot x - b)^{\T} \cdot (A \cdot x - b) = x^{\T} \cdot A^{\T} \cdot A \cdot x - 2 \cdot x^{\T} \cdot A^{\T} \cdot b + b^{\T} \cdot b

Это квадратичная функция по xx. Вектор-градиент равен

f(x)=2ATAx2ATb\nabla f(x) = 2 \cdot A^{\T} \cdot A \cdot x - 2 \cdot A^{\T} \cdot b

В точке минимума градиент равен нулю, поэтому получаем уравнение ATAx=ATbA^{\T} \cdot A \cdot x' = A^{\T} \cdot b. Матрица ATAA^{\T} \cdot A невырождена благодаря полноранговости AA, значит решение единственно.

Геометрическая интерпретация

Рассмотрим линейный оператор A:RnRmA : \RR^n \to \RR^m и его образ ImARm\operatorname{Im} A \subset \RR^m. Вектор AxA \cdot x всегда лежит в ImA\operatorname{Im} A, а bb вообще говоря, не обязан туда попадать.

Решение задачи наименьших квадратов можно описать так: AxA \cdot x' — это ортогональная проекция вектора bb на подпространство ImA\operatorname{Im} A.

Действительно, вектор остатка r=bAxr' = b - A \cdot x' при этом ортогонален всему подпространству ImA\operatorname{Im} A

yImAry=0\forall\, y \in \operatorname{Im} A \quad \langle r' y \rangle = 0

Так как любой yImAy \in \operatorname{Im} A имеет вид y=Azy = A \cdot z, условие ортогональности эквивалентно AT(bAx)=0A^{\T} \cdot (b - A \cdot x') = 0, то есть нормальным уравнениям ATAx=ATbA^{\T} \cdot A \cdot x' = A^{\T} \cdot b.

Связь с псевдообратной матрицей

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

x=A+bx' = A^+ \cdot b

Если rankA=nm\rank A = n \le m, то A+=(ATA)1ATA^+ = (A^{\T} \cdot A)^{-1} \cdot A^{\T} и из нормальных уравнений сразу получаем

x=(ATA)1ATb=A+bx' = (A^{\T} \cdot A)^{-1} \cdot A^{\T} \cdot b = A^+ \cdot b

Это решение единственно и даёт вектор минимальной нормы, если решений задачи наименьших квадратов несколько.

В случае недоопределённой системы Ax=bA \cdot x = b при m<nm < n, когда решений (при совместности) бесконечно много, псевдообратная матрица также выделяет «наиболее короткое» решение:

x=A+bx' = A^+ \cdot b

это единственное решение минимальной нормы x2=min{x2Ax=b}\lVert x' \rVert_2 = \min\limits \{\lVert x \rVert_2 \mid A \cdot x = b\}.

Таким образом, псевдообратная матрица Мура—Пенроуза и метод наименьших квадратов дают единый язык для описания переопределённых и недоопределённых систем: в обоих случаях вектор A+bA^+ \cdot b играет роль «наилучшего» решения — либо в смысле минимизации ошибки Axb\lVert A \cdot x - b \rVert, либо в смысле минимальной нормы среди всех точных решений.

Пример. Пусть

A=(101112),b=(122)A = \begin{pmatrix} 1 & 0\\ 1 & 1\\ 1 & 2 \end{pmatrix}, \quad b = \begin{pmatrix} 1\\ 2\\ 2 \end{pmatrix}

Это переопределённая система в R3\RR^3 с двумя неизвестными. Решим её методом наименьших квадратов.

Составим нормальные уравнения:

ATA=(3335)ATb=(56)A^{\T} \cdot A = \begin{pmatrix} 3 & 3\\ 3 & 5 \end{pmatrix} \quad A^{\T} \cdot b = \begin{pmatrix} 5\\ 6 \end{pmatrix}
ATAx=ATb    {3x1+3x2=53x1+5x2=6A^{\T} \cdot A \cdot x = A^{\T} \cdot b \iff \begin{cases} 3 x_1 + 3 x_2 = 5\\ 3 x_1 + 5 x_2 = 6 \end{cases}

Вычитаем первое уравнение из второго 2x2=12 \cdot x_2 = 1, значит x2=12x_2 = \tfrac{1}{2}. Подставляя обратно, получаем 3x1+312=53 \cdot x_1 + 3 \cdot \tfrac{1}{2} = 5, откуда x1=76x_1 = \tfrac{7}{6}.

x=(7612)x' = \begin{pmatrix} \tfrac{7}{6}\\[2pt] \tfrac{1}{2} \end{pmatrix}

Вектор AxA \cdot x' есть ортогональная проекция bb на \Im} A, а сам xx' совпадает с A+bA^+ \cdot b для данной матрицы AA с полным столбцовым рангом.