本题没有特别困难的内容，但综合了许多初等数论的技巧，步骤较多。
考虑长度为n的序列的最后一格子，若最后一格为数字方格，则将其去除后得一长度为$n-1$的序列，若最后一格为空白方格，则将其去除后得一长度为$n-2$的序列，数字方格共有10种，因此有
$$T(n)=10T(n-1)+T(n-2).$$
令
$$F(n)=\begin{cases} 0 & n\le 0, \\ 1 & n=1, \\ T(n-1) & n\ge 2.\end{cases}$$
则也有
$$F(n)=10F(n-1)+F(n-2).$$
接下来证明这样一个二阶线性递推式的一般性质：
定理1：若$u,v$互素，且$a_n=\begin{cases} 0 & n\le 0, \\ 1 & n=1, \\ ua_{n-1}+va_{n-2} & n\ge 2\end{cases}$, 则对任意$m,n\ge 1$都有$\gcd(a_m,a_n)=a_{\gcd(m,n)}$.
为此需要先证明一些引理，令
$$A=\begin{pmatrix}u & v \\ 1& 0\end{pmatrix},$$
则
$$\begin{pmatrix}a_n \\ a_{n-1}\end{pmatrix}=A\begin{pmatrix}a_{n-1} \\ a_{n-2}\end{pmatrix}=\ldots=A^{n-2}\begin{pmatrix}a_2 \\ a_1 \end{pmatrix}=A^{n-2}\begin{pmatrix}u \\ 1 \end{pmatrix}.$$
从而有
$$A^{n-2}=\begin{pmatrix}a_{n-1} & va_{n-2} \\ a_{n-2} & va_{n-3} \end{pmatrix},$$
即
$$A^{n}=\begin{pmatrix}a_{n+1} & va_{n} \\ a_{n} & va_{n-1} \end{pmatrix}.$$
考虑等式$A^{m+n}=A^mA^n$:
$$\begin{pmatrix}a_{m+n+1} & va_{m+n} \\ a_{m+n} & va_{m+n-1} \end{pmatrix}=A^{m+n}=A^mA^n=\begin{pmatrix}a_{m+1} & va_{m} \\ a_{m} & va_{m-1} \end{pmatrix}\begin{pmatrix}a_{n+1} & va_{n} \\ a_{n} & va_{n-1} \end{pmatrix},$$
比较左下角的元素可得：
引理1： $a_{m+n}=a_ma_{n+1}+va_{m-1}a_n$。
依次取$n=m,2m,\ldots$即可归纳证明：
引理2：对任意$k\in\mathbb N$，$a_m$均可整除$a_{km}$。
另一方面，由递推式$a_n=ua_{n-1}+va_{n-2}$和$\gcd(u,v)=1$很容易归纳检验对任意$n\ge 1$都有$\gcd(a_n,v)=1$，从而再用辗转相除法可得：
$$\gcd(a_{n},a_{n-1})=\gcd(ua_{n-1}+va_{n-2},a_{n-1})=\gcd(va_{n-2},a_{n-1})=\gcd(a_{n-2},a_{n-1})=\ldots=\gcd(a_2,a_1)=1.$$
进一步结合引理2则有
引理3：对任意$k\in\mathbb N$，$\gcd(a_m,a_{km+1})=1$。
最后，为证明定理1，不妨设$n\ge m$，并将其写作$n=km+r=km+1+r-1$，其中$r\in\{0,1,\ldots,m-1\}$, $k\in\mathbb N$，则有
$$\gcd(a_m,a_n)=\gcd(a_m,a_{km+1}a_{r}+va_{km}a_{r-1})=\gcd(a_m,a_{km+1}a_{r})=\gcd(a_m,a_r),$$
其中首个等号系引理1，第二个等号因引理2成立，第三个等号因引理3成立。上式即在下标上作辗转相除法，因此最终会得到$\gcd(a_d,a_{k'd})$的形式，其中$d=\gcd(m,n)$, $k'\in\mathbb N$，由引理2可知$\gcd(a_d,a_{k'd})=a_d$，从而定理1得证。
应用定理1，可得
$$S(L)=\sum_{1\le a,b,c\le L}F_{\gcd(c^a+1,c^b+1)}.$$
要计算$\gcd(c^a+1,c^b+1)$，不妨设$\gcd(a,b)=1$，否则可令$c'=c^{\gcd(a,b)}$。不妨再设$a\ge b$，并将其写作$a=kb+r$，其中$k\in\mathbb N$, $r\in\{0,1,\ldots,b-1\}$。
此时有两种情况：
（1）若$b=1$，则显然a为奇数时$\gcd(c^a+1,c+1)=c+1$，a为偶数时由$\gcd(c,c+1)=1$可得
$$\gcd(c^a+1,c+1)=\gcd(c^a-c,c+1)=\gcd(c^{a-1}-1,c^b+1)=\gcd(c^{a-1}+c,c+1)=\gcd(c^{a-2}+1,c^b+1)=\ldots=\gcd(2,c^b+1)=\begin{cases}2 & c\text{为奇数} \\ 1 & c\text{为偶数}\end{cases}$$
（2）若$b>1$，那么r,b互素，注意到$\gcd(c^b,c^b+1)=1$作辗转相除可得
$$\gcd(c^a+1,c^b+1)=\gcd(c^a-c^b,c^b+1)=\gcd(c^{a-b}-1,c^b+1)=\gcd(c^{a-b}+c^b,c^b+1)$$
若$k\ge2$，那么上式$=\gcd(c^{a-2b}+1,c^b+1)=\ldots=\begin{cases}c^{r} & k\text{为偶数} \\ c^{b+r} & k\text{为奇数} \end{cases}$，否则上式$=\gcd(c^{r}+c^b,c^b+1)=\gcd(c^{b-r}+1,c^b+1)$。
不论何者都又转化为指数上的辗转相除，且指数上始终均为奇数或一奇一偶，因此最终又转化为(1)的情况。
综上所述，若$d=\gcd(a,b)$，则有
$$\gcd(c^a+1,c^b+1)=\begin{cases}c^{d}+1 &\text{若}a/d, b/d\text{均为奇数} \\ 2 & \text{若}a/d, b/d\text{至少有一个偶数，且}c\text{为奇数} \\ 1 & \text{若}a/d, b/d\text{至少有一个偶数，且}c\text{为偶数} \end{cases}$$
因此先遍历a,b，统计每种情况的频次，再遍历c求出对应的$F_{\gcd(c^a+1,c^b+1)}$，此处用线性递推式的特征方程很容易得
$$F_n=\frac{(5+\sqrt{26})^n-(5-\sqrt{26})^n}{2\sqrt{26}},$$
而26是987898789的二次剩余：
$$445663912^2=26\pmod{987898789}$$
2和445663912在模987898789乘法群中的逆分别是493949395和625078636，因此有
$$F_n \bmod 987898789=\left((5+445663912)^n-(5-445663912)^n\right)\cdot493949395 \cdot625078636 \bmod 987898789$$
最后，由于987898789是素数，指数上的n可以模去$\varphi(987898789)=987898788$再用快速幂计算，汇总即得最终结果$970746056$。
注：因使用sympy计算2和445663912的乘法逆，以及便于用内建函数作快速幂，以下代码为Python 3，代码中还打印了进度信息。
import sympy,math
target=2000
mod=987898789
m=445663912
inv2=sympy.mod_inverse(2,mod) #493949395
invM=sympy.mod_inverse(m,mod) #625078636
F=lambda n:((pow(5+m,n,mod)-pow(5-m,n,mod))*inv2*invM)%mod
s=0
q=[0]*(target+1)
for a in range(1,target+1):
for b in range(1,target+1):
d=math.gcd(a,b)
if (a//d)%2==(b//d)%2==1:
q[d]+=1
else:
s+=1
print("initialization completed, start computing....")
res=(s*(F(1)+F(2))*target//2)%mod
for d in range(1,target+1):
if q[d]>0:
res=(res+q[d]*sum(F(pow(c,d,mod-1)+1) for c in range(1,target+1)))%mod
print(f"d={d} checked, current sum: {res}")
print(res)
|
