ZS the Coder has recently found an interesting concept called the Birthday Paradox. It states that given a random set of 23 people, there is around 50% chance that some two of them share the same birthday. ZS the Coder finds this very interesting, and decides to test this with the inhabitants of Udayland.

In Udayland, there are 2^{n} days in a year. ZS the Coder wants to interview k people from Udayland, each of them has birthday in one of 2^{n} days (each day with equal probability). He is interested in the probability of at least two of them have the birthday at the same day. 

ZS the Coder knows that the answer can be written as an irreducible fraction $\frac{A}{B}$. He wants to find the values of A and B (he does not like to deal with floating point numbers). Can you help him?


-----Input-----

The first and only line of the input contains two integers n and k (1 ≤ n ≤ 10^18, 2 ≤ k ≤ 10^18), meaning that there are 2^{n} days in a year and that ZS the Coder wants to interview exactly k people.


-----Output-----

If the probability of at least two k people having the same birthday in 2^{n} days long year equals $\frac{A}{B}$ (A ≥ 0, B ≥ 1, $\operatorname{gcd}(A, B) = 1$), print the A and B in a single line.

Since these numbers may be too large, print them modulo 10^6 + 3. Note that A and B must be coprime before their remainders modulo 10^6 + 3 are taken.


-----Examples-----
Input
3 2

Output
1 8
Input
1 3

Output
1 1
Input
4 3

Output
23 128


-----Note-----

In the first sample case, there are 2^3 = 8 days in Udayland. The probability that 2 people have the same birthday among 2 people is clearly $\frac{1}{8}$, so A = 1, B = 8.

In the second sample case, there are only 2^1 = 2 days in Udayland, but there are 3 people, so it is guaranteed that two of them have the same birthday. Thus, the probability is 1 and A = B = 1.

Difficulty: HARD
Tags: Number theory

Mathematics

Probability

Starter code:

Solutions:
import os
MOD = 1000003
inv2 = pow(2, MOD - 2, MOD)

def logm(n, m):
	ans = 0
	whole = True
	while n >= m:
		whole = whole and n % m == 0
		ans += 1
		n //= m
	if n == 1:
		return (ans, whole)
	return (ans, False)

def fact_exp(n, k):
	ans = 0
	while n != 0:
		n //= k
		ans += n
	return ans

def main():
	(n, k) = map(int, input().split())
	e2 = n + fact_exp(k - 1, 2)
	div = pow(2, n * k - e2, MOD)
	(e, w) = logm(k, 2)
	if e > n or (e == n and (not w)):
		print(1, 1)
		exit()
	num = 1
	Nr = pow(2, n, MOD)
	for t in range(1, k):
		i = (Nr - t) % MOD
		if i == 0:
			num = 0
			break
		p = 0
		while t % 2 == 0:
			p += 1
			t //= 2
		num = num * i * pow(inv2, p, MOD) % MOD
	print((div - num) % MOD, div)
main()


m = 10 ** 6 + 3
(n, k) = map(int, input().split())
p = 1
for i in range(n):
	p *= 2
	if p > k:
		break
if p < k:
	print('1 1')
	exit()
gcd = tmp = k - 1
while tmp:
	gcd -= tmp % 2
	tmp //= 2
b = pow(2, (k - 1) * n - gcd, m)
a = 1
mem = [-1] * 100
for i in range(1, k):
	cnt = 0
	while i % 2 == 0:
		i //= 2
		cnt += 1
	if mem[cnt] == -1:
		mem[cnt] = pow(2, n - cnt, m)
	a = a * (mem[cnt] - i + m) % m
	if a == 0:
		break
print((b - a + m) % m, b)


(n, k) = map(int, input().split())
MOD = 1000003
K = k - 1
max_deg = 0
while K > 0:
	max_deg += K // 2
	K //= 2
den_deg = n * (k - 1) - max_deg
kk = 1
for i in range(n):
	kk *= 2
	if kk >= k:
		break
else:
	print(1, 1)
	return
numerator = 1
two_p_n = pow(2, n, MOD)
for i in range(1, min(k, MOD + 1)):
	numerator *= (two_p_n - i + MOD) % MOD
	if numerator == 0:
		break
	numerator %= MOD
rev = (MOD + 1) // 2
numerator *= pow(rev, max_deg, MOD)
numerator %= MOD
denumerator = pow(2, den_deg, MOD)
numerator = (denumerator + MOD - numerator) % MOD
print(numerator, denumerator)


(n, k) = [int(x) for x in input().split()]
if 2 ** min(n, 62) < k:
	print('1 1')
	exit(0)
mod = int(1000003)
x = n
mult = 2
while mult <= k - 1:
	x += (k - 1) // mult
	mult *= 2

def modpow(x, n, m):
	mult = x % m
	a = 1
	ans = 1
	while a <= n:
		if a & n:
			ans *= mult
			ans %= m
		mult *= mult
		mult %= m
		a *= 2
	return ans
denom = modpow(2, n * k - x, mod)
num = modpow(2, n * k - x, mod)
if k < mod:
	s = 1
	for i in range(1, k):
		val = i
		y = 0
		while not val % 2:
			val //= 2
			y += 1
		s *= modpow(2, (n - y) % (mod - 1), mod) - val
		s %= mod
	num -= s
	num %= mod
	if num < 0:
		s += mod
print(num, denom)


(n, k) = map(int, input().split())
mod = 1000003
if n < 70 and 2 ** n < k:
	print('1 1\n')
	exit(0)

def modpow(a, e):
	ret = 1
	while e > 0:
		if e % 2 == 1:
			ret = ret * a % mod
		a = a * a % mod
		e = e // 2
	return ret

def pw(a, e):
	ret = 1
	while e > 0:
		if e % 2 == 1:
			ret *= a
		a *= a
		e = e // 2
	return ret
par = n
for i in range(1, 100):
	par += (k - 1) // pw(2, i)
mul = 1
cur = modpow(2, n)
for i in range(k):
	mul = cur * mul % mod
	cur -= 1
	if mul == 0:
		break
if mul != 0:
	mul = mul * modpow(modpow(2, par), mod - 2) % mod
up = (modpow(2, n * k - par) - mul) % mod
if up < 0:
	up += mod
print(up, end=' ')
print(modpow(2, n * k - par))


import math
(n, k) = [int(x) for x in input().split()]
if n < 70 and k > 2 ** n:
	print(1, 1)
	exit(0)
mod = int(1000000.0) + 3

def fastpow(a, b):
	(t, ans) = (a, 1)
	while b:
		if b & 1:
			ans = ans * t % mod
		t = t * t % mod
		b >>= 1
	return ans
t = k - 1
cnt = 0
while t:
	cnt += t >> 1
	t >>= 1
x = 0
t = fastpow(2, n)
if k < mod:
	x = 1
	for i in range(1, k):
		x = x * (t - i) % mod
y = fastpow(2, n * (k - 1))
inv = fastpow(2, mod - 2)
inv = fastpow(inv, cnt)
x = (x * inv % mod + mod) % mod
y = (y * inv % mod + mod) % mod
x = (y - x + mod) % mod
print(x, y)


import sys
mod = 10 ** 6 + 3
(n, k) = map(int, input().split())
if n < 100:
	if 2 ** n < k:
		print(1, 1)
		sys.exit()

def factor(n, p):
	if n < p:
		return 0
	return n // p + factor(n // p, p)

def inv(n):
	return pow(n, mod - 2, mod)
two = inv(pow(2, n + factor(k - 1, 2), mod))
v = 1
if k >= mod:
	v = 0
else:
	N = pow(2, n, mod)
	for i in range(k):
		v = v * (N - i) % mod
A = (pow(2, n * k, mod) - v) * two % mod
B = pow(2, n * k, mod) * two % mod
print(A, B)


import math

def euclid_algorithm(a, b):
	(t1, t2) = (abs(a), abs(b))
	(x1, y1, x2, y2) = (int(math.copysign(1, a)), 0, 0, int(math.copysign(1, b)))
	if t1 < t2:
		(t1, t2) = (t2, t1)
		(x1, y1, x2, y2) = (x2, y2, x1, y1)
	while t2 > 0:
		k = int(t1 // t2)
		(t1, t2) = (t2, t1 % t2)
		(x1, y1, x2, y2) = (x2, y2, x1 - k * x2, y1 - k * y2)
	return (t1, x1, y1)

def opposite_element(x, p):
	(gcd, k, l) = euclid_algorithm(x, p)
	if gcd != 1:
		return -1
	return k % p

def solve(n, k):
	if n < 70 and k > 1 << n:
		return (1, 1)
	(s, l) = (0, k - 1)
	while l > 0:
		l >>= 1
		s += l
	p = 10 ** 6 + 3
	x = pow(2, n, p)
	t = pow(opposite_element(2, p), s, p)
	q = pow(2, n * (k - 1), p) * t % p
	r = 1
	if k > p:
		r = 0
	else:
		for i in range(1, k):
			r *= x - i
			r %= p
	return ((q - r * t) % p, q)
(n, k) = map(int, input().split())
(x, y) = solve(n, k)
print(x, y)


(n, k) = map(int, input().split())
if n <= 100 and k > 2 ** n:
	print(1, 1)
	exit(0)
MOD = 1000 * 1000 + 3
INV2 = (MOD + 1) // 2

def add(x, y):
	return (x + y) % MOD

def sub(x, y):
	rez = x - y
	rez %= MOD
	rez += MOD
	rez %= MOD
	return rez

def mult(x, y):
	return x * y % MOD

def binpow(x, y):
	if x == 1 or y == 0:
		return 1
	if x == 0:
		return 0
	rez = binpow(x, y // 2)
	rez = mult(rez, rez)
	if y % 2 == 1:
		rez = mult(rez, x)
	return rez
A = n * k
B = n
temp = k - 1
while temp >= 2:
	B += temp // 2
	temp //= 2
G = min(A, B)
m = binpow(2, n)
P = 1
for i in range(k):
	P = mult(P, sub(m, i))
	if P == 0:
		break
P = mult(P, binpow(INV2, G))
Q = binpow(m, k)
Q = mult(Q, binpow(INV2, G))
P = sub(Q, P)
print(P, Q)