方曉濤，嚴 正，王 晗，徐瀟源，陳 玥，許少倫

（1. 電力傳輸與功率變換控制教育部重點實驗室（上海交通大學），上海市 200240；2. 香港中文大學機械與自動化工程系，香港 999077）

0 引言

電動汽車的應用能夠有效串聯(lián)配電網(wǎng)清潔電力生產(chǎn)和交通網(wǎng)低碳能源消費［1］，是保障“碳達峰、碳中和”目標實現(xiàn)的重要途徑。隨著電動汽車大規(guī)模推廣應用，交通網(wǎng)與配電網(wǎng)的耦合運行特征日益顯著［2］，并受到交通網(wǎng)與配電網(wǎng)多重不確定性的影響。量化分析不確定性因素的影響，對提升交通網(wǎng)運行的可靠性和配電網(wǎng)運行的安全性具有重要意義。

以電動汽車的充電行為作為交互紐帶，交通網(wǎng)與配電網(wǎng)耦合運行特征體現(xiàn)在：1）交通網(wǎng)的擁堵狀況與充電站電價影響電動汽車用戶的駕駛行為和充電站點的選擇，從而改變配電網(wǎng)內(nèi)充電負荷需求；2）大規(guī)模充電負荷的時空分布特性直接影響配電網(wǎng)潮流，同時改變充電站的充電電價，進而影響交通流分布。對上述交通網(wǎng)車流與配電網(wǎng)潮流相互影響的耦合運行方式進行建模并獲得“物質(zhì)-交通流”與“能量-電力流”的聯(lián)合流結果，是目前國內(nèi)外學者關注的熱點問題。文獻［3-5］分別以網(wǎng)絡均衡模型與最優(yōu)潮流模型刻畫交通流與電力流，得到交通-電力耦合系統(tǒng)的聯(lián)合流結果；文獻［6］計及交互的交通流與電力流，構建了交通網(wǎng)與配電網(wǎng)耦合系統(tǒng)的網(wǎng)絡均衡模型并提出了相應求解方法；文獻［7-8］分別通過變分不等式與博弈論分析了耦合系統(tǒng)運行均衡狀態(tài)的存在性與唯一性。以上述交通網(wǎng)與配電網(wǎng)聯(lián)合流模型為基礎，不少學者也開展了耦合系統(tǒng)的風險評估［9］、聯(lián)合規(guī)劃［10］等相關研究。

然而，實際交通-電力耦合系統(tǒng)內(nèi)存在大量不確定性因素，上述文獻所獲得的交通網(wǎng)與配電網(wǎng)聯(lián)合流結果尚未考慮多重不確定性的影響。對于交通網(wǎng)而言，道路通行可用容量易受到車輛事故、道路施工等隨機事件的影響［11］；電動汽車用戶的出行需求、車輛的充電需求受到天氣、環(huán)境溫度等因素影響［12］；多變的道路通行情況（“路”）、隨機的車輛用戶行為（“車”）成為交通網(wǎng)運行的主要不確定性因素［13］。對于配電網(wǎng)而言，間歇性的分布式可再生能源發(fā)電（“源”）、波動性的電力負荷需求（“荷”）直接影響配電網(wǎng)節(jié)點注入功率，并造成配電網(wǎng)運行狀態(tài)的顯著改變。“路-車-源-荷”多重不確定性將伴隨著交通網(wǎng)與配電網(wǎng)的不斷交互傳播不確定性，并影響耦合系統(tǒng)運行。在此情況下，獲得耦合系統(tǒng)的聯(lián)合流計算結果將面臨如下問題：

1）考慮“路-車-源-荷”多重不確定性，如何構建交通網(wǎng)與配電網(wǎng)的耦合運行模型，并實現(xiàn)不確定性環(huán)境下交通-電力聯(lián)合流的量化表征；

2）如何量化評估“路-車-源-荷”多重不確定性因素對交通流、電力流的影響程度，以反映不確定性傳播下交通網(wǎng)與配電網(wǎng)的交互影響。

針對上述問題，本文提出了考慮“路-車-源-荷”多重不確定性的交通網(wǎng)與配電網(wǎng)概率聯(lián)合流分析方法。首先，分別建立了考慮“路-車”不確定性因素的概率交通分配模型和考慮“源-荷”不確定性因素的概率最優(yōu)潮流模型，并實現(xiàn)不確定性環(huán)境下交通網(wǎng)與配電網(wǎng)的概率聯(lián)合流計算；然后，引入基于Sobol’法的全局靈敏度分析方法，量化多重不確定性因素對交通網(wǎng)與配電網(wǎng)耦合運行的影響，辨識關鍵不確定性因素；最后，通過耦合算例仿真計算驗證了所提方法的有效性。

1 考慮“路-車”不確定性因素的概率交通分配模型

1.1 交通網(wǎng)“路-車”不確定性建模

本文采用概率分布函數(shù)表征交通網(wǎng)“路”“車”不確定性因素的隨機特征。對于“路”而言，隨機變量包括道路容量和道路環(huán)境溫度；對于“車”而言，隨機變量包括車輛出行需求和不同出行階段電池荷電狀態(tài)（SOC）。給定交通網(wǎng)絡GT(N，A)，其中N為節(jié)點集合（本文中交通網(wǎng)節(jié)點包括車流交匯處和車輛充電時充電站所在處），A為道路集合。設定S為充電站節(jié)點集合，且S?N；W為所有汽車的起點-終點（O-D）對的集合；KGw和KEw分別為燃油汽車和電動汽車在O-D 對w∈W之間的路徑集合。對交通網(wǎng)“路”“車”不確定性因素建模如下：

1）道路容量

道路容量指道路正常通行下的最大車流量。在交通網(wǎng)運行中，對于道路a=(i，j)∈A（i，j∈N且i和j分別為道路a的起點和終點），其容量的隨機波動ξC，a可采用正態(tài)分布刻畫［14］，記為ξC，a～N(μC，a，σC，a)，其中，μC，a和σC，a分別為ξC，a的均值與標準差。因此，ξC，a的概率密度函數(shù)可表示為：

2）道路環(huán)境溫度

交通網(wǎng)內(nèi)環(huán)境溫度本身具有較大不確定性，文獻［15］指出環(huán)境溫度的隨機特征可利用正態(tài)分布刻畫。由于道路環(huán)境溫度會直接影響電動汽車電池的各項屬性（充、放電效率等），同時也與用戶是否使用空調(diào)等決策行為相關，因此，本文根據(jù)FleetCarma公司在北美地區(qū)收集的用戶出行記錄數(shù)據(jù)，進一步構建了電動汽車能耗率?與環(huán)境溫度ξT之間的函數(shù)關系［16］：

式中：n為冪指數(shù)；bn為給定系數(shù)，其數(shù)值如附錄A 表A1 所示。

3）汽車出行需求

汽車出行需求指O-D 對間的移動需求，燃油汽車和電動汽車在O-D 對w之間出行需求量的隨機波動ξGN，w和ξEN，w可采用正態(tài)分布刻畫［17］。

4）電動汽車電池初始SOC

由于從起點出發(fā)時電動汽車電池的初始SOC各不相同，若統(tǒng)一采用一個概率分布函數(shù)刻畫不同電動汽車初始SOC 的隨機特征會導致較大誤差。為提高建模的準確性，本文將O-D 對w之間起點處的電動汽車分為Tw個類別，各個類別內(nèi)所有電動汽車的電池初始SOC 采用固定數(shù)值近似替代。各個類別內(nèi)的電動汽車數(shù)量占比記為ξS，w，t，t=1，2，…，Tw，所有類別占比總和為1，即

5）電動汽車充電結束時的SOC

電動汽車用戶在充電站的停留時間具有一定的隨機性，在充電結束時電動汽車并不一定能準確完成預計的充電目標。本文假定電動汽車充電結束時的SOC（記為ξF）近似服從正態(tài)分布［18］，且該分布的均值與預計的充電目標相關。

1.2 交通網(wǎng)概率交通分配模型

交通分配指的是將各個O-D 對間的車輛出行需求量按照一定的準則分配到交通網(wǎng)的各個道路上，并得到交通網(wǎng)的流量分布情況。靜態(tài)交通分配模型是常用的交通流計算模型，其針對某一較長時段內(nèi)的交通需求進行分配，所得到的各個道路交通流量為該時段內(nèi)的平均交通流量。在靜態(tài)交通分配模型的基礎上，本文考慮交通網(wǎng)中的“路-車”不確定性因素，提出概率交通分配模型用于獲得交通流的概率統(tǒng)計特征。其中，本文所提概率交通分配模型中同時考慮了燃油汽車與電動汽車的出行成本，并做如下說明：1）不考慮燃油汽車的油耗因素，燃油汽車的出行成本僅包括出行時間成本；2）電動汽車的出行成本不僅包括道路行駛的時間成本，還包括在充電站的充電費用和充電時間成本。

1）隨機道路容量下道路的通行時間

考慮道路容量的隨機波動，對于道路a∈A，采用含隨機道路容量波動的美國聯(lián)邦公路局（BPR）函數(shù)來描述車輛通行時間ta與車流量xa之間的關系［6］：

式中：t為道路a在車流量為0 時的通行時間；ca為道路a的基準容量。

2）擁堵效應下電動汽車的充電花費時間

電動汽車在充電站花費的時間包括排隊時間和充電時間。對于排隊時間，本文考慮擁堵效應的影響［19］，采用如下函數(shù)描述電動汽車在充電站s的排隊時間ts與等待充電的車流量xs之間的關系：

電動汽車在充電站s的充電時間為：

式中：ts，k，w和Fs，k，w分別為O-D 對w之間路徑k上電動汽車在充電站s的充電時間和充電量；Ps為充電站s的充電功率；η為充電效率。

3）基于用戶均衡準則的概率交通分配模型

基于用戶均衡（UE）準則（即Wardrop 第一準則）［20］，考慮交通網(wǎng)“路-車”不確定性的影響，建立概率交通分配模型如下：

式中：cT為交通網(wǎng)所有汽車用戶的出行成本；ω為出行時間等效成本的轉化系數(shù)；λs為充電站s的充電電價；和分別為燃油汽車和電動汽車在O-D 對w之間路徑k上的車流量；為燃油汽車在O-D 對w之間的出行需求基準量；為電動汽車在O-D 對w之間的出行需求基準量；γs，k，w為電動汽車在O-D對w之間路徑k與充電站節(jié)點s的關聯(lián)系數(shù)，當路徑k上的電動汽車需要在充電站s充電時，γs，k，w為1，否則為0；k，w和k，w分別為燃油汽車和電動汽車在O-D 對w之間路徑k與道路a的關聯(lián)系數(shù)，當路徑k∈包含道路a時，，k，w為1，否則為0，當路徑k∈包含道路a時，k，w為1，否則為0。

上述模型的決策變量為各O-D 對汽車在各個路徑上的車流量，目標函數(shù)由4 項構成：第1項為交通網(wǎng)中所有汽車道路通行時間的等效成本；第2、3 項分別為所有電動汽車在充電站排隊時間與充電時間的等效成本；第4 項為所有電動汽車的總充電費用。約束式（8）至式（11）為燃油汽車與電動汽車在各個O-D 對之間的路徑車流量約束；約束式（12）和式（13）分別表示道路、充電站的車流量與路徑車流量之間的關系。求解上述概率交通分配模型即可得到概率交通流結果。

1.3 交通網(wǎng)車輛最優(yōu)路徑生成

為求解上述概率交通分配模型，需要先獲取所有可能的最優(yōu)路徑集合和。實際電動汽車選擇路徑時，除了考慮出行成本外，還需要顧及電池SOC 能否支撐完成全部及后續(xù)行程。因此，電動汽車的路徑規(guī)劃與充電決策是相互關聯(lián)的，本文將其作為一個整體來考慮，在確定電動汽車最優(yōu)路徑時相應考慮充電決策信息，包括其是否在充電站充電、充電時間與充電量等。

本文通過建立并求解優(yōu)化模型來生成可能的最優(yōu)路徑，其原理是尋找一條當前交通流分布下（道路車流為、充電站車流為），燃油汽車或電動汽車出行成本最小的可行路徑。

1）燃油汽車最優(yōu)路徑生成模型

對于O-D 對w之間的燃油汽車，生成路徑k的優(yōu)化模型目標函數(shù)為：

其中

2）電動汽車最優(yōu)路徑生成模型

對于O-D 對w之間的電動汽車，生成路徑k的優(yōu)化模型目標函數(shù)為：

其中

式中：c為該路徑的總出行成本；t為該路徑的總出行時間，包括道路通行時間、充電站排隊時間與充電時間。

模型的約束條件如下：

在傳統(tǒng)的電大遠程教學過程中，是通過內(nèi)容的傳輸來實現(xiàn)遠程教學，需要學生進行下載操作才能對教學視頻、音頻進行觀看。即通過WEB服務器傳輸?shù)綖g覽者的計算機中，下載后存儲到硬盤中再對文件進行觀看，這些文件包含圖片、視頻、音頻等內(nèi)容，一些容量過大的文件會使觀看過程中出現(xiàn)各種問題，存在聲音與視頻不相符的情況，影響觀看的效果。該做法將導致計算機以及WEB的負荷增大，并且使文件的時效性降低。并且傳統(tǒng)的教學過程中無法實現(xiàn)學生與遠程教師的實時互動，不能對問題進行及時解決，大數(shù)據(jù)技術的應用能夠增強多媒體在教學過程中起到的作用，確保教師以及學生之間交流的實時性，同時通過郵件以及平臺留言等方式也可以解決以上問題。

上述約束條件中，式（19）確保了路徑的可行性；式（20）和式（21）為電動汽車在通過道路、充電站后的電池SOC 平衡約束；式（22）為電動汽車在充電站的充電量與充電時間約束，當路徑與充電站節(jié)點的關聯(lián)系數(shù)γi，k，w為0 時，電動汽車的充電量與充電時間也為0；式（23）為輔助變量約束，當路徑與道路的關聯(lián)系數(shù)，w為0 時，輔助變量近似無約束，當，w為1 時，輔助變量為0；式（24）至式（28）對各個節(jié)點處電動汽車電池SOC 進行限制；式（29）確保電動汽車在該路徑上最多只在一個充電站充電。

2 考慮“源-荷”不確定性因素的概率最優(yōu)潮流模型

2.1 配電網(wǎng)“源-荷”不確定性建模

對于含分布式可再生能源的配電網(wǎng)，間歇性可再生能源發(fā)電（風電、光伏發(fā)電）與波動性負荷是影響其運行的主要不確定性因素。給定配電網(wǎng)GD(EN，EL)，其中EN為節(jié)點集合，EL為支路集合。配電網(wǎng)與上級電網(wǎng)的聯(lián)絡節(jié)點為平衡節(jié)點，EG為電源所在節(jié)點集合（包含平衡節(jié)點），且EG?EN。本文將分布式可再生能源等效為功率注入方向相反的負荷，其所在節(jié)點不包含于電源節(jié)點集合EG中。進一步，采用概率分布函數(shù)，對配電網(wǎng)“源”“荷”不確定性因素建模如下：

1）可再生能源不確定性建模

風電機組、光伏電池的輸出功率分別與風速、光照強度密切相關，在長時間尺度下，風速近似服從Weibull 分布［22］：

式中：v為風速；ξk和ξc分別為Weibull 分布的形狀參數(shù)和尺度參數(shù)。

光照強度近似服從Beta 分布［23］：

式中：r和rmax分別為光照強度和最大光照強度；ξa和ξb為Beta 分布的形狀參數(shù)；Γ(·)為伽瑪函數(shù)。

基于風速與光照強度數(shù)據(jù)，依據(jù)文獻［24］的函數(shù)關系能夠得到對應的風電機組與光伏電池有功功率。進一步，通過節(jié)點i∈EN處可再生能源的功率因數(shù)和有功功率PR，i可得到其無功功率QR，i。

2）負荷不確定性建模

配電網(wǎng)中除充電站充電負荷外其他負荷的有功功率PB，i（i∈EN）可認為近似服從正態(tài)分布［24］，記為PB，i～N(μB，i，σB，i)，其中，μB，i和σB，i分別為PB，i的均值和標準差。PB，i的概率密度函數(shù)可表示為：

通過節(jié)點i處除充電站充電負荷外其他負荷的功率因數(shù)和有功功率PB，i可得到其無功功率QB，i。

2.2 配電網(wǎng)概率最優(yōu)潮流模型

相比于輸電網(wǎng)，配電網(wǎng)具有特殊的網(wǎng)絡特性，如線路電阻與電抗的比值較大、輻射性結構等，一般采用交流最優(yōu)潮流算法獲得配電網(wǎng)最優(yōu)潮流。在此基礎上，本文考慮不確定性“源-荷”的影響，建立了配電網(wǎng)概率最優(yōu)潮流模型，以獲得配電網(wǎng)最優(yōu)潮流的概率統(tǒng)計特征：

式中：cD為配電網(wǎng)的運行成本；αi和βi為分布式電源發(fā)電的成本系數(shù)；β0為配電網(wǎng)從上級電網(wǎng)購電的成本系數(shù)；PG，i和QG，i分別為節(jié)點i電源或平衡節(jié)點的有功和無功注入功率；PC，i為節(jié)點i充電站充電負荷有功功率；PL，i和QL，i分別為線路i的有功和無功傳輸功率；Vi為節(jié)點i的電壓幅值；θij為節(jié)點i和節(jié)點j電壓的相角差；Gij和Bij分別為配電網(wǎng)節(jié)點導納矩陣第i行、第j列元素的實部和虛部；G，i、G，i和G，i、G，i分別為節(jié)點i電源或平衡節(jié)點有功和無功注入功率的下限、上限；和分別為節(jié)點i電壓幅值的下限和上限；為線路i的最大傳輸功率。

上述模型的目標函數(shù)為配電網(wǎng)的運行成本。式（34）和式（35）為節(jié)點功率平衡約束；式（36）至式（39）分別為電源與平衡節(jié)點的有功和無功注入約束、節(jié)點電壓幅值約束及線路傳輸功率約束。

本文采用交通網(wǎng)運行時段內(nèi)充電負荷的平均功率作為配電網(wǎng)充電負荷的有功功率，且暫不考慮充電負荷的無功功率。其中，充電站s的充電負荷量Fs可在求解概率交通分配模型后求得：

本文以充電站所在節(jié)點邊際電價作為充電站充電電價［7］，其可通過求解交流最優(yōu)潮流模型中的拉格朗日對偶乘子得到。

3 考慮多重不確定性因素的交通網(wǎng)與配電網(wǎng)概率聯(lián)合流計算

3.1 交通網(wǎng)與配電網(wǎng)耦合運行均衡狀態(tài)分析

交通網(wǎng)與配電網(wǎng)協(xié)同運行中，基于充電站的充電負荷與充電電價耦合互動。在交通網(wǎng)與配電網(wǎng)的相互作用下，耦合系統(tǒng)會達到一個均衡狀態(tài)。該均衡狀態(tài)解的存在性與唯一性的證明并非本文關注的重點，具體可參見文獻［6-8］。

本文采用分散迭代算法，通過迭代求解概率交通分配模型與配電網(wǎng)概率最優(yōu)潮流模型，得到不確定性場景下耦合系統(tǒng)的均衡狀態(tài)解，迭代求解過程如附錄A 圖A1 所示。給定隨機場景(ξ，PD)，其中，ξ為交通網(wǎng)“路-車”不確定性因素，包括道路容量ξC、道路環(huán)境溫度ξT、汽車出行需求ξN、電池初始SOC分類后的各類電動汽車數(shù)量占比ξS、電動汽車充電結束時的SOC 變量ξF；PD為配電網(wǎng)“源-荷”不確定性因素，包括可再生能源有功功率PR與除充電站充電負荷外其他負荷的有功功率PB。配電網(wǎng)側充電站節(jié)點的電價λ會影響交通網(wǎng)汽車用戶的出行成本cT，從而影響用戶的充電決策與路徑規(guī)劃，并體現(xiàn)為交通流的變化；交通流的變化會通過映射η(fE)（即式（40））改變交通網(wǎng)充電站的充電負荷F，進而影響配電網(wǎng)潮流的分布，再經(jīng)由映射π(PG)（即最優(yōu)潮流模型的拉格朗日乘子提?。┓答仦槌潆婋妰rλ的調(diào)整，反過來影響交通網(wǎng)流量分布；通過不斷交互迭代獲得交通網(wǎng)與配電網(wǎng)耦合系統(tǒng)的均衡狀態(tài)解，即隨機場景(ξ，PD)下交通-電力聯(lián)合流結果。

3.2 交通網(wǎng)與配電網(wǎng)概率聯(lián)合流計算流程

為表征“路-車-源-荷”多重不確定性因素影響下交通網(wǎng)與配電網(wǎng)的聯(lián)合流特征，建立交通網(wǎng)與配電網(wǎng)概率聯(lián)合流計算模型：

式中：XT為求解交通分配模型所需的交通網(wǎng)數(shù)據(jù)；XD為求解配電網(wǎng)最優(yōu)潮流模型所需的配電網(wǎng)數(shù)據(jù)；Y為交通網(wǎng)與配電網(wǎng)均衡狀態(tài)下的耦合系統(tǒng)狀態(tài)變量（如交通網(wǎng)的道路車流、出行成本，配電網(wǎng)的節(jié)點電壓、線路功率等）；f(·)為求解交通網(wǎng)與配電網(wǎng)耦合運行均衡狀態(tài)的映射關系（如3.1 節(jié)所示）。

上述模型中，系統(tǒng)的輸入變量ξ、PD具有隨機不確定性，考慮不確定性傳播過程，系統(tǒng)的輸出變量即交通網(wǎng)與配電網(wǎng)耦合系統(tǒng)狀態(tài)變量Y也具有概率特性。本文通過蒙特卡洛模擬（MCS）方法，計算所有不確定性場景下的交通網(wǎng)與配電網(wǎng)聯(lián)合流結果，并以所得交通流和電力流的概率分布量化表征不確定性因素對交通網(wǎng)與配電網(wǎng)耦合系統(tǒng)運行狀態(tài)的影響。交通網(wǎng)與配電網(wǎng)概率聯(lián)合流計算流程如下：

步驟1：獲取交通網(wǎng)、配電網(wǎng)的基本信息，并選取交通網(wǎng)與配電網(wǎng)中的“路-車-源-荷”多重不確定性因素組成系統(tǒng)M維獨立輸入隨機變量(ξ，PD)；

步驟2：基于統(tǒng)計數(shù)據(jù)，采用本文所提概率分布函數(shù)刻畫交通網(wǎng)與配電網(wǎng)中的“路-車-源-荷”多重不確定性因素，建立輸入變量(ξ，PD)的隨機模型；

步驟3：基于所建立的輸入隨機變量模型，通過Sobol’序列采樣或拉丁超立方采樣等低差異序列采樣方法［25］，得到J組輸入隨機變量樣本(ξj，PD，j)（j=1，2，…，J）；

步驟4：基于大規(guī)模的輸入隨機變量樣本，通過所提交通網(wǎng)與配電網(wǎng)耦合運行均衡狀態(tài)求解方法，獲得每組樣本(ξj，PD，j)對應的交通網(wǎng)與配電網(wǎng)聯(lián)合流輸出變量樣本Yj；

步驟5：計算輸出變量樣本Yj的統(tǒng)計特征量，如均值、方差、概率密度函數(shù)（PDF）、累積分布函數(shù)（CDF）等，作為交通網(wǎng)與配電網(wǎng)概率聯(lián)合流計算結果。

4 基于全局靈敏度分析的交通網(wǎng)與配電網(wǎng)多重不確定性因素影響量化

交通網(wǎng)與配電網(wǎng)概率聯(lián)合流計算結果從宏觀上表征了“路-車-源-荷”多重不確定性因素對耦合系統(tǒng)運行狀態(tài)的影響，然而，各個不確定性因素對交通流、電力流的影響程度并不相同，如何量化“路”“車”“源”“荷”各自對交通流、電力流的影響程度，以反映不確定性傳播下耦合系統(tǒng)的交互影響是需要進一步解決的問題。靈敏度分析方法是量化輸入變量影響程度的有效手段［26］，但交通網(wǎng)與配電網(wǎng)耦合系統(tǒng)內(nèi)的不確定性傳播過程并非瞬時完成，且概率聯(lián)合流結果體現(xiàn)的是長時間尺度運行下不確定性因素的影響，這使得傳統(tǒng)以偏導數(shù)為基礎的局部靈敏度方法不能適用。因此，本文引入全局靈敏度分析方法，充分考慮長時間尺度下輸入隨機變量的概率分布特征，以實現(xiàn)交通網(wǎng)與配電網(wǎng)內(nèi)“路-車-源-荷”多重不確定性因素的影響量化。

4.1 基于Sobol’法的全局靈敏度分析

全局靈敏度分析是一種量化輸入隨機變量影響的靈敏度分析手段，目前已用于含隨機源-荷的電力系統(tǒng)潮流分析［27］、小干擾穩(wěn)定分析［28］等。本文采用Sobol’法計算交通網(wǎng)與配電網(wǎng)多重不確定性因素的全局靈敏度指標，通過指標大小獲得不確定性因素重要性排序，從而實現(xiàn)交通網(wǎng)與配電網(wǎng)耦合系統(tǒng)內(nèi)不確定性因素的影響量化。

Sobol’法基于方差分析（ANOVA）理論，采用方差表征系統(tǒng)輸入、輸出變量的不確定性，并將系統(tǒng)輸出變量的方差歸因于各個輸入變量及輸入變量間的交互作用［29］。假設系統(tǒng)模型為Y=f(z)，其中，z=[z1，z2，…，zM]為M維獨立輸入隨機變量?；赟obol 分解理論，可對系統(tǒng)模型做如下分解：

式 中 ：f0為 常 數(shù) ；fi(zi)、fi，j(zi，zj)、f1，2，…，M(z1，z2，…，zM)為各分解子項函數(shù)。

若f(z)平方可積，且展開式中各子項均值為0，則各子項間兩兩相互正交，且各子項可唯一確定。將式（42）兩邊同時取方差可得：

其中

式中：E(Y|zi)為給定zi時Y的條件期望；E(Y|zi，zj)為給定zi與zj時Y的條件期望；D為輸出變量的方差；Di為一階偏方差，表示單一輸入變量zi對輸出變量方差的貢獻程度；Di，j為二階偏方差，表示輸入變量zi與zj之間交互作用對輸出變量方差的貢獻程度。以此類推，其他各階偏方差具有類似含義。輸入變量zi的總偏方差DTi表示變量zi單獨作用以及zi與其他輸入變量交互作用對系統(tǒng)輸出變量方差的共同貢獻程度，如式（45）所示。

式中：Di1，i2，…，is為輸入變量zi1，zi2，…，zis的s階偏方差；Ωi={(i1，i2，…，is)：?k，1 ≤k≤s，ik=i}。

定義標準化的全局靈敏度指標為：

式中：Si（1 ≤i≤M）為輸入變量zi的一階靈敏度指標（FSI），反映了單一輸入隨機變量zi對輸出變量不確定性（方差）的影響；（1 ≤i≤M）為輸入變量zi的總靈敏度指標（TSI），反映了輸入變量zi單獨作用以及zi與其他輸入變量交互作用對系統(tǒng)輸出變量不確定性（方差）的共同影響。

FSI 與TSI 的計算通?；贛CS 方法［26］，然而，大量反復的樣本生成與模型求解使得該方法耗時嚴重。為了提高交通網(wǎng)與配電網(wǎng)多重不確定性因素影響量化過程的計算效率，本文引入稀疏多項式混沌展開（SPCE）代理模型以實現(xiàn)全局靈敏度指標快速計算。

將系統(tǒng)模型Y=f(z)表示為正交多項式的稀疏混沌展開如下：

SPCE 代理模型的構建流程可參見文獻［30］。當SPCE 代理模型確定后，利用展開式的各項系數(shù)即可計算得到系統(tǒng)輸入隨機變量的全局靈敏度指標。對于系統(tǒng)輸入隨機變量ξu（u=1，2，…，M），其FSI 與TSI 計算公式如下［31］：

式中：Au={I∈AS：iu＞0；il=0，?l∈{1，2，…，M}u}；={I∈AS：iu＞0}。

SPCE 代理模型的構建僅需少量系統(tǒng)輸入-輸出變量樣本，因此，應用SPCE 代理模型計算全局靈敏度指標能夠顯著減少計算耗時，提高計算效率。

4.2 “路-車-源-荷”多重不確定性因素的影響量化

結合Sobol’法與SPCE 代理模型，計算交通網(wǎng)與配電網(wǎng)“路-車-源-荷”多重不確定性因素的全局靈敏度指標，以實現(xiàn)不確定性因素的影響量化。需要說明的是，本文側重于量化耦合系統(tǒng)中不確定性因素對耦合系統(tǒng)狀態(tài)變量的單獨影響，因此，后續(xù)計算中僅采用了FSI 作為全局靈敏度指標對輸入隨機變量進行重要性排序，具體計算流程如下：

步驟1：基于MCS 方法獲得少量交通網(wǎng)與配電網(wǎng)耦合系統(tǒng)輸入-輸出樣本，用于SPCE 代理模型構建。該過程即為交通網(wǎng)與配電網(wǎng)概率聯(lián)合流計算過程，具體步驟見3.2 節(jié)。

步驟2：基于Sobol’方法與代理模型，實現(xiàn)全局靈敏度指標的快速計算，具體計算如下。

1）設定步驟1 所生成的輸入-輸出變量樣本規(guī)模為JS，通過累積分布函數(shù)變換建立輸入變量樣本(ξj，PD，j)（j=1，2，…，JS）與標準正態(tài)分布變量樣本zj（j=1，2，…，JS）之間的轉換關系，標準正態(tài)分布變量z與原輸入隨機變量(ξ，PD)的全局靈敏度指標分別對應。

2）選取SPCE 代理模型的基函數(shù)類型為Hermite 正交多項式，將標準正態(tài)分布變量樣本zj作為代理模型輸入變量樣本，輸出變量樣本為Yj（j=1，2，…，JS），構建SPCE 代理模型。

3）利用SPCE 代理模型表達式，基于式（48）計算系統(tǒng)輸入隨機變量的FSI。

4）通過比較輸入隨機變量的全局靈敏度指標數(shù)值大小，獲得交通網(wǎng)與配電網(wǎng)耦合系統(tǒng)不確定性因素的重要性排序，辨識對耦合系統(tǒng)具有顯著影響的關鍵不確定性因素。

上述計算過程的流程圖如附錄A 圖A2 所示。

5 算例仿真

5.1 算例系統(tǒng)

本文用于仿真測試的耦合系統(tǒng)A 由16 節(jié)點環(huán)形交通網(wǎng)與改進的IEEE 33 節(jié)點配電網(wǎng)構成；耦合系統(tǒng)B 由美國Sioux Falls 交通網(wǎng)與改進的IEEE 69節(jié)點配電網(wǎng)構成。耦合系統(tǒng)A、B 的拓撲圖分別如附錄B 圖B1、圖B2 所示。耦合系統(tǒng)A、B 的參數(shù)設置分別見附錄C、附錄D?；谖墨I［24］，本文假定風電機組、光伏電池的輸出功率與風速、光照強度之間的關系如附錄E 所示。本文仿真環(huán)境為Intel Core i7-10 700 8 核CPU、16 GB 內(nèi)存，使用MATLAB R2020b 編譯和測試。

5.2 交通網(wǎng)與配電網(wǎng)概率聯(lián)合流計算

1）耦合系統(tǒng)A

考慮交通網(wǎng)與配電網(wǎng)中的“路-車-源-荷”多重不確定性因素，本文采用MCS 方法計算交通網(wǎng)與配電網(wǎng)的概率聯(lián)合流，采樣規(guī)模J=10 000。耦合系統(tǒng)A中所設置的系統(tǒng)輸入隨機變量為：5 個分布式可再生能源的輸出功率PR，i（i=12，22，18，33，25）；4 個區(qū)域的道路容量隨機波動ξC，i（i=1，2，3，4，分別對應工業(yè)區(qū)、休閑區(qū)、居民區(qū)、商業(yè)區(qū)），區(qū)域i內(nèi)任意一條道路a的容量隨機波動ξC，a均為ξC，i；起點分別為節(jié)點1、3、4 的O-D 對汽車出行需求量隨機波動ξN，i（i=1，3，4），表示起點為i的任意O-D 對w間燃油汽車和電動汽車的出行需求量隨機波動與均為ξN，i；起點分別為節(jié)點1、3、4 的電動汽車電池初始SOC 為0.3、0.6 的數(shù)量占比ξS，i，t（i=1，3，4；t= 1，2）（本算例將O-D 對起點處的電動汽車分為3 個類別，3 個類別內(nèi)電動汽車的電池初始SOC 分別采用0.3、0.6、0.9 代替，起點為i的任意OD 對w間第t個類別的電動汽車數(shù)量占比ξS，w，t均為ξS，i，t，第3 類的電動汽車占比可由總比例1 減去另外2 類的占比之和得到）；電動汽車充電結束時的電池SOC 變量ξF；道路環(huán)境溫度ξT，總的輸入變量維度M=20。在計算交通網(wǎng)與配電網(wǎng)耦合運行均衡狀態(tài)時，靜態(tài)交通分配模型通過調(diào)用IPOPT 求解器求解，配電網(wǎng)交流最優(yōu)潮流模型則借助Matpower 工具包求解。

根據(jù)所得耦合系統(tǒng)狀態(tài)變量的統(tǒng)計結果，分別選取33 節(jié)點配電網(wǎng)各節(jié)點的電壓幅值與16 節(jié)點交通網(wǎng)各道路的車流量作為概率聯(lián)合流分析對象，繪制箱形圖如圖1（a）和（b）所示。由圖1（a）可知，配電網(wǎng)各節(jié)點電壓幅值分布范圍均不一致，其中，節(jié)點17、18 的電壓幅值分布最分散，其波動幅度較大且存在電壓幅值較高的極端場景，相對其他節(jié)點具有較大的電壓越限風險；節(jié)點24、25 的電壓幅值較大概率位于0.95 p.u.以下，是配電網(wǎng)中電壓薄弱的環(huán)節(jié)，也存在一定的安全風險。由圖1（b）可知，交通網(wǎng)各道路車流量的分布范圍差異明顯，道路1→13、3→14 的平均車流量及其波動幅度顯著高于其他道路，極端場景下其車流量最大分別可達2 248 輛/h與2 276 輛/h，在道路可用容量意外減小時（如交通事故、道路維護等）更容易發(fā)生瞬時車流量過大而導致的交通堵塞現(xiàn)象。因此，在道路擴建規(guī)劃時，也需優(yōu)先考慮上述道路。

圖1 交通網(wǎng)與配電網(wǎng)耦合系統(tǒng)A 狀態(tài)變量的箱形圖Fig.1 Box plots of state variables in coupled system A of transportation network and distribution network

對于33 節(jié)點配電網(wǎng)，分別計算節(jié)點8 電壓幅值（V8）、充電站4 充電電價（λ4）的PDF 如圖2（a）、（b）所示；對于16 節(jié)點交通網(wǎng)，分別計算交通網(wǎng)中道路8 →9 車流量（x8→9）、充電站4 充電負荷（F4）的PDF如圖2（c）、（d）所示。對比圖2 中的結果可知，充電站4 節(jié)點電價的PDF 具有雙峰分布的特點，其在交通網(wǎng)與配電網(wǎng)耦合交互過程中，容易受到多重不確定性因素的影響；由圖2（c）可知，道路8→9 的車流量存在超過道路預警容量的風險（本文以道路容量基準值加上隨機波動標準差的3 倍設定為道路預警容量），基于道路8→9 車流量的PDF 能夠計算得到其越限概率（即車流量超過道路預警容量的概率）為28.22%。

圖2 交通網(wǎng)與配電網(wǎng)耦合系統(tǒng)A 狀態(tài)變量的PDFFig.2 PDF of state variables in coupled system A of transportation network and distribution network

2）耦合系統(tǒng)B

對于耦合系統(tǒng)B，分別選取69 節(jié)點配電網(wǎng)中節(jié)點56 的電壓幅值（V56）與Sioux Falls 交通網(wǎng)中道路12 的車流量（x12）作為分析對象，計算其PDF 如附錄F 圖F1 所示。根據(jù)耦合系統(tǒng)B 的概率聯(lián)合流計算結果，得到Sioux Falls 交通網(wǎng)道路車流量均值與越限概率熱力圖如圖3（a）、（b）所示。由圖3（a）可知，僅有少量道路的車流量均值顯著大于其他道路；由圖3（b）可知，部分道路的車流量越限概率較高，容易發(fā)生交通堵塞現(xiàn)象。結合圖3（a）、（b）可知，部分道路由于容量偏小，在車流量不大的情況下也存在較高的交通堵塞風險。為了改善該交通網(wǎng)的通行能力，需要對上述道路進行合理擴建或增強與電網(wǎng)的聯(lián)合調(diào)控?；谏鲜鼋Y果可知，對于較大規(guī)模的交通網(wǎng)與配電網(wǎng)耦合系統(tǒng)，本文所提概率聯(lián)合流計算方法仍具有適用性。

圖3 Sioux Falls 交通網(wǎng)道路車流量均值與越限概率熱力圖Fig.3 Heat maps of traffic flow mean value and limit violation probability of road in Sioux Falls transportation network

綜合耦合系統(tǒng)A、B 的概率聯(lián)合流計算結果可知，由于交通網(wǎng)與配電網(wǎng)內(nèi)“路-車-源-荷”多重不確定性因素的影響，通過不確定性傳播過程，系統(tǒng)的輸出變量（即交通網(wǎng)與配電網(wǎng)耦合系統(tǒng)狀態(tài)變量）也具有概率特性。同時，通過計算交通網(wǎng)與配電網(wǎng)的概率聯(lián)合流結果，系統(tǒng)運行者能夠獲得系統(tǒng)各狀態(tài)變量的概率統(tǒng)計信息，從而為交通網(wǎng)、配電網(wǎng)的聯(lián)合運行分析、運行狀態(tài)感知與預警提供更全面的信息參考。

5.3 交通網(wǎng)與配電網(wǎng)多重不確定性因素影響量化

在交通網(wǎng)與配電網(wǎng)概率聯(lián)合流計算的基礎上，利用基于Sobol’法的全局靈敏度分析方法，量化交通網(wǎng)與配電網(wǎng)中“路-車-源-荷”多重不確定性因素的影響，辨識顯著影響交通網(wǎng)與配電網(wǎng)耦合運行的關鍵不確定性因素。本文利用UQLab［32］構建交通網(wǎng)與配電網(wǎng)聯(lián)合流的SPCE 代理模型，實現(xiàn)全局靈敏度指標的快速計算。

1）耦合系統(tǒng)A

首先，為了驗證基于SPCE 代理模型的全局靈敏度指標計算方法的有效性，將所提方法與傳統(tǒng)的MCS 方法進行對比分析?？紤]耦合系統(tǒng)A 中的20 個系統(tǒng)輸入隨機變量，以V8作為輸出變量，分別計算各個系統(tǒng)輸入隨機變量的FSI。SPCE 代理模型的最高階數(shù)為10，模型構建所用的樣本規(guī)模JS=600；MCS 方法的采樣規(guī)模為10 000。任意選取5 個系統(tǒng)輸入隨機變量：WT1 的輸出功率（PR，12）、區(qū)域1 道路的容量隨機波動（ξC，1）、起點為節(jié)點1 的O-D 對汽車出行需求量隨機波動（ξN，1）、起點為節(jié)點1 的電動汽車電池初始SOC 為0.3 的數(shù)量占比（ξS，1，1）、電動汽車電池充電結束時的SOC（ξF）。通過SPCE 代理模型與MCS 方法分別計算得到FSI結果與計算用時如表1 所示。

表1 SPCE 代理模型與MCS 方法的FSI 計算結果對比Table 1 Comparison of FSI calculation results between SPCE agent model and MCS method

由表1 可知，通過SPCE 代理模型計算得到的FSI 結果與MCS 方法基本一致，其平均數(shù)值誤差為0.004 6。同時，通過SPCE 代理模型計算FSI 能夠大幅減少計算用時（由MCS 方法的191.2 h 減少至2.3 h）。因此，基于SPCE 代理模型的全局靈敏度指標計算方法能夠在保證計算精度的前提下，顯著提高計算效率。

對于33 節(jié)點配電網(wǎng)，以節(jié)點8 電壓幅值（V8）為輸出變量，通過基于SPCE 代理模型的全局靈敏度指標計算方法，得到20 個系統(tǒng)輸入隨機變量的FSI如附錄F 圖F2（a）所示；對于16 節(jié)點交通網(wǎng)，以道路8→9 車流量（x8→9）為輸出變量，計算得到20 個系統(tǒng)輸入隨機變量的FSI 如圖F2（b）所示。

由附錄F 圖F2 可知，僅有少量系統(tǒng)輸入隨機變量的FSI 數(shù)值結果較大，說明并非所有系統(tǒng)輸入變量的隨機性對系統(tǒng)輸出都有顯著影響。以FSI 計算數(shù)值結果較大的系統(tǒng)輸入隨機變量作為關鍵輸入變量，則對于輸出變量V8影響最顯著的關鍵輸入變量為WT1、WT2 的輸出功率（PR，12、PR，22）及起點為節(jié)點1 的電動汽車電池初始SOC 為0.3 的數(shù)量占比（ξS，1，1）；對于輸出變量x8→9影響最顯著的關鍵輸入變量為區(qū)域1 道路的容量波動（ξC，1）及起點為節(jié)點3、4 的O-D 對汽車出行需求量隨機波動（ξN，3、ξN，4）?？梢钥吹?，在本文所設置的參數(shù)條件下，交通網(wǎng)中的“路-車”不確定性因素對配電網(wǎng)狀態(tài)變量的影響相較于配電網(wǎng)中的“源-荷”不確定性因素對交通網(wǎng)狀態(tài)變量的影響更為顯著。

進一步，為了驗證所提交通網(wǎng)與配電網(wǎng)多重不確定性因素影響量化方法在辨識關鍵輸入隨機變量方面的有效性，本文以V8為輸出變量設置如下3 種場景：

場景1：考慮所有20 個輸入變量的隨機性；

場景2：僅考慮3 個關鍵輸入變量的隨機性，其余非關鍵輸入變量取固定值；

場景3：僅考慮其余非關鍵輸入變量的隨機性，3 個關鍵輸入變量取固定值。

基于交通網(wǎng)與配電網(wǎng)概率聯(lián)合流計算，可分別得到3 個場景下V8的CDF 如圖4 所示。

圖4 不同場景下V8的CDFFig.4 CDF of V8 in different cases

由圖4 可知，僅考慮3 個關鍵輸入變量的隨機性時（場景2），輸出變量V8的CDF 與考慮所有20 個輸入變量隨機性時（場景1）基本一致；反之，當3 個關鍵輸入變量取固定值，僅考慮其余輸入變量隨機性時（場景3），輸出變量V8的CDF 與場景1 有很大差異。算例結果表明，全局靈敏度指標所辨識的關鍵輸入變量對輸出變量具有顯著影響，關鍵輸入變量的隨機性是輸出變量隨機性的主要來源，進一步，驗證了本文所提交通網(wǎng)與配電網(wǎng)多重不確定性因素影響量化方法在辨識關鍵不確定性因素方面的有效性。

2）耦合系統(tǒng)B

對于耦合系統(tǒng)B，分別選取69 節(jié)點配電網(wǎng)中節(jié)點56 的電壓幅值（V56）與Sioux Falls 交通網(wǎng)中道路12 的車流量（x12）作為輸出變量，計算系統(tǒng)隨機輸入變量的FSI。取FSI 指標排序前5 的系統(tǒng)輸入隨機變量作為關鍵隨機變量，計算結果如附錄F 表F1 所示。由表F1 可知，對于輸出變量V56影響最顯著的關鍵輸入變量為位于節(jié)點12、51、52、58、59 的可再生能源輸出功率，對于輸出變量x12影響最顯著的關鍵輸入變量為工業(yè)區(qū)2 道路的容量波動及起點為節(jié)點1、2、4、9 的O-D 對汽車出行需求量隨機波動。

6 結語

本文考慮“路-車-源-荷”多重不確定性因素，建立了交通網(wǎng)與配電網(wǎng)概率聯(lián)合流計算模型，并引入了基于SPCE 代理模型的全局靈敏度分析方法，量化多重不確定性因素對交通網(wǎng)與配電網(wǎng)耦合運行的影響。主要結論如下：

1）考慮“路-車-源-荷”多重不確定性因素的概率聯(lián)合流計算模型能夠量化表征交通網(wǎng)與配電網(wǎng)耦合運行的聯(lián)合流概率特性，為交通網(wǎng)、配電網(wǎng)的聯(lián)合運行分析、運行狀態(tài)感知與預警提供更全面的信息參考。

2）基于SPCE 代理模型的全局靈敏度分析方法能夠準確、高效地量化“路-車-源-荷”多重不確定性因素對交通網(wǎng)與配電網(wǎng)耦合運行的影響程度，辨識關鍵不確定性因素，為耦合系統(tǒng)不確定性分析中隨機變量的降維提供方法支撐。

下一步工作是研究汽車用戶成本認知水平差異、用戶行為的隨機-模糊特征對交通網(wǎng)與配電網(wǎng)概率聯(lián)合流的影響。

附錄見本刊網(wǎng)絡版（http：//www.aeps-info.com/aeps/ch/index.aspx），掃英文摘要后二維碼可以閱讀網(wǎng)絡全文。