## Task

Given a positive integer, `n`, return the number of possible ways such that `k` positive integers multiply to `n`. Order matters.

**Examples**
```
n = 24
k = 2
(1, 24), (2, 12), (3, 8), (4, 6), (6, 4), (8, 3), (12, 2), (24, 1) -> 8

n = 100
k = 1
100 -> 1

n = 20
k = 3
(1, 1, 20), (1, 2, 10), (1, 4, 5), (1, 5, 4), (1, 10, 2), (1, 20, 1),
(2, 1, 10), (2, 2, 5), (2, 5, 2), (2, 10, 1), (4, 1, 5), (4, 5, 1),
(5, 1, 4), (5, 2, 2), (5, 4, 1), (10, 1, 2), (10, 2, 1), (20, 1, 1) -> 18
```
**Constraints**
`1 <= n <= 500_000_000`
and `1 <= k <= 1000`

Difficulty: interview

Starter code:

def multiply(n, k):
	

Solutions:
from scipy.special import comb

def multiply(n, k):
    r, d = 1, 2
    while d * d <= n:
        i = 0
        while n % d == 0:
            i += 1
            n //= d
        r *= comb(i + k - 1, k - 1, exact=True)
        d += 1
    if n > 1: 
        r *= k
    return r

from scipy.special import comb
from collections import Counter

def factorize_exponents(n):
    result = Counter()
    while n % 2 == 0:
        result[2] += 1
        n //= 2
    k = 3
    m = int(n ** 0.5)
    while k <= m:
        if n % k == 0:
            result[k] += 1
            n //= k
            m = int(n ** 0.5)
        else:
            k += 2
    if n > 1:
        result[n] += 1
    return list(result.values())

def multiply(n, k):
    factorized = factorize_exponents(n)
    total = 1
    for exp in factorized:
        total *= comb(exp + k - 1, k - 1, exact=True)
    return total


from math import gcd, factorial
from random import randint

def primesfrom2to(n):
    # https://stackoverflow.com/questions/2068372/fastest-way-to-list-all-primes-below-n/3035188#3035188
    """
    Returns an array of primes, 2 <= p < n.
    """
    sieve = [True]*(n//2)
    for i in range(3, int(n**0.5)+1, 2):
        if sieve[i//2]:
            sieve[i*i//2::i] = [False]*((n- i*i -1)//(2*i) + 1)
    return [2] + [2*i + 1 for i in range(1, n//2) if sieve[i]]


bound = 7100
small_primes = primesfrom2to(bound)
small_primes_set = set(small_primes)

def factorint(n):
    """
    Returns a dict containing the prime factors of n as keys and their
    respective multiplicities as values.
    """
    factors = dict()

    limit = int(n**0.5)+1
    for prime in small_primes:
        if prime > limit:
            break
        while n%prime == 0:
            if prime in factors:
                factors[prime] += 1
            else:
                factors[prime] = 1
            n //= prime

    if n >= 2:
        large_factors(n, factors)
    return factors


def large_factors(n, factors):
    while n > 1:
        if isprime(n):
            if n in factors:
                factors[n] += 1
            else:
                factors[n] = 1
            break
        
        factor = pollard_brent(n)
        large_factors(factor, factors)
        n //= factor


def pollard_brent(n):
    # https://comeoncodeon.wordpress.com/2010/09/18/pollard-rho-brent-integer-factorization/
    """
    Returns a factor of n. The returned factor may be a composite number.
    """
    if n%2 == 0:
        return 2

    y, c, m = randint(1, n-1), randint(1, n-1), randint(1, n-1)
    g, r, q = 1, 1, 1
    while g == 1:
        x = y
        for i in range(r):
            y = (pow(y, 2, n)+c)%n

        k = 0
        while k < r and g == 1:
            ys = y
            for i in range(min(m, r-k)):
                y = (pow(y, 2, n)+c)%n
                q = q*abs(x-y)%n
            g = gcd(q, n)
            k += m
        r *= 2

    if g == n:
        while True:
            ys = (pow(ys, 2, n)+c)%n
            g = gcd(abs(x-ys), n)
            if g > 1:
                break

    return g  


def isprime(n):
    # https://en.wikipedia.org/wiki/Miller–Rabin_primality_test#Deterministic_variants
    if n%2 == 0:
        return False
    elif n < bound:
        return n in small_primes_set

    d, r = n-1, 0
    while d%2 == 0:
        d //= 2
        r += 1

    for a in (2, 3, 5, 7):
        x = pow(a, d, n)
        if x == 1 or x == n-1:
            continue

        for i in range(r-1):
            x = pow(x, 2, n)
            if x == n-1:
                break
        else:
            return False

    return True


def partitions(n, k):
    """
    Returns the number of ways to write n as the sum of k nonnegative 
    integers. Order matters.
    """
    result = 1
    for i in range(n+1, n+k):
        result *= i
    return result//factorial(k-1)


def multiply(n, k):
    total, factors = 1, factorint(n)
    for exponent in list(factors.values()):
        total *= partitions(exponent, k)
    return total


def j(n,k):
    i = 1
    for f in range(1,k) : i *= f + n
    for f in range(1,k) : i //= f
    return i
def multiply(n,k):
    i = f = 1
    while f * f < n :
        l,f = 0,-~f
        while n % f < 1 : l,n = -~l,n // f
        if l : i *= j(l,k)
    return i * (n < 2 or j(1,k))

from scipy.special import comb

def multiply(n, k):
    r, d = 1, 2
    while d * d <= n:
        i = 0
        while n % d == 0:
            i += 1
            n //= d
        if i: r *= comb(i + k - 1, k - 1, exact=True)
        d += 1
    if n > 1: 
        r *= k
    return r

from collections import Counter
from functools import reduce
from operator import mul

# The instructions say n<=500_000_000 or 5e8.  However, the tests appear
# to use n<1e12.  We adjust accordingly.

# Use a sieve to pre-compute small primes up to 1e6.  Using trial division
# is sufficent in factoring up to 1e12.  We seed with the first prime to
# simplify is_prime().
small_primes = [2]

# This is suffient for testing primality only up to the square of the
# largest of small_primes.  This is OK to generate the sequence of primes.
def is_prime(n):
    for p in small_primes:
        if p*p > n:
            return True
        if n%p == 0:
            return False

# Sieve to generate the list of small primes.
for n in range(3, 1000000):
    if is_prime(n):
        small_primes.append(n)

# factor an integer up to 1e12
def factor(n):
    if n==1:
        return Counter()
    for p in small_primes:
        if n%p!=0:
            continue
        factors = factor(n//p)
        factors[p] += 1;
        return factors;
    return Counter({ n: 1 })

# count the number of ways of partitioning r objects k ways
def count_partitions(k, r):
    numerator = reduce(mul, list(range(r+1, r+k)), 1)
    denominator = reduce(mul, list(range(1, k)), 1)
    return numerator//denominator

def multiply(n, k):
    factors = factor(n)
    print(n)
    print(factors)
    return reduce(mul, [count_partitions(k, r) for r in list(factors.values())], 1)


from math import sqrt,factorial
def c(n,m):
        return factorial(n)//(factorial(n-m)*factorial(m))

def factor(n):#Amount of prime numbers
        factors={}
        max=int(sqrt(n)+1)
        f=2
        while f<=max:
            if n%f==0:
                factors[f]=0
                while n%f==0:
                    factors[f]+=1
                    n//=f
                max=int(sqrt(n)+1)
            f+=1
        if n!=1:
            factors[n]=1
        return factors
    
def multiply(n, k):
    mul=1
    for m in factor(n).values():
        mul*=c(m+k-1,k-1)
    return mul

import operator as op
from functools import reduce

def factors(n):
    div = 2
    d = dict()

    while n > (div - 1) ** 2:
        if n % div == 0:
            n = n // div
            if div in d:
                d[div] += 1
            else:
                d[div] = 1
        else:
            div += 1    
    if n in d:
        d[n] += 1
    elif n > 1:
        d[n] = 1
    r = []
    for k in d.keys():
        r.append((k, d[k]))
    return r

def get_pascal(r):
    if r == 1:
        return [1]
    elif r == 2:
        return [1,1]
    
    s = [1,1]
    while len(s) < r:
        ns = [1]
        for i in range(len(s) - 1):
            ns.append(s[i] + s[i + 1])
        ns = ns + [1]
        s = ns
    return s

def ncr(n, r):
    if r > n:
        return 0
    r = min(r, n-r)
    numer = reduce(op.mul, range(n, n-r, -1), 1)
    denom = reduce(op.mul, range(1, r+1), 1)
    return numer // denom

def get_cases(r,k):
    p = get_pascal(r)
    sum_v = 0
    for i,v in enumerate(p):        
        sum_v += v * ncr(k,i + 1)
    #print(r, k, p, sum_v)
    return sum_v

def multiply(n,k):
    if n == 1: return 1
    f = factors(n)    
    #print(f)
    s = 1
    for r in f:
        s *= get_cases(r[1], k)
    return s

# I'm slightly sad, I had a simple clean solution, that wasn't quite fast enough, needed 2-3 more seconds
#
# def multiply(n, k):
#     return k if k<2 else sum(multiply(x,k-1) for x in factors(n))
#

from collections import Counter
from functools import lru_cache, reduce
from itertools import combinations
from operator import mul
from math import factorial

def multiply(n, k):
    if k<2: return k
    if k==2: return len(factors(n))
    #get all the combinations that prod to n
    tups = get_factor_tuples(n,k)
    #insert all the extra 1s
    tups = { tuple((1,)*(k-len(x))+x) for x in tups }
    #do factorial math to find # of unique combinations of each tuple (multinomial coefficient)
    f=factorial(k)
    return sum(f // reduce(mul,[ factorial(v) for v in list(Counter(p).values()) if v > 1 ] + [1]) for p in tups)

def get_factor_tuples(n,k):
    pair = {(n,)}
    for kp in range(1,k):
        k_tup = set(x for x in pair if len(x)==kp)
        for grp in k_tup: 
            for i,e in enumerate(grp):
                g = grp[:i]+grp[i+1:]
                if e not in primefactors(n):
                    for a in factors(e):
                        if 1<a<e:
                            b = e//a
                            pair.add(tuple(sorted((a,b)+g)))
    return pair

@lru_cache(maxsize=None)
def primefactors(n):
    i,f = 3,[]
    while n&1 == 0:
        f.append(2)
        n = n >> 1
    while i * i <= n:
        if n % i: i+=2
        else:
            f.append(i)
            n//=i
    if n > 1: f.append(n)
    return f

@lru_cache(maxsize=None)
def factors(n):
    pf = primefactors(n)
    return { 1,n } | { reduce(mul,x) for z in range(1,len(pf)) for x in combinations(pf,z) }