劉 哲，李永明，周 悅，鄒岸新，徐祿文

（1. 重慶大學(xué)電氣工程學(xué)院，重慶 400030；2. 國網(wǎng)重慶市電力公司電力科學(xué)研究院，重慶 401123）

0 引言

我國地緣遼闊，特高壓輸電能夠?qū)崿F(xiàn)大容量、遠距離輸電，降低損耗，節(jié)省線路走廊，有利于構(gòu)建堅強電網(wǎng)[1]，正是高電壓等級輸電的這些優(yōu)勢，我國近些年發(fā)展建設(shè)了很多特高壓線路。隨著高壓輸電線路電壓等級由超高壓升高為特高壓，特高壓輸電系統(tǒng)電壓高、電流大，帶電導(dǎo)體表面及附近空間的電場強度明顯增大，電暈放電產(chǎn)生的可聽噪聲和無線電干擾影響突出[2]。輸變電系統(tǒng)電磁場的工程計算是一項需要追求精度和時間的任務(wù)，為得到更為準(zhǔn)確的結(jié)果，輸電線路的仿真計算趨于三維化，同時考慮的環(huán)境因素也趨于復(fù)雜化，這些所帶來的直接結(jié)果是所要計算的數(shù)據(jù)越來越多，給計算機帶來了前所未有的挑戰(zhàn)，所以怎樣快速并又準(zhǔn)確的計算出周圍的電磁場分布一直是科研工作者追求的目標(biāo)。

計算機技術(shù)日新月異，將計算機的技術(shù)優(yōu)勢最大限度地應(yīng)用于電力系統(tǒng)，是目前國際智能電網(wǎng)的研究重點之一[3]。其中，計算機技術(shù)中的并行計算是解決大數(shù)據(jù)計算與時間復(fù)雜度之間矛盾的有效方法，也只有并行計算機可以解決今天日益增加的計算需要所面臨的困難[4]。并行計算也不是什么新鮮事，在光化學(xué)，以及當(dāng)下熱門的云計算領(lǐng)域都有并行計算的身影[5]-[11]。并行計算就是可以同時對多個任務(wù)，多條指令進行處理，從而縮短計算時間提高計算的效率。目前并行計算主要有 2種方式，MPI和 OpenMP。MPI是通過多機群組合，并在機間通過消息傳遞的方式進行交流，達到并行計算的效果。但這種多機集群環(huán)境，內(nèi)存不共享，在全局共享數(shù)據(jù)讀寫操作中，必須依靠機器間的通信來進行數(shù)據(jù)搬遷[12]，所以往往造成通信開銷遠遠高于計算本身所花的時間。而多核環(huán)境中內(nèi)存是共享的，可以省下通信開銷，鎖操作的效率遠大于通信效率，并且在過去計算機的發(fā)展中，計算機都是單核高頻的追求計算性能，但是在現(xiàn)有工藝下，改善 CPU 性能的傳統(tǒng)方法如提升時鐘速度和指令吞吐量等在摩爾定律限制下已經(jīng)難有大的進展，單核發(fā)展已處于瓶頸狀態(tài)，多核計算將成為未來的發(fā)展趨勢[13]，所以過去的多臺單核計算可以簡化到現(xiàn)在一臺多核計算機上，實現(xiàn)經(jīng)濟優(yōu)化，編程優(yōu)化等優(yōu)勢。

因此，本文采用可應(yīng)用于單機多核的 OpenMP的并行計算方法，結(jié)合模擬電荷法研究特高壓輸電走廊的電場并行計算，仿真分析了離地1.5m工頻電場的分布，對比串行算法驗證其準(zhǔn)確性，并分析并行計算在效率上的優(yōu)勢。

1 OpenMP 并行計算

OpenMP是一個為共享存儲結(jié)構(gòu)計算編寫并行程序的接口，1997年10月，由SGI公司發(fā)起，HP、SUN、IBM和Intel等聯(lián)合開發(fā)和推行一種新的共享主存編程模型的建議標(biāo)準(zhǔn)。其跨平臺，可伸縮的模型給并行程序的設(shè)計者在桌面電腦和巨型機上開發(fā)并行程序提供了簡單靈活的接口，在不降低性能的情況下克服了機器專用編譯指示所帶來的移植性問題，得到了各界的認可和接受[14]。編程模型以線程為基礎(chǔ)[15]，主要應(yīng)用于對稱式共享存儲器多處理機結(jié)構(gòu)的并行計算機（SMP）。SMP是一類常見的共享存儲并行計算機，其內(nèi)部結(jié)構(gòu)圖如下1所示。

圖1 SMP 結(jié)構(gòu)圖Fig.1 Structure chart of SMP

該類結(jié)構(gòu)最重要的兩個特點就是對稱性和共享性，對稱性是指系統(tǒng)中各個處理器的結(jié)構(gòu)相同，地位相等，可以對稱地訪問存儲器以及 I/O設(shè)備；共享性指的是存儲器的地址共享，具有統(tǒng)一且唯一的地址空間，所有的處理器共享該地址空間。

OpenMP最基本的執(zhí)行模式為fork/join并行執(zhí)行模型，它編程簡單并且可以增量開發(fā)。在OpenMP編程中，整個程序被看作一個大的串行執(zhí)行結(jié)構(gòu)，也可以叫做隱式初始任務(wù)區(qū)，并行結(jié)構(gòu)只是其中的一個執(zhí)行單元。并行結(jié)構(gòu)中可以嵌套多個子并行結(jié)構(gòu)，如果環(huán)境變量設(shè)置允許嵌套并行，這些嵌套子并行結(jié)構(gòu)內(nèi)還包含多個子孫并行結(jié)構(gòu)。如圖2所示，執(zhí)行過程大致為，一個OpenMP程序從一個主線程開始，在程序某點需要并行時程序派生出一些額外的線程組成線程組，被派生出來的線程稱為組的從屬線程，并行區(qū)域中的代碼在不同的線程中并行執(zhí)行，當(dāng)并行域中的代碼執(zhí)行完之后，其他并行線程銷毀，控制流程回到單獨的主線程來執(zhí)行[16]。由于是共享內(nèi)存的編程模型，降低了通信開銷，用戶不必處理迭代劃分、數(shù)據(jù)共享、線程調(diào)度及同步等低級別的細節(jié)[17]。

圖2 并行流程Fig.2 Parallel process

該模式的關(guān)鍵是判斷數(shù)據(jù)間的依賴性，如果不存在依賴性，則可以直接插入#pragma omp parallel for 編譯；如果存在依賴性，則需要對數(shù)據(jù)進行分析，查看是否可以通過循環(huán)變換或者分塊等技術(shù)手段消除依賴性，或者進一步發(fā)現(xiàn)其內(nèi)在的并行性。

2 模擬電荷法基本原理

CSM是基于電磁場的唯一性定理和等效原理，其原理簡明，在具備場的基本概念的基礎(chǔ)上，方法應(yīng)用較為簡捷，并可獲得足夠的計算精度[18]。其基本原理就是通過一組離散化的模擬電荷替代介質(zhì)邊界面上連續(xù)分布的束縛電荷，在已知電位的原邊界面上設(shè)定與模擬電荷數(shù)目相當(dāng)?shù)钠ヅ溆嬎泓c，應(yīng)用疊加原理建立各匹配點的電位方程，求解出各模擬點電荷值，并應(yīng)用疊加原理得到場域中其他任意點的電位和電場值。

圖3 單元線電荷Fig.3 Unit of line charge

圖3 為單元線電荷計算示意圖。將導(dǎo)線分為若干個有限長線段，每一線段用線性模擬電荷等效，設(shè)空間中有一直線段。

計算時，首先分段并對大地做鏡像處理，求出各個匹配點產(chǎn)生的電位系數(shù) P，再根據(jù)已知的邊界條件建立如下矩陣方程：

pτφ= （1）

式中，p為輸電線電位系數(shù)矩陣，τ為待求的線電荷列向量；φ為輸電線匹配點電位列向量。

通過求解該方程，解出電荷的值，再根據(jù)電荷與電場之間的關(guān)系，計算出三個方向上的電場值Ex、Ey、Ez，然后再合成，求得空間中任意點P(x,y,z)產(chǎn)生的電場強度為：

3 并行計算

由于OpenMP比較方便的植入性，使原有的串行程序并行化，因此編程的時候首先編寫輸電線下建筑物周圍電場計算的串行程序，然后再充分挖掘程序中可并行部分，對其進行并行程序編寫，實現(xiàn)并行計算、提升計算的效率。

文中程序在裝有 4核主頻 3.30 GHz 的 Xeon CPU的服務(wù)器（內(nèi)存8 G），Windows 10 64位系統(tǒng)上采用 C和 Matlab混合編寫，C語言編譯器Microsoft Visual Studio 2010版，在充分利用Matlab對于矩陣編寫和計算的便利性的基礎(chǔ)上，用VS2010編寫 C語言并加入 OpenMP增加并行性，再用Matlab調(diào)用C語言所得數(shù)據(jù)畫出仿真圖形。

具體的，在Matlab中根據(jù)模型的坐標(biāo)建立相應(yīng)的矩陣，再通過mexFunction函數(shù)調(diào)用Matlab中矩陣的首地址把數(shù)據(jù)傳入到C語言中，在C語言中進行并行程序的編寫，在 Matlab用 mex語句調(diào)用 C語言中所得的結(jié)果繪出仿真圖形。mexFunction函數(shù)的標(biāo)準(zhǔn)結(jié)構(gòu)如下：

void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])

{；//相應(yīng)的語句}

其中整型參數(shù)nlhs[]和nrhs[]分別為輸出參數(shù)個數(shù)和輸入?yún)?shù)個數(shù)，mxArray類型參數(shù)plhs[]和prhs[]分別為輸入?yún)?shù)的指針和輸出參數(shù)的指針。

4 計算實例及仿真結(jié)果

4.1 計算模型

計算模型選取500 kV超高壓輸電走廊，每相采用4分裂導(dǎo)線，導(dǎo)線間距為0.323 m，子導(dǎo)線半徑為0.0148 m；最低點離地14.19 m，檔距為400 m，相間距為11.8 m。兩個桿塔之間輸電線的懸掛形狀呈懸鏈線形狀，當(dāng)輸電線兩端的懸掛點對地高度一致時，由于對稱性，弧垂的最大值則會出現(xiàn)在檔距的中間點處。懸鏈線方程如下：

式中，σ0為導(dǎo)線水平應(yīng)力，MPa；γ為導(dǎo)線比載，N/m·mm2；l檔距，m；h導(dǎo)線最低點對地高度，m。

由于實際的超高壓輸電線路比較復(fù)雜，而有些因素的影響可以忽略不計，所以本文做如下簡化：1、輸電線表面等電位；2、地面是無窮大的導(dǎo)體，電位為零，沿線地面電阻率相同；3、認為電荷分布沿線路分布沒有畸變，不考慮線路上電位的變化；4、只考慮線路主要部分形成的電場，忽略桿塔、金具、絕緣子等附近物體的影響。

簡化后的模型如圖 5所示，以輸電線走向為 y軸，垂直輸電線方向為x軸，垂直向上為z軸，坐標(biāo)原點取在輸電線中心位置。

圖5 輸電走廊Fig.5 Pow er transmission corridor

4.2 計算精度

在精度方面，通過分別用Matlab和C語言并行編程計算，從量化的角度所得每一個點的電場結(jié)果作差比較，兩值的差值最大也只能達到小數(shù)點后 6位；從可視化的角度，圖6是取垂直于輸電線路并通過建筑物的中心平面，通過疊加作圖可以看到兩個數(shù)值作出來的曲線幾乎一模一樣，兩種角度證明了C語言并行計算下精度的可行性，并不會影響計算結(jié)果。

圖6 C 語言和Matlab的結(jié)果對比Fig.6 Comparative results of C and Matlab

4.3 計算效率

仿真分析了輸電走廊(x,y)=(-600~600，-20~20)，離地1.5 m工頻電場的分布首先本次的模型在Matlab下計算時間為 36342 s，而 C語言單核串行下計算時間僅為186.79 s，僅從串行計算的角度，C語言就是Matlab計算速度的194倍。

并行計算中，加速比S(p)以及并行效率E(p)是評價并行效果的重要指標(biāo)，其中S(p)是P個處理器同時計算求解時，單核串行時間與多核并行時間之間的比值，即

式中，Ts是單個處理器在串行模式下計算所花費的時間；TP是P個處理器在并行模式下計算所花費的時間。

并行效率E(p)是加速比與計算時所用到的CPU核數(shù)的比值，如式（5）所示。當(dāng)加速比與核心數(shù)P接近時，并行效率越高。

描述這兩個指標(biāo)的重要定律：Amdahl定律。假設(shè)計算的問題固定規(guī)模不變，串行部分所占整個程序的比例f也不變，我們可以得到加速比S(p)為：

驗證并行計算的有效性，并且為了避免計算的波動性，對每種情況都進行10次計算。

效果圖如下所示：

圖7 加速比Fig.7 Accelerate rate

圖8 加速效率Fig.8 Accelerate efficiency

雙核模式下，加速比基本穩(wěn)定在2.0，說明雙核情況下并行起到了很好的加速效果；同樣的模型在4核模式下，加速比為3.5以上，相比于4核情況下的理想結(jié)果有一點差距，原因之一是點電荷與線電荷計算的差異性，線電荷要進行較為復(fù)雜的積分運算，所以所花的時間相比點電荷較長，造成線程之間的等待；其次根據(jù) Amdahl可知，串行占比也會對并行的效率造成影響。

5 結(jié)論

本文有別與傳統(tǒng)的只采用Matlab進行編程的模式，引入C語言與Matlab混合編程，并在模擬電荷法的基礎(chǔ)上引入OpenMP進行并行編程，既利用了Matlab生成矩陣的便捷性和對數(shù)據(jù)圖像處理的強大功能，又利用了C語言的對于計算的快速編譯性以及對于并行計算的可塑性，計算了500kV高壓輸電走廊下空間電場分布的情況，實現(xiàn)了對輸電線走廊周圍電場環(huán)境的快速求解，為其他類似電磁場數(shù)值計算提供了思路。

[1] 王艷, 王樹歡, 等. 特高壓六相輸電線路工頻電場強度計算[J]. 中國電力, 2015, 48(3): 44-49.WANG Y, WANG S H, et al. Calculation of Power Frequency Electric Field Intensity of UHV Six-phase Transmission Lines[J]. ELECTRIC POWER, 2015, 48(3): 44-49.

[2] 劉振亞. 中國特高壓交流輸電技術(shù)創(chuàng)新[J]. 電網(wǎng)技術(shù),2013, 37(3): T1-T8.LIU Z Y. Innovation of UHVAC Transmission Technology in China[J]. Power System Technology, 2013, 37(3): T1-T8.

[3] 梅生偉, 朱建全. 智能電網(wǎng)中的若干數(shù)學(xué)與控制科學(xué)問題及其展望[J]. 自動化學(xué)報, 2013, 39(2): 119-131.MEI S W, ZHU J Q. Mathematical and Control Scientiflc Issues of Smart Grid and Its Prospects[J]. ACTA AUTOMATICA SINICA, 2013, 39(2): 119-131.

[4] B. Butrylo, F. Musy, L. Nicolas. A survey of parallel solvers for the finite element method in computational electromagnetics[J]. 2004, COMPEL Vol. 23 Iss 2pp. 531-546.

[5] 趙忠獻, 李鴻健, 豆育生. 并行矩陣乘的光化學(xué)反應(yīng)模擬研究[J]. 軟件, 2010, 31(10): 16-20.ZHAO Z X, LI H J, DOU Y S. Research on parallel matrix multiplication Simulation of Photochemical Reactions[J].Software, 2010, 31(10): 16-20.

[6] 劉靜. 基于獨立生成樹的網(wǎng)絡(luò)多路徑傳輸方法研究[J]. 軟件, 2016, 37(4): 25-28.LIU J. Research on the Multi-path Transmission Method Based on Independent Spanning Tree[J]. Software, 2016,37(4): 25-28.

[7] 劉成軍. 基于消息傳遞接口的線性方程組并行計算研究[J].軟件, 2013, 34(1): 119-120.

[8] 唐萬和, 楊海東, 黎展滔, 等. 考慮能耗成本和拖期成本的非同等并行機調(diào)度[J]. 軟件, 2014, 35(3): 52-57.

[9] 曾晉, 呂韜. 一種基于零系數(shù)檢查的通道并行算法[J]. 軟件, 2014, 35(4): 153-156.

[10] 唐曼. 無線環(huán)境下基于多流優(yōu)先級的并行多路傳輸[J]. 軟件, 2014, 35(10): 46-50.

[11] 武濤. 基于云計算的并行動態(tài)路徑搜索算法研究[J]. 軟件,2015, 36(4): 128-132.WU T. Research on Parallel Dynamic Path Search Algorithm Based on Cloud Computing[J]. Software, 2015, 36(4): 128-132.

[12] T. Oyama, T. Kitahara, Y. Serizawa. Parallel processing for power system analysis using band matrix[J]. IEEE TRANSACTION ON POWER SYSTEMS, 1990 Vol. 5 NO. 3 1010-1016.

[13] 周偉明. 多核計算與程序設(shè)計[M]. 武漢: 華中科技大學(xué)出版社, 2009.ZHOU W M. Multi-core Computing and Programming[M].Wuhan, China, Huazhong University of Science and Technology Press, 2009.

[14] The OpenMP API specification for parallel programming[EB/OL]. http://www.openmp. org/wp-content/uploads/openmp-4.5.pdf.2016.11.

[15] Leonardo D, Ramesh M. OpenMP: An Industry Standard API for Shared-Memory Programming[J]. IEEE Computational Science and Engineering, 1998, 5(1): 46-55.

[16] 鄒競, 馬華. 一種基于OpenMP的并行混合PVS算法[J]. 計算機應(yīng)用研究, 2016, 33(1): 56-59, 91.ZOU J, MA H. Parallel hybrid PVS algorithm based on OpenMP[J]. Application Research of Computers, 2016, 33(1):56-59, 91.

[17] Q uinn M J. Parallel programming in C with MPI and OpenMP[M]. 北京. 清華大學(xué)出版社, 2005.

[18] 倪光正, 錢秀英. 電磁場數(shù)值計算(第2版)[M]. 北京: 機械工業(yè)出版社, 2010. 1.