РЕШЕНИЕ СИСТЕМ ЛИНЕЙНЫХ ДИОФАНТОВЫХ УРАВНЕНИЙ ^{∗}

\begin{matrix}c \\ ◯ \end{matrix} Г.И. Малашонок

Malaschonok G.I. Solving systems of linear Diophantine equations. New methods to solve linear systems of Diophantine equations are suggested. These problems are considered from two points of view – modular and padic. Both methods allow to get solutions of a linear system with the size n × m with the complexity O(n^{β+1}m). For quasi-square system, the p-adic method allows to get solutions with the complexity O(n^{3}). Both estimations have the accuracy up to the logarithmic multipliers, β – is a power in the estimation of the time of matrix multiplication.

∗Вестник ТГУ, Т.5, N. 5, 2000, 620–627

ВВЕДЕНИЕ

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

Основной метод для решения линейных систем диофантовых уравнений основан на вычислении формы Смита матрицы коэффициентов системы [1, 2]. В недавних работах были получены вероятностные методы для нахождения решений линейной системы диофантовых уравнений [3, 4]. Веpоятностные методы являются самыми быстpыми из известных сегодня методов. Однако они не гаpантиpуют нахождение pешения. Если pешение не найдено веpоятностным методом, то всегда существует веpоятность того, что система совместна, но ее pешение не найдено случайно.

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

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

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

Показано, что для квазиквадратных систем сложность p-адического метода O(n^{3}), а модулярного O(n^{β+1}). Для существенно прямоугольных систем, когда число неизвестных m более чем вдвое пpевышает число уpавнений n, оба метода имеют одинаковый порядок сложности – O(n^{β+1}m). Здесь оценки пpиводятся с точностью до логарифмического множителя и пpедполагается пpименение стандаpтного умножения. β – показатель степени в оценке сложности матричного умножения.

Пpиводится алгоpитм со сложностью O(n^{β−1}m) для pешения систем над кольцами вычетов.

Как обычно, все pезультаты пеpеносятся с кольца целых чисел на кольцо полиномов над полем.

КВАЗИЕВКЛИДОВЫ КОЛЬЦА

Пусть R – коммутативное кольцо с единицей. Будем называть диофантовой матрицей для ненулевой пары α из R^{2} унимодулярную матрицу E_{α} ∈ R^{2×2} такую, что

E_{α}α = ( \begin{matrix}g \\ 0 \end{matrix} ), g ≠ 0, det E_{α} = 1.

Для нулевой пары диофантовой матрицей будем называть единичную.

Определение 1 Кольцо R назовем квазиевклидовым, если в кольце определен алгоритм, вычисляющий для любой пары диофантову матрицу.

Из обратимости диофантовой матрицы следует, что идеал (a, b) является главным идеалом, порожденным g. Поэтому всякое квазиевклидово кольцо является кольцом главных идеалов.

Наоборот. Пусть R кольцо главных идеалов, в которой существует алгоритм GI для вычисления генератора g идеала (p, q), порожденного двумя элементами p, q ∈ R:

{g, a, b} = GI(p, q), (g) = (p, q), ap + bq = g

Такое кольцо будем называть GI кольцом. Диофантовой матрицей для ненулевой пары

α = ( \begin{matrix}p \\ q \end{matrix} ) будет E_{α} = ( \begin{matrix}a & b \\ −q/g & p/g \end{matrix} ).

Заметим, что числа a, b и g определены неоднозначно. Поэтому и диофантова матрица определена неоднозначно. Если в R ввести порядок, то этот выбор можно сделать однозначным.

Если R евклидово кольцо, то в качестве алгоритма GI можно взять расширенный алгоритм Евклида.

Предложение 1 Псевдоевклидовость кольца с отношением порядка наследуется его факторкольцом.

Доказательство.

Пусть R′ – факторкольцо псевдоевклидова кольца R c некоторым отношением порядка и ϕ : R → R′ канонический гомоморфизм. Определим функцию ψ : R′ → R, выбирающую представителя класса смежности для p′ ∈ R′ в прообразе ϕ^{−1}(p′):

ψ(p′) = min(ϕ^{−1}(p′)).

Пусть ψ(p′) = p, ψ(q′) = q, α = ( \begin{matrix}p \\ q \end{matrix} ), E_{α} –

диофантова  матрица  для  α  и

E_{α}α = ( \begin{matrix}g \\ 0 \end{matrix} ). Определим для ненулевой пары

α′ = ( \begin{matrix}p′ \\ q′ \end{matrix} ) ∈ R′^{2} диофантову матрицу E_{α′}

следующим образом E_{α′} = ϕ(E_{α}). Очевидно, такое определение корректно, так как

E_{α′}α′ = ( \begin{matrix}g′ \\ 0 \end{matrix} ), det E_{α′} = 1, g′ ≠ 0.

Покажем, что g′ = ϕ(g) ≠ 0.

Из-за обратимости диофантовых матриц p и q делятся на g. Поэтому если ϕ(g) = 0, то ϕ(p) = ϕ(q) = 0, но по условию хоть одно из чисел p′ = ϕ(p), q′ = ϕ(q) не должно быть нулем.

Следовательно, множество квазиевклидовых колец содержит все GI кольца и их факторкольца, в том числе евклидовы кольца и их факторкольца.

Будем называть матрицу A = (a_{ij}), размера n × m, верхней треугольной, если все элементы под главной диагональю нулевые: (i, j) = 0 для всех n ≥ j > i ≥ 1. Заметим, что элементы на главной диагонали и над ней могут тоже быть нулевыми.

Будем говорить, что матрица A разложима, если она является произведением унимодулярной матрицы Ã на верхнюю треугольную матрицу Ā: A = ÃĀ.

Предложение 2 В квазиевклидовом кольце все матрицы разложимы.

Доказательство.

Построим такое разложение для матрицы A = (a_{ij}), размера n × m. Пусть

α_{ij} = ( \begin{matrix}a_{jj} \\ a_{ij} \end{matrix} ), i > j, E_{α_{i,j}} – диофантова

матрица для α_{ij}. Обозначим через E_{α_{i,j}}^{n} унимодулярную матрицу, которая получается после замены в единичной матрице порядка n четырх элементов (jj), (ji), (ij), (ii) соответственно на четыре элемента (11), (12), (21),

(22) матрицы E_{α_{i,j}}. Пусть E_{j} = E_{α_{n,j}} ⋯ E_{α_{j+1,j}}, U = E_{n−1} ⋯ E_{1}. Тогда U – унимодулярная матрица, UA – верхняя треугольная матрица. Эти вычисления приводят к последовательному исключению ненулевых элементов сначала в первом столбце, затем во втором столбце, и так далее до n − 1-го столбца.

Обозначим C_{D} – асимптотическую оценку сложности вычисления диофантовой матрицы, M(n) – асимптотическую оценку сложности умножения матриц порядка n, M(n) = αn^{β}. Известное сегодня значение для β меньше 2,356 [5].

Лемма 1 Пусть R – квазиевклидово кольцо, A, C ∈ R^{n×n}, A – верхняя треугольная матри-

ца. Тогда матрицу ( \begin{matrix}A \\ C \end{matrix} ) можно разложить

( \begin{matrix}A \\ C \end{matrix} ) = Q( Ā ; 0 ) с асимптотической оцен-

кой сложности, не превосходящей

C_{s}(n) = n^{2}C_{D} + αλ_{β}n^{β},  (1)

где λ_{β} = 7(2^{β−2} − 1)^{−1}.

Доказательство. Для разложения матрицы размера 2 × 1 достаточно умножить на ее диофантовую матрицу. Пусть теперь n = 2^{p}, p > 0.

Разобьем A и C на квадратные блоки

A = ( \begin{matrix}a & b \\ 0 & d \end{matrix} ), C = ( \begin{matrix}e & f \\ g & h \end{matrix} ).

Пусть имеются четыре разложения:

( \begin{matrix}a \\ e \end{matrix} ) = Q_{1}^{−1}( \begin{matrix}a_{1} \\ 0 \end{matrix} ), ( \begin{matrix}b_{1} \\ f_{1} \end{matrix} ) = Q_{1}( \begin{matrix}b \\ f \end{matrix} ),

( \begin{matrix}a_{1} \\ g \end{matrix} ) = Q_{2}^{−1}( \begin{matrix}a_{2} \\ 0 \end{matrix} ), ( \begin{matrix}b_{2} \\ h_{1} \end{matrix} ) = Q_{2}( \begin{matrix}b_{1} \\ h \end{matrix} ),

( \begin{matrix}d \\ f_{1} \end{matrix} ) = Q_{3}^{−1}( \begin{matrix}d_{1} \\ 0 \end{matrix} ), ( \begin{matrix}d_{1} \\ h_{1} \end{matrix} ) = Q_{4}^{−1}( \begin{matrix}d_{2} \\ 0 \end{matrix} ).

Матрицы Q_{i} = ( \begin{matrix}α_{i} & β_{i} \\ γ_{i} & δ_{i} \end{matrix} ), i = 1, .., 4,

унимодулярные, a, a_{1}, a_{2}, d, d_{1}, d_{2} – верхние тре-

угольные. Тогда получим Ā = ( \begin{matrix}a_{2} & b_{2} \\ 0 & d_{2} \end{matrix} ) и

Q = ( \begin{matrix}1 & 0 & 0 & 0 \\ 0 & α_{4} & 0 & β_{4} \\ 0 & 0 & 1 & 0 \\ 0 & γ_{4} & 0 & δ_{4} \end{matrix} ) ( \begin{matrix}1 & 0 & 0 & 0 \\ 0 & α_{3} & β_{2} & 0 \\ 0 & γ_{2} & δ_{2} & 0 \\ 0 & 0 & 0 & 1 \end{matrix} ) ×

× ( \begin{matrix}α_{2} & 0 & 0 & β_{2} \\ 0 & 1 & 0 & 0 \\ 0 & 0 & 1 & 0 \\ γ_{3} & 0 & 0 & δ_{3} \end{matrix} ) ( \begin{matrix}α_{1} & 0 & β_{1} & 0 \\ 0 & 1 & 0 & 0 \\ γ_{1} & 0 & δ_{1} & 0 \\ 0 & 0 & 0 & 1 \end{matrix} ) =

( \begin{matrix}α_{1}α_{2} & 0 & β_{1}α_{2} & β_{2} \\ η_{1} & α_{4}α_{3} & η_{2} & β_{4}δ_{2} \\ δ_{3}γ_{1} & γ_{3} & δ_{3}δ_{1} & 0 \\ η_{3} & γ_{4}α_{3} & η_{4} & δ_{4}δ_{2} \end{matrix} ) .

Используются обозначения η_{1} = α_{4}β_{3}γ_{1} + β_{4}γ_{2}α_{1}, η_{2} = α_{4}β_{3}δ_{1} + β_{4}γ_{2}β_{1}, η_{3} = γ_{4}β_{3}γ_{1} + δ_{4}γ_{2}α_{1}, η_{4} = γ_{4}β_{3}δ_{1} + δ_{4}γ_{2}β_{1}.

Следовательно, достаточно выполнить 28 умножений блоков порядка n/2 и четыре разложения для блоков размера n × n/2.

Продолжая выполнять разбиения на все более мелкие блоки, вплоть до блоков 2 × 1, получим асимптотическую оценку сложности вычислений

C_{s}(n) = n^{2}C_{D} + 7 ∑_{i=1}^{p} 4^{i}M(\frac{n}{2^{i}}) < n^{2}C_{D} + αλ_{β}n^{β}.

Здесь число λ_{β} = 7(2^{β−2}−1)^{−1} находится в пределах от λ_{3} = 7 до λ_{log 6} = 14.

Следствие 1 Пусть R – квазиевклидово

кольцо. Тогда матрицу ( \begin{matrix}A \\ C \end{matrix} ) над кольцом R,

размера n×m, с верхним треугольным блоком A ∈ R^{m×m} можно разложить с асимптотической оценкой сложности, не превосходящей C_{s}^{1}(n, m) = = (n − m)mC_{D} + αm^{β−1}(λ_{β}(n − m) + 4n − 8m), при n ≥ 2m и C_{s}^{2}(n, m) = (n − m)mC_{D}+ +α(n − m)^{β−2}(λ_{β}m(n − m) + 10mn − 4n^{2} − 4m^{2}), при m ≤ n ≤ 2m.

Доказательство.

a. Пусть n ≥ 2m и m делит n. Cначала выполним разложение для матрицы, образованной блоком A и верхним квадратным блоком C_{1} матрицы C: ( \begin{matrix}A \\ C_{1} \end{matrix} ) = Q_{1}( Ā ; 0 ). Затем

(n−2m)/m раз проведем разложение матрицы, образованной блоком Ā на предыдущем шаге и i-тым сверху квадратным блоком C_{i} матрицы C, i = 2, … k. Перемножение унимодулярных матриц, получаемых на каждом шаге, требует еще четырех умножений матриц размера m×m. Всего получим

C_{s}^{1}(n, m) = C_{s}(m)\frac{(n − m)}{m} + 4M(m)\frac{(n − 2m)}{m}.

b. Пусть m ≤ n ≤ 2m и n − m делит m. Разобьем матрицы A и C на блоки размера (n − m) × (n − m). Матрица C содержит одну блочную строку. Обозначим C = (C_{1}, C̄_{1}), где C_{1} – левый блок матрицы C. Матрица A содержит k = m/(n − m) блочных строк и столбцов.

Cначала разложим матрицу, образованную верхним левым диагональным блоком A_{1} матрицы A и левым блоком C_{1} матрицы C:

( \begin{matrix}A_{1} \\ C_{1} \end{matrix} ) = Q_{1}( Ā_{1} ; 0 ), Q_{1} – унимодулярная ма-

трица. И вычислим произведение ( Ā′_{1} ; C̄′_{2} ) =

Q_{1}( Ā_{1} ; C̄_{1} ), где Ā_{1} - первая блочная строка ма-

трицы A начиная со второго блока. Обозначим C̄′_{1}=(C_{1}, C̄_{1}), где C_{1} – левый квадратный блок. Всего необходимо 4(k −1) блочных умножений.

Затем последовательно раскладываем матрицы, образованные i-тым (i=2,…, k) диагональным блоком A_{i} матрицы A и левым блоком

C_{i} матрицы C̄′_{i}: ( \begin{matrix}A_{i} \\ C_{i} \end{matrix} ) = Q_{i}( Ā_{i} ; 0 ). И вычис-

ляем произведение ( Ā′_{i} ; C̄_{i+1} ) = Q_{1}( Ā_{i} ; C̄_{i} ), где

Ā_{i} - i-я блочная строка матрицы A, начиная с i + 1-го блока, C̄′_{i}=(C_{i}, C̄_{i}).

Кроме того, для вычисления унимодулярной матрицы Q требуется перемножить унимодулярные матрицы Q′_{i}, размера n × n, соответствующие Q_{i}. Матрица Q′_{i} получается из единичной матрицы I_{n} заменой четырех блоков (i,i), (i,k+1), (k+1,i), (k+1,k+1) соответственно на четыре блока (1,1), (1,2), (2,1), (2,2) матрицы Q_{i}. Умножение на Q′_{i} требует четырех блочных умножений, а всего на i-том шаге нужно 4(k − i) + 4 блочных умножений, i = 2, …, k.

Суммируя по всем шагам, получим C_{s}^{2}(n, m) : kC_{s}(n − m)+ (∑_{i=1}^{k} 4(k − i + 1) − 4) M(n − m) = = kC_{s}(n − m) + 2(k^{2} + k − 2)M(n − m) = = m(n − m)C_{D} + M(n − m)(λ_{β}m/(n − m)+ +2(5mn − 2n^{2} − 2m^{2})/(n − m)^{2}).

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

Теорема 1 В квазиевклидовом кольце всякую квадратную матрицу A порядка n можно разложить с асимптотической оценкой сложности C_{t}(n), не превосходящей

μ_{β}M(n) + \frac{n^{2} − n}{2}C_{D}, μ_{β} = \frac{9 · 2^{β−2} − 2}{(2^{β} − 2)(2^{β−2} − 1)}.

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

Не уменьшая общности, будем полагать, что n = 2^{p}. Разобьем матрицу A на четыре ква-

дратных блока A = ( \begin{matrix}A & B \\ C & D \end{matrix} ).

(1). Пусть имеется разложение блока A :

A = Ã^{−1}Ā.

Обозначим A_{1} = ÃA, B_{1} = ÃB.

(2). Пусть разложение блока ( \begin{matrix}A_{1} \\ C \end{matrix} ) имеет вид

( \begin{matrix}A_{1} \\ C \end{matrix} ) = Q^{−1}( \begin{matrix}A_{2} \\ 0 \end{matrix} ),

Q – унимодулярная, A_{2} – верхняя треугольная.

Обозначим ( \begin{matrix}B_{2} \\ D_{1} \end{matrix} ) = Q( \begin{matrix}B_{1} \\ D \end{matrix} ).

(3). Пусть имеется разложение блока D_{1} :

D_{1} = D̃_{1}^{−1}D̄_{1}.

Обозначим D_{2} = D̄_{1}.

Тогда разложима матрица A:

A = Ã^{−1}Ā,

Ã = diag(1, D̃_{1}) Q diag(Ã, 1), Ā = ( \begin{matrix}A_{2} & B_{2} \\ 0 & D_{2} \end{matrix} ).

Всего, при таком разложении, требуется 9 блочных умножений, два блочных разложения и одно разложение для блока размера n × n/2. Обозначим C_{t}(n) сложность разложения матрицы A. Получим

C_{t}(n) = 2C_{t}(n/2) + 9M(n/2) + C_{s}(n/2). (2)

По лемме 1 получим

C_{t}(n) = 2C_{t}(n/2) + (n/2)^{2}C_{D} + α(9 + λ_{β})(n/2)^{β}.

Будем разбивать блоки A и D далее на квадратные блоки. И так будем поступать вплоть до блоков второго порядка. При этом асимптотическая оценка сложности вычислений будет C_{t}(n) = ∑_{i=1}^{p} 2^{i−1}((\frac{n}{2^{i}} )^{2}C_{D} + α(9 + λ_{β})(n/2^{i})^{β}) < \frac{n^{2}−n}{2}C_{D} + αμ_{β}n^{β}, где число μ_{β} = \frac{1}{2^{β}−2}(9 + \frac{7}{2^{β−2}−1}) находится в пределах от μ_{3} = 8/3 до μ_{log 6} = 23/4.

Следствие 2 В квазиевклидовом кольце разложение прямоугольной матрицы A размера n × m можно выполнить с асимптотической оценкой сложности, не превосходящей C_{t}^{0}(n, m) = C_{D}(n^{2} − n)/2 + αn^{β}(μ_{β} + (m − n)/n) при n < m, C_{t}^{1}(n, m) = C_{D}(nm − (m^{2} + m)/2) + αμ_{β}m^{β}+ α(n − m)^{β−1}(λ_{β}m(n − m) + 10mn − 4n^{2} − 4m^{2}), при m ≤ n ≤ 2m и C_{t}^{2}(n, m) = C_{D}(nm − (m^{2} + m)/2)+ αm^{β−1}(μ_{β}m+λ_{β}(n−m)+4n−8m), при n ≥ 2m.

Доказательство.

Доказательство сразу следует из оценок для C_{t}(n) и C_{s}(n).

При n < m достаточно разложить левый квадратный блок, а затем умножить правый блок из m − n столбцов матрицы A на Ã.

При n > m достаточно разложить верхний квадратный блок и воспользоваться следствием Леммы 1.

Теорема 2 В квазиевклидовом кольце всякая система линейных уравнений Ay = c с верхней треугольной матрицей коэффициентов A ∈ R^{n×n} может быть решена с асимптотической оценкой сложности

C_{n} = C_{D}(n^{2} − n)/2 + αϕ_{β}n^{β},

ϕ_{β} = \frac{1}{2^{β} − 2}(5 + \frac{7}{2^{β−2} − 1}).

Доказательство.

Решение системы Ay = c можно свести к решению двух систем вдвое меньшего размера и одному матричному разложению.

Пусть ( \begin{matrix}A & B \\ 0 & D \end{matrix} )( \begin{matrix}y_{1} \\ y_{0} \end{matrix} ) = ( \begin{matrix}f \\ h \end{matrix} ) – блоч-

ная запись системы Ay = c. Все блоки матрицы A квадратные, блоки A, D – верхние треугольные.

(1). Пусть решением системы Dy_{0} = h будет y_{0} = V s_{0} + v, V ∈ R^{n/2×n/2}, v, s_{0} ∈ R^{n/2}, s_{0} – вектор параметров. Подставим его в систему Ay_{1} + By_{0} = f. Получим Ay_{1} + Gs_{0} = g, где G = BV, g = f − Bv.

(2). Разложим матрицу (A, G): (A, G)P = (Ā, 0). P - унимодулярная матрица, Ā - верхняя треугольная. По лемме 1 и замечанию к лемме оценка сложности такого разложения матрицы (A, G) составляет C_{s}(n/2).

(3). Пусть решением системы Āz_{1} = g будет z_{1} = W t_{1} + w, W ∈ R^{n/2×n/2}, w, t_{1} ∈ R^{n/2}, t_{1} – вектор параметров.

Тогда ( \begin{matrix}z_{1} \\ z_{0} \end{matrix} ) = diag(W, I)( \begin{matrix}t_{1} \\ t_{0} \end{matrix} ) + ( \begin{matrix}w \\ 0 \end{matrix} )

будет решением системы (Ā, 0)( \begin{matrix}z_{1} \\ z_{0} \end{matrix} ) = g.

Пусть P = ( \begin{matrix}P_{a} & P_{b} \\ P_{c} & P_{d} \end{matrix} ) – блочная запись

матрицы P. Так как ( \begin{matrix}y_{1} \\ s_{0} \end{matrix} ) = P ( \begin{matrix}z_{1} \\ z_{0} \end{matrix} ) и

y_{0} = V s_{0} + v, то решение исходной системы Āy = c получим в виде

y = ( \begin{matrix}P_{a}W & P_{b} \\ V P_{c}W & V P_{d} \end{matrix} )( \begin{matrix}t_{0} \\ t_{1} \end{matrix} ) + ( \begin{matrix}P_{a}w \\ V P_{c}w + v \end{matrix} ).

Всего требуется выполнить пять блочных умножений. Обозначим C_{n} сложность решения системы порядка n. Тогда C_{n} = 2C_{n/2} + 5M(n/2) + C_{s}(n/2) = = 2C_{n/2} + α(5 + λ_{β})n^{n/2} + (n/2)^{2}C_{D}. Отсюда получаем асимптотическую оценку

C_{n} = C_{D}(n^{2} − n)/2 + αϕ_{β}n^{β},

ϕ_{β} = \frac{1}{2^{β} − 2}(5 + \frac{7}{2^{β−2} − 1}).

Теорема 3 В квазиевклидовом кольце всякая система линейных уравнений Ax = c, A ∈ R^{n×m} может быть решена со следующими асимптотическими оценками сложности

C_{nm} = (nm − l)C_{D} + αψ_{β}l^{β},

ψ_{β} = μ_{β} + ϕ_{β} + λ_{β} \frac{n−m}{m} + \frac{4n}{m} − 8 при 2m ≤ n; ψ_{β} = μ_{β} + ϕ_{β} + λ_{β} \frac{m−n}{n} + \frac{6m}{n} − 9 при 2n ≤ m;

C_{nm} = (nm − l)C_{D} + αψ_{β}l^{β}+

+α|n − m|^{β−2}(λ_{β}l|m − n| + 10nm − 4m^{2} − 4n^{2}),

ψ_{β} = μ_{β} + ϕ_{β} при m ≤ n ≤ 2m, ψ_{β} = μ_{β} + ϕ_{β} + 2\frac{m}{n} − 1 при n < m ≤ 2n, где l = min(m, n).

Доказательство. Матрица (A, c) раскладывается в произведение (A, c) = P (A′, c′) унимодулярной и верхней треугольной за C_{t}^{∗}(n, m) операций по следствию 2 теоремы 1. Следовательно, достаточно решить систему

A′x = c′  (3)

с верхней треугольной матрицей коэффициентов.

a. Пусть m > n, A′ = (A, B), A – квадратный блок, имеющий верхний треугольный вид. По замечанию к лемме 1 можно построить унимодулярную матрицу Q, такую, что (A, B)Q = (Ā, 0), где Ā – верхняя треугольная, за C_{s}^{∗}(n, m) операций. Решения системы (3) и системы Āy = c′ связаны матрицей Q:

x = Q( \begin{matrix}y \\ t_{0} \end{matrix} ), где t_{0} – вектор параметров.

Пусть решением системы Āy = c′ будет Ut_{1} +u, U ∈ R^{n×n}, t_{1}, u, y ∈ R^{n}, t_{1} – вектор параметров, t_{1} ∈ R^{n}. Тогда решением исходной си-

стемы будет x = Qdiag(U, I)t + Q( \begin{matrix}u \\ 0 \end{matrix} ), где

t ∈ R^{m} – вектор параметров. Для вычисления матричного произведения Qdiag(U, I) нужно αn^{β−1}m операций. По теореме 2 решение системы Āy = c′ можно получить за C_{n} операций.

Cуммируя все операции, получим оценку сложности такого решения

C_{nm}^{1i} = C_{t}^{0}(n, m) + C_{s}^{i}(m, n) + C_{n} + M(n)(m/n),

i = 1, 2, где i = 1 соответствует случаю m ≥ 2n, i = 2 соответствует случаю 2n ≥ m > n.

b. Пусть n ≥ m, A′ = ( \begin{matrix}A \\ B \end{matrix} ), A – ква-

дратный блок, имеющий верхний треугольный вид. По теореме 2 решение системы Ax = c′ можно получить за C_{m} операций, согласно теореме 2. Следовательно, общее число операций для этого решения будет

C_{nm} = C_{t}^{i}(n, m) + C_{m} i = 1, 2,  (4)

где i = 1 соответствует случаю 2m ≥ n ≥ m, i = 2 соответствует случаю n ≥ 2m.

Утверждение теоремы получается подстановкой всех оценок в правую часть.

В частности для случая квадратной матрицы коэффициентов системы получим

C_{nn} = n^{2}C_{D} +αψ_{β}n^{β}, ψ_{β} = \frac{7 · 2^{β−1}}{(2^{β} − 2)(2^{β−1} − 2)}.

Для β = 3: ψ_{β} < 4, 67 . Для β = log 7: ψ_{β} < 6, 54 . Для β = log 6: ψ_{β} = 10, 5.

ПОДЪЕМ ПО ГЕНЗЕЛЮ

Пусть R – коммутативная область с единицей, K – поле частных области R, A ∈ R^{n×m}, c ∈ R^{n}

Ax = c  (5)

– система линейных уравнений над R.

Будем полагать, что n < m, rankA = n, k = m − n + 1 – число линейно независимых решений системы (5). Пусть A = (A_{0}, A_{1}), блок A_{0} квадратный, det A_{0} ≠ 0 и последняя строка матрицы A_{0}^{−1}(A_{1}, c) не содержит нулей. Общий случай сводится к данному, например, так, как предлагается в [4]. Пусть e_{j}– столбцы единичной матрицы I_{k} = (e_{1}, …, e_{k}) порядка k.

Если искать решение системы (5) в ви-

де ( \begin{matrix}x_{j} \\ ξ_{j}e_{j} \end{matrix} ), x_{j} ∈ K^{n−1}, ξ_{j} ∈ K\\{0\},

(j = 1, …, k), то придем к определенной системе порядка k.

Решение такой системы в поле K с помощью p-адического подъема по Гензелю можно получить за

D_{Z} = O(n^{3}(log n + log ∥A∥ + log p)^{2})

бит-операций, если R – кольцо целых чисел Z. Если же R = F [x] – кольцо полиномов над полем F, то достаточно

D_{F [x]} = O(n^{3}(deg A + deg p)^{2})

операций в поле F [3, 6]. Для получения k решений потребуется соответственно kD_{Z} и kD_{F [x]} операций.

Так как ξ_{j} ≠ 0, j = 1, …, k, то все k решений будут линейно независимы. Приводя каждое из решений к общему знаменателю, который является делителем некоторого минора матрицы A, запишем эти решения в виде

x_{j} = ( \begin{matrix}x_{j} \\ α_{j}e_{j} \end{matrix} )χ_{j}^{−1}, j = 1, …, k,

x_{j} ∈ R^{n−1}, α_{j}, χ_{j} ∈ R.

Обозначим δ = (α_{1}, …, α_{k}) и отметим, что δ – делитель минора матрицы (A, c), образованного первыми n − 1 столбцами и последним столбцом.

Домножим числитель и знаменатель каждого решения j на δα_{j}^{−1}. Обозначим x_{j}δα_{j}^{−1} = v_{j}, χ_{j}δα_{j}^{−1} = δ_{j}. Получим базисное множество ре-

шений системы (5) в виде x_{j} = ( \begin{matrix}v_{j} \\ δe_{j} \end{matrix} )δ_{j}^{−1}.

Обозначим u = (δ_{1}, …, δ_{k}) – вектоp и V = (v_{1}, …, v_{k}) – матpицу, u ∈ R^{k}, V ∈ R^{(n−1)×k}.

Теорема 4 Для того чтобы была совместна система (5) в кольце R, необходимо и достаточно, чтобы была совместна в R система

\{ \begin{matrix}V y & = & 0 mod δ \\ uy & = & δ. \end{matrix}  (6)

Если Y – множество решений ситемы (6), то

X = \{( \begin{matrix}y \\ V yδ^{−1} \end{matrix} ) | y ∈ Y \} – множество ре-

шений системы (5).

Доказательство. Множество pешений системы – это гипеpплоскость pазмеpности m − n в пpостpанстве K^{m}. Она опpеделяется m − n + 1 линейно независимым pешением. Следовательно, для всякого решения z ∈ K^{m} системы (5) найдутся числа σ_{i} ∈ K, i = 1, …, k, такие что ∑_{i=1}^{k} x_{i}σ_{i} = z и ∑_{i=1}^{k} σ_{i} = 1. Обозначим σ_{i}δ_{i}^{−1} = τ_{i} ∈ K, z = (y, w), y ∈ R^{k}, w ∈ R^{m−k}, t = (τ_{1}, …, τ_{k}). Тогда

∑_{i=1}^{k} δ_{i}τ_{i}x_{i} = z, ut = 1.

Если рассматривать (t, y, w) как вектор неизвестных, то получим систему уравнений

( \begin{matrix}δI_{k} & −I_{k} & 0 \\ V & 0 & −I_{n−1} \\ u & 0 & 0 \end{matrix} ) ( \begin{matrix}t \\ y \\ w \end{matrix} ) = ( \begin{matrix}0 \\ 0 \\ 1 \end{matrix} )

Домножая каждую часть на обратимую в K матрицу

( \begin{matrix}I_{k} & 0 & 0 \\ −V & δI_{n−1} & 0 \\ −u & 0 & δ \end{matrix} ) ,

приходим к системе

( \begin{matrix}δI_{k} & −I_{k} & 0 \\ 0 & V & −δI_{n−1} \\ 0 & u & 0 \end{matrix} ) ( \begin{matrix}t \\ y \\ w \end{matrix} ) = ( \begin{matrix}0 \\ 0 \\ δ \end{matrix} )

(7)

Пусть система (5) имеет диофантово решение z = (y, w) ∈ R^{m}. Тогда совместна система (7). Следовательно совместна и система (6). Наоборот, пусть совместна система (6) и y ∈ R^{k} – ее решение, тогда (δ^{−1}y, y, δ^{−1}V y) решение системы (7). Следовательно z = (y, δ^{−1}V y) – диофантово решение системы (5).

Оценим сложность решения системы (6) в кольце целых чисел и кольце многочленов.

Матрица V имеет размеры (n − 1) × k, а модуль δ является делителем минора порядка n. Следовательно, для R = Z

log δ ≤ n(log ∥A∥ + 0.5 log n),  (8)

где ∥A∥ – модуль наибольшего элемента матрицы A. Для R = F [z]

deg δ ≤ n deg ∥A∥,  (9)

где deg ∥A∥ – степень наибольшего по степени элемента матрицы A.

Пусть применяется стандартный алгоритм умножения целых чисел и многочленов: для умножения N- битных чисел требуется N^{2} битопераций, а для умножения двух многочленов степени N требуется N^{2} операций умножения в поле F.

При k ≥ n по формуле (3) получим асимптотические оценки числа операций соответственно для Z и F [x] :

η_{β}knM(n)(log ∥A∥ + 0.5 log n)^{2}
η_{β}knM(n) deg^{2} ∥A∥.

При k ≤ n по формуле (4) получим асимптотические оценки числа операций соответственно для кольца целых чисел Z : δk^{β−1}(log ∥A∥ + \frac{log n}{2})^{2}(n^{3}(4 + λ_{β}) +n^{2}k(ν_{β} + μ_{β} − λ_{β} − 4)) и для кольца многочленов над полем F [x]: δk^{β−1} deg^{2} ∥A∥(n^{3}(4 + λ_{β}) +n^{2}k(ν_{β} + μ_{β} − λ_{β} − 4)).

Отметим, что в общем случае сложность вычисления k решений систем в поле K не влияет на эти асимптотические оценки.

Однако, для квазиквадратных систем это не так, ибо величина k мала по сравнению с n, и ей можно пренебречь. В этом случае для Z и F [x] получим, соответственно, следующие оценки числа операций:

O(n^{3}(log n + log ∥A∥ + log p)^{2})
O(n^{3}(deg ∥A∥ + deg p)^{2}).

Если применять быстрые методы умножения целых чисел и многочленов, то эти оценки можно еще улучшить. Так, метод Карацубы требует O(N^{log 3}) операций, а сегодняшний рекорд составляет – O(N log N log log N).

МОДУЛЯРНЫЙ МЕТОД

Будем решать систему (5), полагая что A = (A_{0}, A_{1}), блок A_{0} квадратный, det A_{0} = δ, δ ≠ 0. Рекурсивный метод [7] позволяет привести систему (5) к диагональному виду

(δI_{n}, A′_{1})x = c′  (10)

за O(n^{β−1}m) операций в кольце R.

Все коэффициенты системы (10) – это миноры порядка n матрицы (A, c). Следовательно, на них распространяются оценки (8) и (9), соответственно, для кольца целых чисел и кольца многочленов над полем F.

Если применять стандартное умножение, то количество операций, которое потребуется для вычисления коэффициентов системы (10) для колец Z и F [x], соответственно будет

O(mn^{β+1}(log n+log ∥A∥)^{2}) и O(mn^{β+1} deg^{2} ∥A∥).

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

O(mn^{β}(log n + log ∥A∥)^{2})иO(mn^{β} deg^{2} ∥A∥). (11)

Решение системы (10) сводится к pешению системы

A′_{1}y = c′ mod δ.  (12)

в фактоpкольце R/δR. Пусть Y – множеcтво решений системы (12), тогда

X = \{( \begin{matrix}(c′ − A′_{1}y)δ^{−1} \\ y \end{matrix} )| y ∈ Y \}

– множество решений системы (10).

Оценим сложность решения системы (12). Матрица A′ имеет размеры n × (m − n), а модуль δ является минором порядка n. По теореме 3 получим следующие оценки.

При m − n ≥ n получим соответственно для Z и F [x]:

O(n^{β+1}m(log ∥A∥ + log n)^{2}),
O(n^{β+1}m deg^{2} ∥A∥).  (13)

При m − n ≤ n получим соответственно для Z и F [x]:

O(n^{3}(m − n + 1)^{β−1}(log ∥A∥ + log n)^{2}),

O(n^{3}(m − n + 1)^{β−1} deg^{2} ∥A∥).  (14)

В общем случае оценки (13), (14) являются оценками сложности решения системы модулярным методом, так как оценка (11) имеет меньший порядок по n и m−n. Однако для квазиквадратных систем величиной m − n можно пренебречь и оценка (11) будет давать основной вклад. В этом случае для Z и F [x] получим соответственно следующие оценки числа операций:

O(n^{β+1}(log n + log ∥A∥)^{2}) и O(n^{β+1} deg^{2} ∥A∥).

ЗАКЛЮЧЕНИЕ

Сравнение p-адического и модулярного метода показывает, что для квазиквадратных систем p-адический метод имеет существенное преимущество, он быстрее в O(n^{β−2}) раз.

Другим достоинством p-адического метода для систем прозвольного вида является тот факт, что общее решение системы может быть получено с большой вероятностью за время O((m−n+1)n^{3}) [6]. Ожидаемое количество итераций при этом будет

O(\frac{log n + log log ∥A∥}{m − n + 1}).

Каждая новая итерация начинается со случайного выбора ведущего минора порядка n.

Для существенно прямоугольных систем, когда m > 2n, оба метода имеют одинаковый порядок сложности O(mn^{β+1}).

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

ЛИТЕРАТУРА

1. Gregory R.T., Krishnamurthy E.V. Methods and Applications of Error-Free Computation. Berlin, 1984.

2. Hurt M.F.,Waid C.A. A generalized inverse wich give all the integral solutions to a system of linear equations// SIAM J. Appl. Math. V. 10. 1970. P. 547–550.

3. Mulders T. and Storjohann A. Diophantine Linear System Solving// Proceedings of ISSAC’99: ACM International Symposium on Symbolic and Algebraic Computation. July 1999, Vancouver, Canada. P. 181–188.

4. Malaschonok G. Effective Matrix Methods in Commutative Domains// Formal Power Series and Algebraic Combinatorics/ Ed. by D. Krob, A.A. Mikhalev, A.V. Mikhalev. Springer, 2000. P. 506–517.

5. Coppersmith D. and Winograd S. Matrix multiplication via arithmetic progressions// Journal of Symbolic Computation. 1990. V.

9. P. 251–280.

6. Dixon J. Exact solution of linear equations using p-adic expansions// Numer. Math. 1982. V. 40. P. 137–141.

7. Malaschonok G. Recursive Method for the Solution of systems of Linear Equations// Computational Mathematics/ Ed. by A. Sydow. Proceedings of the 15th IMACS World Congress. Berlin, 1997. V. 1. P. 475– 480.

Поступила в редакцию 1 октября 2000 г.