Problem:
 Score : 400 points    Problem Statement You are given positive integers N and M .  How many sequences a of length N consisting of positive integers satisfy a_1 \times a_2 \times ... \times a_N = M ? Find the count modulo 10^9+7 .  Here, two sequences a' and a'' are considered different when there exists some i such that a_i' \neq a_i'' .      Constraints  All values in input are integers.  1 \leq N \leq 10^5  1 \leq M \leq 10^9         Input Input is given from Standard Input in the following format:  N  M       Output Print the number of the sequences consisting of positive integers that satisfy the condition, modulo 10^9 + 7 .        Sample Input 1 2 6      Sample Output 1 4  Four sequences satisfy the condition: \{a_1, a_2\} = \{1, 6\}, \{2, 3\}, \{3, 2\} and \{6, 1\} .       Sample Input 2 3 12      Sample Output 2 18       Sample Input 3 100000 1000000000      Sample Output 3 957870001

Solution (Python):
from collections import Counter


def main():
    """
    1 <= N <= 10^5
    1 <= M <= 10^9

    正整数 N, M
    a1 * a2 * ... * aN = M となる長さNの数列a
    何通りあるか
    数列A1-ai != 数列A2-ai
    """
    N, M = map(int, input().split())

    ans = editorial(N, M)
    print(ans)


def factorize(N):
    """素因数分解

    TODO: もっと速い方法
    素因数がある場合は範囲をsqrtで狭められるが、ない場合愚直に探してしまう
    N以下の素数を見つけて、それだけで対応すれば速いはず
    """
    import math
    factors = []
    while True:
        if N % 2 == 0:
            factors.append(2)
            N //= 2
        else:
            break
    m = int(math.sqrt(N))
    i = 3
    while i <= m:
        add = False
        while True:
            if N % i == 0:
                factors.append(i)
                N //= i
                add = True
            else:
                break
        if add:
            m = int(math.sqrt(N))
        i += 2
    if N > 1:
        factors.append(N)
    return factors


def f(N, M):
    """
    M を素因数分解すると約数とそれらの個数がわかる
    それらの組合せを制約を満たすようにする
    1 * M, M * 1 をわすれないこと:
        1 * 1 * 1 * M もOK

    3 12
        1,1,12: 3P3 / 2!
        2,2,3:  3P3 / 2!
        3,4,1:  3P3
        6,2,1:  3P3

    http://sucrose.hatenablog.com/entry/2014/10/10/235805
    """
    ans = 0
    m = 10 ** 9 + 7
    factors = factorize(M)
    c = Counter(factors)
    print(c)
    return ans


def editorial(N, M):
    i = 2
    rest_M = M
    ans = 1
    mod = 10 ** 9 + 7
    while i ** 2 <= rest_M:
        div_count = 0
        while rest_M % i == 0:
            div_count += 1
            rest_M //= i

        if div_count > 0:
            
            ans *= calc_comb(N + div_count - 1, div_count)
            ans %= mod

        i += 1

    if rest_M != 1:
        ans *= calc_comb(N + 1 - 1, N - 1)
        ans %= mod

    return ans


def calc_comb(n, r, mod=10 ** 9 + 7):
    if r > n - r:
        return calc_comb(n, n - r)

    ans_mul = 1
    ans_div = 1

    for i in range(r):
        
        ans_mul = ans_mul * (n - i) % mod
        
        ans_div = ans_div * (i + 1) % mod

    
    return ans_mul * mod_pow(ans_div, mod - 2) % mod


def mod_pow(a, p, mod=10 ** 9 + 7):
    if p == 0:
        return 1

    if p % 2 == 0:
        half_p = p // 2
        half = mod_pow(a, half_p)
        return half ** 2 % mod
    else:
        return a * mod_pow(a, p - 1) % mod


def mod_pow_bit(a, p, mod=10 ** 9 + 7):
    """
    http://keita-matsushita.hatenablog.com/entry/2016/12/05/184011
    """
    ans = 1
    while p > 0:
        if (p & 1) == 1:
            ans = ans * a % mod

        a = a ** 2 % mod
        p >>= 1

    return ans


def ref(N, M):
    """
    https://mathtrain.jp/tyohukuc
        n 種類のものから重複を許して r 個選ぶ場合の数
        n H r = (n+r−1) C r

    3 12 の例で示せているのに、重複組合せという言葉が頭に浮かばなかった

    解説動画より
    modの世界における逆元?
        重複組合せ
            フェルマーの小定理
                ref: http://keita-matsushita.hatenablog.com/entry/2016/12/05/184011
            階乗と階乗逆元
        拡張ユークリッド互除法
            今回は出番なし
    """
    def create_tables(N, MOD):
        """
        http://drken1215.hatenablog.com/entry/2018/06/08/210000
        """
        fact_mod = [0] * (N + 1)
        simod = [0] * (N + 1)
        fact_inv_mod = [0] * (N + 1)

        fact_mod[0] = fact_mod[1] = 1
        simod[0] = simod[1] = 1
        fact_inv_mod[0] = fact_inv_mod[1] = 1

        for i in range(2, N + 1):
            fact_mod[i] = (fact_mod[i - 1] * i) % MOD
            simod[i] = MOD - simod[MOD % i] * (MOD // i) % MOD
            fact_inv_mod[i] = (fact_inv_mod[i - 1] * simod[i]) % MOD
            
            
            

        return fact_mod, fact_inv_mod

    def create_fact_mod(N, MOD):
        fact_mod = [0] * (N + 1)
        fact_mod[0] = 1
        fact_mod[1] = 1

        
        for i in range(2, N + 1):
            fact_mod[i] = (fact_mod[i - 1] * i) % MOD

        return fact_mod

    def create_seq_inv_mod(N, MOD):
        simod = [0] * (N + 1)
        simod[0] = 1  
        simod[1] = 1

        
        for i in range(2, N + 1):
            simod[i] = (MOD - MOD // i) * simod[MOD % i] % MOD
            
            

        return simod

    def create_fact_inv_mod(N, MOD, simod):
        fact_inv_mod = [0] * (N + 1)
        fact_inv_mod[0] = 1
        fact_inv_mod[1] = 1

        
        for i in range(2, N + 1):
            fact_inv_mod[i] = (fact_inv_mod[i - 1] * simod[i]) % MOD

        return fact_inv_mod

    def comb_mod(N, R, fact_mod, fact_inv_mod, MOD):
        """
        組合せ nCr: n! / ((n-r)! * r!)
        上記のmodを逆元を使って表現

        f[n] * f[r] % MOD * f[n-r] % MOD 
        これのために事前にテーブルを作っていた
        """
        return fact_mod[N] * fact_inv_mod[N - R] * fact_inv_mod[R] % MOD

    ans = 1
    m = 10 ** 9 + 7
    factors = factorize(M)
    c = Counter(factors)

    
    
    
    fact_mod = create_fact_mod(N + 100, m)
    simod = create_seq_inv_mod(N + 100, m)
    fact_inv_mod = create_fact_inv_mod(N + 100, m, simod)

    
    
    for i in c.values():
        ans *= comb_mod(N - 1 + i, i, fact_mod, fact_inv_mod, m)
        ans %= m

    return ans


if __name__ == '__main__':
    main()

Solution (Python):
import math
import functools


def mmul(a, b, mod):
    return (a % mod) * (b % mod) % mod


def mfact(a, mod):
    fact = 1
    for m in range(a + 1):
        fact *= m % mod
        fact %= mod
    return fact


def mpow(a, b, mod):
    if b == 0:
        return 1
    if b == 1:
        return a % mod
    if b % 2 == 0:
        return mpow(a, b / 2, mod) ** 2 % mod
    return mpow(a, b // 2, mod) ** 2 * a % mod


def mdiv(a, b, mod):
    return mmul(a, mpow(b, mod - 2, mod), mod)



def mnpaca(n, a, mod):
    res = 1
    for num in range(n + 1, n + a + 1):
        res *= num % mod
        res %= mod
    return mdiv(res, math.factorial(a), mod)


mod = 10 ** 9 + 7

N, M = map(int, input().split())

prime_counts = []
i = 2
while i * i <= M:
    count = 0
    while M % i == 0:
        M /= i
        count += 1
    if count:
        prime_counts.append(count)
    i += 1
if M > 1:
    prime_counts.append(1)

res = functools.reduce(lambda x, y: x * y, [mnpaca(N - 1, c, mod) for c in prime_counts], 1)

print(res % mod)

Solution (Python):
N, M = [int(i) for i in input().split()]
mod = 10 ** 9 + 7
def func(M):
  res = []
  i = 2
  while i * i <= M:
    c=0
    while M % i == 0:
      M /= i
      c += 1
    if c > 0:
      res.append(c)
    i += 1
  if M > 1:
    res.append(1)
  return res

def conb(n, r):
  N, R = n, r
  for i in range(1, r):
    N *= n-i
    R *= r-i
  return N // R

res = func(M)
ans = 1
for i in res:
  ans = (ans * conb(i+N-1, i)) % mod
print(ans)

Solution (Python):

from math import factorial
from functools import reduce


def mmul(a, b, mod):
    return (a % mod) * (b % mod) % mod


def mfact(a, mod):
    fact = 1
    for m in range(2, a + 1):
        fact *= m % mod
        fact %= mod
    return fact


def mpow(a, b, mod):
    if b == 0:
        return 1
    if b == 1:
        return a % mod
    if b % 2 == 0:
        return mpow(a, b / 2, mod) ** 2 % mod
    return mpow(a, b // 2, mod) ** 2 * a % mod


def mdiv(a, b, mod):
    return mmul(a, mpow(b, mod - 2, mod), mod)



def mnpaca(n, a, mod):
    res = 1
    for num in range(n + 1, n + a + 1):
        res *= num % mod
        res %= mod
    return mdiv(res, factorial(a), mod)

mod = 1000000007
N, M = [int(s) for s in input().split()]

prime_counts = []
i = 2
while i * i <= M:
    count = 0
    while M % i == 0:
        M /= i
        count += 1
    if count:
        prime_counts.append(count)
    i += 1
if M > 1:
    prime_counts.append(1)

res = reduce(lambda x, y: x * y, [mnpaca(N - 1, c, mod) for c in prime_counts],
             1)

print(res % mod)

Solution (Python):
import sys
from string import ascii_lowercase

lim = 10**5 + 100
mod = 10**9 + 7


def main():
    n, m = map(int, input().split())

    f = [1]*(lim + 1)
    for i in range(2, lim + 1):
        f[i] = (i * f[i-1]) % mod
    rf = [1]*(lim + 1)
    rf[lim] = pow(f[lim], mod - 2, mod)
    for i in range(lim - 1, 1, -1):
        rf[i] = ((i+1) * rf[i+1]) % mod

    def comb(n, k):
        return f[n] * rf[k] * rf[n-k] % mod

    divs = factorize(m)

    ans = 1

    for p, c in divs:
        ans *= comb(n - 1 + c, n - 1)
        ans %= mod

    print(ans)


def factorize(n):
    res = []
    for i in range(2, n + 1):
        if i*i > n:
            break
        cnt = 0
        while n % i == 0:
            n //= i
            cnt += 1
        res.append((i, cnt))

    if n > 1:
        res.append((n, 1))

    return res


if __name__ == '__main__':
    main()

Solution (Python):
N,M=map(int,input().split())
MAX_N=2*10**5
P=10**9+7
def egcd(a, b):
    (x, lastx) = (0, 1)
    (y, lasty) = (1, 0)
    while b != 0:
        q = a // b
        (a, b) = (b, a % b)
        (x, lastx) = (lastx - q * x, x)
        (y, lasty) = (lasty - q * y, y)
    return (lastx, lasty, a)
def inv(x):
    return egcd(x,P)[0]
Fact=[0 for i in range(MAX_N+1)]
Finv=[0 for i in range(MAX_N+1)]
Fact[0]=1
Finv[0]=1
for i in range(MAX_N):
    Fact[i+1]=(Fact[i]*(i+1))%P
    Finv[i+1]=inv(Fact[i+1])%P
def C(n,k):
    return (Fact[n]*(Finv[k]*Finv[n-k])%P)%P
def H(n,r):
    return C(n+r-1,r-1)
ans=1
def factorize(n):
    fct = []  
    b, e = 2, 0  
    while b * b <= n:
        while n % b == 0:
            n = n // b
            e = e + 1
        if e > 0:
            fct.append((b, e))
        b, e = b + 1, 0
    if n > 1:
        fct.append((n, 1))
    return fct
D=factorize(M)
for seq in D:
    k=seq[1]
    ans=ans*C(k+N-1,N-1)
    ans=ans%P
print(ans)

Solution (Python):


mod = 10 ** 9 + 7


def mod_pow(a: int, p: int):
    '''a ** p % mod.'''

    if p == 0:
        return 1

    if p % 2 == 0:
        half_p = p // 2
        half = mod_pow(a, half_p)
        return (half ** 2) % mod
    else:
        return (a * mod_pow(a, p - 1)) % mod


def calc_combination(a: int, b: int):
    '''aCb % mod.'''

    if b > (a - b):
        return calc_combination(a, a - b)

    multiple = 1
    division = 1

    for i in range(b):
        multiple *= a - i
        division *= i + 1
        multiple %= mod
        division %= mod

    return multiple * mod_pow(division, mod - 2) % mod


def solve(n: int, m: int):
    from math import sqrt

    ans = 1
    remain = m

    for j in range(2, int(sqrt(m)) + 1):
        if remain % j == 0:
            count = 0

            while remain % j == 0:
                count += 1
                remain //= j

            ans *= calc_combination(n + count - 1, n - 1)
            ans %= mod

    if remain != 1:
        ans *= calc_combination(n, 1)
        ans %= mod

    return ans


def main():
    n, m = map(int, input().split())

    
    
    print(solve(n, m))


if __name__ == '__main__':
    main()

Solution (Python):
from collections import Counter
N, M = map(int, input().split())
mod = 10 ** 9 + 7



def prime_factorization(n):
    table = []
    for x in range(2, int(n**0.5)+2):
        while n % x == 0:
            table.append(x)
            n //= x
    if n > 1:
        table.append(n)
    return table



def nHk(n, k):
    value = 1
    nk = n + k - 1
    k = min(nk-k, k)
    for i in range(k):
        value = (value * nk) % mod
        nk -= 1
    for x in range(1, k + 1):
        value = (value * pow(x, mod-2, mod)) % mod
    return value


X = Counter(prime_factorization(M)).values()
ans = 1
for x in X:
    ans *= nHk(N, x) % mod

print(ans % mod)

Solution (Python):
def prime_decomposition(n):
  i = 2
  d = {}
  while i * i <= n:
    while n % i == 0:
      n /= i
      if i not in d:
          d[i] = 0
      d[i] += 1
    i += 1
  if n > 1:
    if n not in d:
        d[n] = 1
  return d

def egcd(a, b):
    (x, lastx) = (0, 1)
    (y, lasty) = (1, 0)
    while b != 0:
        q = a // b
        (a, b) = (b, a % b)
        (x, lastx) = (lastx - q * x, x)
        (y, lasty) = (lasty - q * y, y)
    return (lastx, lasty, a)


def modinv(a, m):
    (inv, q, gcd_val) = egcd(a, m)
    return inv % m

def modfact(n,N):
    r = 1
    for i in range(1,n+1):
        r = r % N * i % N
    return r


def modconb(n,r,N):
    p = modfact(n,N)
    q1 = modinv(modfact(r,N),N)
    q2 = modinv(modfact(n-r,N),N)
    return p*q1%N*q2%N


N,M = map(int,input().split())
R = 10**9+7

d = prime_decomposition(M)
ans = 1
for p in d:
    ans *= modconb(d[p]+N-1,N-1,R)
    ans %= R
print(ans)

Solution (Python):


import math, collections
N, M = map(int,input().split())
MOD = 7+10**9

if M == 1:
    print(1)
    exit(0)


Q = []
k = 2
while k < M:
    if M%k == 0:
        Q.append(k)
        M //= k
    else:
        k += 1
Q.append(k)

R = collections.Counter(Q)

def nCr(n,r):
    return (math.factorial(n)//(math.factorial(r)*math.factorial(n-r)))%MOD

ans = 1
for f in R:
    ans *= nCr(N+R[f]-1,R[f])
    ans %= MOD
print(ans)

Solution (Python):
N,M = map(int,raw_input().split())

def factorize(n):
	fct = []  
	b, e = 2, 0  
	while b * b <= n:
		while n % b == 0:
			n = n / b
			e = e + 1
		if e > 0:
			fct.append((b, e))
		b, e = b + 1, 0
	if n > 1:
		fct.append((n, 1))
	return fct

def comb(n,r):
	ans = 1
	for i in range(r):
		ans=ans*(n-i)
		ans=ans/(i+1)
	return ans%(10**9 + 7)

fct= factorize(M)

ans = 1

for f in fct:
	ans=(ans *comb(f[1]+N-1,f[1]))%(10**9+7)

print ans

Solution (Python):
from collections import Counter as c
from scipy.misc import comb
from functools import reduce

def f(n):
  d,i=c(),2
  while i*i<=n:
    while n%i==0:
      n//=i
      d[i]+=1
    i+=1
  if n>1:d[n]+=1
  return d.values()

n,m=map(int,input().split())
mod,v=10**9+7,1
for x in f(m):
  v*=comb(x+n-1,x,exact=True)%mod
print(v%mod)

Solution (Python):
N,M = map(int,input().split())

m = M
f = {}
while m%2==0:
    f[2] = f.get(2,0)+1
    m = m//2

l = int(m**0.5)+1
k = 3
while m != 1:
    if k>l:
        f[m] = f.get(m,0)+1
        break
    while m%k==0:
        f[k] = f.get(k,0)+1
        m = m//k
    k+=2
ans = 1
for k in f:
    for i in range(N-1+f[k],N-1,-1):
        ans = (ans * i) 
    for i in range(f[k]):
        ans = ans // (i+1)
    ans = ans % ((10**9)+7)
print(ans)

Solution (Python):
from collections import defaultdict
import math
N,M=map(int,input().split())
dd=defaultdict(int)
m=M
mod=10**9+7
for i in range(2,10**5):
    s=0
    while m%i==0 and m>1:
        m=m//i
        s+=1
    if s>0:
        dd[i]=s
if m>1:
    dd[10**9]=1
def f(n,k):
    a=1
    for i in range(1,k+1):
        a=a*(n-i+1)
        a%=mod
        a*=pow(i,mod-2,mod)
        a%=mod
    return a
ans=1
for l,t in dd.items():
    ans*=f(N-1+t,t)
    ans%=mod
print(ans)

Solution (Python):

import math
N, M = map(int, input().split())

ips = []
M_tmp = M
for p in range(2, int(math.sqrt(M)) + 1):
    if M_tmp % p > 0: continue
    for i_p in range(1, 33):
        if M_tmp % (p ** (i_p + 1)) > 0: break
    ips.append(i_p)
    M_tmp = M_tmp // (p ** i_p)
    if M_tmp == 1: break
else:
    ips.append(1)




if len(ips) > 0: ip_max = max(ips)
else: ip_max = 0
memo = [[None] * N for i in range(ip_max + 1)] 
for j in range(N): memo[0][j] = 1
for i in range(1, ip_max + 1): memo[i][N - 1] = 1

for i in range(1, ip_max + 1):
    sum_tmp = 0
    vals = list(memo[i - 1])[::-1]
    ptr_memo = memo[i]
    for j, val in enumerate(vals):
        sum_tmp += val
        ptr_memo[N - 1 - j] = sum_tmp

        

ans = 1
for i_p in ips:
    ans = (ans * memo[i_p][0]) % (10 ** 9 + 7)

if M == 1: print(1)
elif ip_max < 1: print(N)
else: print(ans)

Solution (Python):
def C(x, y):
    z = 1
    for i in range(x, x-y, -1):
        z *= i
    for i in range(y, 1, -1):
        z //= i
    return z
def H(x, y):
    return C(x+y-1, y) % 1000000007
n, m = list(map(int, input().split()))
p = 2
prime_factorization = []
while p * p <= m:
    count = 0
    while m % p == 0:
        m //= p
        count += 1
    if count > 0:
        prime_factorization.append(count)
    p += 1
if m > 1:
    prime_factorization.append(1)
ans = 1
for r in prime_factorization:
    ans *= H(n, r)
print(ans % 1000000007)

Solution (Python):
def C(x, y):
    z = 1
    for i in range(x, x-y, -1):
        z *= i
    for i in range(y, 1, -1):
        z //= i
    return z
def H(x, y):
    return C(x+y-1, y)
n, m = list(map(int, input().split()))
p = 2
prime_factorization = []
while p * p <= m:
    count = 0
    while m % p == 0:
        m //= p
        count += 1
    if count > 0:
        prime_factorization.append(count)
    p += 1
if m > 1:
    prime_factorization.append(1)
ans = 1
for r in prime_factorization:
    ans *= H(n, r)
print(ans % 1000000007)

Solution (Python):
n,m=map(int,input().split())
x=10**9+7
def factoring(m):
  i=2
  a=[]
  while i*i<=m:
    count=0
    while m%i==0:
      count+=1
      m=m//i
    if count>0:
      a.append(count)
    i+=1
  if m>1:
    a.append(1)
  return a
def cmb(n,r):
  if r>n-r:
    r=n-r
  if r==0:
    return 1
  if r==1:
    return n
  result=1
  denominator=[i for i in range(1,r+1)]
  numerator=[i for i in range(n-r+1,n+1)]
  for n in numerator:
    result*=n
  for d in denominator:
    result//=d
  return result
fac=factoring(m)
res=1
for b in fac:
  res*=cmb(n+b-1,b)
  res%=x
print(res)

Solution (Python):





"""
RC110 D
"""


n,m = map(int,input().split())

from math import sqrt
from collections import deque
def divisors(x):
    divs = deque()
    sx = int(sqrt(x))
    if x == sx**2:
        divs.append(sx)
        sx = sx-1
    for i in range(sx,0,-1):
        if x%i==0:
            divs.append(i)
            divs.appendleft(x//i)
    return(divs)

mdivs = list(divisors(m))
from collections import defaultdict
divsdict = defaultdict(list)
for d in mdivs:
    divsdict[d] = [dd for dd in mdivs if d%dd==0]

nbin = list(map(int,format(n,'b')))
nnbin = len(nbin)
dpdicts = [defaultdict(int) for i in range(nnbin)]
for d in mdivs:
    dpdicts[0][d]=1

cut = 10**9+7

for i in range(1,nnbin):
    for d in mdivs:
        dpdicts[i][d] = sum(map(lambda x:dpdicts[i-1][x]*dpdicts[i-1][d//x],divsdict[d]))%(10**9+7)

dpdictbin = defaultdict(int)
for d in mdivs:
    dpdictbin[d]=0
dpdictbin[1]=1

nbin.reverse()
for i in filter(lambda x:nbin[x]==1,range(nnbin)):
    for d in mdivs:
        dpdictbin[d] = sum(map(lambda x:dpdictbin[x]*dpdicts[i][d//x],divsdict[d]))%(10**9+7)



print(dpdictbin[m])

Solution (Python):
from collections import Counter


def main():
    """
    1 <= N <= 10^5
    1 <= M <= 10^9

    正整数 N, M
    a1 * a2 * ... * aN = M となる長さNの数列a
    何通りあるか
    数列A1-ai != 数列A2-ai
    """
    N, M = map(int, input().split())

    ans = editorial(N, M)
    print(ans)


def factorize(N):
    """素因数分解

    TODO: もっと速い方法
    素因数がある場合は範囲をsqrtで狭められるが、ない場合愚直に探してしまう
    N以下の素数を見つけて、それだけで対応すれば速いはず
    """
    import math
    factors = []
    while True:
        if N % 2 == 0:
            factors.append(2)
            N //= 2
        else:
            break
    m = int(math.sqrt(N))
    i = 3
    while i <= m:
        add = False
        while True:
            if N % i == 0:
                factors.append(i)
                N //= i
                add = True
            else:
                break
        if add:
            m = int(math.sqrt(N))
        i += 2
    if N > 1:
        factors.append(N)
    return factors


def f(N, M):
    """
    M を素因数分解すると約数とそれらの個数がわかる
    それらの組合せを制約を満たすようにする
    1 * M, M * 1 をわすれないこと:
        1 * 1 * 1 * M もOK

    3 12
        1,1,12: 3P3 / 2!
        2,2,3:  3P3 / 2!
        3,4,1:  3P3
        6,2,1:  3P3

    http://sucrose.hatenablog.com/entry/2014/10/10/235805
    """
    ans = 0
    m = 10 ** 9 + 7
    factors = factorize(M)
    c = Counter(factors)
    print(c)
    return ans


def editorial(N, M):
    """
    https://mathtrain.jp/tyohukuc

    3 12 の例で示せているのに、重複組合せという言葉が頭に浮かばなかった

    解説動画より
    modの世界における逆元?
        重複組合せ
            フェルマーの小定理
            階乗と階乗逆元
        拡張ユークリッド互除法
            今回は出番なし
    """
    def create_tables(N, MOD):
        """
        http://drken1215.hatenablog.com/entry/2018/06/08/210000
        """
        fact_mod = [0] * (N + 1)
        simod = [0] * (N + 1)
        fact_inv_mod = [0] * (N + 1)

        fact_mod[0] = fact_mod[1] = 1
        simod[0] = simod[1] = 1
        fact_inv_mod[0] = fact_inv_mod[1] = 1

        for i in range(2, N + 1):
            fact_mod[i] = (fact_mod[i - 1] * i) % MOD
            simod[i] = MOD - simod[MOD % i] * (MOD // i) % MOD
            fact_inv_mod[i] = (fact_inv_mod[i - 1] * simod[i]) % MOD
            
            
            

        return fact_mod, fact_inv_mod

    def create_fact_mod(N, MOD):
        fact_mod = [0] * (N + 1)
        fact_mod[0] = 1
        fact_mod[1] = 1

        
        for i in range(2, N + 1):
            fact_mod[i] = (fact_mod[i - 1] * i) % MOD

        return fact_mod

    def create_seq_inv_mod(N, MOD):
        simod = [0] * (N + 1)
        simod[0] = 1  
        simod[1] = 1

        
        for i in range(2, N + 1):
            simod[i] = (MOD - MOD // i) * simod[MOD % i] % MOD
            
            

        return simod

    def create_fact_inv_mod(N, MOD, simod):
        fact_inv_mod = [0] * (N + 1)
        fact_inv_mod[0] = 1
        fact_inv_mod[1] = 1

        
        for i in range(2, N + 1):
            fact_inv_mod[i] = (fact_inv_mod[i - 1] * simod[i]) % MOD

        return fact_inv_mod

    def comb_mod(N, R, fact_mod, fact_inv_mod, MOD):
        """
        組合せ nCr: n! / ((n-r)! * r!)
        上記のmodを逆元を使って表現

        f[n] * f[r] % MOD * f[n-r] % MOD 
        これのために事前にテーブルを作っていた
        """
        return fact_mod[N] * fact_inv_mod[N - R] * fact_inv_mod[R] % MOD

    ans = 1
    m = 10 ** 9 + 7
    factors = factorize(M)
    c = Counter(factors)

    
    
    
    fact_mod = create_fact_mod(N + 100, m)
    simod = create_seq_inv_mod(N + 100, m)
    fact_inv_mod = create_fact_inv_mod(N + 100, m, simod)

    
    
    for i in c.values():
        ans *= comb_mod(N - 1 + i, i, fact_mod, fact_inv_mod, m)
        ans %= m

    return ans


if __name__ == '__main__':
    main()

Solution (Python):
import math
mod = 10**9 + 7
n, m = map(int, input().split())


def factorize(n):
    fct = []
    b, e = 2, 0
    while b*b <= n:
        while n%b == 0:
            n //= b
            e += 1
        if e > 0:
            fct.append((b, e))
        b += 1
        e = 0
    if n > 1:
        fct.append((n, 1))
    return fct

num = factorize(m)
ans = 1

for p in num:
    for i in range(p[1]):
        ans *= n+i
    for i in range(p[1]):
        ans //= i+1
print(ans%mod)

Solution (Python):

n,m=map(int,input().split())
l={}
p=2

while m>=p*p:
    if m%p==0:
        if not p in l.keys():
            l[p]=0
        l[p]+=1
        m//=p
    else:
        p+=1
if not m in l.keys():
    l[m]=0
l[m]+=1

ans=1
for p,q in l.items():
    if p==1:
        continue
    for i in range(q):
        ans*=n+i
    for i in range(q):
        ans//=i+1

print(ans%1000000007)

Solution (Python):
N, M = map(int, input().split())
MAX_N = 10 ** 5 + 2000
factors, MOD = {}, 10 ** 9 + 7
fact, fact_inv = [1] * (MAX_N+1), [1] * (MAX_N+1)


def prime_factorization(n):
    i = 2
    while i * i <= n:
        while not n % i:
            factors[i] = factors.get(i, 0) + 1
            n //= i
        i += 1

    if n != 1:
        factors[n] = 1


def combination(n, r):
    return fact[n] * fact_inv[r] * fact_inv[n-r] % MOD


def main():
    prime_factorization(M)

    for i in range(1, N+2001):
        fact[i] = fact[i-1] * i % MOD

    fact_inv[N+2000] = pow(fact[N+2000], MOD-2, MOD)
    for i in range(N+1999, -1, -1):
        fact_inv[i] = fact_inv[i+1] * (i+1) % MOD

    ans = 1
    for i, j in factors.items():
        ans = ans * combination(j+N-1, j) % MOD

    print(ans)


if __name__ == '__main__':
    main()

Solution (Python):
import math
def fact(a,b):
    ans = 1
    while a != b:
        ans *= a
        a -= 1
    return ans
def nCr(n,r):
    return (fact(n,r)) // (math.factorial(n-r))
def nHr(n,r):
    return nCr(n+r-1, r-1)
def prime(n): 
    import math
    num_list = [i + 1 for i in range(2,n,2)]
    list_prime = [2]
    limit = math.sqrt(n)
    if n == 2:
        return list_prime
    else:
        while True:
            p = num_list[0]
            if p >= limit:
                return list_prime + num_list
            list_prime.append(p)
            num_list = [num for num in num_list if num % p != 0]

def primeFactorization(n):
    import math
    list_prime = prime(int(math.sqrt(n)))
    i = 0
    dict_pF = {}
    dict_primeFactorization = {}
    for pri in list_prime:
        dict_pF[pri] = 0
    while True:
        if n == 1:
            for key, value in dict_pF.items():
                if value != 0:
                    dict_primeFactorization[key] = value
            return dict_primeFactorization
        elif i >= len(list_prime):
            dict_pF[n] = 1
            for key, value in dict_pF.items():
                if value != 0:
                    dict_primeFactorization[key] = value
            return dict_primeFactorization
        p = list_prime[i]
        if n % p == 0:
            n //= p
            dict_pF[p] += 1
            continue
        else:
            i += 1

ans = 1
n,m = map(int, input().split())
if m == 1:
    ans = 1
else:
    for factor,degree in primeFactorization(m).items():
        ans *= (nHr(degree,n)) % (10**9+7)
        
print(ans % (10**9+7))

Solution (Python):
import time
import collections
import scipy.misc

def prime_decomposition(n):
  i = 2
  table = []
  while i * i <= n:
    while n % i == 0:
      n /= i
      table.append(int(i))
    i += 1
  if n > 1:
    table.append(int(n))
  return table   

N, M = map(int, input().split())
t_strt = time.time() 

factors = collections.Counter(prime_decomposition(M))
t_elps = time.time()

num_pattern = 1;

for exponent in factors.values():
    tmp_num_pattern = scipy.misc.comb(exponent + N -1, N-1, exact = True) 
    num_pattern = num_pattern * int(tmp_num_pattern)
    
print(num_pattern%(10**9+7))

Solution (Python):
import math




def get_primes(n):
    primes = []
    checks = [True for _ in range(int(math.sqrt(n)) + 1)]
    checks[0] = checks[1] = False
    for i in range(2, len(checks)):
        if checks[i]:
            primes.append(i)
            tmp = 0
            while i + tmp * i < len(checks):
                checks[i + tmp * i] = False
                tmp += 1
    return primes


def combination(a, b):
    i = 1
    for _i in range(a - b + 1, a + 1):
        i *= _i
    return i // math.factorial(b)


n, m = map(int, input().split())

primes = get_primes(m)
elems = []
m_copy = m
divides = dict()
for _p in primes:
    while m_copy % _p == 0:
        m_copy = m_copy // _p
        if str(_p) in divides.keys():
            divides[str(_p)] += 1
        else:
            divides[str(_p)] = 1
if m_copy > 1:
    divides[str(m_copy)] = 1



ans = 1
for i in divides.values():
    ans *= combination(i + n - 1, i)
print(ans % (1000000007))

Solution (Python):
N,M = map(int, input().split())
mod = 10 ** 9 + 7

import math
class PrimeUtils:
    @staticmethod
    def primes_sieve(n):
        sieve = [True] * (n + 1)
        for i in range(2, n + 1):
            if sieve[i]:
                for j in range(i + i, n + 1, i):
                    sieve[j] = False
        return [i for i, x in enumerate(sieve[2:], 2) if x]
    @staticmethod
    
    def pfd(n):
        rn = math.ceil(math.sqrt(n))
        primes = PrimeUtils.primes_sieve(rn)

        factors = []
        cur = n
        for p in primes:
            cnt = 0
            while True:
                div, rem = divmod(cur, p)
                if rem == 0:
                    cnt += 1
                    cur = div
                else:
                    break
            if cnt != 0:
                factors.append((p, cnt))
        if cur != 1 or n == 1:
            factors.append((cur, 1))
        return factors

import itertools, operator, math
class Combination:
    def __init__(self, fact=math.factorial, div=operator.itruediv):
        self.fact = fact
        self.div = div
    def C(self, n, m):
        if n < m: return 0
        if m == 0 or n == m: return 1
        return self.div(
                self.fact(n),
                self.fact(n - m) * self.fact(m)
            )
    def H(self, n, r):
        return self.C(n - 1 + r, r)

class GFp:
    def __init__(self, p):
        self.p = p
    def inv(self, v):
        return pow(v, self.p - 2, self.p)

class CachedFactorial:
    def __init__(self, n, mod):
        self.mod = mod
        self.fact_cache = list(
                itertools.accumulate(
                        itertools.chain(
                                [1],
                                range(1, n+1)
                        ),
                        lambda a,x: (a * x) % mod
                )
        )
    def factorial(self, n):
        return self.fact_cache[n]

class CachedComb(Combination):
    def __init__(self, n, mod):
        self.mod = mod
        self.gfp = GFp(mod)

        self.n = n
        self.fact = CachedFactorial(n, mod)
        super().__init__(
             self.fact.factorial,
             lambda x, y: (x * self.gfp.inv(y)) % self.mod
        )

ps = PrimeUtils.pfd(M)
cc = CachedComb(N + max(x[1] for x in ps), mod)
ans = 1
for p, c in ps:
    if p == 1: continue
    cmb = cc.H(N, c)
    ans = ans * cmb % mod
print(ans)