### smax's blog

By smax, 6 weeks ago,

• +220

 » 6 weeks ago, # | ← Rev. 2 →   +91 Thank you for the blog! The text is very easy to read in my opinion. I would like to add something.You provided one algorithm of finding the $k$-th term of a linear recurrence, which works in $O(M(n)\log{k})$, where $M(n)$ is the time complexity of polynomial multiplication. It is actually not the fastest algorithm (anymore), check out this paper: what you described seems to be Fiduccia's algorithm, and it takes $3M(n)\log{k} + O(n\log{k})$ time to evaluate (according to this paper). Its authors propose another algorithm, which runs in $2M(n)\log(k+1) + M(n)$. It works as follows:Let $Q(x)$ be the characteristic polynomial of our recurrence, and $F(x) = \sum_{i=0}^{\infty}a_ix^i$ be the generating formal power series of our sequence. Then it can be seen that all nonzero terms of $F(x)Q(x)$ are of at most $(n-1)$-st power. This means that $F(x) = P(x)/Q(x)$ for some polynomial $P(x)$. Moreover, we know what $P(x)$ is: it is basically the first $n$ terms of $F(x)Q(x)$, that is, can be found in one multiplication of $a_0 + \ldots + a_{n-1}x^{n-1}$ and $Q(x)$, and then trimming to the proper degree.So what remains is to find $[x^k]\dfrac{P(x)}{Q(x)}$ (recall that $[x^k]G(x)$ is the $k$-th coefficient of a formal power series or a polynomial $G(x)$). We do this recursively: $[x^k]\dfrac{P(x)}{Q(x)} = [x^k]\dfrac{P(x)Q(-x)}{Q(x)Q(-x)} = [x^k]\dfrac{S(x)}{T(x^2)} = [x^k]\dfrac{S_0(x^2) + xS_1(x^2)}{T(x^2)}.$Here $S(x)$ is $P(x)Q(-x)$, $S_0(x)$ and $S_1(x)$ are polynomials composed of only even/odd-indexed terms of $S(x)$, and we notice that $Q(x)Q(-x)$ is an even function, thus its polynomial representation contains only even-indexed terms, which represent a polynomial $T(x)$.Now, if $k$ is even, then $[x^k]\dfrac{S_0(x^2) + xS_1(x^2)}{T(x^2)} = [x^k]\dfrac{S_0(x^2)}{T(x^2)} = [x^{k/2}]\dfrac{S_0(x)}{T(x)},$otherwise $[x^k]\dfrac{S_0(x^2) + xS_1(x^2)}{T(x^2)} = [x^k]\dfrac{xS_1(x^2)}{T(x^2)} = [x^{(k-1)/2}]\dfrac{S_1(x)}{T(x)}.$Thus we divided $k$ by $2$ with the cost of only two multiplications. Moreover, one can probably even use fewer FFTs to compute $Q(x)Q(-x)$ and maybe apply other tricks to speed this up.
•  » » 6 weeks ago, # ^ |   0 Woah that's a neat trick, thanks for sharing!
 » 6 weeks ago, # |   +23 In the first code shown you write that if $n$ is $1$ (which means $f(x) = x - c_0$) then $a = [c[0]]$ without mentioning it. For some of us it might not be trivial that $x^k \equiv c_0^k \; mod(x - c_0)$. If someone else didn't get that, notice that $x^k - r^k = (x-r)(x^{k-1} + x^{k-2}r + x^{k-3}r^2 + \dots + xr^{k-2} + r^{k-1})$. It can be extended to arbitrary polynomials to get $f(x) \equiv f(k) \; mod (x-k)$, which is known as polynomial remainder theorem.
•  » » 6 weeks ago, # ^ |   +14 Good catch, I did forget to explain that part of the code. You can also think of it as $x^k \bmod (x - c_0) = (x \bmod (x - c_0))^k \bmod (x - c_0) = c_0^k \bmod (x - c_0)$. $x \bmod (x - c_0) = c_0$ by doing one round of long division.