Метод Видемана-Копперсмита

{Автор - Участник:GKormakov}

Свойства и структура алгоритма

Общее описание алгоритма

Алгоритм Видемана-Копперсмита был предложен Дугласом Видеманом в статье 1986 года . Алгоритм позволял получить решение системы линейных уравнений [math]\displaystyle{ Ax = b }[/math] над полем GF(2) (в общем случае, над конечным полем GF(q)) в случае, когда матрица A являлась разреженной. В частности, совместный алгоритм Видемана-Копперсмита решал задачу поиска нетривиальных векторов из ядра матрицы [math]\displaystyle{ A }[/math], т.е. нетривиальные решения уравнения [math]\displaystyle{ Ax = 0 }[/math]

Стоит отметить, что алгоритм описывался для квадратной матрицы [math]\displaystyle{ A }[/math], однако его легко можно обобщить на случай прямоугольной матрицы, переходя к матрице [math]\displaystyle{ A^T A }[/math]

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

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

Алгоритм получил широкое применение в алгоритмах криптографии, в том числе, в алгоритме факторизации, использующего факторные базы

Использование стандартных методов Крылова (например, метода Ланцоша) не могли гарантировать условий на матрицу в GF(2), которые гарантировали сходимость, в связи с этим и были предложены перечисленные выше методы.

Математическое описание алгоритма

Метод Видемана

Невырожденный случай

Пусть [math]\displaystyle{ Ax = b, A\in\mathbb{F}_2^{n \times n}, b\in \mathbb{F}_2^n }[/math] - согласованная СЛАУ с квадратной невырожденной матрицей.

[math]\displaystyle{ \begin{align} f(\lambda) = 1 + \sum_{i=1}^n f_i \lambda^i \end{align} }[/math] - характеристический полином матрицы [math]\displaystyle{ A }[/math] с коэффициентами [math]\displaystyle{ f_i\in\mathbb{F}_2 }[/math]

По теореме Гамильтона-Кэли [math]\displaystyle{ f(A) = 0 }[/math], следовательно, [math]\displaystyle{ f(A)b = 0 }[/math]. И, учитывая вид характеристической функции:

[math]\displaystyle{ \begin{align} b = \sum_{i=1}^n f_i A^i b = A(\sum_{i=1}^{n-1} f_{i+1}A^i b ) \end{align} }[/math]

Таким образом, если коэффициенты характеристического полинома известны, то единственное решение исследуемой системы можно найти как линейную комбинацию векторов Крылова: [math]\displaystyle{ \begin{align} x = \sum_{i=1}^{n-1} f_{i+1}A^i b \end{align} }[/math]

Но как искать эти коэффициенты? Для случая СЛАУ с разреженной матрицей Видеман предложил следующую идею: Выбрать случайный вектор [math]\displaystyle{ y \in \mathbb{F}_2^n }[/math], тогда для [math]\displaystyle{ j = 0...n }[/math] верны следующие соотношения:

[math]\displaystyle{ \begin{align} y^Tf(A)A^jb = y^T(\sum_{i=0}^{n} f_iA^i)A^jb = \sum_{i=0}^{n}(y^TA^{i+j}b)f_i = 0 \end{align} }[/math]

Пусть [math]\displaystyle{ \alpha_i = y^T A^ib \in \mathbb{F}_2, \forall i=0...n }[/math], тогда верхняя система преобразуется в однородную СЛАУ:

[math]\displaystyle{ \begin{align} \sum_{i=0}^{n}\alpha_{i+j}f_i = 0 , j = 0...n \end{align} }[/math]

с матрицей [math]\displaystyle{ H_{ij}=\alpha_{i+j} }[/math] порядка n + 1, имеющую ганкелеву структуру, что позволяет проще (за [math]\displaystyle{ O(n^2) }[/math] - алгоритм Тренча) решать СЛАУ, поскольку матрица задаётся 2n + 1 параметрами. Если переставить строки матрицы [math]\displaystyle{ H }[/math] в обратном порядке, то матрица станет теплицевой. Над GF(2) также известны быстрые алгоритмы решения СЛАУ (алгоритм Берлекампа - Масси)

Вырожденный случай

Для случая нахождения произвольного вектора из ядра матрицы [math]\displaystyle{ A }[/math] идея Видемана применяется в следующем виде:

Поскольку [math]\displaystyle{ A }[/math] теперь вырожденная, то характеристический полином будет выглядеть иначе (занулятся младшие коэффициенты):

[math]\displaystyle{ \begin{align} f(\lambda) = \lambda^s(\sum_{i=0}^{n-s} f_{i+s} \lambda^i), причём f_s \neq 0 \end{align} }[/math]

В этом случае выбираем произвольный вектор [math]\displaystyle{ x\in\mathbb{F}_2^n }[/math] и находим вектор [math]\displaystyle{ \omega \in \mathbb{F}_2^n }[/math]: [math]\displaystyle{ \begin{align} \omega = (\sum_{i=0}^{n-s} f_{i+s} A^i)x \end{align} }[/math]

Т.к. размерность ядра матрицы [math]\displaystyle{ \sum_{i=0}^{n-s} f_{i+s} A^i }[/math] не может быть больше n - 1, а размерность пространства [math]\displaystyle{ \mathbb{F} }[/math] равна n, то с вероятностью 1/2 [math]\displaystyle{ \omega \neq 0 }[/math]. Если же он равен 0, то выбираем другой x и ищем степень матрицы [math]\displaystyle{ A }[/math]: [math]\displaystyle{ A^\delta \omega \neq 0, а A^{\delta + 1} \omega = 0 }[/math], тогда [math]\displaystyle{ A^\delta \omega }[/math] - искомый вектор из ядра.

Перепишем итоговый алгоритм Видемана для нахождения вектора из ядра A:

Шаг 1. Выбираем случайные векторы [math]\displaystyle{ x, y \in \mathbb{F}_2^n }[/math]

Шаг 2. [math]\displaystyle{ \forall s = 0, ... 2n }[/math] вычисляем [math]\displaystyle{ \alpha_s = y^TA^sx }[/math] и строим полином [math]\displaystyle{ a(x) = \sum\alpha_ix^i \textit{формальной} }[/math] степени 2n

Шаг 3. Характеристическому полиному [math]\displaystyle{ f(x) }[/math] ставим в соответствие полином [math]\displaystyle{ \widetilde{f(x)} = x^n f(\frac{1}{x}), тогда f_i = \widetilde{f(x)}_{n-i} }[/math] и

[math]\displaystyle{ \sum_{i=0}^n\alpha_{i+j}f_i = \sum_{i=0}^n\alpha_{i+j}\widetilde{f(x)}_{n-i} = 0, j=0...n }[/math], т.е. все коэффициенты в [math]\displaystyle{ a(x)\widetilde{f(x)} }[/math] при [math]\displaystyle{ x^n, x^{n+1}, ..., x^{2n} }[/math] будут равны

Т.е. [math]\displaystyle{ a(x)\widetilde{f(x)} = g(x) mod \,x^{2n + 1}, где\, deg(g) \lt n }[/math]

Таким образом построили полином [math]\displaystyle{ \widetilde{f(x)} }[/math]

Шаг 4. Получаем решение [math]\displaystyle{ \omega = A^\delta \sum_if_iA^ix : A\omega = 0 }[/math] алгоритмом, описанным выше

Алгоритм Видемана-Копперсмита

Алгоритм Видемана по поиску вектора из ядра матрицы A сводился к решению СЛАУ с генкелевой или теплицевой матрицей, что ускоряло алгоритм, однако алгоритм Видемана является последовательным алгоритмом (коэффициенты полинома [math]\displaystyle{ a(x) }[/math] вычисляются последовательно), т.е. он не имеет выигрыша по скорости перед методом Монтгомери при распараллеливании на несколько узлов. Поэтому Копперсмит добавил в алгоритм Видемана работу с блочными матрицами.

Пусть также задана однородная СЛАУ [math]\displaystyle{ Ax = 0\, , A\in\mathbb{F}_2^{n\times n} }[/math], так же требуется найти нетривиальное решение. Опишем предложения Копперсмита:

Шаг 1. Рассмотрим два случайных блока [math]\displaystyle{ Z\in \mathbb{F}_2^{n\times m_b} }[/math] и [math]\displaystyle{ Y\in \mathbb{F}_2^{n\times n_b} }[/math], [math]\displaystyle{ m_b \neq n_b }[/math]. Тогда столбцы матрицы [math]\displaystyle{ AZ = X\in \mathbb{F}_2^{n\times m_b} }[/math] будут принадлежать линейной оболочке столбцов матрицы [math]\displaystyle{ A }[/math]. На практике числа [math]\displaystyle{ m_b, n_b }[/math] кратны размеру машинного слова, это будет использовано при параллельной реализации.

Далее вычисляются матрицы [math]\displaystyle{ a_i = Y^TA^iX \in \mathbb{F}_2^{n_b\times m_b} \forall\, i=0...L }[/math], где [math]\displaystyle{ L = \frac{n}{n_b} + \frac{n}{m_b} + \Delta }[/math], где [math]\displaystyle{ \Delta }[/math] - небольшое натуральное число, например, 10. С помощью этих матриц определяем полином [math]\displaystyle{ a(x) = \sum_{i=0}^La_ix^i }[/math]

Определение Степенью матричного полинома назовём степень полинома максимальной степени.

Однако на практике оказывается более важным смотреть на формальную степень. Работая с ними можно использовать следующие правила: 1) Формальная степень суммы полиномов равна максимуму из формальных степеней каждого слагаемого; 2) Умножение полинома на одночлен увеличивает его формальную степень на 1. Таким образом, построенный на шаге 1 полином является формальной степени L.

Шаг 2. Построим матричный полином [math]\displaystyle{ f(x) = \sum_{i=0}^Df_ix^i }[/math] формальной степени D с матричными коэффициентами [math]\displaystyle{ f_i\in\mathbb{F}_2^{m_b\times(m_b + n_b)} }[/math], для которого выполняются следующие условия:

1. Справедливо полиномиальное соотношение [math]\displaystyle{ a(x)f(x) = g(x) + e(x)x^M }[/math], где [math]\displaystyle{ g(x) }[/math] и [math]\displaystyle{ e(x) }[/math] - матричные полиномы;
2. Степень матричного полинома [math]\displaystyle{ g }[/math] не превосходит [math]\displaystyle{ D }[/math];
3. [math]\displaystyle{ M }[/math] и [math]\displaystyle{ D }[/math] удовлетворяют соотношению [math]\displaystyle{ M - D - 1 \gt \frac{n}{n_b} }[/math]

(Метод построения матричного полинома [math]\displaystyle{ f(x) }[/math] изложим ниже)

Т.к. [math]\displaystyle{ deg(g) \leq D }[/math], то матричные коэффициенты при степенях, больших D, равны 0. Следовательно, матричные коэффициенты a(x) и f(x) удовлетворяют системе:

[math]\displaystyle{ \sum\limits_{k=0}^D a_{i + k + 1}f_{D-k} = 0, \; \forall i=\overline{0, M-D-2} }[/math], что эквивалентно, при подстановке [math]\displaystyle{ a_{i+k+1} }[/math], следующему уравнению:

[math]\displaystyle{ Y^T A^i A\underbrace{\left(\sum\limits_{k=0}^D (A^{k} X) f_{D - k}\right)}_{W} = 0, \; \forall i = \overline{0, M - D - 2} }[/math]

Таким образом, W удовлетворяет условиям А-ортагональности со всеми блоками векторов вида [math]\displaystyle{ Y, A^TY, (A^T)^2Y... (A^T)^{M - D - 2}Y }[/math]. Далее возможно два варианта:

1) Если линейная оболочка блоков вышеперечисленных векторов совпала со всем пространством [math]\displaystyle{ \mathbb{F}_2^n }[/math], то выбором n линейно независимых столбцов и составлением невырожденной матрицы [math]\displaystyle{ \hat{Y} }[/math] получили бы следующее равенство: [math]\displaystyle{ \hat{Y}AW = 0 }[/math]. Домножив на обратную, получим [math]\displaystyle{ AW = 0 }[/math] и если [math]\displaystyle{ W \neq 0 }[/math], то задача решена.

2)Оболочка блоков не охватила всего пространства, тогда можно вероятностно оценить, что: Т.к. [math]\displaystyle{ M - D - 1 \gt \frac{n}{n_b} }[/math], то почти наверняка размерность линейной оболочки [math]\displaystyle{ (A^T)^iY, \; \forall i=\overline{0,M-D-2} }[/math] близка к n. Этого достаточно, чтобы утверждать, что задача решена. (Для деталей см. )

Метод построения матричного полинома

Копперсмит предложил итерационную процедуру:

1. Инициализация:[math]\displaystyle{ M := 0,\; f(x) := [I_{m_b}\; 0],\; e(x) := [a(x)\; I_{n_b}],\; \delta_j = 0 }[/math], где j -- номер столбца [math]\displaystyle{ f(x), \delta_j }[/math] -- формальная степень;
2. Итерация:
   1. Вычислить [math]\displaystyle{ e_0 = e(0) }[/math];
   2. Найти P: [math]\displaystyle{ e_0P = [0\; F_{n_b}] }[/math], причём преобразование P должно быть допустимым;
   3. Вычислить [math]\displaystyle{ e:=\frac{eP(x)}{x}, \; f:=fP(x) }[/math], где [math]\displaystyle{ P(x) = P \begin{pmatrix}I_{m_b} & 0\\ 0 & xI_{n_b}\end{pmatrix} }[/math];
   4. [math]\displaystyle{ M:=M+1 }[/math]
   5. [math]\displaystyle{ \delta_{j_k}:= \delta_{j_k} + 1 }[/math] для тех столбцов, которые были умножены на x;
   6. Вычислить формальную степень D матричного полинома f(x);
3. Останов Проверить критерий останова: [math]\displaystyle{ M - D \gt \frac{n}{n_b} + \Delta }[/math]. Если выполнен, то остановиться, нет - продолжить.

Сколько всего потребуется итераций? Если [math]\displaystyle{ m_b = n_b }[/math], то будет выполнена оценка [math]\displaystyle{ s \approx 2\left(\frac{n}{n_b} + 1\right) }[/math]. Если же брать [math]\displaystyle{ m_b = 2n_b }[/math], как советовал Копперсмит, то [math]\displaystyle{ s \approx \frac{3}{2}\left(\frac{n}{n_b} + 1\right) }[/math] (Более подробный анализ см. в ). С точки зрения вероятностной устойчивости и скорости логичнее использовать второй случай, однако, как показывал опыт, и исходя из указаний Копперсмита (глава 5 ), для практического применения подходит первый случай.

Литература

<references \>