Preparation: adjunction of \(\sqrt{5}\)

As already known, the explicit formula of \(a_n\) can be represented using some constants \(c_1\) and \(c_2\) as:

\[a_n = c_1 \left(\frac{1 + \sqrt{5}}{2}\right)^n + c_2 \left(\frac{1 - \sqrt{5}}{2}\right)^n.\]

Although we want to use this beautiful formula to solve the problem, the issue arises when it comes to representing \(\sqrt{5} \bmod 998244353\). If there is \(\sqrt{5}\) in the modulo-\(998244353\) world, i.e. there exists an integer \(n\) such that \(n^2 = 5 \pmod{998244353}\), then we can replace \(\sqrt{5}\) with \(n\), but unfortunately, such an \(n\) does not exist.

Stalemate? No. It can be resolved by a mathematical concept called adjunction. Roughly speaking, we can manage the value in the following form: (modint998244353) + (modint998244353) \(\times \sqrt{5}\), which can actually be handled just as the ordinary modint998244353. Such an algebraic structure is called a “field obtained by adjoining \(\sqrt{5}\).”

• A bit more formally, for \(F = \lbrace a + b \sqrt{5} \vert a, b \in F_{998244353} \rbrace\) with arithmetic operations properly defined, \(F\) can be proven to form a field, so any algorithm against a field can be applied against \(F\) too.

Adjunction occasionally occur in solutions of advanced problems (especially, adjunction of the \(k\)-th root of \(1\) is typical), so it is worth learning.

Solution

From now on, consider that a value is an element of \(F = \lbrace a + b \sqrt{5} \vert a, b \in F_{998244353} \rbrace\).

\(a_n\) can be represented as

\[a_n = c_1 A^n + c_2 B^n\]

for some constants \(A, B, c_1\), and \(c_2\). (\(c_1\) and \(c_2\) can be found by evaluating the equation for \(n=1, 2\).)

Let \(\mathrm{mod} = 998244353\). Nontrivally, we can assert that

\[A^{2(\mathrm{mod} + 1)} = B^{2(\mathrm{mod} + 1)} = 1,\]

so the period of the sequence \(a_n\) is a divisor of \(2(\mathrm{mod} + 1)\). Thus, if \(N \gt 2(\mathrm{mod} + 1)\), the answer the \(q\)-th power of the answer for \(N = 2 (\mathrm{mod} + 1)\), multiplied by the answer for \(N = r\), where \(N = q \times 2 (\mathrm{mod} + 1) + r\). Thus, it is sufficient if we can solve the problem for \(N \leq 2 (\mathrm{mod} + 1)\). From no on, we assume \(N \leq 2 (\mathrm{mod} + 1)\).

Let \(D = \frac{A}{B}\); then

\[\frac{a_n}{B^n} = c_1 D^n + c_2,\]

so it is sufficient to compute \(\prod_{n=1}^N (c_1 D^n + c_2)\) to evaluate the sought expression; let us consider how to evaluate it.

Let \(M = \lfloor \sqrt{N} \rfloor\). If we can solve the case for \(N = M^2\), all that left is to find \(\mathrm{O}(\sqrt{N})\) more elements \(a_n\) to find the answer, so hereinafter we consider \(N = M^2\) case.

\[ \begin{aligned} &\prod_{n=1}^{M^2} (c_1 D^n + c_2) \\ &= \prod_{i=0}^{M-1} \prod_{j=1}^{M} (c_1 D^{iM + j} + c_2) \\ &= \prod_{i=0}^{M-1} \left\vert \prod_{j=1}^M (c_1 D^j X + c_2) \right \vert_{X=D^{iM}}. \end{aligned} \]

Here, \(\vert f(X) \vert_{X=n}\) means assigning \(X=n\) to \(f(X)\). Let \(F(X) = \prod_{j=1}^M (c_1 D^j X + c_2)\), then the following two steps are required to find the sought value:

1. expand \(F(X)\);
2. evaluate \(F(X)\) at \(X = D^0, D^M, D^{2M}, \dots, D^{M(M-1)}\).

We explain each step.

The former step can be accomplished in a total of \(\mathrm{O}(M \log M)\) in the manner of binary lifting.

Let \(F_n(X) = \prod_{j=1}^n (c_1 D^j X + c_2)\). What we want to find is \(F_M(X)\).
If \(M\) is even, we have the following relation:

\[ \begin{aligned} F_M(X) &= \prod_{j=1}^M (c_1 D^j X + c_2) \\ &= \prod_{j=1}^{M/2} (c_1 D^j X + c_2) \prod_{j=M/2+1}^M (c_1 D^j X + c_2) \\ &= F_{M/2}(X) F_{M/2} (D^{M/2} X). \end{aligned} \]

This allows us to find \(F_M(X)\) based on \(F_{M/2}(X)\) in \(\mathrm{O}(M \log M)\) time. If \(M\) is odd, we can use the relation \(F_M(X) = F_{M-1}(X) (c_1 D^M X + c_2)\). By implementing them, \(F_M(X)\) can be expanded in a total of \(\mathrm{O}(M \log M)\) time.

Now we consider the latter step. It essentially asks to do multi-point evaluation of a polynomial on a geometric sequence. While an ordinary multi-point evaluation costs \(\mathrm{O}(n \log^2 n)\) time, where \(n\) is the degree of the polynomial and the number of points to evaluate, but that on a geometric sequence can be done in \(\mathrm{O}(M \log M)\) time with an algorithm called chirp z-transform. Article in Kyopro Encyclopedia of Algorithms (Japanese), Article in CodeForces (English).

By implementing it appropriately, the original problem can be solved in a total of \(\mathrm{O}(\sqrt{\mathrm{mod}} \log \mathrm{mod})\) time, which is fast enough.