換熱網絡合成問題的并行BB/SQP混合算法

姜楠，劉永忠,2，朱天鴻

?

換熱網絡合成問題的并行BB/SQP混合算法

姜楠1，劉永忠1,2，朱天鴻1

（1西安交通大學化學工程與技術學院，陜西西安 710049；2熱流科學與工程教育部重點實驗室，陜西西安 710049）

換熱網絡合成問題通常可用非凸、非線性、不可微的混合整數非線性規劃模型描述。基于GPU的并行計算技術為求解大規模模型提供了高效支撐。針對已有并行SQP算法求解換熱網絡合成問題中存在二元變量組合數過多、并行SQP算法求解結果嚴重依賴初值等問題，提出了BB/SQP混合并行算法。該算法采用BB算法代替枚舉法，不但大大減少了模型求解中可能的二元變量組合，而且為SQP算法選出了可行的初值，從而提高了算法的求解質量。研究表明，所提出的混合并行算法能夠有效求解換熱網絡合成問題，且并行計算相比串行計算的求解速度顯著增加，加速比可達39。

換熱網絡合成；混合整數非線性規劃；GPU；并行算法

引 言

為了獲得最小公用工程費用、最少換熱單元數、最小換熱面積和總投資費用的換熱網絡[1-4]，其問題可抽象為復雜的混合整數非線性規劃（MINLP）模型[5-6]。該模型具有非線性和非凸性的特點。采用啟發式算法難以得到全局最優解，而確定性算法的求解時間過長。因此，目前的趨勢是通過并行加速策略來求解大規模優化模型。目前所采用的主要方法為：直接使用易于并行化的啟發式算法對模型進行求解；或通過模型分割的方法將原問題拆解為可并行求解的一系列子問題，之后對子問題采取并行化處理。

在使用并行的啟發式算法時，由于圖形處理器（GPU）表現出優越的并行計算性能，被引入到通用計算的領域中。一些學者提出了基于GPU并行加速的粒子群算法求解多目標優化問題[7-10]，提高了問題的求解效率和算法的尋優概率。Wong[11]實現了基于GPU的并行遺傳算法，且提出了一種基于GPU的多目標進化算法。Munawar等[12]基于GPU使用改進的遺傳算法求解MINLP問題，使用熵測量的方法調整基因搜索的方向，達到了雙精度下20和單精度下40的加速比。康麗霞等[13]提出了一個基于OpenMP系統的并行遺傳算法，在基本遺傳算法的基礎上引入了一系列調節和控制策略，用于改善算法的收斂性，提高了算法獲得最優解的概率。然而，上述并行求解優化問題的算法僅限于進化類算法，這是由進化類算法天然并行的特點決定的，并未克服進化類算法易于陷入局部最優和易于發散的缺陷。

采用模型分割方法實現并行化，首先需要將MINLP問題拆解成眾多可并行求解的子問題，主要途徑是通過分解模型中的約束條件來達到目的。針對較大規模的MINLP問題，Wah等[14]提出了約束分割法，通過將模型中的約束按一定規則分割，從而得到一定數量的子問題，之后即可對子問題進行并行求解。Bogataj等[15]提出了一種通過將非凸MINLP問題轉化為求解一類單凸MINLP問題的方法，并證明了該方法在模型規模較小時能夠保證全局最優解。Bjorkqvist等[16]提出了基于Internet的分布式并行求解方法，其目的還是需要整合更多的計算資源來輔助并行解決MINLP問題中存在的大量計算密集型任務。由于GPU更擅長處理計算密集型任務，當問題規模不大或者算法處理的任務數目或迭代的次數較少時，GPU并不能顯著提升計算性能。而隨著問題規模的增大，處理數據迭代規模增大時，GPU的加速效果越發明顯。

目前對使用確定性算法求解換熱網絡MINLP問題的研究相對較少，針對確定性算法求解大規模MINLP問題所需時間過長甚至無法在有限時間內求解的問題，康麗霞等[17]提出了基于GPU加速求解MINLP問題的SQP并行算法，采用CPU/GPU異構并行的策略，以枚舉的方式保證了整數變量解空間的完備性，構建了異步并行的內核函數實現了NLP子問題的求解。然而，該算法占用內存空間較大，且初值依賴性強。

針對啟發式算法求解質量差、確定性算法求解時間長等問題，通過分析分支定界算法（BB）的可并行性和序貫二次規劃算法（SQP）求解約束優化問題的全局收斂性及其超線性的收斂速度，本文提出了基于GPU加速求解換熱網絡合成問題的BB/SQP混合并行算法，目標是提高求解精度、加快計算速度。本文將分別介紹串行和并行的BB/SQP混合算法的設計，然后采用串行和并行算法分別對換熱網絡優化問題進行求解，通過對比，闡明本文算法的準確性和求解效率。

1 BB與SQP混合算法流程

分支定界算法（BB）是一種求解純整數或混合整數規劃的有效算法，主要通過分支、定界、比較與剪枝3個關鍵的步驟進行尋優，通常其計算量大大少于窮舉法。序貫二次規劃算法（SQP）是求解有約束的非線性規劃問題（NLP）的有效算法之一，是一種循環迭代算法，具有良好的全局收斂性和超線性的收斂速度。

換熱網絡合成問題可描述為

模型中包含連續型變量和二元變量兩類變量。換熱量或流股溫度等參數用連續型變量表示，是否存在換熱單元用二元變量表示。

根據分支定界法的思想，構造原問題的松弛問題0，表示為

求解該松弛問題的結果可分為3種情況：

（1）該松弛問題無解，則原問題也無解，計算結束；

（2）該松弛問題有解，且最優解滿足原問題的約束，則得到了原問題的最優解，計算結束；

（3）該松弛問題有解，但最優解不滿足原問題的約束，則需要進行下一步的分支。

求解1和2子問題，如果1得到了滿足原問題約束的可行解，則將其作為當前得到的最優解，其相應的目標函數值1作為當前的上界。接著對結點2進行檢查，如果其有解，再分為3種情況：

（1）若2的解不滿足原問題的約束，且2不小于當前的上界，則將2剪枝，1的解即為所求解；

（2）若2的解滿足原問題的約束，且2小于當前的上界，則2的解即為所求的全局最優解；

（3）若2的解不滿足原問題的約束，且2小于當前的上界，則其仍有望得到最優解，以2作為父結點繼續進行分支。

對于第3種情況，選取2的解中某一不滿足整數約束的解，重復分支步驟。如此不斷重復上述的分支、定界和剪枝的過程，如果出現更好的上界，則隨時更新上界為更小的目標函數值。

對換熱網絡合成問題的混合整數非線性規劃（MINLP）模型，本文將其分解為一系列的混合整數規劃（MIP）和非線性規劃（NLP）子問題，其中的MIP問題可以使用BB算法對混合整數組合進行有選擇的取舍，選取其中有可能通往最優解的結點進行分支，同時每一個結點都是一個NLP問題，可以使用SQP算法對NLP問題進行有效的求解。如此一來，可以設計出一個主子結構，主結構使用BB算法作為框架，子結構使用SQP算法求解每個結點中的NLP子問題，隨著BB算法的進行，分支、定界、剪枝的程度愈加深入，不斷增加的約束條件使NLP問題逐漸滿足原MINLP問題的約束，在遍歷了所有的活結點后，最終可求得原MINLP問題的最優解。

2 串行和并行BB/SQP混合算法設計

2.1 串行BB/SQP混合算法的設計

本文所需解決問題的目標是最小化換熱網絡的年度費用，采取優先隊列式分支定界法。將各結點按照優先級進行排序，目標函數值越小的結點優先級越高。因此，本文以最小堆的形式組織活結點表，將活結點按照各自的目標函數值排序。如此可以保證每次取出的表首結點都是活結點表中目標函數值最小的結點，結點一旦被取出后則從活結點表中刪除。

按照優先隊列式分支定界法的原則，混合算法的流程如圖1所示。

圖1 串行BB/SQP混合算法流程

由圖1可以看出，串行BB/SQP混合算法的基本流程如下。

（1）構造換熱網絡優化MINLP問題的松弛問題，記為NLP0。

（2）調用SQP算法模塊求解NLP0，結果分為3種情況：若NLP0無解，則原MINLP問題無解，計算結束；若NLP0有解且解滿足原MINLP問題的所有約束，則得到原問題的最優解，其目標函數值即為最優目標函數值，計算結束；若NLP0有解但解不滿足原問題的約束條件，則轉步驟（3）。

（3）將該結點加入活結點表。

（4）若活結點表不為空，取出活結點表中的優先級最高的結點，將其分支產生兩個子問題NLP1和NLP2，之后將該結點從活結點表中刪除，調用SQP求解模塊對子問題分別求解。

（5）對NLP1和NLP2的求解結果分別進行檢查，若其目標函數值大于當前的上界，將該結點舍棄；否則，轉步驟（6）。

（6）進一步判斷子問題的解是否滿足原問題的約束，若不滿足，轉步驟（3）；若滿足，則得到原問題的一組可行解，將其記錄下來，同時得到一個更優的上界，將上界更新為該結點的目標函數值，轉步驟（4）。

（7）重復步驟（3）～步驟（6），直到活結點表為空時，計算終止。

在分支步驟中，新產生的子結點復制其父結點的最優解作為子結點中SQP算法的迭代初值。相比一般的初值，從父結點繼承的最優解通常情況下能夠更接近子結點的最優解，為SQP算法提供更好的初值，加速問題的求解。

2.2 并行BB/SQP混合算法的設計

由于并行化程序的要求，需要將算法設計為若干適于并行獨立求解的部分，因此如何將任務進行合理分解是本算法的重點。

分支定界算法的搜索起始于一個根結點，這一根結點是將原混合整數規劃中整數約束的部分去掉，將整數變量松弛為連續型變量進行求解的，其解空間會相應增大。通常情況下，根結點的解難以滿足原問題的約束，但是其目標函數值可作為該優化問題最優解的一個下界，可先通過并行SQP算法對根結點求解，如果根結點的解不是原問題的解，則將原問題分解為若干子問題，再利用并行混合算法進行求解。

在換熱網絡合成問題的MINLP模型中，其中的二元變量數目為已知，假設為，則個二元變量的所有組合形式即為完整的換熱網絡結構序列，可通過窮舉得到，其總數為2，或者可以看做是深度為1的滿二叉樹的所有葉子結點所構成的集合，如圖2所示。

圖2中，結點內的序列表示已確定的二元變量組合即部分網絡結構。可以看出，深度不同，結點數目不同，結點的深度與其已確定的不完整的結構序列位數一一對應，且左子結點是在其父結點結構序列的基礎上增加一位，以“0”占位；右子結點是在其父結點結構序列的基礎上增加一位，以“1”占位，即下層的結點均是由上層的結點分支得到，且不會有遺漏。這與需要獲得某一位的組合序列的思路是一致的。為了利用GPU并行計算的特性，可根據問題的規模，選取分支樹的某一深度，列舉出深度為時的所有結點的組合序列共2-1種，它們是最終網絡結構的前-1位的所有組合形式。將這些結點加入初始活結點表中，之后分別初始化于每個block中作為并行深度優先分支定界法搜索的根結點，深度的選取需要同時考慮所要求解問題的整型變量數目以及硬件自身的能力。

圖2 分支樹的深度與產生組合數目的關系

經并行初始化之后，每個block擁有各自預先確定的部分網絡結構。因此，每個block都可看作是一個深度優先分支定界法的根結點，可通過核函數（kernel）執行并行的深度優先搜索（DFS）[18]，如圖3所示。

圖3 根結點與block的關系

圖3中，初始活結點表中的每個結點均作為深度優先搜索的根結點，同時分布于各個block中，之后以其為根結點的深度優先搜索DFS在每個block中執行，除了搜索方式不同，其他執行過程與串行混合算法相似，如此每一個子空間都可以同時進行求解，且每個block中的解均為原問題解空間中的子集。這里所執行的深度優先搜索是一個回溯的過程，其本質是先序遍歷一棵狀態樹，同時將狀態樹的建立隱含在遍歷過程中，使用它可以避免窮舉式搜索，這種方式適用于求解一些組合數相當大的問題[19]。本文所設計的結合并行SQP算法的深度優先BB/SQP搜索流程如圖4所示。

從圖4可以看出，在執行分支時每次只產生一個子結點，子結點的迭代初值可設置為其父結點的最優解，如果該子結點經過比較之后被剪枝，則將搜索返回其父結點，分支產生另一子結點。

圖4 深度優先搜索算法流程

下面給出本文所采用的深度優先BB/SQP搜索的基本步驟：

（1）設置各根結點的初始信息，包括對初始結構序列和變量賦初值，讀入已知數據等；

（2）選取其中某一分支進行試探，即在缺省的網絡結構序列中添加一位“0”或“1”；

（3）使用并行SQP算法對該結點進行求解，之后判斷此分支是否可行，若不可行則轉步驟（6）；

（4）若分支可行則前進一步繼續試探；

（5）找到一組可行解則進行記錄，若未找到則轉步驟（2）；

（6）選取另一分支進行搜索，轉步驟（3），若該結點的分支全部試探完則轉步驟（7）；

（7）退回一步返回父結點，若未退到頭則轉步驟（6）；

（8）已退到頭則結束或輸出無解。

當所有block完成搜索之后，會產生大量的解，此時需要執行同步步驟來確保所有block均已計算完畢，并且將所有block的解匯總比較，進行最優解的篩選。

由于DFS的過程發生在device端，而GPU的流處理器相對簡單，若在device端使用復雜的分支與定界程序會大大影響GPU的執行效率。因此，本文通過使用有序的分支策略和一個簡單的定界函數來替代串行部分復雜的分支定界過程，以此緩解GPU的這一局限。線程間的同步和問題的規模是影響程序執行時間的兩個重要因素，因此應該設置更緊的上界以剪掉更多的分支，創建更小的分支定界樹，避免上界的多次更新同步。因此初始上界的設置就顯得尤為重要，通常可通過啟發式算法或多項式近似的方法來得到相對合適的初始上界。如果問題的規模很大，所需求解的并行的子問題過多時，可按子問題數目選擇合適的block數目，通過多輪實現的形式，每輪計算一部分子問題，將問題分批次處理。圖5給出了并行混合算法的計算流程。

圖5 并行BB/SQP混合算法流程

3 算例分析

本文算例計算所用設備為一臺安裝了Nvidia Tesla K20c和Quadro K600的工作站，計算環境為：CPU為Pentium(R) Dual-Core CPU E5800，3.20 GHz。軟件環境為：Windows 7操作系統，CUDA toolkit版本為6.0。

GPU為Nvidia Tesla K20c顯卡一塊，流處理器內核數量104，CUDA核心數量2496，計算能力3.5，單精度浮點性能3.52 Tflop·s-1，容量4 GB，時鐘頻率0.71 GHz；Quadro K600顯卡一塊，流處理器內核數量8，CUDA核心數量192，計算能力3.0，單精度浮點性能194.84 Gflop·s-1，顯存容量1 GB，時鐘頻率0.88 GHz。

3.1 對小型換熱網絡模型的計算測試

本節采用一個規模較小的換熱網絡案例來進行測試，算例來源于文獻[20]。該換熱網絡包括2股熱流和2股冷流，基礎數據如表1所示。

表1 算例1的流股數據

本例換熱網絡模型規模為2×2×3，無分流，最小傳熱溫差為10 K，模型中共包含68個變量，其中有16個整型變量。采用本文算法得到的換熱網絡結構見圖6。表2中列出了串行程序與并行程序計算結果的對比，并與文獻中的結果進行比較。

由表可見，串行計算和并行計算均得到了比文獻中更優的目標函數值，表明本文算法能夠很好地求解小型換熱網絡優化問題，且并行程序相比串行程序執行時間大大縮短，達到了約39的加速比。

圖6 算例1的換熱網絡結構

表2 算例1的計算結果對比

3.2 對中型換熱網絡模型的計算測試

本節將對一個規模較大換熱網絡的MINLP模型進行測試，該算例來源于文獻[21]，基礎數據見表3。

表3 算例2的流股數據

本算例的規模為4×3×2，有分流，最小傳熱溫差為1 K，模型中共有整型變量31個，連續變量95個。由于串行程序的計算時間過長，本節僅給出采用并行計算進行測試的結果，所得到的換熱網絡結構如圖7所示，結果對比如表4所示。

由表可見，本文設計的并行混合算法可求得算例2的最優整型變量組合，其相應的目標函數值不僅優于文獻中的最優目標函數值，且得到了切實可行的換熱網絡設計方案，證明了本文設計算法可以對換熱網絡的MINLP問題進行有效求解。

圖7 算例2的換熱網絡結構

4 結 論

針對換熱網絡的合成問題，本文分析并設計了將分支定界算法與序貫二次規劃算法相結合的混合算法，同時實現了混合算法的串行計算和基于GPU的并行計算。通過對兩個規模不同的算例進行測試，驗證了本文提出的BB/SQP并行混合算法對于求解換熱網絡優化問題的有效性。對比串行計算和并行計算，本文提出的并行混合算法具有顯著的加速效果。

[1] BAKAR S H A, HAMID M K A, ALWI S R W. Selection of minimum temperature difference (Δmin) for heat exchanger network synthesis based on trade-off plot [J]. Applied Energy, 2016, 162 (15): 1259-1271.

[2] NOVAK P Z, KRAVANJA Z. A methodology for the synthesis of heat exchanger networks having large numbers of uncertain parameters [J]. Energy, 2015, 92 (3): 373-382.

[3] GROSSMANN I E, APAP R M, CALFA B A. Recent advances in mathematical programming techniques for the optimization of process systems under uncertainty [J]. Computers & Chemical Engineering, 2015, 37 (1): 1-4.

[4] LEE S, GROSSMANN I E. A global optimization algorithm for nonconvex generalized disjunctive programming and applications to process systems [J]. Computers & Chemical Engineering, 2001, 25 (11): 1675-1697.

[5] YEE T F, GROSSMANN I E, KRAVANJA Z. Simultaneous optimization models for heat integration (Ⅰ): Area and energy targeting and modeling of multi-stream exchangers [J]. Computers & Chemical Engineering, 1990, 14 (10): 1151-1164.

[6] YEE T F, GROSSMANN I E. Simultaneous optimization models for heat integration (Ⅱ): Heat exchanger network synthesis [J]. Computers & Chemical Engineering, 1990, 14 (10): 1165-1184.

[7] ZHOU Y, TAN Y. GPU-based parallel multi-objective particle swarm optimization [J]. International Journal of Artificial Intelligence, 2011, 7 (11): 125-141.

[8] WANG J, WU Z. Opposition-based particle swarm optimization for solving large scale optimization problems on GPU [J]. Journal of Wuhan University(Natural Science), 2011, 57 (2): 148-154.

[9] MUSSI L, DAOLIO F, CAGNONI S. Evaluation of parallel particle swarm optimization algorithms within the CUDA architecture [J]. Information Sciences, 2011, 181 (20): 4642-4657.

[10] MUSSI L, NASHED Y S G, CAGNONI S. GPU-based asynchronous particle swarm optimization [C]//Proceedings of the 13th Annual Conference on Genetic and Evolutionary Computation. Dublin: ACM. 2011: 1555-1562.

[11] WONG M L. Parallel multi-objective evolutionary algorithms on graphics processing units [C]//Proceedings of the Conference on Genetic and Evolutionary Computation. Canada: Companion Material, 2009: 2515-2522.

[12] MUNAWAR A, WAHIB M, MUNETOMO M. Advanced genetic algorithm to solve MINLP problems over GPU [C]//Congress of Evolutionary Computation. LA: IEEE, 2011: 318-325.

[13] 康麗霞, 姜楠, 夏明星, 等. 基于OpenMP的并行GA加速求解換熱網絡設計 [J]. 高校化學工程學報, 2016, 30 (2): 423-438. KANG L X , JIANG N, XIA M X,Paralleled genetic algorithm for design of HEN based on OpenMP system [J]. Journal of Chemical Engineering of Chinese Universities, 2016, 30 (2): 423-438.

[14] WAH B W, CHEN Y X. Constraint partitioning in penalty formulations for solving temporal planning problems [J]. Artificial Intelligence, 2006, 170 (3): 187-231.

[15] BOGATAJ M, KRAVANJA Z. An alternative strategy for global optimization of heat exchanger networks [J]. Applied Thermal Engineering, 2012, 43 (43): 75-90.

[16] BJORKQVIST J, WESTERLUND T. Parallel solution of disjunctive MINLP problems [J]. Chemical Engineering Communications, 2002, 185 (1): 115-124.

[17] 康麗霞, 張燕蓉, 唐亞哲, 等.基于GPU加速求解MINLP問題的SQP并行算法 [J]. 化工學報, 2012, 63 (11): 3597-3601. KANG L X, ZHANG Y R, TANG Y Z,. Paralleled SQP algorithm for solution of MINLP problems based on GPU acceleration [J].CIESC Journal, 2012, 63 (11): 3597-3601.

[18] CHAKROUN I, MELAB N. Operator-level GPU-accelerated branch and bound algorithms [J]. Procedia Computer Science, 2013, 18 (1): 280-289.

[19] CARNEIRO T, MURITIBA A E, NEGREIROS M. A new parallel schema for branch-and-bound algorithms using GPGPU [C]//Computer Architecture and High Performance Computing. Espirito Santo: IEEE, 2011: 41-47.

[20] ZAMORA J M, GROSSMANN I E. A global MINLP optimization algorithm for the synthesis of heat exchanger networks with no stream splits [J]. Computers & Chemical Engineering, 1998, 22 (3): 367-384.

[21] VIDYASHANKAR K, NARASIMHAN S. Comparison of heat exchanger network synthesis using Floudas and Yee superstructures [J]. Indian Chemical Engineer, 2010, 52 (1): 1-22.

Parallel algorithm of hybrid BB/SQP for heat exchanger network synthesis

JIANG Nan1, LIU Yongzhong1,2, ZHU Tianhong1

(1School of Chemical Engineering and Technology, Xi’an Jiaotong University, Xi’an 710049, Shaanxi, China; 2Key Laboratory of Thermo-Fluid Science and Engineering, Ministry of Education, Xi’an 710049, Shaanxi, China)

Heat exchanger network synthesis can be described by a mixed integer non-linear programming (MINLP) model, which features non-convex, non-linear and non-differentiable optimization. The parallel computing technology based on GPU provides an efficient support for solving large scale models. In this work, a hybrid algorithm combined branch and bound method (BB) with sequential quadratic programming (SQP) is proposed to overcome difficulties in the existing parallel SQP algorithm, such as too many combinations of integer variables, dependency of initial values and. The BB method is adopted in the hybrid algorithm instead of the exhaustive method. It can not only reduce the combinations of integer variables, but also select feasible initial values for the SQP algorithm. The solution quality is improved. The results of the examples show that the proposed parallel hybrid algorithm can solve the heat exchanger network synthesis problems efficiently. Compared to the serial algorithm, the proposed parallel algorithm has much higher executive speed with the speedup ratio of 39.

heat exchanger network synthesis; mixed integer non-linear programming; GPU; hybrid algorithm

date: 2016-09-13.

Prof. LIU Yongzhong, yzliu@mail.xjtu.edu.cn

10.11949/j.issn.0438-1157.20161287

TQ 021.8

A

0438—1157（2016）12—5169—07

國家自然科學基金項目（21376188，21676211）。

supported by the National Natural Science Foundation of China (21376188, 21676211).

2016-09-13收到初稿，2016-09-22收到修改稿。

聯系人：劉永忠。第一作者：姜楠（1990—），男，碩士研究生。