优化 MATLAB 代码

Question

此代码需要很长时间才能运行（超过 10 分钟）。有什么方法可以优化它，使其在不到一分钟的时间内完成？

clear all;
for i = 1:1000000
    harmonicsum = 0;
    lhs = 0;
    for j = 1:i
        % compute harmonic sum
        harmonicsum = harmonicsum + 1/j;
        % find sum of factors
        if (mod(i,j)==0)
            lhs = lhs + j;
        end
    end
    %define right hand side (rhs) of Riemann Hypothesis
    rhs = harmonicsum + log(harmonicsum) * exp(harmonicsum);

    if lhs > rhs
        disp('Hypothesis violated')
    end
end

Answer 1

内部循环执行大约 1000000*(1000000+1)/2 = 500000500000 次！难怪它很慢。也许你应该尝试不同的近似方法。

Answer 2

向量化您的算法，我可以将执行时间略微减少到约 8.5 分钟。在一个语句中计算所有谐波和：

Harmonicsum = cumsum(1 ./ (1:1e6));

您现在可以在一个语句中计算右侧：

rhs =harmonicsum + log(harmonicsum) .* exp(harmonicsum);

我无法向量化这些因素的确定，所以这是我能想到的对它们求和的最快方法。 MATLAB 的FACTOR 命令允许您为每次迭代生成所有主要因子。然后，我们使用UNIQUE 和NCHOOSEK 计算所有可能组合的唯一乘积集。这避免了将每个整数作为一个因素进行测试。

lhs = 零（1e6，1）； 对于 ii = 1:1e6 primeFactor = 因子（ii）； numFactor = 长度（primeFactor）； allFactor = []; 对于 jj = 1:numFactor-1 allFactor = [allFactor;独特的（产品（nchoosek（primeFactor，jj），2））]； 结尾 lhs(ii) = sum(allFactor) + 1 + ii; 结尾 lhs(1) = 1;

找出违反黎曼假设的指数：

isViolated = find(lhs > rhs);

Answer 3

@b3 has a great vectorization of rhs.

不过有一个错字，需要使用times 而不是mtimes：

harmonicsum = cumsum(1 ./ (1:1e6));
rhs = harmonicsum + log(harmonicsum) .* exp(harmonicsum);

对于lhs，我大致基于 Eratosthenes 的筛子提出以下建议：

lhs = 1 + [1:1e6];
lhs(1) = 1;
for iii = 2:numel(lhs)/2
    lhs(2*iii:iii:end) = lhs(2*iii:iii:end) + iii;
end;

执行时间仅为 2.45 秒（对于这半个问题）。包括 rhs 和 find 的计算在内的总计不到 3 秒。

我目前正在运行其他版本以确保结果相同。

---

编辑：发现lhs(1) 的错误并对其进行特殊处理（这是一种特殊情况，唯一的自然数，其中 1 和 N 不是不同的因数）