You are given a sequence whose $n^{\text{th}}$ term is 

$T_n=n^K\times R^n$
You have to evaluate the series
$S_n=T_1+T_2+T_3+\cdots+T_n$
Find $S_n\:\text{mod}\:(10^9+7)$.

Input Format  

The first line of input contains $\mathbf{T}$, the number of test cases. 

Each test case consists of three lines, each containing $\mbox{K}$, $n$ and $\mbox{R}$ respectively.  

Output Format 

For each test case, print the required answer in a line.  

Constraints  

$1\leq T\leq10$ 

$1\leq K\leq10^3$ 

$1\leq n\leq10^{16}$ 

$2\leq R\leq10^{16}$ 

$R\ \text{mod}\ (10^9+7)\neq1$  

Sample Input  

2
2
5
2
3
4
3

Sample Output  

1146
5988

Explanation  

Case 1: $1146=1^2\times2^1+2^2\times2^2+3^2\times2^3+4^2\times2^4+5^2\times2^5$ 

Case 2: $5988=1^3\times3^1+2^3\times3^2+3^3\times3^3+4^3\times3^4$

Difficulty: HARD
Tags: Matrices

Starter code:

Solutions:
mod = 1000000007
factlist = [1]
for i in range(2005):
	factlist.append(factlist[-1] * (i + 1) % mod)
factinvlist = [0] * len(factlist)
factinvlist[-1] = pow(factlist[-1], mod - 2, mod)
for i in range(len(factlist) - 1, 0, -1):
	factinvlist[i - 1] = factinvlist[i] * i % mod

def choose(n, r):
	if r < 0 or n < r:
		return 0
	else:
		return factlist[n] * factinvlist[r] * factinvlist[n - r] % mod
t = int(input())
for _ in range(t):
	k = int(input())
	n = int(input())
	r = int(input())
	rmininv = pow(r - 1, mod - 2, mod)
	rinv = pow(r, mod - 2, mod)
	oneminrinv = pow(1 - rinv, mod - 2, mod)
	arr = []
	for i in range(k + 1):
		arr.append([0] * (k + 1))
		arr[i][i] = 1 - rinv
		for j in range(i - 1, -1, -1):
			arr[j][i] = -pow(-1, i - j) * choose(i, j) * rinv
	outputvec = [0] * k + [1]
	inputvec = [0] * (k + 1)
	sumzero = 0
	toprint = 0
	for i in range(k, -1, -1):
		inputvec[i] = oneminrinv * outputvec[i] % mod
		for j in range(i):
			outputvec[j] -= arr[j][i] * inputvec[i] % mod
		toprint += inputvec[i] * pow(n, i, mod) * pow(r, n, mod)
	sumzero = inputvec[0]
	print((toprint - sumzero) % mod)


from itertools import chain

def formula(base, exp, n, mod=10 ** 9 + 7):
	b_inv = pow(base - 1, -1, mod)
	val_by_exp = [base * (pow(base, n, mod) - 1) * b_inv % mod]
	over = [1]
	for e in range(1, exp + 1):
		over = [(a + b) % mod for (a, b) in zip(chain([0], over), chain(over, [0]))]
		_ = pow(1 + n, e, mod) * pow(base, n, mod) - 1
		_ -= sum((over[i] * val_by_exp[i] for i in range(e)))
		val_by_exp.append(base * _ * b_inv % mod)
	return val_by_exp

def solve(k, n, r):
	return formula(r, k, n)[-1]
for _ in range(int(input())):
	k = int(input())
	n = int(input())
	r = int(input())
	print(solve(k, n, r))


mod = 1000000007
factlist = [1]
for i in range(2005):
	factlist.append(factlist[-1] * (i + 1) % mod)
factinvlist = [0] * len(factlist)
factinvlist[-1] = pow(factlist[-1], mod - 2, mod)
for i in range(len(factlist) - 1, 0, -1):
	factinvlist[i - 1] = factinvlist[i] * i % mod

def choose(n, r):
	if r < 0 or n < r:
		return 0
	else:
		return factlist[n] * factinvlist[r] * factinvlist[n - r] % mod
t = int(input())
for _ in range(t):
	k = int(input())
	n = int(input())
	r = int(input())
	rmininv = pow(r - 1, mod - 2, mod)
	rinv = pow(r, mod - 2, mod)
	oneminrinv = pow(1 - rinv, mod - 2, mod)
	arr = []
	for i in range(k + 1):
		arr.append([0] * (k + 1))
		arr[i][i] = 1 - rinv
		for j in range(i - 1, -1, -1):
			arr[j][i] = -pow(-1, i - j) * choose(i, j) * rinv
	outputvec = [0] * k + [1]
	inputvec = [0] * (k + 1)
	sumzero = 0
	toprint = 0
	for i in range(k, -1, -1):
		inputvec[i] = oneminrinv * outputvec[i] % mod
		for j in range(i):
			outputvec[j] -= arr[j][i] * inputvec[i] % mod
		toprint += inputvec[i] * pow(n, i, mod) * pow(r, n, mod)
	sumzero = inputvec[0]
	print((toprint - sumzero) % mod)


mod = 1000000007
factlist = [1]
for i in range(2005):
	factlist.append(factlist[-1] * (i + 1) % mod)
factinvlist = [0] * len(factlist)
factinvlist[-1] = pow(factlist[-1], mod - 2, mod)
for i in range(len(factlist) - 1, 0, -1):
	factinvlist[i - 1] = factinvlist[i] * i % mod

def choose(n, r):
	if r < 0 or n < r:
		return 0
	else:
		return factlist[n] * factinvlist[r] * factinvlist[n - r] % mod
t = int(input())
for _ in range(t):
	k = int(input())
	n = int(input())
	r = int(input())
	rmininv = pow(r - 1, mod - 2, mod)
	rinv = pow(r, mod - 2, mod)
	oneminrinv = pow(1 - rinv, mod - 2, mod)
	arr = []
	for i in range(k + 1):
		arr.append([0] * (k + 1))
		arr[i][i] = 1 - rinv
		for j in range(i - 1, -1, -1):
			arr[j][i] = -pow(-1, i - j) * choose(i, j) * rinv
	outputvec = [0] * k + [1]
	inputvec = [0] * (k + 1)
	sumzero = 0
	toprint = 0
	for i in range(k, -1, -1):
		inputvec[i] = oneminrinv * outputvec[i] % mod
		for j in range(i):
			outputvec[j] -= arr[j][i] * inputvec[i] % mod
		toprint += inputvec[i] * pow(n, i, mod) * pow(r, n, mod)
	sumzero = inputvec[0]
	print((toprint - sumzero) % mod)


mod = 1000000007
factlist = [1]
for i in range(2005):
	factlist.append(factlist[-1] * (i + 1) % mod)
factinvlist = [0] * len(factlist)
factinvlist[-1] = pow(factlist[-1], mod - 2, mod)
for i in range(len(factlist) - 1, 0, -1):
	factinvlist[i - 1] = factinvlist[i] * i % mod

def choose(n, r):
	if r < 0 or n < r:
		return 0
	else:
		return factlist[n] * factinvlist[r] * factinvlist[n - r] % mod
t = int(input())
for _ in range(t):
	k = int(input())
	n = int(input())
	r = int(input())
	rmininv = pow(r - 1, mod - 2, mod)
	rinv = pow(r, mod - 2, mod)
	oneminrinv = pow(1 - rinv, mod - 2, mod)
	arr = []
	for i in range(k + 1):
		arr.append([0] * (k + 1))
		arr[i][i] = 1 - rinv
		for j in range(i - 1, -1, -1):
			arr[j][i] = -pow(-1, i - j) * choose(i, j) * rinv
	outputvec = [0] * k + [1]
	inputvec = [0] * (k + 1)
	sumzero = 0
	toprint = 0
	for i in range(k, -1, -1):
		inputvec[i] = oneminrinv * outputvec[i] % mod
		for j in range(i):
			outputvec[j] -= arr[j][i] * inputvec[i] % mod
		toprint += inputvec[i] * pow(n, i, mod) * pow(r, n, mod)
	sumzero = inputvec[0]
	print((toprint - sumzero) % mod)


mod = 1000000007
factlist = [1]
for i in range(2005):
	factlist.append(factlist[-1] * (i + 1) % mod)
factinvlist = [0] * len(factlist)
factinvlist[-1] = pow(factlist[-1], mod - 2, mod)
for i in range(len(factlist) - 1, 0, -1):
	factinvlist[i - 1] = factinvlist[i] * i % mod

def choose(n, r):
	if r < 0 or n < r:
		return 0
	else:
		return factlist[n] * factinvlist[r] * factinvlist[n - r] % mod
t = int(input())
for _ in range(t):
	k = int(input())
	n = int(input())
	r = int(input())
	rmininv = pow(r - 1, mod - 2, mod)
	rinv = pow(r, mod - 2, mod)
	oneminrinv = pow(1 - rinv, mod - 2, mod)
	arr = []
	for i in range(k + 1):
		arr.append([0] * (k + 1))
		arr[i][i] = 1 - rinv
		for j in range(i - 1, -1, -1):
			arr[j][i] = -pow(-1, i - j) * choose(i, j) * rinv
	outputvec = [0] * k + [1]
	inputvec = [0] * (k + 1)
	sumzero = 0
	toprint = 0
	for i in range(k, -1, -1):
		inputvec[i] = oneminrinv * outputvec[i] % mod
		for j in range(i):
			outputvec[j] -= arr[j][i] * inputvec[i] % mod
		toprint += inputvec[i] * pow(n, i, mod) * pow(r, n, mod)
	sumzero = inputvec[0]
	print((toprint - sumzero) % mod)