April 21, 2017
При решении задач с помощью производящих функций (и вообще при работе с формальными степенными рядами) может возникнуть такая задача: по известным степенным рядам двух функций построить степенной ряд их частного.
$f$ и $h$ известны, найти $g$.
С алгоритмической точки зрения есть две естественные постановки:
Мы займёмся второй задачей.
Из определения понятно, что $ h_{k} = \sum_{i=0}^{k} f_{i} g_{k-i} $. Тогда можно найти $g_{0} = \frac{h_{0}}{f_{0}}$, потом при известном $g_{0}$ вычислить $g_{1} = \frac{h_{1} - f_{1} g_{0}}{f_{0}}$, затем вычислить $g_{2}$ и так далее.
Если все $g_{0}$, $g_{1}$, $\ldots$, $g_{k-1}$ уже найдены, то $g_{k}$ можно посчитать за время $O(k)$. Общее время работы составит $O(n^{2})$.
Заметим также, что $g_{k}$ зависит только от $f_{i}$ и $h_{i}$ при $i \le k$.
Мы будем решать задачу поиска обратного степенного ряда, то есть по $f_{i}$ найти $g_{i}$ такие, что
Поймём сначала, что этого достаточно для решения исходной задачи. А именно: пусть $G(x) = \sum_{i=0}^{n-1} g_{i} x^{i}$ — многочлен, составленный из первых $n$ коэффициентов степенного ряда, обратного к $F$, $H(x) = \sum_{i=0}^{n-1} h_{i} x^{i}$ — многочлен, полученный отбрасыванием старших степеней степенного ряда $\sum_{i=0}^{\infty} h_{i} x^{i}$. Тогда первые коэффициенты произведения $G(x) \cdot H(x)$ и будут искомыми коэффициентами в исходной задаче:
Последний переход верен, так как $ \sum_{i=0}^{k-j} f_{i} g_{(k-j)-i} $ — это $(k-j)$-й коэффициент многочлена, равного $1$.
Умножить $G(x)$ на $H(x)$ можно за время $O(n \log n)$ с помощью FFT.
Мы тоже будем строить коэффициенты $G(x)$ постепенно, используя предыдущие, однако будем делать это не по одному, а группами по $2^{k}$. Основной шаг алгоритма — по первым $2^{k}$ коэффициентам получить следующие $2^{k}$ коэффициентов за время $O \left( k 2^{k} \right)$.
Но сначала нужно найти хотя бы один коэффициент. Это сделать легко: $g_{0} = \frac{1}{f_{0}}$.
Обозначим $m=2^{k}$, $F(x) = \sum_{i=0}^{2m-1} f_{i} x^{i}$, $G(x) = \sum_{i=0}^{2m-1} g_{i} x^{i}$. Какое условие на $G(x)$ покажет, что это правильные $2m$ коэффициентов обратного степенного ряда для $\sum_{i=0}^{\infty} f_{i} x^{i}$ ? Мы знаем, что $G(x)$ однозначно определяется первыми $2m$ коэффициентами $F(x)$ и $H(x)$; следовательно, $G(x)$ должен быть таким, чтобы первые $2m$ коэффициентов $F(x) \cdot G(x)$ совпали с первыми $2m$ коэффициентами $H(x)$, т.е. $F(x) \cdot G(x) = 1 + x^{2m} \cdot P(x)$, где $P(x)$ — некоторый многочлен.
Введём обозначения:
На предыдущих шагах алгоритма мы по $F_{1}(x)$ построили $G_{1}(x)$ такой, что $F_{1}(x) \cdot G_{1}(x) = 1 + x^{m} \cdot P(x)$. Сейчас мы хотим построить подходящий $G_{2}(x)$.
Заметим, что $F(x) = F_{1}(x) + x^{m} \cdot F_{2}(x)$, $G(x) = G_{1}(x) + x^{m} \cdot G_{2}(x)$. Тогда
$F_{2}(x)G_{2}(x)$ — это некоторый многочлен, поэтому на последнее слагаемое можно не обращать внимания. А ещё мы помним, что
Таким образом, нужно найти $G_{2}(x)$ такой, что
Поделим на $x^{m}$ и домножим на $G_{1}(x)$:
Но $F_{1}(x)G_{1}(x) = 1 + x^{m} \cdot P(x)$, поэтому
Таким образом, для получения $G_{2}(x)$ нужно выполнить 3 умножения многочленов размера $2m$ (после умножения старшие коэффициенты можно отбросить). При помощи FFT это можно сделать за $O(m \log m)$.
Общее время работы алгоритма составит
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
vector<double> multiply(vector<double> A, vector<double> B);
vector<double> inverse(vector<double> F)
{
int oldSz = (int)F.size();
int n = oldSz;
while(n & (n - 1))
{
F.push_back(0);
n++;
}
vector<double> G;
G.push_back(1. / F[0]);
while((int)G.size() < (int)F.size())
{
n = G.size();
G.resize(2 * n);
for (int i = n; i < 2 * n; i++)
G[i] = 0;
vector<double> F1, F2;
F1.resize(2 * n);
F2.resize(2 * n);
for (int i = 0; i < n; i++)
{
F1[i] = F[i];
F2[i] = F[n + i];
F1[n + i] = F2[n + i] = 0;
}
vector<double> P = multiply(G, F1);
for (int i = 0; i < n; i++)
{
P[i] = P[n + i];
P[n + i] = 0;
}
F2 = multiply(F2, G);
for (int i = 0; i < n; i++)
P[i] += F2[i];
P = multiply(P, G);
for (int i = 0; i < n; i++)
G[n + i] = -P[i];
}
G.resize(oldSz);
return G;
}
Все вышеизложенное верно и для модульной арифметики. В коде нужно использовать правильное умножение многочленов (смотри статью про FFT), сложение и вычитание производить по модулю.
Несмотря на то, что сложность алгоритма составляет $O(n \log n)$, константа под $O$ весьма велика.