結果
| 問題 |
No.2305 [Cherry 5th Tune N] Until That Day...
|
| コンテスト | |
| ユーザー |
👑 Kazun
|
| 提出日時 | 2022-10-27 20:02:14 |
| 言語 | PyPy3 (7.3.15) |
| 結果 |
AC
|
| 実行時間 | 3,285 ms / 10,000 ms |
| コード長 | 23,679 bytes |
| コンパイル時間 | 529 ms |
| コンパイル使用メモリ | 82,176 KB |
| 実行使用メモリ | 130,768 KB |
| 最終ジャッジ日時 | 2024-12-15 19:19:35 |
| 合計ジャッジ時間 | 14,530 ms |
|
ジャッジサーバーID (参考情報) |
judge2 / judge5 |
(要ログイン)
| ファイルパターン | 結果 |
|---|---|
| other | AC * 21 |
ソースコード
class Modulo_Polynomial():
__slots__=("Poly", "max_degree")
def __init__(self, Poly=[], max_degree=2*10**5):
""" 多項式の定義
P: 係数のリスト
max_degree
※Mod: 法はグローバル変数から指定
"""
if Poly:
self.Poly=[p%Mod for p in Poly[:max_degree]]
else:
self.Poly=[0]
self.max_degree=max_degree
def __str__(self):
return str(self.Poly)
def __repr__(self):
return self.__str__()
def __iter__(self):
yield from self.Poly
#=
def __eq__(self,other):
from itertools import zip_longest
return all([a==b for a,b in zip_longest(self.Poly,other.Poly,fillvalue=0)])
#+,-
def __pos__(self):
return self
def __neg__(self):
return self.scale(-1)
#items
def __getitem__(self, index):
if isinstance(index, slice):
return Modulo_Polynomial(self.Poly[index], self.max_degree)
else:
if index<0:
raise IndexError("index is negative (index: {})".format(index))
elif index>=len(self.Poly):
return 0
else:
return self.Poly[index]
def __setitem__(self, index, value):
if index<0:
raise IndexError("index is negative (index: {})".format(index))
elif index>=self.max_degree:
return
if index>=len(self.Poly):
self.Poly+=[0]*(index-len(self.Poly)+1)
self.Poly[index]=value%Mod
#Boole
def __bool__(self):
return any(self.Poly)
#簡略化
def reduce(self):
""" 高次の 0 を切り捨て
"""
P=self.Poly
for d in range(len(P)-1,-1,-1):
if P[d]:
break
self.resize(d+1)
return
#シフト
def __lshift__(self,other):
if other<0:
return self>>(-other)
if other>self.max_degree:
return Modulo_Polynomial([0],self.max_degree)
G=[0]*other+self.Poly
return Modulo_Polynomial(G,self.max_degree)
def __rshift__(self,other):
if other<0:
return self<<(-other)
return Modulo_Polynomial(self.Poly[other:],self.max_degree)
#次数
def degree(self):
P=self.Poly
for d in range(len(self.Poly)-1,-1,-1):
if P[d]:
return d
return -float("inf")
#加法
def __add__(self,other):
P=self; Q=other
if Q.__class__==Modulo_Polynomial:
N=min(P.max_degree,Q.max_degree)
A=P.Poly; B=Q.Poly
else:
N=P.max_degree
A=P.Poly; B=Q
return Modulo_Polynomial(Calc.Add(A,B),N)
def __radd__(self,other):
return self+other
#減法
def __sub__(self,other):
P=self; Q=other
if Q.__class__==Modulo_Polynomial:
N=min(P.max_degree,Q.max_degree)
A=P.Poly; B=Q.Poly
else:
N=P.max_degree
A=P.Poly; B=Q
return Modulo_Polynomial(Calc.Sub(A,B),N)
def __rsub__(self,other):
return (-self)+other
#乗法
def __mul__(self,other):
P=self
Q=other
if Q.__class__==Modulo_Polynomial:
a=b=0
for x in P.Poly:
if x:
a+=1
for y in Q.Poly:
if y:
b+=1
if a>b:
P,Q=Q,P
P.reduce();Q.reduce()
U,V=P.Poly,Q.Poly
M=min(P.max_degree,Q.max_degree)
if a<2*P.max_degree.bit_length():
B=[0]*(len(U)+len(V)-1)
for i in range(len(U)):
if U[i]:
for j in range(len(V)):
B[i+j]+=U[i]*V[j]
if B[i+j]>Mod:
B[i+j]-=Mod
else:
B=Calc.Convolution(U,V)[:M]
B=Modulo_Polynomial(B,M)
B.reduce()
return B
else:
return self.scale(other)
def __rmul__(self,other):
return self.scale(other)
#除法
def __floordiv__(self,other):
if not other:
raise ZeroDivisionError
if isinstance(other,int):
return self/other
self.reduce()
other.reduce()
return Modulo_Polynomial(Calc.Floor_Div(self.Poly, other.Poly),
max(self.max_degree, other.max_degree))
def __rfloordiv__(self,other):
if not self:
raise ZeroDivisionError
if isinstance(other,int):
return Modulo_Polynomial([],self.max_degree)
#剰余
def __mod__(self,other):
if not other:
return ZeroDivisionError
self.reduce(); other.reduce()
r=Modulo_Polynomial(Calc.Mod(self.Poly, other.Poly),
min(self.max_degree, other.max_degree))
r.reduce()
return r
def __rmod__(self,other):
if not self:
raise ZeroDivisionError
r=other-(other//self)*self
r.reduce()
return r
def __divmod__(self,other):
q=self//other
r=self-q*other; r.reduce()
return (q,r)
#累乗
def __pow__(self,other):
if other.__class__==int:
n=other
m=abs(n)
Q=self
A=Modulo_Polynomial([1],self.max_degree)
while m>0:
if m&1:
A*=Q
m>>=1
Q*=Q
if n>=0:
return A
else:
return A.inverse()
else:
P=Log(self)
return Exp(P*other)
#逆元
def inverse(self, deg=None):
assert self.Poly[0], "定数項が0"
if deg==None:
deg=self.max_degree
return Modulo_Polynomial(Calc.Inverse(self.Poly, deg), self.max_degree)
#除法
def __truediv__(self,other):
if isinstance(other, Modulo_Polynomial):
if Calc.is_sparse(other.Poly):
d,f=Calc.coefficients_list(other.Poly)
K=len(d)
H=[0]*self.max_degree
alpha=pow(other[0], Mod-2, Mod)
H[0]=alpha*self[0]%Mod
for i in range(1, self.max_degree):
c=0
for j in range(1, K):
if d[j]<=i:
c+=f[j]*H[i-d[j]]%Mod
else:
break
c%=Mod
H[i]=alpha*(self[i]-c)%Mod
H=Modulo_Polynomial(H, min(self.max_degree, other.max_degree))
return H
else:
return self*other.inverse()
else:
return pow(other,Mod-2,Mod)*self
def __rtruediv__(self,other):
return other*self.inverse()
#スカラー倍
def scale(self, s):
return Modulo_Polynomial(Calc.Times(self.Poly,s),self.max_degree)
#最高次の係数
def leading_coefficient(self):
for x in self.Poly[::-1]:
if x:
return x
return 0
def censor(self, N=-1, Return=False):
""" N 次以上の係数をカット
"""
if N==-1:
N=self.max_degree
N=min(N, self.max_degree)
if Return:
return Modulo_Polynomial(self.Poly[:N],self.max_degree)
else:
self.Poly=self.Poly[:N]
def resize(self, N, Return=False):
""" 強制的に Poly の配列の長さを N にする.
"""
N=min(N, self.max_degree)
P=self
if Return:
if len(P.Poly)>N:
E=P.Poly[:N]
else:
E=P.Poly+[0]*(N-len(P.Poly))
return Modulo_Polynomial(E,P.max_degree)
else:
if len(P.Poly)>N:
del P.Poly[N:]
else:
P.Poly+=[0]*(N-len(P.Poly))
#代入
def substitution(self, a):
""" a を (形式的に) 代入した値を求める.
a: int
"""
y=0
t=1
for p in self.Poly:
y=(y+p*t)%Mod
t=(t*a)%Mod
return y
#=================================================
class Calculator:
def __init__(self):
self.primitive=self.__primitive_root()
self.__build_up()
def __primitive_root(self):
p=Mod
if p==2:
return 1
if p==998244353:
return 3
if p==10**9+7:
return 5
if p==163577857:
return 23
if p==167772161:
return 3
if p==469762049:
return 3
fac=[]
q=2
v=p-1
while v>=q*q:
e=0
while v%q==0:
e+=1
v//=q
if e>0:
fac.append(q)
q+=1
if v>1:
fac.append(v)
g=2
while g<p:
if pow(g,p-1,p)!=1:
return None
flag=True
for q in fac:
if pow(g,(p-1)//q,p)==1:
flag=False
break
if flag:
return g
g+=1
#参考元: https://judge.yosupo.jp/submission/72676
def __build_up(self):
rank2=(~(Mod-1)&(Mod-2)).bit_length()
root=[0]*(rank2+1); iroot=[0]*(rank2+1)
rate2=[0]*max(0, rank2-1); irate2=[0]*max(0, rank2-1)
rate3=[0]*max(0, rank2-2); irate3=[0]*max(0, rank2-2)
root[-1]=pow(self.primitive, (Mod-1)>>rank2, Mod)
iroot[-1]=pow(root[-1], Mod-2, Mod)
for i in range(rank2)[::-1]:
root[i]=root[i+1]*root[i+1]%Mod
iroot[i]=iroot[i+1]*iroot[i+1]%Mod
prod=iprod=1
for i in range(rank2-1):
rate2[i]=root[i+2]*prod%Mod
irate2[i]=iroot[i+2]*prod%Mod
prod*=iroot[i+2]; prod%=Mod
iprod*=root[i+2]; iprod%=Mod
prod=iprod = 1
for i in range(rank2-2):
rate3[i]=root[i + 3]*prod%Mod
irate3[i]=iroot[i + 3]*iprod%Mod
prod*=iroot[i + 3]; prod%=Mod
iprod*=root[i + 3]; iprod%=Mod
self.root=root; self.iroot=iroot
self.rate2=rate2; self.irate2=irate2
self.rate3=rate3; self.irate3=irate3
def Add(self, A, B):
""" 必要ならば末尾に元を追加して, [A[i]+B[i]] を求める.
"""
if type(A)==int:
A=[A]
if type(B)==int:
B=[B]
m=min(len(A), len(B))
C=[(A[i]+B[i])%Mod for i in range(m)]
C.extend(A[m:])
C.extend(B[m:])
return C
def Sub(self, A, B):
""" 必要ならば末尾に元を追加して, [A[i]-B[i]] を求める.
"""
if type(A)==int:
A=[A]
if type(B)==int:
B=[B]
m=min(len(A), len(B))
C=[0]*m
C=[(A[i]-B[i])%Mod for i in range(m)]
C.extend(A[m:])
C.extend([-b%Mod for b in B[m:]])
return C
def Times(self,A, k):
""" [k*A[i]] を求める.
"""
return [k*a%Mod for a in A]
#参考元 https://judge.yosupo.jp/submission/72676
def NTT(self, A):
""" A に Mod を法とする数論変換を施す
※ Mod はグローバル変数から指定
References:
https://github.com/atcoder/ac-library/blob/master/atcoder/convolution.hpp
https://judge.yosupo.jp/submission/72676
"""
N=len(A)
H=(N-1).bit_length()
l=0
I=self.root[2]
rate2=self.rate2; rate3=self.rate3
while l<H:
if H-l==1:
p=1<<(H-l-1)
rot=1
for s in range(1<<l):
offset=s<<(H-l)
for i in range(p):
x=A[i+offset]; y=A[i+offset+p]*rot%Mod
A[i+offset]=(x+y)%Mod
A[i+offset+p]=(x-y)%Mod
if s+1!=1<<l:
rot*=rate2[(~s&-~s).bit_length()-1]
rot%=Mod
l+=1
else:
p=1<<(H-l-2)
rot=1
for s in range(1<<l):
rot2=rot*rot%Mod
rot3=rot2*rot%Mod
offset=s<<(H-l)
for i in range(p):
a0=A[i+offset]
a1=A[i+offset+p]*rot
a2=A[i+offset+2*p]*rot2
a3=A[i+offset+3*p]*rot3
alpha=(a1-a3)%Mod*I
A[i+offset]=(a0+a2+a1+a3)%Mod
A[i+offset+p]=(a0+a2-a1-a3)%Mod
A[i+offset+2*p]=(a0-a2+alpha)%Mod
A[i+offset+3*p]=(a0-a2-alpha)%Mod
if s+1!=1<<l:
rot*=rate3[(~s&-~s).bit_length()-1]
rot%=Mod
l+=2
#参考元 https://judge.yosupo.jp/submission/72676
def Inverse_NTT(self, A):
""" A を Mod を法とする逆数論変換を施す
※ Mod はグローバル変数から指定
References:
https://github.com/atcoder/ac-library/blob/master/atcoder/convolution.hpp
https://judge.yosupo.jp/submission/72676
"""
N=len(A)
H=(N-1).bit_length()
l=H
J=self.iroot[2]
irate2=self.rate2; irate3=self.irate3
while l:
if l==1:
p=1<<(H-l)
irot=1
for s in range(1<<(l-1)):
offset=s<<(H-l+1)
for i in range(p):
x=A[i+offset]; y=A[i+offset+p]
A[i+offset]=(x+y)%Mod
A[i+offset+p]=(x-y)*irot%Mod
if s+1!=1<<(l-1):
irot*=irate2[(~s&-~s).bit_length()-1]
irot%=Mod
l-=1
else:
p=1<<(H-l)
irot=1
for s in range(1<<(l-2)):
irot2=irot*irot%Mod
irot3=irot2*irot%Mod
offset=s<<(H-l+2)
for i in range(p):
a0=A[i+offset]
a1=A[i+offset+p]
a2=A[i+offset+2*p]
a3=A[i+offset+3*p]
beta=(a2-a3)*J%Mod
A[i+offset]=(a0+a1+a2+a3)%Mod
A[i+offset+p]=(a0-a1+beta)*irot%Mod
A[i+offset+2*p]=(a0+a1-a2-a3)*irot2%Mod
A[i+offset+3*p]=(a0-a1-beta)*irot3%Mod
if s+1!=1<<(l-2):
irot*=irate3[(~s&-~s).bit_length()-1]
irot%=Mod
l-=2
N_inv=pow(N,Mod-2,Mod)
for i in range(N):
A[i]=N_inv*A[i]%Mod
def non_zero_count(self, A):
""" A にある非零の数を求める. """
return len(A)-A.count(0)
def is_sparse(self, A, K=None):
""" A が疎かどうかを判定する. """
if K==None:
K=25
return self.non_zero_count(A)<=K
def coefficients_list(self, A):
""" A にある非零のリストを求める.
output: ( [d[0], ..., d[k-1] ], [f[0], ..., f[k-1] ]) : a[d[j]]=f[j] であることを表している.
"""
f=[]; d=[]
for i in range(len(A)):
if A[i]:
d.append(i)
f.append(A[i])
return d,f
def Convolution(self, A, B):
""" A, B で Mod を法とする畳み込みを求める.
※ Mod はグローバル変数から指定
"""
if not A or not B:
return []
N=len(A)
M=len(B)
L=M+N-1
if min(N,M)<=50:
if N<M:
N,M=M,N
A,B=B,A
C=[0]*L
for i in range(N):
for j in range(M):
C[i+j]+=A[i]*B[j]
C[i+j]%=Mod
return C
H=L.bit_length()
K=1<<H
A=A+[0]*(K-N)
B=B+[0]*(K-M)
self.NTT(A)
self.NTT(B)
for i in range(K):
A[i]=A[i]*B[i]%Mod
self.Inverse_NTT(A)
return A[:L]
def Autocorrelation(self, A):
""" A 自身に対して,Mod を法とする畳み込みを求める.
※ Mod はグローバル変数から指定
"""
N=len(A)
L=2*N-1
if N<=50:
C=[0]*L
for i in range(N):
for j in range(N):
C[i+j]+=A[i]*A[j]
C[i+j]%=Mod
return C
H=L.bit_length()
K=1<<H
A=A+[0]*(K-N)
self.NTT(A)
for i in range(K):
A[i]=A[i]*A[i]%Mod
self.Inverse_NTT(A)
return A[:L]
def Multiple_Convolution(self, *A):
""" A=(A[0], A[1], ..., A[d-1]) で Mod を法とする畳み込みを行う.
"""
from collections import deque
if not A:
return [1]
Q=deque(list(range(len(A))))
A=list(A)
while len(Q)>=2:
i=Q.popleft(); j=Q.popleft()
A[i]=self.Convolution(A[i], A[j])
A[j]=None
Q.append(i)
i=Q.popleft()
return A[i]
def Inverse(self, F, length=None):
if length==None:
M=len(F)
else:
M=length
if M<=0:
return []
if self.is_sparse(F):
"""
愚直に漸化式を用いて求める.
計算量: F にある係数が非零の項の個数を K, 求める最大次数を N として, O(NK) 時間
"""
d,f=self.coefficients_list(F)
G=[0]*M
alpha=pow(F[0], Mod-2, Mod)
G[0]=alpha
for i in range(1, M):
for j in range(1, len(d)):
if d[j]<=i:
G[i]+=f[j]*G[i-d[j]]%Mod
else:
break
G[i]%=Mod
G[i]=(-alpha*G[i])%Mod
del G[M:]
else:
"""
FFTの理論を応用して求める.
計算量: 求めたい項の個数をNとして, O(N log N)
Reference: https://judge.yosupo.jp/submission/42413
"""
N=len(F)
r=pow(F[0],Mod-2,Mod)
m=1
G=[r]
while m<M:
A=F[:min(N, 2*m)]; A+=[0]*(2*m-len(A))
B=G.copy(); B+=[0]*(2*m-len(B))
Calc.NTT(A); Calc.NTT(B)
for i in range(2*m):
A[i]=A[i]*B[i]%Mod
Calc.Inverse_NTT(A)
A=A[m:]+[0]*m
Calc.NTT(A)
for i in range(2*m):
A[i]=-A[i]*B[i]%Mod
Calc.Inverse_NTT(A)
G.extend(A[:m])
m<<=1
G=G[:M]
return G
def Floor_Div(self, F, G):
assert F[-1]
assert G[-1]
F_deg=len(F)-1
G_deg=len(G)-1
if F_deg<G_deg:
return []
m=F_deg-G_deg+1
return self.Convolution(F[::-1], Calc.Inverse(G[::-1],m))[m-1::-1]
def Mod(self, F, G):
while F and F[-1]==0:
F.pop()
while G and G[-1]==0:
G.pop()
if not F:
return []
return Calc.Sub(F, Calc.Convolution(Calc.Floor_Div(F,G),G))
def Polynominal_Coefficient(P,Q,N):
""" [X^N] P/Q を求める.
References:
http://q.c.titech.ac.jp/docs/progs/polynomial_division.html
https://arxiv.org/abs/2008.08822
https://arxiv.org/pdf/2008.08822.pdf
"""
P=P.Poly.copy(); Q=Q.Poly.copy()
m=1<<((len(Q)-1).bit_length())
P.extend([0]*(2*m-len(P)))
Q.extend([0]*(2*m-len(Q)))
while N:
R=[Q[i] if i&1==0 else -Q[i] for i in range(2*m)]
Calc.NTT(P); Calc.NTT(Q); Calc.NTT(R)
for i in range(2*m):
P[i]*=R[i]; P[i]%=Mod
Q[i]*=R[i]; Q[i]%=Mod
Calc.Inverse_NTT(P); Calc.Inverse_NTT(Q)
if N&1==0:
for i in range(m):
P[i]=P[2*i]
else:
for i in range(m):
P[i]=P[2*i+1]
for i in range(m):
Q[i]=Q[2*i]
for i in range(m,2*m):
P[i]=Q[i]=0
N>>=1
if Q[0]==1:
return P[0]
else:
return P[0]*pow(Q[0],Mod-2,Mod)%Mod
def Nth_Term_of_Linearly_Recurrent_Sequence(A, C, N, offset=0):
""" A[i]=C[0]*A[i-1]+C[1]*A[i-2]+...+C[d-1]*A[i-d] で表される数列 (A[i]) の第 N 項を求める.
A=(A[0], ..., A[d-1]): 最初の d 項
C=(C[0], ..., C[d-1]): 線形漸化式
N: 求める項数
offset: ずらす項数 (初項が第 offset 項になる)
"""
assert len(A)==len(C)
d=len(A)
N-=offset
if N<0:
return 0
elif N<d:
return A[N]%Mod
A=Modulo_Polynomial(A,d+1)
Q=Modulo_Polynomial([-C[i-1] if i else 1 for i in range(d+1)], d+1)
P=A*Q; P[d]=0
return Polynominal_Coefficient(P,Q,N)
#==================================================
def solve():
global Mod,Calc
Mod=998244353; Calc=Calculator()
N=int(input())
P=[-1]+list(map(int,input().split()))
W=[-1]+list(map(int,input().split()))
W_sum=[0]*(N+1)
for i in range(1,N+1):
W_sum[P[i]]+=W[i]
B=[0]*(N+1); B[0]=1
dep=[0]*(N+1)
leaf=[True]*(N+1)
for i in range(1, N+1):
coef=W[i]*pow(W_sum[P[i]], Mod-2, Mod)
B[i]=B[P[i]]*coef%Mod
dep[i]=dep[P[i]]+1
leaf[P[i]]=False
D=max(dep)
C=[0]*(D+2)
for i in range(N+1):
if leaf[i]:
C[dep[i]]+=B[i]
for d in range(D+1,-1,-1):
if d:
C[d]-=C[d-1]
else:
C[d]+=1
X=[1 if t==0 else 0 for t in range(-(D+1),0+1)]
# 本計算
Q=int(input())
Ans=[0]*Q
for q in range(Q):
A,K=map(int,input().split())
Ans[q]=Nth_Term_of_Linearly_Recurrent_Sequence(X, C, K-dep[A], -(D+1))*B[A]%Mod
if A==0:
Ans[q]=(Ans[q]-1)%Mod
return Ans
#==================================================
import sys
input=sys.stdin.readline
write=sys.stdout.write
print(*solve(), sep="\n")
Kazun