Матрицы и системы линейных уравнений: операции, определитель и методы решения

Матрица — это прямоугольная таблица чисел из m строк и n столбцов; над матрицами определены сложение, умножение на число и матричное умножение. Основной прикладной сценарий — через матрицы компактно записать и решить систему линейных уравнений (СЛУ), где несколько уравнений вида a1x1 + a2x2 + … = b собраны в одну запись A·X = B. Ниже разберем операции, которые для этого нужны, определитель и обратную матрицу, а затем три рабочих метода решения (Гаусс, Крамер, матричный) с кодом на Python.

Сразу разведем два близких термина, чтобы примеры не казались противоречивыми. Система линейных уравнений — это набор уравнений относительно неизвестных чисел x1, x2, …, xn; ее решение — набор чисел. Матричное уравнение — это уравнение вида A·X = B, где A, B, X сами являются матрицами, и искать надо неизвестную матрицу X. Классическая СЛУ — частный случай, где X и B являются столбцами.

Операции с матрицами: минимум для решения СЛУ

Чтобы решать системы, достаточно четырех операций. Сложение и вычитание выполняются поэлементно и только для матриц одинакового размера m×n; результат — матрица того же размера. Умножение на число k умножает каждый элемент и допустимо для матрицы любого размера.

Матричное умножение A·B определено, только если число столбцов A равно числу строк B. Элемент результата в строке i и столбце j — это сумма произведений элементов i-й строки A на соответствующие элементы j-го столбца B. Важно не путать его с поэлементным умножением: это разные операции, и в numpy у них разные операторы.

import numpy as np

A = np.array([[1, 2], [3, 4]])
B = np.array([[5, 6], [7, 8]])

print("A + B =\n", A + B)
print("A @ B (матричное) =\n", A @ B)
print("A * B (поэлементно) =\n", A * B)
print("A.T (транспонирование) =\n", A.T)

Вывод:

A + B =
 [[ 6  8]
 [10 12]]
A @ B (матричное) =
 [[19 22]
 [43 50]]
A * B (поэлементно) =
 [[ 5 12]
 [21 32]]
A.T (транспонирование) =
 [[1 3]
 [2 4]]

Оператор @ дает матричное произведение, а * — поэлементное. Это частая ошибка новичков: A * B для квадратных матриц не падает с ошибкой, а молча возвращает не то, что нужно. Транспонирование A.T меняет строки и столбцы местами и обозначается на письме как A с верхним индексом T.

Определитель и обратная матрица

Определитель (детерминант, обозначается |A| или det A) — это число, которое можно вычислить только для квадратной матрицы. Для матрицы 2×2 он равен a11·a22 — a12·a21. Для 3×3 работает правило треугольников (правило Саррюса), а для больших размеров определитель на практике считают не по формуле, а численно — например, через ту же процедуру исключения, что и в методе Гаусса.

Определитель важен по одной причине: он показывает, есть ли у системы единственное решение. Если |A| не равен нулю, матрица называется невырожденной, у нее существует обратная матрица A^(-1), и система A·X = B имеет ровно одно решение. Если |A| = 0, матрица вырождена, обратной у нее нет, а у системы либо нет решений, либо их бесконечно много.

Обратная матрица A^(-1) — это такая матрица, что A·A^(-1) = A^(-1)·A = E, где E — единичная матрица (единицы на главной диагонали, нули вне ее). Это матричный аналог обратного числа: для числа a обратным является 1/a, и a·(1/a) = 1. Условие существования тоже аналогично: как 1/a существует при a не равном нулю, так A^(-1) существует при |A| не равном нулю.

import numpy as np

A = np.array([[2, 1, -1],
              [-3, -1, 2],
              [-2, 1, 2]], dtype=float)

print("det A =", round(np.linalg.det(A), 6))
print("A^-1 =\n", np.linalg.inv(A))

Вывод:

det A = -1.0
A^-1 =
 [[ 4.  3. -1.]
 [-2. -2.  1.]
 [ 5.  4. -1.]]

Определитель равен -1 (не нулю), значит обратная матрица существует и система с этой A будет решаться однозначно.

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

Возьмем сквозной пример через всю статью — систему из трех уравнений:

2x +  y -  z =  8
-3x -  y + 2z = -11
-2x +  y + 2z = -3

В матричной форме это A·X = B, где A — матрица коэффициентов (та самая, что выше), X — столбец неизвестных (x, y, z), B — столбец свободных членов (8, -11, -3). Ниже три способа найти X. Какой выбрать — зависит от задачи; сравнение в таблице.

Метод Гаусса

Метод Гаусса — это последовательное исключение неизвестных: элементарными преобразованиями строк систему приводят к треугольному виду (прямой ход), а затем находят неизвестные снизу вверх (обратная подстановка). Это основной численный метод; именно он лежит в основе numpy.linalg.solve. Реализуем его вручную, с частичным выбором ведущего элемента — это защищает от деления на очень малое число.

def gauss(A, b):
    n = len(A)
    # расширенная матрица: коэффициенты + столбец свободных членов
    M = [row[:] + [b[i]] for i, row in enumerate(A)]
    for col in range(n):
        # частичный выбор: наибольший по модулю элемент столбца в ведущие
        pivot = max(range(col, n), key=lambda r: abs(M[r][col]))
        if abs(M[pivot][col]) < 1e-12:
            raise ValueError("матрица вырождена: единственного решения нет")
        M[col], M[pivot] = M[pivot], M[col]
        # обнуляем столбец ниже ведущего элемента
        for r in range(col + 1, n):
            factor = M[r][col] / M[col][col]
            for c in range(col, n + 1):
                M[r][c] -= factor * M[col][c]
    # обратная подстановка снизу вверх
    x = [0.0] * n
    for i in range(n - 1, -1, -1):
        s = sum(M[i][j] * x[j] for j in range(i + 1, n))
        x[i] = (M[i][n] - s) / M[i][i]
    return x

A = [[2, 1, -1], [-3, -1, 2], [-2, 1, 2]]
b = [8, -11, -3]
print(gauss(A, b))

Вывод (числа с плавающей точкой, поэтому вместо ровно 3 и -1 появляются близкие значения):

[2.0, 3.0000000000000004, -0.9999999999999999]

То есть x = 2, y = 3, z = -1. Ветка с ValueError не даст выдать мусор за ответ, если матрица окажется вырожденной: тогда ведущий элемент близок к нулю, и функция честно сообщает, что единственного решения нет.

Матричный метод

Если A квадратная и невырожденная, обе части уравнения A·X = B умножают слева на A^(-1). Получается A^(-1)·A·X = A^(-1)·B, а поскольку A^(-1)·A = E и E·X = X, ответ равен X = A^(-1)·B. Порядок умножения важен: для уравнения X·A = B умножать надо справа, и ответ будет X = B·A^(-1); для A·X·B = C ответ X = A^(-1)·C·B^(-1).

import numpy as np

A = np.array([[2, 1, -1], [-3, -1, 2], [-2, 1, 2]], dtype=float)
B = np.array([8, -11, -3], dtype=float)

X = np.linalg.inv(A) @ B
print("X =", X)
# на практике solve точнее и быстрее, чем явная inv:
print("solve =", np.linalg.solve(A, B))

Вывод:

X = [ 2.  3. -1.]
solve = [ 2.  3. -1.]

На практике для решения системы предпочтителен np.linalg.solve(A, B), а не явное вычисление обратной матрицы: он точнее численно и быстрее. Явная A^(-1) оправдана, когда обратная матрица нужна сама по себе или систему решают многократно с одной A и разными B.

Метод Крамера

Метод Крамера выражает каждое неизвестное через определители: xi = det(Ai) / det(A), где Ai получается из A заменой i-го столбца на столбец свободных членов B. Метод нагляден для вывода формул и удобен для систем 2×2 и 3×3, но для больших размеров невыгоден: вычислять много определителей дорого, поэтому в реальном коде берут Гаусса.

import numpy as np

A = np.array([[2, 1, -1], [-3, -1, 2], [-2, 1, 2]], dtype=float)
B = np.array([8, -11, -3], dtype=float)

det = np.linalg.det(A)
x = []
for i in range(A.shape[1]):
    Ai = A.copy()
    Ai[:, i] = B          # заменяем i-й столбец на свободные члены
    x.append(np.linalg.det(Ai) / det)
print([round(v, 6) for v in x])

Вывод:

[2.0, 3.0, -1.0]

Все три метода дали один ответ (2, 3, -1) — это нормальная взаимная проверка. Если бы det был равен нулю, деление в методе Крамера потеряло бы смысл, и это тот же сигнал вырожденности, что и в методах выше.

Что делать, если система не решается

Единственное решение есть не всегда. Типовые симптомы и что это значит:

• LinAlgError: Singular matrix от np.linalg.solve или det, равный нулю. Матрица вырождена: строки/столбцы линейно зависимы. Решений либо нет, либо бесконечно много — однозначного ответа не существует.
• Матрица не квадратная (уравнений и неизвестных разное число). Обратной матрицы и определителя для нее нет; тут нужен ранг матрицы и метод Гаусса на прямоугольной системе, а для приближенного решения переопределенной системы — метод наименьших квадратов (np.linalg.lstsq).
• Числа-ответы выглядят как огромные или почти совпадают между строками. Матрица близка к вырожденной (плохо обусловлена); мелкая погрешность входных данных сильно искажает ответ.

Где применяется / связь с практикой

Матрицы и СЛУ — это не только упражнение из линейной алгебры. В компьютерной графике и игровых движках повороты, масштабирование и перенос объектов задаются умножением координат на матрицы преобразований (обычно 3×3 и 4×4 с однородными координатами). В машинном обучении данные, веса и градиенты хранятся как матрицы/тензоры, а обучение — это по сути потоки матричных умножений; линейная регрессия сводится к решению системы. В инженерных и физических расчетах, экономических моделях и оптимизации СЛУ возникают постоянно, а численные методы вроде Гаусса — это то, что реально исполняется под капотом библиотек.

Выводы

• Матрица — таблица чисел m×n; операции сложения и вычитания требуют равного размера, матричное умножение — совпадения «столбцы первой = строки второй»; в numpy матричное произведение это @, а * — поэлементное.
• Определитель считается только для квадратной матрицы; при |A| не равном нулю существует обратная матрица и система A·X = B имеет единственное решение, при |A| = 0 матрица вырождена.
• Три метода решения СЛУ: Гаусс (универсальный, основа библиотек), матричный (X = A^(-1)·B), Крамер (через определители, для малых систем).
• В коде для решения системы берут np.linalg.solve(A, B) — он точнее и быстрее явной обратной матрицы.
• Если det = 0 или матрица не квадратная, единственного решения нет; нужен ранг, метод Гаусса на прямоугольной системе или наименьшие квадраты.

FAQ

Чем матричное умножение отличается от поэлементного в numpy?
A @ B (или np.dot) дает матричное произведение по правилу «строка на столбец», A * B умножает элементы на одинаковых позициях. Для квадратных матриц оба не выдают ошибку, но результаты разные, поэтому оператор надо выбирать осознанно.

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

Можно ли найти обратную матрицу у неквадратной матрицы?
Обычной обратной — нет, она определена только для квадратных невырожденных матриц. Для прямоугольных используют псевдообратную матрицу (np.linalg.pinv), например при решении переопределенных систем методом наименьших квадратов.