10

在回答另一个问题时,我偶然发现了一个问题,即在没有Symbolic Math Toolbox的情况下,我实际上如何能找到整数的所有因数。

例如:

factor(60)

返回:

 2     2     3     5

unique(factor(60))

因此将返回所有素因子,“1”缺失。

 2     3     5

而且我正在寻找一个可以返回所有因素的函数(1数字本身并不重要,但它们会很好)

预期输出x = 60

 1     2     3     4     5     6    10    12    15    20    30    60     

我想出了那个相当庞大的解决方案,除了它可能可以被矢量化,没有任何优雅的解决方案吗?

x = 60;

P = perms(factor(x));
[n,m] = size(P);
Q = zeros(n,m);
for ii = 1:n
    for jj = 1:m
        Q(ii,jj) = prod(P(ii,1:jj));
    end
end

factors = unique(Q(:))'

另外我认为,对于某些大数字,此解决方案将失败,因为 perms 需要矢量长度<11。

4

3 回答 3

13

n您可以通过将数字除以包含整数 1 到 的向量来找到数字的所有因数n,然后找到除以 1 后的余数正好为零(即整数结果):

>> n = 60;
>> find(rem(n./(1:n), 1) == 0)

ans =

     1     2     3     4     5     6    10    12    15    20    30    60
于 2014-01-09T19:11:52.763 回答
11

以下是查找整数因子的六种不同实现的比较:

function [t,v] = testFactors()
    % integer to factor
    %{45, 60, 2059, 3135, 223092870, 3491888400};
    n = 2*2*2*2*3*3*3*5*5*7*11*13*17*19;

    % functions to compare
    fcns = {
        @() factors1(n);
        @() factors2(n);
        @() factors3(n);
        @() factors4(n);
        %@() factors5(n);
        @() factors6(n);
    };

    % timeit
    t = cellfun(@timeit, fcns);

    % check results
    v = cellfun(@feval, fcns, 'UniformOutput',false);
    assert(isequal(v{:}));
end

function f = factors1(n)
    % vectorized implementation of factors2()
    f = find(rem(n, 1:floor(sqrt(n))) == 0);
    f = unique([1, n, f, fix(n./f)]);
end

function f = factors2(n)
    % factors come in pairs, the smaller of which is no bigger than sqrt(n)
    f = [1, n];
    for k=2:floor(sqrt(n))
        if rem(n,k) == 0
            f(end+1) = k;
            f(end+1) = fix(n/k);
        end
    end
    f = unique(f);
end

function f = factors3(n)
    % Get prime factors, and compute products of all possible subsets of size>1
    pf = factor(n);
    f = arrayfun(@(k) prod(nchoosek(pf,k),2), 2:numel(pf), ...
        'UniformOutput',false);
    f = unique([1; pf(:); vertcat(f{:})])'; %'
end

function f = factors4(n)
    % http://rosettacode.org/wiki/Factors_of_an_integer#MATLAB_.2F_Octave
    pf = factor(n);                    % prime decomposition
    K = dec2bin(0:2^length(pf)-1)-'0'; % all possible permutations
    f = ones(1,2^length(pf));
    for k=1:size(K)
      f(k) = prod(pf(~K(k,:)));        % compute products 
    end; 
    f = unique(f);                     % eliminate duplicates
end

function f = factors5(n)
    % @LuisMendo: brute-force implementation
    f = find(rem(n, 1:n) == 0);
end

function f = factors6(n)
    % Symbolic Math Toolbox
    f = double(evalin(symengine, sprintf('numlib::divisors(%d)',n)));
end

结果:

>> [t,v] = testFactors();
>> t
t =
    0.0019        % factors1()
    0.0055        % factors2()
    0.0102        % factors3()
    0.0756        % factors4()
    0.1314        % factors6()

>> numel(v{1})
ans =
        1920

尽管第一个矢量化版本最快,但基于循环的等效实现 ( factors2) 也紧随其后,这要归功于自动 JIT 优化。

请注意,我必须禁用蛮力实现 ( factors5()),因为它会引发内存不足错误(1:3491888400以双精度存储向量需要超过 26GB 的内存!)。这种方法对于大整数显然是不可行的,无论是空间还是时间。

结论:使用以下矢量化实现:)

n = 3491888400;
f = find(rem(n, 1:floor(sqrt(n))) == 0);
f = unique([1, n, f, fix(n./f)]);
于 2014-01-10T14:50:32.293 回答
10

对@gnovice 答案的改进是跳过除法操作:rem单独就足够了:

n = 60;
find(rem(n, 1:n)==0)
于 2014-01-09T20:27:42.190 回答