演算法設計與分析期末專案

1. 題目 固定容量設施分配問題

1.1 題目描述

​假設有\(n\)個設施和\(m\)位客戶，希望從中選擇：

1. 在\(n\)個設施中選擇開啟哪些設施
2. 將客戶分配給開啟的設施

要求開啟設施的總費用和客戶分配到設施的總費用之和最小，且設施的容量不能超過設施中客戶分配費用的總和。

1.2 題目要求

1. 至少兩種演算法解答題目

2. 獲得71種測例的執行結果(具體格式略，在文章中體現)

2. 數學模型分析

2.1 變數定義

在本題目中，為方便分析，定義如下的變數:

1. \(I = \lbrace 0, 1, \cdots, n-1 \rbrace 表示設施 0, 1, \cdots, n-1\)

2. $J = \lbrace 0, 1, \cdots, m-1 \rbrace 表示客戶0, 1, \cdots, m-1 $

3. \(y_{i}\)表示設施\(i\)的開啟狀況，具體如下所示：\[y_{i}=\begin{cases} 1, 設施i開啟 \\ 0, 設施i未開啟 \\ \end{cases}\]

4. \(x_{ij}\)表示客戶\(j\)的分配情況，具體如下所示：\[x_{ij}=\begin{cases} 1, 表示客戶j被分配到設施i \\ 0, 表示客戶j未被分配到設施i \\ \end{cases}\]

5. \(V_{i}\)表示設施\(i\)的容量

6. \(C_{i}\)

   表示設施\(i\)的開啟費用

7. \(D_{ij}\)表示客戶\(j\)分配到設施\(i\)的需求容量

8. \(G_{ij}\)表示客戶\(j\)分到到設施\(i\)的分配費用

9. \(L\)表示分配的總費用

2.2 目標函式

題目的目標函式可以寫成如下數學形式，即求所以開啟設施所需要的費用和所用客戶分配到相應工廠所需的費用之和：
\[L=min \lbrace \sum_{i=0}^{n-1}C_{i}y_{i}+\sum_{i=0}^{n-1}\sum_{j=0}^{m-1}G_{ij}x_{ij} \rbrace\]

2.3 約束條件

目標函式的約束條件如下所示，即每個設施中客戶需求量\(D_{ij}\)

3. 演算法設計與分析

3.1 貪心演算法

解決這個問題的第一種解法是採用貪心演算法，貪心的方式是遍歷每一位客戶\(j \in J\),在設施集合\(I\)中尋找最佳的設施\(i\)，使得到分配客戶\(j\)為止，所得到的總費用\(F_{j}\)最小，則：\[i=\underset{i \in I}{\operatorname{argmin}}\,\lbrace F_{j-1}+y_{i}C_{i}+G_{ij}\rbrace [V_{i}-D_{j} \geq 0] \\ F_{-1} = 0,G = F_{m-1}\]

3.1.1 演算法分析

固定容量的設施問題在常數時間下是不可解的，因此可以使用一些近似的演算法來求解，貪心演算法就是其中的一種。在實驗中，使用顧客設施的分配序列\(S=[s_1,s_2,\cdots,s_m],s_i \in I\)，作為問題的解，使用\(E=[e_1,e_2,\cdots,e_n],n \in \lbrace0,1\rbrace\)表示設施的開啟狀態，總在貪心的過程中，總共需要經歷兩層迴圈，因此演算法的時間複雜度為\(O(n^2)\)。但是貪心演算法畢竟是不是求解固定容量設施的準確解法，因此不可避免會存在誤差，隨著資料規\(n\)模的增大，誤差會以\(lgn\)增大，因此在大資料下，貪心演算法的準確度不高，跟最優解之間的差距可以達到40%，但是在規模較小或者中等的時候，貪心演算法的速度和解的質量還是非常好的，能夠在10%以內。
演算法在C++程式碼中的關鍵部分如下所示：

    // Greedy Algorithm.
    // For every Costomer t.
    for (int t = 0; t < cosNum; t++) {
        int minCos = INT_MAX;
        int bestServe = -1;
        // Costomers t choose facility j.
        // Assume all the facility are open at first.
        for (int j = 0; j < facNum; j++) {
            // fct[j].volumn-csm[t].demand >= 0 means 
            // the fct[j] have enough space to satisfy the demand.
            if ((fct[j].volumn - csm[t].demand >= 0) 
            && (minCos >= G+csm[t].assignment_cost[j])) {
                bestServe = j;
                minCos = G+csm[t].assignment_cost[j];
            }
        }

        serveList[t] = bestServe;
        // Update the current volumns of fct.
        fct[serveList[t]].volumn -= csm[t].demand; 
        // Update the specific fct to open.
        isOpen[serveList[t]] = 1;    
        // Update the current total cost.              
        G = minCos;                                
    }

3.1.2 測例分析

3.2 多鄰域搜尋與模擬退火演算法

多鄰域搜尋演算法常與模擬退火演算法相結合，是一種啟發式演算法。模擬退火的出發點是基於物理中固體物質的退火過程與一般組合優化問題的相似性，其物理退火過程由加溫過程、等溫過程和冷卻過程兩部分組成。多鄰域搜尋是指利用多種鄰域搜尋演算法對解進行遍歷搜尋，這種方法容易陷入區域性最優解中，因此利用模擬退火演算法可以在一定的概率上接受差解，從而跳出區域性最優解，最終收斂到全域性最優解。

3.2.1 多鄰域搜尋

3.2.1.1 更改分配設施

第一種鄰域搜尋的方式是隨機選擇一個客戶\(j \in J\)，在設施集\(I\)中隨機選擇一個不同於當前分配設施的設施\(i\)，如果\(i\)的剩餘容量能夠滿足\(j\)的需求\(D_{ij}\)，則將\(x_{ij}\)的值置為1，原來的設施置為0。實際操作中發現，如果在選擇設施的時候能夠加入貪心，效果會更好，更好後的相關程式碼如下所示：

    if (p <= 0.5) {
        int cosIndex = random(0, cosNum);
        //int newServe = random(0, facNum);


        // New serve should have enough space
        /*
        while ((newServe == serveList[cosIndex]) 
        || (fctTemp[newServe].volumn < csm[cosIndex].demand)) {
            // If do not have enough space.
            // create a new index.
            newServe = random(0, facNum);
        }*/

        // Find the best new serve.
        int newServe = serveList[cosIndex];
        int minCost = currentBestSolution;
        for (int i = 0; i < facNum; i++) {
            // The i-th facility.
            // If the facility has enough space and not currently facility.
            if (serveList[cosIndex] != i 
            && fctTemp[i].volumn >= csm[cosIndex].demand) {
                // Check whether current serve facility is close after csm[cosIndex] is not assigned.
                int currentOpen = (fctTemp[serveList[cosIndex]].volumn
                                +csm[cosIndex].demand == fctTemp[serveList[cosIndex]].Max) ? -1 : 0;
                int newCost = minCost 
                            - csm[cosIndex].assignment_cost[serveList[cosIndex]] 
                            + csm[cosIndex].assignment_cost[i] 
                            + currentOpen * fctTemp[serveList[cosIndex]].cost;
                if (newCost < minCost) {
                    // Accept the new solution.
                    minCost = newCost;
                    newServe = i;
                }
            }
        }

3.2.1.2 交換設施

第二種鄰域搜尋的方式，是隨機選取兩個不相同的客戶\(j_1\)和\(j_2\)，如果原先兩個客戶所服務的設施\(i_1\)和\(i_2\)不相同，而且容量能夠分別滿足\(j_2\)和\(j_1\)的需求，則將兩個客戶所服務的設施交換，在程式碼中的部分如下所示：

    // Case2: swap.
    else if (p > 0.33 && p <= 0.66) {
        // Random select two costomer position.
        int cosNumOne = random(0, cosNum);
        int cosNumTwo = random(0, cosNum);

        // Prevent the two position to be same.
        // Comfiren that the two position have enough space.
        while ((cosNumOne == cosNumTwo)
        || (fctTemp[serveList[cosNumOne]].volumn + csm[cosNumOne].demand < csm[cosNumTwo].demand)
        || (fctTemp[serveList[cosNumTwo]].volumn + csm[cosNumTwo].demand < csm[cosNumOne].demand)) {
            cosNumOne = random(0, cosNum);
            cosNumTwo = random(0, cosNum);
        }
        fctTemp[serveList[cosNumOne]].volumn = fctTemp[serveList[cosNumOne]].volumn + csm[cosNumOne].demand - csm[cosNumTwo].demand;
        fctTemp[serveList[cosNumTwo]].volumn = fctTemp[serveList[cosNumTwo]].volumn + csm[cosNumTwo].demand - csm[cosNumOne].demand;

        newList[cosNumOne] = serveList[cosNumTwo];
        newList[cosNumTwo] = serveList[cosNumOne];
    }

3.2.2 模擬退火

由於多鄰域搜尋容易陷入區域性最優值，因此可以引入模擬退火來跳出區域性最優解，並最終收斂到全域性最優解。
模擬退火可以按照如下步驟進行描述：

給定初溫T=T0，隨機產生出事狀態s=s0,令k=0;
    Repeat:
        Repeat:
            產生新狀態sj=Generate(s)
                if min{1, exp[-(C(sj)-C(s))/tk]}>=random[0,1] s=sj;
        Until 抽樣穩定準則滿足;
        退溫tk+1=update(tk)，並令k=k+1;
    Until 演算法終止準則滿足
輸出演算法搜尋結果

結合三種鄰域搜尋，在程式碼中模擬退火的過程程式碼如下所示：

    while (T_begin > T_end) {
        vector<facility> oldTemp = fct;
        vector<int> oldOpen = isOpen;
        vector<int> serveTemp = serveList;
        for (int i = 0; i < iterator; i++) {
            // According to serveList.

            fctTemp = fct;
            vector<int> newSolution = createNewSolution(fctTemp); // A new serveList.
            
            for (int i = 0; i < fctTemp.size(); i++) {
                openTemp[i] = (fctTemp[i].volumn == fctTemp[i].Max) ? 0 : 1;
            }

            int newSolutionLength = getCurrentCost(newSolution, openTemp);
            int delta = newSolutionLength - currentBestSolution;

            // Annealing.
            if (delta >= 0) {
                double r = (double)rand() / (RAND_MAX);
                // Keep the original solution.
                if (exp(-delta / T_begin) <= r) {
                    // Accept the bad solution.
                    serveList = newSolution;
                    fct = fctTemp;
                    isOpen = openTemp;
                    currentBestSolution = newSolutionLength;
                }
            }
            // Else Accept the solution.
            else {
                serveList = newSolution;
                fct = fctTemp;
                isOpen = openTemp;
                currentBestSolution = newSolutionLength;
            }
        }
        
        if (globalBest > currentBestSolution) {
            globalBest = currentBestSolution;
            fct = oldTemp;
            isOpen = oldOpen;
            serveList = serveTemp;
        }

        cout << T_begin << " 當前全域性最優結果為: " << globalBest << endl;

        // Update the temperature.
        T_begin *= q_rate;
    }

演算法分析

相對於貪心演算法，模擬退火演算法和多鄰域演算法在某些測例中的準度都要高很多，但是程式碼執行的速度會比較慢，而且經驗表明，初溫通常設定為50000，退火率為0.99，內迴圈次數為5000，因此程式碼執行的時間往往會比較長，但是通常如果多鄰域演算法比較好的話，會收斂到一個比較好的結果。

測例分析