2013-12-12 149 views
4

我真的很努力在MATLAB上優化微積分代碼。嵌套循環優化

獲取非線性材料的材料屬性是一項繁重的計算。

此計算需要超過240萬步。 它本身相當簡單,因爲它包含了大量的sum。 唯一的問題是數字存儲在各種數組和列表中,這有點令人困惑。

下面是原來的代碼:

Tensor=zeros(3,3,3,3); 
for m=1:3 
    for n=1:3 
     for o=1:3 
      for p=1:3 
       for x=1:16 
        for y=1:16 
         for z=1:16 
          for i=1:3 
           for j=1:3 
            for k=1:3 
             for l=1:3 
              for r=1:3 
               for s=1:3 
                sum=sum+(1/(8*(pi()^2))*P{x,y,z}(i,m)*P{x,y,z}(j,n)*P{x,y,z}(k,o)*P{x,y,z}(l,p)*(TensorC(i,j,k,l)-TensorC0(i,j,r,s))*TensorA(r,s,k,l)*sin(omega(x))*p_omega(x)*p_phi(y)*p_beta(z); 
               end 
              end 
             end              
            end 
           end 
          end 
         end 
        end 
       end 
       Tensor(m,n,o,p)=sum; 
      end 
     end 
    end 
end 

P{x,y,z}(i,m)是依據式I的變化(同爲其他):i and m確定式的類型並且將結果與xyz變量計算。

總和中的所有其他數字都是在數組和張量中拾取的實數。

我試圖以儘快計算出它們的減少操作次數,提取從最後的循環儘可能多的變量:

Tensor=zeros(3,3,3,3); 
CO1=1/(8*(pi()^2)); 
for m=1:3 
    for n=1:3 
     for o=1:3 
      for p=1:3 
      sum=C0_tensor(m,n,o,p); 
       for x=1:16 
        CO7=sin(omega(x)); 
        CO8=p_omega(x); 
        for y=1:16 
         CO9=p_phi(y); 
         for z=1:16 
          CO10=p_beta(z); 
          for i=1:3 
           CO2=P{x,y,z}(i,m); 
           for j=1:3 
            CO3=P{x,y,z}(j,n); 
            for k=1:3 
             CO4=P{x,y,z}(k,o); 
             for l=1:3 
              CO5=P{x,y,z}(l,p); 
              CO6=TensorC(i,j,k,l); 
              for r=1:3 
               for s=1:3 
                CO11=TensorC0(i,j,r,s); 
                CO12=TensorA(r,s,k,l); 
                sum=sum+CO1*CO2*CO3*CO4*CO5*(CO6-CO11)*CO12*CO7*CO8*CO9*CO10; 
               end 
              end 
             end              
            end 
           end 
          end 
         end 
        end 
       end 
       Tensor(m,n,o,p)=sum; 
      end 
     end 
    end 
end 

儘管如此,計算實在是太長了。

我看不出並行或矢量化計算的任何方式......

從一個陣列或一個矩陣檢索一個特定值的操作似乎很慢...

待辦事項你認爲我應該建立一個包含所有值而不是使用多個值的巨大張量?

+1

如果真的很難避免循環,你應該開始考慮使用其他有利於循環的編程語言,比如'C++' – Ray

+1

8o ...好親切......並祝你好運。你對matlab代碼有很大的依賴性嗎? – Brian

+0

量化*太長*。 –

回答

1

由於您正在用相同的名稱覆蓋有用的函數,因此不應將sum用作變量名稱。

以剛內環這裏,你計算每個rs單個值被添加到您的輸出值:

for r=1:3 
     for s=1:3 
       CO11=TensorC0(i,j,r,s); 
       CO12=TensorA(r,s,k,l); 
       sum=sum+CO1*CO2*CO3*CO4*CO5*(CO6-CO11)*CO12*CO7*CO8*CO9*CO10; 
     end 
    end 

但是,你可以立刻得到C011/C012爲3 ×3個矩陣,做一個總和,然後添加到您的輸出: (改變了你的sumout這裏,請注意*,而不是*在適當的地方。):

C011 = squeeze(TensorCO(i,j,:,:)); 
C012 = squeeze(TensorCO(:,:,k,l)); 
s = CO1*CO2*CO3*CO4*CO5*(CO6-CO11).*CO12*CO7*CO8*CO9*CO10; 
out = out + sum(s(:)); 

而且,當你這樣做:

CO7=sin(omega(x)); 
CO8=p_omega(x); 

(這只是那麼C07 * C08在以後的公式) - 你不需要((x)的歐米茄),每n,m,p.重新計算的罪,因此可以取出完全的循環。

預計算sin和由p_omega(外側環)相乘:

omega78 = sin(omega).*p_omega; 

然後,只需檢索在x環路和使用而不是C07*C08C78 = omega78(x)

+0

非常感謝您的幫助! 我運用了你的想法,並調整了一些數學和數組/張量結構。 基本上,主循環速度快了10倍。 我會繼續努力尋找新的想法! – Empirehell

+1

我優化了張量和數組結構。 現在主循環運行9分鐘而不是2小時! – Empirehell