Given two integers, $m$ and ${k}}$, Alice loves to calculate their power sum using the following formula:
$\textit{PowerSum}(m,k)\equiv\sum\limits_{i=1}^m i^k$

Bob has a set, ${s}$, of $n$ distinct pairwise coprime integers. Bob hates multiples of these integers, so he subtracts ${i}^{k}$ from Alice's power sum for each ${i} \in [1,m]$ whenever there exists at least one $j\in[1,n]$ such that $i\:{mod}\:s_j\equiv0$.

Alice and Bob are now confused about the final value of the power sum and decide to turn to Eve for help. Can you write a program that helps Eve solve this problem? Given ${q}$ queries consisting of $n$, $m$, and ${k}}$, print the value of the power sum modulo $10^9+7$ on a new line for each query.

Input Format

The first line contains an integer, ${q}$, denoting the number of queries. The $2\cdot q$ lines describe each query over two lines:

The first line contains three space-separated integers denoting the respective values of $n$ (the number of integers in Bob's set), ${k}$ (the exponent variable in the power sum formula), and $m$ (the upper range bound in the power sum formula). 
The second line contains $n$ distinct space-separated integers describing the respective elements in set ${s}$. 

Constraints

$1\leq q\leq2$  
$1\leq n\leq50$  
$0\leq k\leq10$  
$1\leq m\leq10^{12}$  
$1\leq s_j\leq10^{12}$     
$s_i\neq s_j,\text{where}\ i\neq j$     
$gcd(s_i,s_j)\equiv1,\text{where}\ i\neq j$  

Output Format

For each query, print the resulting value of the power sum after Bob's subtraction, modulo $10^9+7$.

Sample Input
2
2 1 10
2 3
3 2 18
4 13 9

Sample Output
13
1055

Explanation

We perform the following $q=2$ queries:

Alice first calculates the sum $1^1+2^1+...+10^1=55$. Bob's set contains $2$ and $3$ only, so he subtracts the power of all numbers that are multiples of $2$ and/or $3$ from Alice's sum to get: $55-2^1-3^1-4^1-6^1-8^1-9^1-10^1=13$. We then print the result of $13\:\text{mod}\: (10^9+7)=13$ on a new line.
Alice first calculates the sum $1^2+2^2+...+18^2=2109$. Bob then subtracts multiples of ${4}$, ${9}$, and $13$ from Alice's sum to get: $2109-4^2-8^2-9^2-12^2-13^2-16^2-18^2=1055$. We then print the result of $1055\: \text{mod}\:(10^9+7)=1055$ on a new line.

Source: hackerrank | Difficulty: HARD
Tags: Number theory


Solutions:
from collections import defaultdict
p = 1000000007

def inv(i):
	return pow(i, p - 2, p)
mat = []
matinv = []
numpows = 11
for i in range(numpows):
	row = [0] * numpows
	mat.append(row)
	matinv.append(row.copy())
	matinv[i][i] = 1
for i in range(numpows):
	mat[i][0] = 1
	mat[i][i] = i + 1
	for j in range(1, i):
		mat[i][j] = mat[i - 1][j] + mat[i - 1][j - 1]
for i in range(numpows):
	for j in range(i + 1):
		if (i + j) % 2 == 1:
			mat[i][j] = p - mat[i][j]
for i in range(numpows):
	pivnum = mat[i][i]
	for j in range(numpows):
		mat[i][j] = mat[i][j] * inv(pivnum) % p
		matinv[i][j] = matinv[i][j] * inv(pivnum) % p
	for j in range(i + 1, numpows):
		subnum = p - mat[j][i]
		for k in range(numpows):
			mat[j][k] = (mat[j][k] + subnum * mat[i][k]) % p
			matinv[j][k] = (matinv[j][k] + subnum * matinv[i][k]) % p
powsums = dict()

def sumpow(i, k):
	try:
		return powsums[i, k]
	except:
		try:
			powsums[i, k] = (powsums[i - 1, k] + pow(i, k, p)) % p
			return powsums[i, k]
		except:
			toret = 0
			for a in range(k + 1):
				toret = (toret + matinv[k][a] * pow(i, a + 1, p)) % p
			powsums[i, k] = toret
			return toret

def noreps(a):
	toret = []
	currval = -1
	for x in a:
		if x != currval:
			toret.append(x)
		currval = x
	return toret

def merge(a, b):
	bind = 0
	if len(b) == 0:
		return a
	left = -1
	right = len(a)
	while right - left > 1:
		mid = (left + right) // 2
		if a[mid] > b[0]:
			left = mid
		else:
			right = mid
	toret = a[:left]
	aind = left
	while aind < len(a) and bind < len(b):
		if a[aind] == b[bind]:
			toret.append(a[aind])
			aind += 1
			bind += 1
		elif a[aind] > b[bind]:
			toret.append(a[aind])
			aind += 1
		elif b[bind] > a[aind]:
			toret.append(b[bind])
			bind += 1
	while aind < len(a):
		toret.append(a[aind])
		aind += 1
	while bind < len(b):
		toret.append(b[bind])
		bind += 1
	return toret
t = int(input())
for _ in range(t):
	(n, k, m) = list(map(int, input().split(' ')))
	s = list(map(int, input().split(' ')))
	s.sort()
	s.reverse()
	s = s[:n]
	tosum = defaultdict(lambda : 0)
	tosum[m] = 1
	tosumkeys = [m]
	for x in s:
		adddict = defaultdict(lambda : 0)
		addkeys = []
		for y in tosumkeys:
			if y < x:
				break
			adddict[y // x] += tosum[y]
			addkeys.append(y // x)
		addkeys = noreps(addkeys)
		mult = pow(x, k, p)
		for y in addkeys:
			tosum[y] = (tosum[y] - adddict[y] * mult) % p
		tosumkeys = merge(tosumkeys, addkeys)
	toprint = 0
	for x in reversed(tosumkeys):
		toprint += sumpow(x, k) * tosum[x]
	print(toprint % p)