#	User	Rating
1	tourist	4009
2	jiangly	3823
3	Benq	3738
4	Radewoosh	3633
5	jqdai0815	3620
6	orzdevinwang	3529
7	ecnerwala	3446
8	Um_nik	3396
9	ksun48	3390
10	gamegame	3386

#	User	Contrib.
1	cry	167
2	Um_nik	163
3	maomao90	162
3	atcoder_official	162
5	adamant	159
6	-is-this-fft-	158
7	awoo	157
8	TheScrasse	154
9	Dominater069	153
9	nor	153

Hi everyone, today I'd like to elaborate something more on the problem I on the recently ended div 1 contest.

As described in the official solution, a bottleneck of the problem is the following: Given a sequence $$${c_j}_{-n \leq j \leq n}$$$, compute

$$$ g_k = \sum_{-k\leq j\leq k} ([x^j] (x + 1/x)^k) c_j $$$

for all $$$0 \leq k \leq n$$$.

In the official solution, an $$$O(n\log^2 n)$$$ time algorithm is presented, and we also claimed that it is possible to compute $$$g_k$$$ in $$$O(n\log n)$$$ time. In this post, I'd like to elaborate on how to do that.

By the way, another bottleneck of the problem is the relaxed convolution problem, where the usual implementation also takes $$$O(n\log^2 n)$$$ time. I mentioned in comment that the relaxed convolution has a somehow practical solution that has time complexity

$$$ O\left(\frac{n\log^2 n}{\log \log n}\right), $$$

and there are more intricated optimization that have theoretical complexity $$$O(n(\log n)^{1+o(1)})$$$, we have theoretically faster algorithms for problem I.

Well, let's now focus on the previous problem. If you know that methodology called transposition principle (aka. Tellegen's principle), you will immediately see that the problem is equivalent to the following: Given a polynomial $$$A(x)$$$ of degree $$$n$$$, compute the coefficients of $$$A(x + 1/x)$$$.

Well. If you don't know that, hope there will be someone writing a blog post about that. And just forget the previous problem, and focus on the problem of computing the coefficients of $$$A(x + 1/x)$$$. :)

So here comes our protagonist of this post: the composition problem.

Composition problem

The problem is described as follows. Given a polynomial $$$A(x)$$$ of degree $$$n$$$, and a rational function $$$P(x)/Q(x)$$$, where $$$P, Q$$$ has a constant degree, compute the coefficients of the polynomial $$$A(P(x)/Q(x))$$$. Here we mean that, since $$$Q(x)^n A(P(x)/Q(x))$$$ must be a polynomial, the goal is to compute its coefficients.

Another intuitive interpretation is that, $$$A(P(x)/Q(x))$$$ is a rational function, and we want to compute the coefficients of its numerator. This naturally leads to the following general algorithm: Consider a divide-and-conquer algorithm, let $$$A(x) = A_0(x) + x^{n/2} A_1(x)$$$, we compute $$$A_0(P/Q)$$$ and $$$A_1(P/Q)$$$ separately, and merge them via

$$$ A(P/Q) = A_0(P/Q) + (P/Q)^{n/2}A_1(P/Q). $$$

The time complexity of this algorithm is given by the recurrence relation $$$T(n) = 2T(n/2) + O(n\log n)$$$, which has solution $$$O(n\log^2 n)$$$.

Several years ago I started to think about this problem. Unfortunately, I don't know how to do better than this in general. However, I do discover some special cases that have $$$O(n\log n)$$$ algorithms, where indeed handle the case $$$x+1/x = (1 + x^2)/x$$$.

You might suspect that there is some magic theory behind this, but actually I don't have. The only thing I did was to find some algebraic tricks, and applying them just works.

Building Blocks

We first start with several building blocks that are useful for the composition problem.

The first case is $$$A(x+c)$$$. Actually this is almost the only nontrivial tool we can use to speed up the composition problem. I guess most of the audience already knows how to do this in $$$O(n\log n)$$$, but for completeness I will still present the idea here. Just look at the expansion

$$$ b_k = \sum_j a_j \binom j k = \frac 1{k!}\sum_j a_j j! \cdot \frac{c^{j-k}}{(j-k)!}, $$$

which is a convolution.

The second case is $$$A(x^k)$$$. This only requires $$$O(n)$$$ time, since it essentially needs no arithmetic operations.

The third case is $$$A(\alpha x)$$$, where $$$\alpha$$$ is a constant. This can be done by scaling the coefficients of $$$A(x)$$$, which also takes $$$O(n)$$$ time.

Basic Tricks

With the three building blocks above, we can already handle something that seems to be much more.

The first problem is $$$A(x^2 + ax + b)$$$, i.e., the composition of $$$A(x)$$$ with a quadratic polynomial. This can be done by an algebraic trick:

$$$ x^2 + ax + b = \left(x + \frac a 2\right)^2 + \left(b - \frac{a^2}{4}\right). $$$

What does this mean? This means that we can compute $$$A(x^2 + ax + b)$$$ by computing the compositions

$$$ A(x^2 + ax + b) = A(x) \circ \left(x +\left( b - \frac{a^2}{4}\right)\right) \circ x^2 \circ \left(x + \frac a 2\right) $$$

in the order from left to right.

The second problem is $$$A((ax+b)/(cx+d))$$$, i.e., the Mobius transformation. Just look at

$$$ \frac{ax+b}{cx+d} = \frac ac + \frac{bc-ad}{c(cx+d)}, $$$

and we can compute $$$A((ax+b)/(cx+d))$$$ by computing the compositions

$$$ A\left(\frac{ax+b}{cx+d}\right) = A(x) \circ \left(x + \frac ac\right) \circ \frac{bc-ad}{c}x^{-1} \circ (cx+d). $$$

Well, there are some technical details on how to compute the compositions after composing $$$x^{-1}$$$, but these are just tedious and I think you clever readers can figure them out by yourself. :)

The Trick for $$$A(x + 1/x)$$$

Now we are going to solve $$$A(x + 1/x)$$$, and here is the trick:

$$$ x + 1/x = 2\frac{x+1}{x-1} \circ x^2 \circ \frac{x+1}{x-1}. $$$

Using this special case, we actually can handle general $$$A(P(x)/Q(x))$$$ in $$$O(n\log n)$$$ time, where $$$P, Q$$$ are quadratic polynomials.

For $$$(ax^2 + bx + c) / (dx^2 + ex + f)$$$, we can first move out a constant term, then reduce to the cases $$$(b' x + c') / (dx^2 + ex + f)$$$.

Then substitute by $$$x - c'/b'$$$, we reduce to the case $$$b'x / (d'x^2 + e'x + f')$$$.

Then take reciprocal, we reduce to the case $$$(d'x^2 + e'x + f') / b'x = (d'x^2 + f') / b'x + e'/b'$$$.

Seems that $$$d'x/b' + f'x^{-1}/b'$$$ is different with $$$x+x^{-1}$$$, but after scaling, we can again reduce to the case $$$x+x^{-1}$$$.

Well, for some special cases, we cannot do division, but then it can also be reduced to other cases, such as composing a quadratic polynomial, which is also already settled.

(Acknowledgement: This blog (Chinese) reminded me that it's possible to settle all $$$(ax^2 + bx + c) / (dx^2 + ex + f)$$$).

The Trick for $$$A(x^3 + ax^2 + bx + c)$$$

Actually, we can also handle $$$A(x^3 + ax^2 + bx + c)$$$ in $$$O(n\log n)$$$ time.

The first two tricks are just applying $$$x+c$$$, then we can figure out that the only thing we need to handle is to compute $$$A(x^3 + cx)$$$. Since we can also do scaling, we only need to solve some specific $$$c$$$.

The game ends with the following trick:

$$$ (x^3 - 3x) \circ \frac{x+x^{-1}}{2} = \frac{x^3+x^{-3}}{2}. $$$

This means that

$$$ x^3-3x = \frac{x+x^{-1}}{2} \circ x^3 \circ \left(\frac{x+x^{-1}}{2}\right)^{\langle -1\rangle}. $$$

You may ask, how can I composite $$$\left(\frac{x+x^{-1}}{2}\right)^{\langle -1\rangle}$$$? Well, the problem is that it will occur the composition $$$x^{1/2}$$$ when you expand the composition chain of $$$x+x^{-1}$$$, but at that time, it can be proved that all terms of the polynomial have the form $$$x^{2k}$$$, then the composition $$$x^{1/2}$$$ has clear meanings.

Well, there are many technical details that I omitted, for example, it will occur in some essential cases that we need to do scaling with a parameter is a square root of something. This can be handled by introducing a formal element $$$R = \sqrt \alpha$$$, and there might be some other technical details that I omitted.

Whatever, the long chain of composition and the technical details make the $$$O(n\log n)$$$ algorithm not so practical according to my experimental test, but I think it is still an interesting result.

Can we go further?

Is it possible to solve the composition problem in $$$O(n\log n)$$$ time for general rational functions? At least for the current techniques we have, the answer is no.

Here is an argument. What we have done before is to reduce the composition problem to the composition chain with the Mobius transform $$$(ax+b)/(cx+d)$$$, and $$$x^k$$$, and $$$x^{1/k}$$$. Note that if we want to solve the equation $$$P(x) = 0$$$. Suppose $$$P(x)$$$ has such a composition chain, then we can solve the equation $$$P(x) = 0$$$ by solving the composition chain step by step, this will result in an expression of the roots of $$$P(x)$$$ in radicals. However, by the celebrated Abel–Ruffini theorem, we know that there are polynomials (actually general polynomials of degree $$$\geq 5$$$) that cannot be solved by radicals, and this means that these polynomials definitely cannot be represented by the composition chain we want!

So there might be some interesting questions: Try to find some characterization of the polynomials that can be represented by the composition chain. This seems to let me recall the Galois theory, or a characterization of the composition chain of a polynomial by polynomials, known as the "Ritt decomposition theorem". However, they seem to not be directly related to the composition chain we are talking about, I wonder if there might or might not be some well-developed mathematical theory behind this problem.

Another question is, can we find some other building blocks that can be used to solve the composition problem, thus extending our ability to solve the composition problem for more cases? Or more ambitiously, can we find an entirely different algorithm that can solve the composition problem in $$$O(n\log n)$$$ time for general rational functions?

I'll be happy to see if any of these questions have positive or negative answers.

Anyway, have a nice day!

Comments (8)

Write comment?

adamant

9 months ago, # |

+17

you will immediately see that the problem is equivalent to the following
If you don't know that, hope there will be someone writing a blog post about that

Very briefly, transposition principle simply says that if there is a scheme that computes $$$Ax$$$ for any $$$x$$$ in $$$T(n)$$$, then there is also a scheme that computes $$$A^\top x$$$ for any $$$x$$$ in $$$T(n) + O(n)$$$, or something like that. This is done by representing $$$A = A_1 \dots A_n$$$, a sequence of elementary linear transforms, and then representing the transpose as $$$A^\top = A_n^\top \dots A_1^\top$$$.

Finding $$$A(B(x))$$$ for any $$$A(x)$$$ with a fixed $$$B(x)$$$ means finding the sum of $$$a_i [x^j] B^i(x)$$$ over $$$i$$$ for every $$$j$$$. Its transpose would be finding the sum of $$$b_j [x^j] B^i(x)$$$ over $$$j$$$ for every $$$i$$$, which is exactly what happens in this article's example.

I have a draft which is waiting its time for 21 months already, but I really struggle to wrap it up.

I guess it would really help if someone could provide anything of the following:

Some examples where it is useful. There are standard examples, like that polynomial interpolation/evaluation are transpose of partial fraction decomposition, or how Bostan-Mori is the transpose of $$$x^n \bmod q$$$. But those are really more like "fun facts" rather than actual examples of problems where it significantly simplifies or improves something.
How do people actually use it? As in, do you implement a transposed multiplication and use it as a building block for other transposed algorithms by meticulously transposing every elementary operation, like Tellegen's principle suggest? Or do you somehow just keep the idea in mind while solving the transposed problem directly?
It seems that $$$(A^\top)^{-1} = (A_1^\top)^{-1} \dots (A_n^\top)^{-1}$$$ tends to "resemble" $$$A$$$ in some way, as it even uses the operations in the same order. So, you could guess that e.g. if applying $$$A$$$ involves a particular technique, applying $$$(A^\top)^{-1}$$$ should be doable with a "similar" technique. I think on practice it means that if the "direct" computation of $$$A$$$ is cheap and $$$A^{-1}$$$ is not, it is likely that computing $$$A^{-1}$$$ as the transpose of $$$(A^\top)^{-1}$$$ would be faster than the direct computation. Are there more examples of this?

→ Reply

Kapt

9 months ago, # ^ |

← Rev. 2 →

+25

I've only seen one problem with Tellegen's principle in the wild: 102978D - Do Use FFT (though you probably have seen that one as well, seeing as your team solved it).

I implement mulT ($$$mul^t$$$ from the article) and build stuff on it. I imagine polynomials having infinite trailing zeroes and transpose elementary operations with polynomials (i.e. addition, multiplication by a constant polynomial), don't know if that counts as transposing every elementary operation or not.

From what I've seen in people's codes for multipoint, everyone implements mulT as well, though a lot of very strong coders for some reason implement mulT using standart multiplication instead of transposing algorithm for it, which strictly loses efficiency.

Everyone I know who knows this principle learnt it either from the problem I linked or from me and I couldn't find any resources explaining it other than the original article and I don't know how one is supposed to find out about it without seeing editorial for the linked problem. At some point I even considered writing a blog about it myself, so even if you think that you don't have enough examples or there is still some parts you haven't figured out yet, publishing a blog about it would still be very useful as it'd be the only easy to find resource and it'll likely attract attention of people who might understand the subject better

+10

Oh, we solved it without Tellegen's principle and I didn't look through the tutorial. Thanks for pointing it out!

My solution was like this. Consider

$$$ S_k(t) = \sum\limits_{i=1}^n C_i A_i^{n-t} \prod\limits_{j=1}^k (A_i + B_j), $$$

the problem asks to find $$$S_k(n)$$$ for all $$$k$$$. If we define $$$F_k(x) = \sum\limits_{t=0}^n S_k(t) x^t$$$, then

$$$ F_k(x) = (x + B_k) F_{k-1}(x). $$$

Now, we need to find $$$[x^n] F_k(x)$$$ for all $$$k$$$ from $$$1$$$ to $$$n$$$. First of all, we notice that

$$$ F_0(x) = \sum\limits_{i=1}^n C_i \sum\limits_{t=0}^n A_i^{n-t} x^t = \left(\sum\limits_{i=1}^n \frac{C_i}{1-A_i x} \bmod{x^{n+1}}\right)^R, $$$

which is computable in $$$O(n \log^2 n)$$$ with D&C as a sum of rational functions. After this, we notice that $$$[x^n] F_i(x)$$$ only depends on $$$[x^{(n-t)..n}] F_{i-t}(x)$$$, rather than on the full $$$F_{i-t}(x)$$$, and this allows us to find all $$$[x^n] F_i(x)$$$ recursively in $$$O(n \log^2 n)$$$ if we know $$$F_0(x)$$$ and sub-products of $$$(x+B_1)\dots(x+B_n)$$$.

Now that I compare it with the intended solution, operations seem quite similar, but I don't really see where Tellegen's principle come into play.

As for how people find out about it, one other way is by seeing that Chinese people (well, mostly Elegia actually) mention it in the comments every now and then xD (e.g. here).

nor

About examples, ecnerwala's comment here explains how the transposition principle can be used to motivate the construction used while computing the prefix sums of the Dirichlet convolution of two not-necessarily-multiplicative functions. Hopefully this is an example of what you're looking for.

tranquility

By the way, if we forget about the polynomial $$$A$$$ and solve a simpler problem — just calculate the free coefficient of the polynomial (ax^2+bx+c+d/x+e/x^2)^n for all n, then it can be done in linear time, because they satisfy a polynomial recurrence. Maybe this can help somehow?

Wielomian

Very interesting and educational blog, thank you! :) One small thing: in the section for $$$A(x + c)$$$ you dropped the coefficient $$$c^{k-j}$$$ under the summation :)

Elegia

+11

Thanks for the correction!