Линейные уравнения

Некоторые линейные функции обратимы: например, поворот на угол. Другие — необратимы: например, проекция. Понятие обратимости можно продолжить и на матрицы.

Определение. Матрица $A$ является обратимой, если существует матрица $A^{-1}$ такая, что

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

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

Понятие обратимости матриц важно для решения линейных уравнений:

В противном случае будет либо ноль, либо бесконечное количество решений, в общем случае живущих на каком-то многомерном пространстве. Об этом удобно думать так: каждое условие (строчка-уравнение) задает какую-то гиперплоскость в $n$-мерном пространстве, и пересечения этих гиперплоскостей либо пустые, либо совпадают, либо дают какое-то подпространство меньшей размерности.

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

#Метод Гаусса

Метод Гаусса (англ. Gaussian elimination) позволяет решать системы линейных уравнений, а также находить обратные к матрицам и решать смежные задачи.

Он основывается на том факте, что с системами уравнений мы можем свободно делать следующие операции:

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

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

Цель алгоритма — используя эти три операции привести систему к виду

и тогда решением будет просто $x = b’$.

#Алгоритм

Будем действовать в предположении, что решение существует и единственно.

Припишем вектор $b$ справа от матрицы $A$, и будем итеративно приводить получившуюся $n \times (n+1)$ матрицу к требуемому виду. На шаге $i$, мы хотим сделать так, чтобы в ячейке $a_{ii}$ стояла единица, а во всех остальных ячейках $i$-того столбца стояли нули.

Чтобы выполнить $i$-тый шаг:

1. найдем какой-нибудь ненулевой элемент на $i$-том столбце и поменяем его строку местами с $i$-той;
2. поделим всю $i$-тую строку на элемент $a_{ii}$ (он ненулевой: мы его специально искали на предыдущем шаге);
3. пройдемся по всем остальным строкам $j \ne i$ и вычтем $i$-тую строку, помноженную на коэффициент $a_{ji}$.

В результате после $i$-того шага элемент $a_{ii}$ будет единичным (потому что мы его поделили на самого себя), а во всех остальных ячейках $i$-того столбца стали нули (потому что мы в этих позициях вычли $a_{ij} \cdot a_{ii} = a_{ij}$).

Таким образом, через $n$ шагов мы приведем матрицу к нужному виду, и ответом будет последний, $(n+1)$-ый приписанный столбец.

#Реализация

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

typedef vector< vector<double> > matrix;

vector<double> gauss(matrix a) {
    int n = (int) a.size();
    for (int i = 0; i < n; i++) {
        // находим опорный элемент
        int best = i;
        for (int j = i + 1; j < n; j++)
            if (abs(a[j][i]) > abs(a[best][i]))
                best = j;
        // меняем строчки местами
        swap(a[best], a[i]);
        // нормализуем i-тую строчку (слева нули, с ними ничего делать не надо)
        for (int j = n; j >= i; j--)
            a[i][j] /= a[i][i];
        // зануляем все остальные ячейки в i-том столбце
        for (int j = 0; j < n; j++)
            if (j != i)
                // слева только нули, так что можно пройтись только по правой части
                for (int k = n; k >= i; k--)
                    a[j][k] -= a[i][k] * a[j][i];
    }
    vector<double> ans(n);
    for (int i = 0; i < n; i++)
        ans[i] = a[i][n];
    return ans;
}

Алгоритм работает $O(n^3)$.

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

#Булевы матрицы

В контексте олимпиад, в большинстве задач, где требуется решать системы линейных уравнений, на самом деле это нужно делать над полем $\mathbb{Z}_2$, то есть по модулю 2.

Задача. Есть $n$ лампочек и $n$ переключателей: каждый активированный переключатель меняет состояние (включает или выключает) какого-то подмножества лампочек. Известно состояние всех лампочек, нужно восстановить состояние переключателей.

Нас по сути просят решить следующую систему:

Здесь $x$ — состояния переключателей, $b$ — состояния лампочек, $A$ — информация о том, влияет ли переключатель на лампочку.

В таком случае можно значительно ускорить и упростить обычный метод Гаусса через битсет:

typedef bitset<maxn> t;
typedef array<t, maxn> matrix;

t gauss(matrix a) {
    for (int i = 0; i < n; i++) {
        int nonzero = i;
        for (int j = i + 1; j < n; j++)
            if (a[j][i])
                nonzero = j;
        swap(a[nonzero], a[i]);
        for (int j = 0; j < n; j++)
            if (j != i && a[j][i])
                a[j] ^= a[i];
    }
    t x;
    for (int i = 0; i < n; i++)
        x[i] = a[i][n] ^ a[i][i];
    return x;
}

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

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

В случае входных данных из «реального мира» можно добавить небольшой шум $\epsilon_{ij}$ к коэффициентам: это гарантирует, что решение, причем уникальное, всегда будет существовать.

#Базисы

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

Базисы есть не только в линейной алгебре. Например, ${1, x, x^2}$ является базисом всех квадратных трёхчленов. Или ${\neg, \land, \lor}$ является базисом всех логических выражений (то есть всё можно выразить через И, ИЛИ и НЕ). В произвольном языке программирования можно выделить какой-то набор команд, через который можно будет написать всё что угодно, и он тоже в этом смысле будет базисом.

Определение. Система векторов ${a_1, a_2, \dots, a_n}$ называется линейно зависимой, если существует такой набор действительных коэффициентов $(\lambda_1, \lambda_2, \dots, \lambda_n)$, что $\lambda_i \neq 0$ хотя бы для одного $i$ и

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

Утверждение. Базис $n$-мерного пространства — это линейно независимый набор из $n$ векторов.

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

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

Упражнение. Дан массив из $10^5$ целых чисел от $0$ до $(2^{30}-1)$. Требуется найти количество различных подпоследовательностей этого массива, xor-сумма которых равна заданному числу $x$.