Read problems statements in Mandarin Chinese  and Russian. 

This problem's statement is really a short one.

You are given an integer S. Consider an infinite sequence S, 2S, 3S, ... . Find the first number in this sequence that can be represented as Q^{3}, where Q is some positive integer number. As the sought number could be very large, please print modulo (10^{9} + 7).

The number S will be given to you as a product of N positive integer numbers A_{1}, A_{2}, ..., A_{N}, namely S = A_{1} * A_{2} * ... * A_{N}

------ Input ------ 

The first line of the input contains an integer T denoting the number of test cases. The description of T test cases follows.

The first line of each test case contains an integer N.

Then there is a line, containing N space separated integers, denoting the numbers A_{1}, A_{2}, ..., A_{N}.

------ Output ------ 

For each test case, output a single line containing the first term of the sequence which is the perfect cube, modulo 10^{9}+7.

------ Constraints ------ 

$1 ≤ T ≤ 10$
$(Subtask 1): N = 1, 1 ≤ S ≤ 10^{9} - 15 points.$
$(Subtask 2): N = 1, 1 ≤ S ≤ 10^{18} - 31 point.$
$(Subtask 3): 1 ≤ N ≤ 100, 1 ≤ A_{i} ≤ 10^{18} - 54 points.$

----- Sample Input 1 ------ 
2

2

2 2

2

2 3
----- Sample Output 1 ------ 
8

216
----- explanation 1 ------ 
Example case 1. First few numbers in the infinite sequence 4, 8, 12, 16, 20, , etc. In this sequence, 8 is the first number which is also a cube (as 23 = 8).
Example case 2. First few numbers in the infinite sequence 6, 12, 18, 24, , etc. In this sequence, 216 is the first number which is also a cube (as 63 = 216).

Difficulty: VERY_HARD
Tags: Number theory

Mathematics

Starter code:

Solutions:
import math as m9

def prime_sieve(n):
	if n % 2 == 1:
		n += 1
	sieve = [True] * int(n / 2)
	n2 = sqrt_int(n)
	for i in range(3, n2 + 1, 2):
		if sieve[int(i / 2)]:
			sieve[int(i * i / 2)::i] = [False] * int((n - i * i - 1) / (2 * i) + 1)
	return [2] + [2 * i + 1 for i in range(1, int(n / 2)) if sieve[i]]

def sqrt_int(n):
	n2 = int(n ** 0.5)
	if (n2 + 1) ** 2 == n:
		print('         sqrt_int: rounding issues for n = 8d, n2 = %9.4f' % (n, n2))
	return n2

def issq(n):
	if n < 0:
		print('issq-error - CANNOT calc the square root of n = %8.3f' % n)
	valid = n == verifytruncatednthroot(n, 2) ** 2
	return valid

def verifytruncatednthroot(n, root):
	n1 = n ** (1 / root)
	n2 = int(n1)
	check1 = (n2 + 1) ** root
	if check1 <= n:
		n2 += 1
	check2 = (n2 + 0) ** root
	if check2 > n:
		print('         verifytruncatednthroot: ROUNDING issue(2)', n, root, n1, n2, check2)
	return n2

def primefactorsuptop_typ2(n, primes):
	v = []
	k = n
	for q in primes:
		while k % q == 0:
			v.append(q)
			k //= q
	if k > 1:
		v.append(k)
	return v

def fcube(A, primes, modx=0):
	W = []
	R0 = []
	mp = max(primes)
	for x in A:
		if x == 1:
			x = x
		else:
			pf = primefactorsuptop_typ2(x, primes)
			k9 = len([q for q in pf if q > 100])
			if pf[-1] <= mp:
				W.extend(pf)
			else:
				W.extend(pf[:-1])
				R0.append(pf[-1])
	G0 = []
	if len(R0) == 1:
		G0 = list(R0)
	else:
		for (i, r1) in enumerate(R0):
			for (j, r2) in enumerate(R0):
				if i < j:
					g = m9.gcd(r1, r2)
					if g > 1:
						if issq(g):
							g2 = sqrt_int(g)
						else:
							g2 = g
						if not g2 in G0:
							G0.append(g2)
	loops = 0
	G00 = list(sorted(G0, reverse=True))
	if len(G0) <= 1:
		if len(G0) == 1:
			t1 = G0[0]
			if issq(t1):
				t2 = sqrt_int(t1)
				G1 = [t2]
			else:
				G1 = list(G0)
		else:
			G1 = list(G0)
	else:
		while True:
			loops += 1
			G1 = []
			for (i, r1) in enumerate(G00):
				if issq(r1):
					t2 = sqrt_int(t1)
					if not t2 in G1:
						G1.append(t2)
				for (j, r2) in enumerate(G00):
					if i < j:
						g = m9.gcd(r1, r2)
						T0 = [g, r1 // g, r2 // g]
						T1 = [q for q in T0 if q > 1]
						for t in T1:
							if issq(t):
								t2 = sqrt_int(t)
							else:
								t2 = t
							if not t2 in G1:
								G1.append(t2)
			if sorted(G00) == sorted(G1):
				break
			G00 = sorted(G1)
	G11 = list(G1)
	for (i, r1) in enumerate(G1):
		for (j, r2) in enumerate(G1):
			p = max(r1, r2)
			q = min(r1, r2)
			g = m9.gcd(p, q)
			if g > 1:
				f1 = p // g
				if f1 in G11 and p in G11:
					G11.remove(p)
	G = list(G11)
	R1 = []
	for r in R0:
		x = r
		for g in G:
			while x % g == 0:
				W.append(g)
				x //= g
		if x > 1:
			if issq(x):
				x2 = sqrt_int(x)
			else:
				x2 = x
			if not x2 in R1:
				R1.append(x2)
	W2 = []
	s = 1
	for x in sorted(list(set(W))):
		y = W.count(x)
		z = 3 * m9.ceil(y / 3)
		s *= x ** z
		if modx > 0:
			s %= modx
	for x in list(set(R1)):
		y = R1.count(x)
		z = 3 * m9.ceil(y / 3)
		s *= x ** z
		if modx > 0:
			s %= modx
	return s
m = 0
modfcube = 10 ** 9 + 7
m = int(input())
maxprimex = 20
maxprimex = 1000 * 10 ** 3
primes = prime_sieve(maxprimex)
for i in range(1, m + 1):
	n = int(input())
	s = input()
	u = s.split(' ')
	A = [int(x) for x in u]
	y = fcube(A, primes, modx=modfcube)
	print(y)


from random import randrange
from fractions import gcd
from collections import defaultdict, Counter
K = 30

def miller_rabin(n):
	if n == 2:
		return True
	if n == 1 or n % 2 == 0:
		return False
	d = n - 1
	while d % 2 == 0:
		d //= 2
	for i in range(K):
		x = randrange(1, n)
		t = d
		y = pow(x, t, n)
		while t != n - 1 and y != 1 and (y != n - 1):
			y = y * y % n
			t *= 2
		if y != n - 1 and t % 2 == 0:
			return False
	return True

def f(x, m):
	return (x * x + 11) % m

def pollard_rho(n):
	if n % 2 == 0:
		return 2
	x = 2
	y = 2
	d = 1
	while d == 1:
		x = f(x, n)
		y = f(f(y, n), n)
		d = gcd(abs(x - y), n)
	return d

def factor(n):
	if n == 1:
		return defaultdict(int)
	temp = [n]
	factors = defaultdict(int)
	while temp:
		x = temp.pop()
		if x <= M:
			while x > 1:
				factors[ds[x]] += 1
				x //= ds[x]
			continue
		if miller_rabin(x):
			factors[x] += 1
			continue
		y = pollard_rho(x)
		if x == y:
			temp.append(x)
		else:
			temp.append(y)
			temp.append(x // y)
	return factors
MOD = 1000000007
M = 1000000
ds = [1] * (M + 1)
for i in range(2, M + 1):
	if ds[i] == 1:
		for j in range(i, M + 1, i):
			ds[j] = i
for _ in range(int(input())):
	n = int(input())
	ct = defaultdict(int)
	a = [int(i) for i in input().split()]
	r = 1
	for v in a:
		for (k, v) in factor(v).items():
			ct[k] += v
			r *= pow(k, v, MOD)
			r %= MOD
	for (k, v) in ct.items():
		r *= pow(k, (3 - v % 3) % 3, MOD)
	print(r % MOD)