from random import shuffle,randrange rd=randrange class primes(): def __init__(self, n): self.prime_num = n self.min_prime = [-1] * (self.prime_num + 1) # 2以上の自然数に対して最小の素因数を表す self.min_prime[0] = 0 self.min_prime[1] = 1 i = 2 self.prime = [] self.memo_prifac = {} self.memodiv=dict() while i <= self.prime_num: if self.min_prime[i] == -1: self.min_prime[i] = i self.prime.append(i) for j in self.prime: if i * j > self.prime_num or j > self.min_prime[i]: break self.min_prime[j * i] = j i += 1 def prifac(self, n): # 素因数分解した結果を返す if n in self.memo_prifac: return self.memo_prifac[n] res = {} x = n while x > 1: p = self.min_prime[x] if p in res: res[p] += 1 else: res[p] = 1 x //= p # self.memo_prifac[n] = res #場合によってはこの行を消すと高速化 return res def divisors(self, n): # 約数列挙 メモした方がいいかも if n in self.memodiv:return self.memodiv[n] if n== 1: return [1] prf = self.prifac(n) keys = [key for key in prf] def divsearch(i): if i == len(keys) - 1: return [keys[i] ** j for j in range(prf[keys[i]] + 1)] else: res = [] subres = divsearch(i + 1) p = keys[i] for j in range(prf[p] + 1): for node in subres: res.append(node * p ** j) return res self.memodiv[n]=divsearch(0) return self.memodiv[n] pri=primes(10**6+100) u=[0]*(10**6+100) u[1]=1 for x in range(2,10**6+10): u[x]=1 d=pri.prifac(x) for p in d: if d[p]>=2:u[x]=0 u[x]*=-1 from math import gcd mod=998244353 def naive(n,p): ans=0 for bit in range(1,2**n): x=[] for i in range(n): if (bit>>i)&1:x.append(p[i]) k=len(x) flag=1 for i in range(k-1): if gcd(x[i],x[i+1])==1:flag=0 ans+=flag return ans%mod def sol1(n,P): p=[0]+P[:] ma=max(p) dp=[0]*(ma+1) for i in range(1,n+1): res=0 for j in range(1,ma+1): if gcd(p[i],j)==1:res+=dp[j] dp[p[i]]+=sum(dp)+1-res dp[p[i]]%=mod return sum(dp)%mod def sol2(n,P): p=[0]+P[:] for i in range(1,n+1): res=1 for q in pri.prifac(p[i]):res*=q p[i]=res ma=max(p) dp=[0]*(ma+1) g=[0]*(ma+1) sumdp=0 for i in range(1,n+1): f1=0 div=[] for m in pri.divisors(p[i]):div.append(m) for m in div: f1+=u[m]*g[m] delta=sumdp+1-f1 dp[p[i]]+=delta dp[p[i]]%=mod sumdp+=delta sumdp%=mod for m in div: g[m]+=delta g[m]%=mod return sum(dp)%mod n=int(input()) a=list(map(int,input().split())) print(sol2(n,a)) cnt=0 while 0: cnt+=1 print(cnt) n=randrange(1,19) p=[rd(1,100) for i in range(n)] ansn=naive(n,p) ans1=sol1(n,p) ans2=sol2(n,p) if ans1!=ansn: print(n) print(*p) exit()