显然能表示为某个$A_n$的连续子列需要以1开头。考察以1开头，长度为k的连续子列，假设这样的子列共有$L_k$种，则根据OEIS A005942中的信息可得
$$L_{2k}=L_k+L_{k+1}, \quad L_{2k+1}=2L_{k+1}$$
因此给定$n$，首先需要根据上述信息确定$A_n$二进制表示的位数$k_n$，然后找出在所有以1开头且长度为$k_n$连续子列中它是第几个数。
最后，根据OEIS A010060中的信息可知，若将Thue-Morse序列分割为长度为$4$的子列的拼接，那么只有0110和1001两种子列。更一般地，若将其分割为长度为$2^t$的子列的拼接，那么同样也只有两种子列不断重复。
利用这一性质就可以快速定位到$A_n$，不过具体实现仍很繁琐，以下代码改编自官方论坛，最终结果是$178476944$。
from math import sqrt
MOD=10**9
powers=[2]
chunks=[(0,1)]
ptable=[0,1]
for i in range(1, 64):
a,b=chunks[-1]
chunks.append(((a*powers[-1]+b)%MOD,(b*powers[-1]+a)%MOD))
powers.append((powers[-1]*powers[-1])%MOD)
for x in range(1, 128):
ptable.append(ptable[x])
ptable.append(1-ptable[x])
def getPow(n):
x=1
while x<=n:
x<<=1
return x>>2
def getLen(n):
if n<4:
return (0,2,2,3)[n]
p=getPow(n-1)
return n-1+2*p if n-1>3*p else 2*n-2-p
def getSum(n):
if n<5:
return (0,0,2,4,7)[n]
p=getPow(n-1)
return 4+(19*p//2+5)*(p-1)//3+((n-1-2*p)*(n+p-2) if n-1<3*p else (n-3*p-1)*(n-3*p-2)//2+p*(5*n-11*p-6))
def invPos(n):
a,b=1,int(sqrt(n*1.5))+2
while a<b-1:
c=(a+b)//2
if (getSum(c)>n):
b=c
else:
a=c
return (a,n-getSum(a))
def findSeq(n,offset,bit):
if n<4:
return (((-1,-1)),((1,0),(0,1)),((0,2),(5,1)),((0,4),(3,2),(5,1)))[n][offset][bit]
return 2*findSeq((n+1)//2,offset,bit) if offset<getLen((n+1)//2) else (2*findSeq(n//2+1,getLen(n)-offset-1,1-bit)+1)
def parity(n):
while n>=(1<<32):
n=(n>>32)^(n&((1<<32)-1))
n=(n>>16)^(n&65535)
n=(n>>8)^(n&255)
return ptable[n]
def value(n):
n,offset=invPos(n)
pos=findSeq(n,offset,1)
x=0
for p in range(pos,pos+n):
x<<=1
if parity(p):
x|=1
return x
res=0
for p in range(1,19):
n,offset=invPos(10**p)
pos=findSeq(n,offset,1)
end=pos+n
x=0
shift=-1
while True:
shift+=1
a=1<<shift
if a+pos>end:
break
if a&pos:
x*=powers[shift]
x+=chunks[shift][parity(pos)]
x%=MOD
pos+=a
t=pos^end
while shift:
shift-=1
a=1<<shift
if a+pos<=end and a&t:
x*=powers[shift]
x+=chunks[shift][parity(pos)]
x%=MOD
pos+=a
t^=a
res+=x
print res%MOD
|
