matlab 降低演算法時間複雜度的方法？

01-29

第一次來知乎求教大神們！懇請多多指教！！我在用matlab編寫一個小演算法，這個演算法裡面可能多次循環的嵌套，導致得到最終結果（輸入Reader=800,Tag=1000，r=30，範圍為[1,900]的時候），花費了將近800多秒！！！勞駕各方神聖給我指點迷津，降低我這個演算法的時間複雜度，有什麼好點子好方法么？
更新：原題是RFID網路冗餘閱讀器去除演算法，即要去除掉系統網路中冗餘的閱讀器，就是圖中的紅色圈圈，下圖是已經去除後的效果。
【演算法執行結果圖如下】

演算法設計：
① definition:
cc(covered count) ——每個tag被reader覆蓋的數量（黑點到紅圈圓心的距離：D&<=r）
nc(neighbor count)——每個reader鄰居的數量（兩個圓圈原點距離：L&<=2*r）
ncc(neighbor cover count)——每個reader的鄰居所覆蓋的所有tag數量的總和
② steps:
a.找出所有cc=1的tag，由於這些tag只被獨有的reader覆蓋，所以這些對應的reader為非冗餘reader，其下所覆蓋的所有tag均被該reader鎖定(hold)；
b.當所有cc=1對應的reader被找到後，即剩下的tag的cc均&>1，剩下的reader根據nc值由低到高依次循環執行c,d,e操作（鄰居越多的閱讀器在實際中越容易對其他閱讀器產生干擾）；
c.當有多個nc值相等的reader，判斷相同nc值的reader的ncc值大小，根據ncc的值由高到低循環進行d,e操作；

d.該reader中所有的tag的cc值減1，該reader為冗餘reader；
e.找出c操作中cc=1的tag，重複a操作，假如依舊沒有cc=1的tag，重複b操作；
f.直到所有tag均被reader鎖定後，去除結束。
【附代碼：】
%x1,y1,x2,y2是隨機數，在其他函數里產生的已知變數，Reader, Tag, r均為GUI上輸入的變數，輸出為冗餘閱讀器數量和整個演算法執行的時間 function [tcba_redundant,TCBAtime]=TCBA(x1,y1,x2,y2,Reader,Tag,r) %TCBA演算法： tic; %%%%%%%%%%%%%【計算tagcc的值】%%%%%%%%%%%%%%%% tagcc=zeros(1,Tag); for j=1:Reader for i=1:Tag if showdistance(x1(j),y1(j),x2(i),y2(i))&<=r tagcc(i)=tagcc(i)+1; end end end % disp(strcat(tagcc: ,num2str(tagcc))); %%%%%%%%%%%%%%%%%%%%%【找出cc=1的tag對應的reader（判斷reader是否冗餘），對於非冗餘的reader所覆蓋的tag進行鎖定】%%%%%%%%%%%%%%%%%%%%%% tagholder=zeros(1,Tag); readerclosed=zeros(1,Reader); readerredundant=ones(1,Reader); for j=1:Reader for i=1:Tag if showdistance(x1(j),y1(j),x2(i),y2(i))&<=r tagcc(i)==1 tagholder(i)==0 tagholder(i)=j; for k=1:Tag if showdistance(x1(j),y1(j),x2(k),y2(k))&<=r tagholder(k)==0 tagholder(k)=j; end end readerclosed(j)=1; readerredundant(j)=0; end end end % disp(strcat(tagholder: ,num2str(tagholder))); % disp(strcat(readerclosed： ,num2str(readerclosed))); % disp(strcat(readerredundant： ,num2str(readerredundant))); %%%%%%%%%%%%%【NC】%%%%%%%%%%%%%%%% for a=1:Reader reader.nc(a)=0; reader.ncc(a)=0; for b=1:Reader if showdistance(x1(a),y1(a),x1(b),y1(b))&<=2*r showdistance(x1(a),y1(a),x1(b),y1(b))&>0 readerclosed(a)==0 readerredundant(a)==1 reader.nc(a)=reader.nc(a)+1; %%%%%%%%%%%%%【NCC】%%%%%%%%%%%%%%%% for i=1:Tag if showdistance(x1(b),y1(b),x2(i),y2(i))&<=r reader.ncc(a)=reader.ncc(a)+1; end end end end end % disp(strcat(reader.nc: ,num2str(reader.nc))); % disp(strcat(reader.ncc: ,num2str(reader.ncc))); %%%%%%%%%%%%%【對於冗餘的reader進行操作】%%%%%%%%%%%%%%%% %%%%%%%%%%%%%【根據reader的nc和ncc值由高到低進行篩選】%%%%%%%%%%%%%%%% while (length(find(reader.nc==0))~=Reader) reader.ncmax=find(reader.nc==max(reader.nc)); reader.ncc2=reader.ncc(reader.ncmax); reader.nccmax=reader.ncmax(find(reader.ncc2==max(reader.ncc2))); ncc=reader.nccmax(1); % disp(strcat(reader.ncmax對應rid: ,num2str(reader.ncmax))); % disp(strcat(reader.ncmax對應ncc值(reader.ncc2): ,num2str(reader.ncc2))); % disp(strcat(reader.nccmax對應rid(ncc)，即要去除的reader為: ,num2str(ncc))); for k=1:Tag if showdistance(x1(ncc),y1(ncc),x2(k),y2(k))&<=r tagcc(k)=tagcc(k)-1; end end reader.nc(ncc)=0; readerclosed(ncc)=1; %對剩下的reader重複最開始的操作進行判斷是否為冗餘reader for j=1:Reader if readerclosed(j)==0 readerredundant(j)==1 for i=1:Tag if showdistance(x1(j),y1(j),x2(i),y2(i))&<=r tagcc(i)==1 tagholder(i)==0 tagholder(i)=j; for h=1:Tag if showdistance(x1(j),y1(j),x2(h),y2(h))&<=r tagholder(h)==0 tagholder(h)=j; end end readerclosed(j)=1; readerredundant(j)=0; end end end end end % disp(strcat(~tagcc: ,num2str(tagcc))); % disp(strcat(~readerclosed： ,num2str(readerclosed))); % disp(strcat(~readerredundant： ,num2str(readerredundant))); % disp(strcat(~reader.nc: ,num2str(reader.nc))); % disp(strcat(~reader.ncc: ,num2str(reader.ncc))); % disp(strcat(~tagholder: ,num2str(tagholder))); %%%%%%%%%%%%%%%%%%%%%【去除沒有hold任何標籤的reader】%%%%%%%%%%%%%%%%%%%%%% a=unique(tagholder); holderid=a(find(a~=0)); disp(strcat(覆蓋有標籤的reader的rid為：,num2str(holderid))); cla(gca); for k=holderid sita=0:pi/20:2*pi;%角度[0,2*pi] xr=x1(k)+r*cos(sita); yr=y1(k)+r*sin(sita); readers=plot(xr,yr,-r); hold on; end %對坐標軸進行相關設置 [tags,TO,tid]=tag(x2,y2,Tag); h=legend([readers,tags],reader,tag); set(h,color,[1 1 0]); %%%%%%%%%%%%%【計算冗餘的reader數量】%%%%%%%%%%%%%%%% tcba_redundant=Reader-length(holderid); fprintf(TCBA演算法得到冗餘閱讀器個數：%d ,tcba_redundant); %%%%%%%%%%%%%%%%%%%%%【TCBA演算法執行所需時間】%%%%%%%%%%%%%%%%%%%%%% time=toc; TCBAtime=num2str(time); fprintf(TCBA演算法執行所需時間為：%s秒 ,TCBAtime); fprintf(------------------------------------------------ );
代碼寫的不漂亮望大神們見諒！！
------------華麗的分割線----------（根據知友@小耿的幫助進行修改結果圖）---------------↓↓↓↓↓↓↓

【具體修改代碼見下面我的回答】

（如果大神來答過就不會有人看我的答案了。所以我要趕在大神之前把答案發出來，這樣就會有人看了來指正，這樣我就會進步。

以上是我的小心思，恩。）

首先抱歉，我不熟悉matlab。以下不能保證語法都正確。

題主你代碼的主要問題在於執行了太多次showdistance()函數。我一塊一塊說。

假設reader數是m，tag數是n。

一，計算tagcc部分

如果是我的話，我會用一個n*m的矩陣（比如A = zeros[n,m]）把每一個reader是否覆蓋某一個tag保存起來（如果覆蓋，就將矩陣中對應的值設為1）。這樣，我之後再也無需調用showdistance()函數來尋找reader覆蓋的tag，而是直接從矩陣中讀取結果：

首先，第j個reader是否覆蓋第i個tag就只要看是否A(i,j)&>0；

其次，第i個tag的CC值是sum(A(i,:))，第j個reader覆蓋量為sum(A(:,j))；

第三，覆蓋第i個tag的reader為find(A(i,:)&>0)，第j個reader覆蓋的tag是find(A(:,j)&>0)。

全部轉化為了矩陣操作，不僅減少了循環次數而且效率高，空間換取時間。

二，判斷reader是否冗餘部分

上面的矩陣就派上用場了：

for j=1:Reader for i = find(A(:,j)&>0) if tagcc(i) == 0 %別忘了先把不覆蓋任何tag的reader去掉 elseif tagcc(i) == 1 %鎖定非冗餘reader ...

三，計算reader.nc部分

第43行：

for b=1:Reader

這裡沒必要從1開始循環，可以從a+1開始。但要注意修改：

1，除了reader.nc(a)，reader.nc(b)也要+1。

2，第44行無需showdistance(x1(a),y1(a),x1(b),y1(b))&>0這個判斷條件了。

另外第44行可以把readerclosed(a)==0這個條件放在showdistance(x1(a),y1(a),x1(b),y1(b))&<=2*r前面。這樣當前者為假時，可以短路跳過後者，減少showdistance()執行次數。

四，演算法部分

後半部分判斷冗餘reader和前面一樣可以優化。我是想說一下前半部分「尋找nc最大值的reader中ncc最大值的那個」。這部分尋找一次的複雜度是O(m)，加上外層while循環就是O(m^2)。

我不太熟悉matlab。如果是別的語言，我會用priority queue，這樣複雜度可以降到O(m·logm)。

你需要兩個priority queue，第一個出隊nc最大值（可以叫ncPq），第二個出隊ncc最大值（可以叫nccPq）。我寫個偽代碼：

所有未鎖定reader加入ncPq while nccPq或ncPq非空 if nccPq為空 ncPq最大值出隊，加入nccPq while ncPq.max == nccPq.max ncPq最大值出隊，加入nccPq end end nccPq最大值出隊，這就是要刪除的reader end

matlab也可以實現priority queue，可參考http://stackoverflow.com/questions/12469502/how-to-implement-priority-queue-in-matlab

分析得比較匆忙，可能有些地方有問題或方法值得改進。還請各位大神在評論區賜教。

============================================

@小耿其實已經把很多需要優化的點都提到了。

準確來說你的計算中性能消耗比較多的部分可能出現在：

判斷 reader 是否冗餘的過程中：對每個 reader 的進行了一次「全局 tag 的遍歷」，再對滿足要求的 tag 與 reader 進行操作（其中「對所有 tag 又進行了一次全局遍歷」）
在尋找 reader 臨近 reader 的過程中：對每個 reader 進行了一次「全局 reader 遍歷」，再對「滿足要求」的 reader 進行了一次「全局 tag 遍歷」
篩選過程中：這部分相對計算較多，我們忽略掉一些小的操作，那麼剩下的是，每一輪篩選過程中，對於滿足要求的 reader，對「全局的 tag 進行遍歷」，又對「滿足要求」的 reader 「在該嵌套內」進行了一次「對全局的 tag 遍歷」

注意我標黑體部分的提示。

1. 對全局判斷的優化

-------------------------------------------------------------------

所謂的「全局」意味著每次遍歷是不加選擇的對整體目標進行一次操作，實際上這個例子中所有需要比較的部分只有相鄰的 reader 與 reader ， reader 與 tag 間。解決這個問題的方案我暫時想到一種（不一定是最好的）：對整體區域進行劃分。

將整個區域劃分成 $a*b$ 個子區域，每個區域大小 $W*H$ 。每個 reader （的圓心）與 tag 都屬於一個區域（判定方法很簡單，只需要 $(x, y)$ 坐標分別對 $(W,H)$ 相除取整），這個區域結構用一個二維Cell Array（Cell Arrays - MATLAB Simulink）來表示。它其中的每個元素包含兩個數組，分別是該區域內 reader 的 id 數組與 tag 的 id 數組。

這樣，每次遍歷就只需要對reader所處區域與它的鄰域內的 reader 或 tag 進行對比就行了。

當然， W 與 H 的選取有一些講究，它要求至少滿足鄰域外的 reader 都不滿足相交， tag 都不滿足覆蓋要求，具體值至少是多少題主可以自行計算一下。

生成這個二維區域數組只需要遍歷一次 reader ，再遍歷一次 tag 。

2. 對覆蓋與相交關係判斷的優化

-------------------------------------------------------------------

同樣「滿足要求」意味著每次對距離進行計算並與某些值比較大小。

這個方面 @小耿提到用 $n*m$ 矩陣來存儲 reader 與 tag 間的覆蓋關係。這是一種用空間換時間的方法，而且在這個例子中 $n*m$ 也不算太大，所以肯定是可行的。

我的方法想法是類似的。利用兩個 Cell Array 來存儲關係。一個存儲對特定的 reader ，它所覆蓋的 tag 的 id 數組；另一個存儲對特定的 reader，它與附近 reader 相交的 reader 的 id 數組。

利用之前對全局判斷的優化，這個工作可以更快完成。

3. 對嵌套的優化

-------------------------------------------------------------------

注意這段代碼

%對剩下的reader重複最開始的操作進行判斷是否為冗餘reader for j=1:Reader if readerclosed(j)==0 readerredundant(j)==1 for i=1:Tag if showdistance(x1(j),y1(j),x2(i),y2(i))&<=r tagcc(i)==1 tagholder(i)==0 tagholder(i)=j; for h=1:Tag if showdistance(x1(j),y1(j),x2(h),y2(h))&<=r tagholder(h)==0 tagholder(h)=j; end end readerclosed(j)=1; readerredundant(j)=0; end end end end

最裡層的 loop 並沒有利用中間層 loop 的結果，因此可以把它從裡面提出來與中間層並列。

============================================

具體 implementation 就交給聰明的題主了（逃

感謝 @小耿知友的幫助！修改了部分代碼後的TCBA演算法執行時間瞬間下降到了14秒，如果沒有其他大程序佔用CPU的話估計能到8，9秒這樣！十分不錯！

在演算法的代碼上還沒做太大修改，還在學慣用profile和priority queue這個東東，邊學邊試著再修改修改。

【修改後結果圖】（這次跑的時候可能電腦在同時運行其他大程序因此比之前耗時都長一點）

【修改後的代碼】

①tagcc的修改：

%%%%%%%%%%%%%【計算tagcc, readercc的值】%%%%%%%%%%%%%%%% rht=zeros(Tag,Reader); for j=1:Reader for i=1:Tag if showdistance(x1(j),y1(j),x2(i),y2(i))&<=r rht(i,j)=1; end end end % disp(num2str(rht)); tagcc=zeros(1,Tag); readercc=zeros(1,Reader); %tagcc：某個標籤被多少個reader覆蓋 for i=1:Tag tagcc(i)=sum(rht(i,:)); end %readercc：某個reader覆蓋多少個tag for j=1:Reader readercc(j)=sum(rht(:,j)); end disp(strcat(readercc: ,num2str(readercc))); disp(strcat(tagcc: ,num2str(tagcc)));

②篩選tagcc=1的代碼修改：

%%%%%%%%%%%%%【找出cc=1的tag對應的reader（判斷reader是否冗餘），對於非冗餘的reader所覆蓋的tag進行鎖定，同時對沒有鎖定任何tag的reader進行去除】%%%%%%%%%%%%% tagholder=zeros(1,Tag); readerclosed=zeros(1,Reader); readerredundant=ones(1,Reader); %找到並去除沒有覆蓋任何tag的reader allzero=find(all(rht==0)); for j=allzero readerclosed(j)=1; end

%鎖定僅被一個reader覆蓋的tag對應的reader only=find(tagcc==1); for i=only oc=find(rht(i,:)==1); readerclosed(oc)=1; readerredundant(oc)=0; tagholder(i)=oc; end for j=oc for k=find(rht(:,j)&>0) if tagholder(k)==0 tagholder(k)=j; end end end

③nc和ncc值的操作：

%%%%%%%%%%%%%【NC】%%%%%%%%%%%%%%%% reader.nc=zeros(1,Reader); reader.ncc=zeros(1,Reader); for a=1:Reader %從a+1開始循環，不用從1開始了 for b=a+1:Reader %把readerclosed(a)==0 readerredundant(a)==1提前優先進行判斷 if readerclosed(a)==0 readerredundant(a)==1 showdistance(x1(a),y1(a),x1(b),y1(b))&<=2*r %找到鄰居後兩個reader.nc同時+1 reader.nc(a)=reader.nc(a)+1; reader.nc(b)=reader.nc(b)+1; %%%%%%%%%%%%%【NCC】%%%%%%%%%%%%%%%% for i=1:Tag %showdistance()函數用rht矩陣代替 if rht(i,b)&>0 reader.ncc(a)=reader.ncc(a)+1; end end end end end

（之後的演算法部分也是都把showdistance()函數利用矩陣rht進行替換，就不貼代碼了。）

（修改待續）

不如先做下profit？看下到底哪裡慢？