優化MATLAB代碼

此代碼需要很長時間才能運行（超過10分鐘）。有什麼方法可以優化它，以便在不到一分鐘內完成？優化MATLAB代碼

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

來源

2011-10-03 icobes

關鍵是不要問「 N「的因素是什麼，但是」什麼數字有'j作爲因素「，這更容易。 –

@b3 has a great vectorization of rhs.

一個錯字，就更加需要使用times，而不是mtimes：

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

對於lhs，我提出以下建議，鬆散的基礎上埃拉托色尼篩：

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是不顯着的因素）

來源

2011-10-03 21:29:56

感謝您的錯字。在'lhs'計算上幹得不錯。我認爲你需要的只是給整個向量加1。 –

對於諧波總和，我需要把它放在循環內嗎？我不遵循你的代碼如何爲我的每個值做一個和諧總和？ – icobes

@icobes：這是'cumsum'的魔力。請注意'harmonicsum（i）= harmonicsum（i-1）+ 1/i'。 –

內循環執行大約1000000 *（1000000 + 1）/ 2 = 500000500000次！難怪它很慢。也許你應該嘗試一種不同的近似方法。

來源

2011-10-03 20:08:25

我不知道如何測試因素，並測試1,000,000數字中的每一個，而不使用雙循環。 – icobes

向量化你的算法，我可以將執行時間稍微縮短到〜8.5分鐘。計算所有諧波款項在一個聲明：

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

現在可以計算出右側在一個聲明：

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

我無法向量化的因素決定，所以這是我可以用最快的方式來總結它們。 MATLAB的FACTOR命令允許您爲每次迭代生成所有素數因子。然後我們使用UNIQUE和NCHOOSEK來計算所有可能組合的唯一產品組。這避免了測試每個整數作爲一個因素。

 
lhs = zeros(1e6, 1); 
for ii = 1:1e6 
    primeFactor = factor(ii); 
    numFactor = length(primeFactor); 
    allFactor = []; 
    for jj = 1:numFactor-1 
     allFactor = [allFactor; unique(prod(nchoosek(primeFactor, jj), 2))]; 
    end 
    lhs(ii) = sum(allFactor) + 1 + ii; 
end 
lhs(1) = 1;

尋找在這黎曼假設是違反了指數：

 
isViolated = find(lhs > rhs);

來源

2011-10-03 21:05:23

據我所知，你會得到'lhs（1）'，'lhs（2）'，'lhs（3）'的錯誤結果，也可能是所有其他的結果。 –

冉從問題的原始代碼找到'lhs（1:10）'，我現在匹配，你的不是，不是一個單一的元素:( –

我認爲這個問題（除了'ii == 1' case）是'jj == numFactor'產生'ii'本身，這是你分別計算的。將循環改爲'jj = 1：numFactor'並且應該更好。 –

優化MATLAB代碼

回答

相關問題