找出一个数的所有两个因子的分解。

Question

给定一个 N，我试图找到每个正整数 a 和 b，使得 N = a*b。

我开始使用 sympy.ntheory.factorint 分解为质因数，它给我一个字典因子 -> 指数。

我已经有了这段代码，但我不想获得重复项（a 和 b 扮演相同的角色）：

import itertools

from sympy.ntheory import factorint


def find_decompositions(n):
    prime_factors = factorint(n)
    cut_points = {f: [i for i in range(1+e)] for f, e in prime_factors.items()}
    cuts = itertools.product(*cut_points.values())
    decompositions = [((a := np.prod([f**e for f, e in zip(prime_factors, cut)])), n//a) for cut in cuts]
    return decompositions

例子：

In [235]: find_decompositions(12)
Out[235]: [(1, 12), (3, 4), (2, 6), (6, 2), (4, 3), (12, 1)]

我想要得到的是：

Out[235]: [(1, 12), (3, 4), (2, 6)]

我试图通过使用范围扩展，如e//2，1 + e//2，(1+e)//2，1 + (1+e)//2来减少cut_points中的范围一半。但是这些方法都没有奏效。

一个简单的解决方案显然是计算相同的结果并返回：

decompositions[:(len(decompositions)+1)//2]

但我正在寻找最终能减少计算次数的解决方案。

- Nick Skywalker

---

在你的算法中，“sympy.ntheory.factorint”这个因式分解操作是比较耗费资源的。在此之后，其他的操作都可以被视为“廉价”的。（至少对于具有多个“大”质数因子的数字而言） - MrSmith42

Answer 1

你正在使用模块 `sympy`，该模块已经有一个名为 `divisors` 函数 的函数：

from sympy import divisors

print(divisors(12))
# [1, 2, 3, 4, 6, 12]

def find_decompositions(n):
    divs = divisors(n)
    half = (len(divs) + 1) // 2  # +1 because perfect squares have odd number of divisors
    return list(zip(divs[:half], divs[::-1]))

print(find_decompositions(12))
# [(1, 12), (2, 6), (3, 4)]

print(find_decompositions(100))
# [(1, 100), (2, 50), (4, 25), (5, 20), (10, 10)]

- Stef

---

太棒了！我对sympy不太了解，因为这是我第一次使用。谢谢！ - Nick Skywalker

需要修复以被接受：half = (len(divs)+1)//2（例如，find_decompositions(100) 错过了 (10, 10)） - Nick Skywalker

@NickSkywalker 当然！已修复 - Stef

Answer 2

您可以利用生成“cuts”的方式的优势，这意味着在位置i处生成的每个“cut”都有一个在位置-i处的互补部分。如果 n 是一个完全平方数，您需要一些额外的逻辑。

def find_decompositions_paired(n):
    prime_factors = factorint(n)
    primes = np.array(list(prime_factors.keys()))
    prime_exponents = np.array(list(prime_factors.values()))
    cut_points = [range(1+e) for e in prime_factors.values()]

    cuts = [c for c in itertools.product(*cut_points)]
    decompositions = [
        (a, b) for a, b in zip(cuts[:len(cuts)//2], cuts[::-1])
    ]
    # if n is a square there's going to be an orphan
    if len(cuts) % 2:
        decompositions.append((cuts[len(cuts)//2], cuts[len(cuts)//2]))

    # calculate the value of the primes from the factors
    decompositions = [
        (np.power(primes, d[0]).prod(), np.power(primes,d[1]).prod())
        for d in decompositions
    ]
    return decompositions

测试36的结果是正确的分解

In [156]: find_decompositions_paired(36)
[((0, 0), (2, 2)), ((0, 1), (2, 1)), ((0, 2), (2, 0)), ((1, 0), (1, 2)), ((1, 1), (1, 1))]
Out[156]: [(1, 36), (3, 12), (9, 4), (2, 18), (6, 6)]

大约有1/3的加速效果。

In [177]: %timeit find_decompositions(np.prod(range(1,20))) 
474 ms ± 12.3 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)

In [178]: %timeit find_decompositions_paired(np.prod(range(1,20)))
303 ms ± 7.47 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)

为了完整起见，我们可以看一下原生的sympy divisors函数的表现（根据Stef的回答）。

def find_decompositions_sympy(n):
    div = divisors(n)
    half = len(div) // 2 + 1
    return [(div[i], div[-i-1]) for i in range(half)]

毫不奇怪，它表现得好得多。

In [179]: %timeit find_decompositions_sympy(np.prod(range(1,20)))
13.1 ms ± 443 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)

- Elliot

---

唉...我真是太愚蠢了。谢谢你的指出，我会立即修正答案。 - Elliot

请注意，标准Python库中还有math.prod函数，它似乎更适合用于range，并且返回的是Python整数而不是np.int64。 - Stef

感谢您详细的回答。使用 divisors 函数似乎是最佳方法。 - Nick Skywalker

Answer 3

除非你的N非常大，否则你可以使用一个简单的列表推导式，它迭代到N的平方根，并选择因子及其对应的数值。

def factorPairs(N):
    return [ (f,N//f) for f in range(1,int(N**0.5)+1) if N%f == 0 ]

print(factorPairs(12))
[(1, 12), (2, 6), (3, 4)]

print(factorPairs(36))
[(1, 36), (2, 18), (3, 12), (4, 9), (6, 6)]

- Alain T.

---

我建议使用from math import isqrt，然后使用1 + isqrt(N-1)代替int(N**0.5)+1。 - Stef