5

Julia から始めて、MATLAB コードを Julia に (基本的には行ごとに) 変換しました。Julia コードはかなり遅い (50x など) ことに気付きました。元の問題は、値関数を補間する動的プログラミングの問題です。その補間は、コードがほとんどの時間を費やす場所です。そこで、パフォーマンスの違いを示す最小限のサンプル コードを作成しようとしました。注意すべき重要な点は、これが補間のスプライン近似であり、グリッドが好ましくは不規則である、つまり等間隔ではないということです。MATLAB コード:

tic
spacing=1.5;
Nxx = 300; 
Naa = 350;
Nalal = 200; 
sigma = 10;
NoIter = 500;

xx=NaN(Nxx,1);
xmin = 0.01;
xmax = 400;
xx(1) = xmin;
for i=2:Nxx
    xx(i) = xx(i-1) + (xmax-xx(i-1))/((Nxx-i+1)^spacing);
end

f_U = @(c) c.^(1-sigma)/(1-sigma);  
W=NaN(Nxx,1);
W(:,1) = f_U(xx);

xprime = ones(Nalal,Naa);
for i=1:NoIter
     W_temp = interp1(xx,W(:,1),xprime,'spline');
end
toc

Elapsed time is 0.242288 seconds.

ジュリアコード:

using Dierckx
function performance()

const spacing=1.5 
const Nxx = 300 
const Naa = 350
const Nalal = 200 
const sigma = 10 
const NoIter = 500 

xx=Array(Float64,Nxx)
xmin = 0.01
xmax = 400
xx[1] = xmin
for i=2:Nxx
    xx[i] = xx[i-1] + (xmax-xx[i-1])/((Nxx-i+1)^spacing)
end

f_U(c) =  c.^(1-sigma)/(1-sigma)
W=Array(Float64,Nxx,1)
W[:,1] = f_U(xx)

W_temp = Array(Float64,Nalal,Naa)

interp_spline = Spline1D(xx,W[:,1])
xprime = ones(Nalal,Naa)
for i=1:NoIter
    for j=1:Naa
         for k=1:Nalal
         W_temp[k,j] = evaluate(interp_spline, xprime[k,j])
         end
    end
end

end

@time(performance())
4.200732 seconds (70.02 M allocations: 5.217 GB, 4.35% gc time)

元の問題では行列なので、W を 2 次元配列にしました。さまざまな補間パッケージについて調査しましたが、不規則なグリッドとスプラインのオプションはあまりありませんでした。MATLAB の interp1 は明らかに利用できません。

私の問題は、Julia コードを記述して MATLAB と同じ結果が得られる場合、Julia の方が一般的に高速であると考えていたことです。しかし、明らかにそうではないので、コーディングに注意を払う必要があります。私はプログラマーではありません。もちろん最善を尽くしていますが、ここで簡単に修正できる明らかな間違いを犯しているのか、それとも（あまりにも）頻繁に発生するのか、注意を払う必要があるのか​​ を知りたいです。私の Julia コーディング - それを学ぶ価値がないかもしれないからです。同様に、ここで Julia を高速化できれば (たとえば、割り当てが少し大きく見えるなど)、おそらく MATLAB も高速化できます。Julia に対する私の希望は、同様のコードの場合、MATLAB よりも高速に実行されることでした。

タイミングについていくつかコメントした後、次のコードも実行しました。

using Dierckx

tic()
const spacing=1.5 
const Nxx = 300 
const Naa = 350
const Nalal = 200 
const sigma = 10 
const NoIter = 500 

xx=Array(Float64,Nxx)
xmin = 0.01
xmax = 400
xx[1] = xmin
for i=2:Nxx
     xx[i] = xx[i-1] + (xmax-xx[i-1])/((Nxx-i+1)^spacing)
end

f_U(c) =  c.^(1-sigma)/(1-sigma)
W=Array(Float64,Nxx,1)
W[:,1] = f_U(xx)

W_temp = Array(Float64,Nalal,Naa)

interp_spline = Spline1D(xx,W[:,1])
xprime = ones(Nalal,Naa)
for i=1:NoIter
    for j=1:Naa
         for k=1:Nalal
         W_temp[k,j] = evaluate(interp_spline, xprime[k,j])
         end
     end
 end
 toc()

elapsed time: 
7.336371495 seconds

さらに遅い、うーん...

別の編集: この場合、ループを 1 つ削除すると、実際には高速になりますが、それでも MATLAB に匹敵するものではありません。コード：

function performance2()

const spacing=1.5 
const Nxx = 300 
const Naa = 350
const Nalal = 200 
const sigma = 10 
const NoIter = 500 

xx=Array(Float64,Nxx)
xmin = 0.01
xmax = 400
xx[1] = xmin
for i=2:Nxx
    xx[i] = xx[i-1] + (xmax-xx[i-1])/((Nxx-i+1)^spacing)
end

f_U(c) =  c.^(1-sigma)/(1-sigma)
W=Array(Float64,Nxx,1)
W[:,1] = f_U(xx)

W_temp = Array(Float64,Nalal,Naa)

interp_spline = Spline1D(xx,W[:,1])
xprime = ones(Nalal,Naa)
for i=1:NoIter
    for j=1:Naa
         W_temp[:,j] = evaluate(interp_spline, xprime[:,j])
    end
 end

 end

@time(performance2())
1.573347 seconds (700.04 k allocations: 620.643 MB, 1.08% gc time)

別の編集で、ループを同じ回数繰り返します。

function performance3()

const spacing=1.5 
const Nxx = 300 
const Naa = 350
const Nalal = 200 
const sigma = 10 
const NoIter = 500 

xx=Array(Float64,Nxx)
xmin = 0.01
xmax = 400
xx[1] = xmin
for i=2:Nxx
    xx[i] = xx[i-1] + (xmax-xx[i-1])/((Nxx-i+1)^spacing)
end

f_U(c) =  c.^(1-sigma)/(1-sigma)
W=Array(Float64,Nxx,1)
W[:,1] = f_U(xx)

W_temp = Array(Float64,Nalal,Naa)
W_tempr = Array(Float64, Nalal*Naa)

interp_spline = Spline1D(xx,W[:,1])
xprime = ones(Nalal,Naa)
xprimer = reshape(xprime, Nalal*Naa)

for i=1:NoIter
        W_tempr = evaluate(interp_spline, xprimer)
end

W_temp = reshape(W_tempr, Nalal, Naa)
end

tic()
performance3()
toc()

elapsed time: 
1.480349334 seconds

それでも、MATLAB とはまったく比較になりません。ちなみに、私の実際の問題では、ループを簡単に50k回実行しており、より大きなxprime行列にアクセスしていますが、その部分が違いを生むかどうかはわかりません.

4

1 に答える 1