厲文秀

（河海大學(xué)能源與電氣學(xué)院，江蘇南京 210098）

在自然界，系統(tǒng)狀態(tài)的未來發(fā)展趨勢往往既取決于當(dāng)前運(yùn)行狀態(tài)，也與過去的狀態(tài)密切相關(guān)，這類現(xiàn)象稱為時(shí)滯現(xiàn)象。時(shí)滯現(xiàn)象在電力系統(tǒng)的控制回路中非常普遍，在過去分析電力系統(tǒng)時(shí)，往往忽略時(shí)滯的影響，這主要是因?yàn)榭刂破饔杀镜亓繕?gòu)成，此時(shí)的通信延時(shí)在10 ms以下[1-4]。但是，隨著互聯(lián)電力系統(tǒng)規(guī)模和復(fù)雜度的不斷增加，大量功率需要遠(yuǎn)距離傳輸，區(qū)域間功率振蕩的可能性和危害性日益增加，僅靠局部反饋信號(hào)設(shè)計(jì)的控制器越來越難保證互聯(lián)電力系統(tǒng)的穩(wěn)定性，必須充分利用廣域量測系統(tǒng)所能提供的遠(yuǎn)方設(shè)備信息，進(jìn)行系統(tǒng)協(xié)同控制，方能確保系統(tǒng)的安全穩(wěn)定運(yùn)行。而廣域量測信息存在明顯的延時(shí)，不能完全忽略，因此研究時(shí)滯環(huán)節(jié)對系統(tǒng)穩(wěn)定性的影響，具有十分重要的現(xiàn)實(shí)意義[5-8]。

現(xiàn)實(shí)生活中電力系統(tǒng)總是存在著各種各樣隨機(jī)干擾，這些干擾有的是由擾動(dòng)引起的系統(tǒng)內(nèi)部結(jié)構(gòu)、參數(shù)變化，有的是由可再生能源、電動(dòng)汽車等接入引起的外部隨機(jī)激勵(lì)。當(dāng)前隨著可再生能源的大規(guī)模并網(wǎng)，其隨機(jī)性對電力系統(tǒng)的影響不容忽略，有必要研究隨機(jī)因素作用下電力系統(tǒng)的穩(wěn)定性問題[9-10]。

隨機(jī)電力系統(tǒng)中加入時(shí)滯，其研究的難度大大增加。近年來隨著隨機(jī)H∞控制問題的解決，隨機(jī)時(shí)滯系統(tǒng)的穩(wěn)定性和魯棒控制問題的研究取得了一定進(jìn)展，出現(xiàn)了一些新的研究成果。文獻(xiàn)[11]研究了一類帶Markovian切換的非線性時(shí)滯隨機(jī)系統(tǒng)的穩(wěn)定性問題；文獻(xiàn)[12-13]通過使用模型轉(zhuǎn)換和線性矩陣不等式的方法，研究了一類時(shí)滯不確定隨機(jī)系統(tǒng)的平均時(shí)間延遲和指數(shù)穩(wěn)定性的魯棒控制；文獻(xiàn)[14-15]研究了不確定隨機(jī)系統(tǒng)的魯棒穩(wěn)定性條件。

目前在電力系統(tǒng)領(lǐng)域，隨機(jī)時(shí)滯系統(tǒng)的穩(wěn)定性研究還未見報(bào)道。本文以勵(lì)磁系統(tǒng)輸入信號(hào)延時(shí)作為時(shí)滯環(huán)節(jié)，并通過引入負(fù)荷擾動(dòng)項(xiàng)，建立單機(jī)無窮大隨機(jī)時(shí)滯系統(tǒng)的數(shù)學(xué)模型；仿真獲得隨機(jī)時(shí)滯電力系統(tǒng)的受擾軌跡，分析不同時(shí)滯、不同擾動(dòng)強(qiáng)度下隨機(jī)時(shí)滯電力系統(tǒng)的穩(wěn)定性；負(fù)荷模型分別采用恒阻抗、恒功率、恒電流模型，分析比較采用不同負(fù)荷模型對系統(tǒng)穩(wěn)定性的影響。仿真表明，3種負(fù)荷模型所能承受的時(shí)滯上限能力從大到小依次為：恒阻抗模型、恒電流模型、恒功率模型。

1 隨機(jī)時(shí)滯電力系統(tǒng)模型

1.1 負(fù)荷模型

用于描述負(fù)荷特性的數(shù)學(xué)模型稱之為負(fù)荷模型[16]。如果所采用的負(fù)荷模型不夠準(zhǔn)確，在系統(tǒng)仿真過程中，會(huì)影響電力系統(tǒng)的動(dòng)態(tài)響應(yīng)以及系統(tǒng)安全穩(wěn)定問題，甚至仿真得到的結(jié)果與系統(tǒng)的真實(shí)情況存在很大的差異[17-20]。

一個(gè)靜態(tài)負(fù)荷模型表示任意瞬時(shí)的負(fù)荷特性是該瞬時(shí)的母線電壓幅值和頻率的代數(shù)函數(shù)。有功功率分量P和無功功率分量Q分別予以考慮。本文只考慮與電壓相關(guān)的負(fù)荷模型。

傳統(tǒng)上負(fù)荷對電壓的依賴特性用指數(shù)模型來表示

式（1）中，P和Q為當(dāng)母線電壓幅值為V時(shí)的負(fù)荷有功和無功分量；下標(biāo)“0”表示初始運(yùn)行條件時(shí)相關(guān)變量的值。

這個(gè)模型的參數(shù)是指數(shù)a和b。當(dāng)這些指數(shù)等于0、1或2，該模型分別表示恒功率、恒電流或恒阻抗特性。對于合成負(fù)荷，它們的值取決于負(fù)荷分量的綜合特性。指數(shù)a（或b）接近等于當(dāng)U=U0時(shí)的斜率dP/dU（或dQ/dU）。對于合成的系統(tǒng)負(fù)荷，指數(shù)a通常在0.5～1.8的范圍內(nèi)；指數(shù)b的典型值在1.5～6.0之間。指數(shù)b的一個(gè)突出特點(diǎn)是，它的變化是電壓的非線性函數(shù)。這是配電變壓器和電動(dòng)機(jī)的磁飽和所引起的。電壓越高，Q也就明顯的越大。

廣泛用來表示負(fù)荷的電壓依賴性的另一種模型是多項(xiàng)式模型

這個(gè)模型通常稱為ZIP模型，因?yàn)樗怯珊阕杩梗╖）、恒電流（I）和恒功率（P）分量組成的。該模型的參數(shù)是系數(shù)p1～p3和q1～q3，它們定義了每個(gè)分量的比率。

1.2 三節(jié)點(diǎn)電力系統(tǒng)隨機(jī)時(shí)滯模型

以圖1所示3節(jié)點(diǎn)電力系統(tǒng)，其中發(fā)電機(jī)采用三階實(shí)用模型，計(jì)及勵(lì)磁系統(tǒng)動(dòng)態(tài)。

圖1 典型3節(jié)點(diǎn)系統(tǒng)Fig. 1 A classical three-bus power system

系統(tǒng)中負(fù)荷模型采用的是集結(jié)負(fù)荷模型，包括商業(yè)居民負(fù)荷和工業(yè)負(fù)荷。商業(yè)負(fù)荷和居民負(fù)荷均采用恒功率靜態(tài)模型；工業(yè)負(fù)荷包括靜態(tài)工業(yè)負(fù)荷和動(dòng)態(tài)工業(yè)負(fù)荷兩部分，其中靜態(tài)工業(yè)負(fù)荷采用（ZIP）模型（即恒阻抗、恒電流和恒功率負(fù)荷模型），以恒功率模型為例，動(dòng)態(tài)工業(yè)負(fù)荷則可以表示為如下形式[21-22]：

式（3）中，P0，Q0分別為商業(yè)居民負(fù)荷的有功、無功負(fù)荷；P1，Q1分別為工業(yè)負(fù)荷中靜態(tài)部分的恒功率有功、無功負(fù)荷；U、θ為負(fù)荷母線處電壓及相角。Kpω，Kqω，Kpv，Kqv，Kqv2和T為負(fù)荷特性常數(shù)。

在式（3）的基礎(chǔ)上，建立數(shù)學(xué)模型

式（4）中，Xd和Xq分別為dq軸同步電抗；E′q為q軸暫態(tài)電勢；X′d為d軸暫態(tài)電抗；δ和ω分別為發(fā)電機(jī)轉(zhuǎn)子的功角和角速度；ω0為角速度的初始穩(wěn)態(tài)值；Tj為發(fā)電機(jī)的慣性時(shí)間常數(shù)；Pm為機(jī)械功率；Pe為電磁功率；D為阻尼系數(shù)；T′d0為d軸的暫態(tài)時(shí)間常數(shù)；Ef為勵(lì)磁電勢。PL和QL是負(fù)荷母線從網(wǎng)絡(luò)中吸收的有功功率和無功功率，是U，θ，δ的函數(shù)。其表達(dá)式為

式（5）中，E′0、Y′0和θ′0是根據(jù)戴維南等效定理，將負(fù)荷點(diǎn)上并聯(lián)的電容C并入無窮大母線后系統(tǒng)的等值參數(shù)

其中Y=1/X0，X0=Xl/2+XT2，XT2為降壓變壓器T2的電抗。

式（4）中的Pe、Id、Vg可以根據(jù)發(fā)電機(jī)的電壓向量圖，得到關(guān)系式

式（7）中，EQ為假想電動(dòng)勢；Eq為空載電動(dòng)勢；V為無窮大母線電壓；Vg為勵(lì)磁系統(tǒng)端電壓。

一般認(rèn)為功角穩(wěn)定性問題與有功負(fù)荷的變化有關(guān)，而電壓穩(wěn)定性與無功負(fù)荷有關(guān)，所以本文只考慮無功負(fù)荷發(fā)生隨機(jī)變動(dòng)。在系統(tǒng)已知負(fù)荷模型中加入隨機(jī)擾動(dòng)項(xiàng)，這樣模型系統(tǒng)中無功負(fù)荷形式為

1.3 隨機(jī)微分方程的數(shù)值計(jì)算方法

對于隨機(jī)微分方程來說，只有滿足一些特殊表達(dá)的簡單系統(tǒng)才能夠解析求解，更多的情況需通過數(shù)值計(jì)算方法獲得解過程的軌跡。EM數(shù)值方法是隨機(jī)微分方程最簡單的數(shù)值求解的方法[23-24]。設(shè)隨機(jī)過程X（t）是隨機(jī)微分方程（6）的解過程，對于某個(gè)正整數(shù)N，記Δt=（T-t0）/N，τj=jΔt，Xj=X（τj），j=0，1，2，…，N，則EM數(shù)值方法的差分迭代公式為:

2 系統(tǒng)動(dòng)態(tài)響應(yīng)分析

系統(tǒng)參數(shù)設(shè)置如下：

Tj=10.0 s，Xd=0.982 pu，Xq=0.982 pu，X′d=0.344 pu，Xe=0.604 pu，D=0.5，T′d0=5.0 s，Vg0=1.05 pu，Pm=0.45 pu，ω0=314.15 rad/s，Ka=5，Ta=0.5 s，C=12，Kpω=-0.4，Kqω=-0.03，Kpv=0.3，Kqv，=-2.8，Kqv2=2.1，T=8.5，P0=0.6，Q0=1.3，P1=0，Q0=10.946。

采用EM數(shù)值方法仿真獲得系統(tǒng)電壓變化曲線，仿真1 000次并平均，根據(jù)平均響應(yīng)曲線分析時(shí)滯對穩(wěn)定性的影響。

1）σ0=0.01

從圖2中，我們可以看出，當(dāng)激勵(lì)強(qiáng)度σ0=0.01時(shí)，電壓振蕩幅度比較小，當(dāng)時(shí)滯比較小時(shí)，負(fù)荷母線電壓變化比較平衡，系統(tǒng)穩(wěn)定。而當(dāng)隨著時(shí)滯τ逐漸增大到τ=0.139 s時(shí)，母線電壓大幅度降低，系統(tǒng)失穩(wěn)。

2）σ0=0.5

3）σ0=1

從圖3中可以看出，當(dāng)激勵(lì)強(qiáng)度σ0=0.5時(shí)，電壓振蕩幅度比較小，當(dāng)時(shí)滯比較小時(shí)，負(fù)荷母線電壓變化比較平衡，系統(tǒng)穩(wěn)定。而當(dāng)隨著時(shí)滯τ逐漸增大到τ=0.121 s時(shí)，母線電壓大幅度降低，系統(tǒng)失穩(wěn)。圖4中，當(dāng)激勵(lì)強(qiáng)度σ0=1時(shí)，無論有沒有時(shí)滯，或者時(shí)滯數(shù)值為多大，系統(tǒng)都失穩(wěn)。比較圖2、圖3可以發(fā)現(xiàn)，隨著激勵(lì)強(qiáng)度的增加，負(fù)荷母線電壓呈現(xiàn)明顯的隨機(jī)振蕩，這主要是由于隨機(jī)擾動(dòng)添加到了Q1中，負(fù)荷點(diǎn)的電壓則對無功波動(dòng)反應(yīng)迅速，而且隨機(jī)激勵(lì)強(qiáng)度越大，系統(tǒng)振蕩頻率越快。

圖4 σ0=1時(shí)電壓響應(yīng)曲線Fig. 4 The voltage curve under σ0=1

表1列出了不同隨機(jī)擾動(dòng)強(qiáng)度，不同時(shí)滯作用下系統(tǒng)的穩(wěn)定域與不穩(wěn)定域。

表1 恒功率模型時(shí)不同時(shí)滯、不同激勵(lì)強(qiáng)度下的電壓穩(wěn)定性Tab. 1 Stability under different excitation intensities and time delays with constant power load

從表1中，可分析出：隨機(jī)時(shí)滯電力系統(tǒng)的穩(wěn)定性不僅取決于激勵(lì)的強(qiáng)度，還與時(shí)滯的大小密切有關(guān)。且隨著激勵(lì)強(qiáng)度的增加，系統(tǒng)可承受的時(shí)滯上限逐漸降低，穩(wěn)定域縮小。

3 負(fù)荷模型對系統(tǒng)穩(wěn)定性的影響

當(dāng)負(fù)荷采用恒阻抗模型時(shí)，即

分析在這種模型條件下，不同時(shí)滯、不同激勵(lì)強(qiáng)度下系統(tǒng)的穩(wěn)定性。仿真結(jié)果如表2所示。

表2 恒阻抗模型時(shí)不同時(shí)滯、不同激勵(lì)強(qiáng)度下的電壓穩(wěn)定性Tab. 2 Stability under different excitation intensities and time delays with constant impedance load

當(dāng)負(fù)荷采用恒阻抗模型時(shí)，即

分析在這種模型條件下，不同時(shí)滯、不同激勵(lì)強(qiáng)度下系統(tǒng)的穩(wěn)定性。仿真結(jié)果如表3所示。

表3 恒電流模型時(shí)不同時(shí)滯、不同激勵(lì)強(qiáng)度下的電壓穩(wěn)定性Tab. 3 Stability under different excitation intensities and time delays with constant current load

綜合表1—表3可知，采用不同的負(fù)荷模型，對時(shí)滯隨機(jī)電力系統(tǒng)的穩(wěn)定性會(huì)有不同的影響。在相同的隨機(jī)激勵(lì)強(qiáng)度下，不同的負(fù)荷模型所能承受的時(shí)滯上限有所差異。3種負(fù)荷模型所能承受的時(shí)滯上限的能力由大到小的順序?yàn)椋汉阕杩鼓Ｐ?恒電流模型>恒功率模型。其原因在于：系統(tǒng)中負(fù)荷母線處安裝的電容，其作用是為了支撐負(fù)荷點(diǎn)的電壓，但負(fù)荷模型中隨機(jī)擾動(dòng)的存在會(huì)導(dǎo)致系統(tǒng)電壓產(chǎn)生隨機(jī)波動(dòng)并略有下降。由于恒阻抗模型所需無功功率與電壓平方值成正比，恒電流模型所需無功功率與電壓值成正比，因此這兩種模型所需要的無功功率都會(huì)隨著系統(tǒng)電壓的降低而減少，且恒阻抗模型減少的程度多于恒電流模型。而恒功率模型所需的無功功率并不隨系統(tǒng)電壓變化而變化，導(dǎo)致系統(tǒng)電壓進(jìn)一步下降，因此會(huì)更加不利于系統(tǒng)的電壓穩(wěn)定。

4 結(jié)語

本文基于一個(gè)經(jīng)典3節(jié)點(diǎn)電力系統(tǒng)模型，以勵(lì)磁系統(tǒng)輸入信號(hào)延時(shí)作為時(shí)滯環(huán)節(jié)，并通過引入負(fù)荷擾動(dòng)項(xiàng)，建立隨機(jī)時(shí)滯系統(tǒng)的數(shù)學(xué)模型，并運(yùn)用隨機(jī)歐拉法對該微分方程組進(jìn)行數(shù)值求解，仿真分析不同時(shí)滯、不同擾動(dòng)強(qiáng)度下隨機(jī)時(shí)滯電力系統(tǒng)的穩(wěn)定性；同時(shí)負(fù)荷模型分別采用恒阻抗、恒功率、恒電流模型，分析比較不同負(fù)荷模型對系統(tǒng)穩(wěn)定性的影響。仿真表明，3種負(fù)荷模型所能承受的時(shí)滯上限能力從大到小依次為：恒阻抗模型、恒電流模型、恒功率模型。

[1] 江全元，鄒振宇，曹一家，等. 考慮時(shí)滯影響的電力系統(tǒng)分析和時(shí)滯控制研究進(jìn)展[J]. 電力系統(tǒng)自動(dòng)化，2005，29（3）: 2-7.JIANG Quanyuan，ZOU Zhenyu，CAO Yijia，et al.Overview of power system stability analysis and wide-area control in consideration of time delay[J]. Automation of Electric Power System，2005，29（3）: 2-7（in Chinese）.

[2] 關(guān)海平，焦曉亮. 600 MW汽輪發(fā)電機(jī)組運(yùn)行狀況的調(diào)研分析[J]. 節(jié)能技術(shù)，2011，29（1）: 57－65.GUAN Haiping，JIAO Xiaoliang. Study and analysis on operating conditions of 600 MW turbo-generator sets[J].Energy Conservation Technology，2011，29（1）: 57－65（in Chinese）.

[3] 劉春梅，李勇，辛楊，等. 電力系統(tǒng)低頻振蕩的影響因素研究[J]. 陜西電力，2011，39（6）: 12-15.LIU Chunmei，LI Yong，XIN Yang，et al. Study on influencing factors of low frequency oscillation of power system[J]. Shaanxi Electric Power，2011，39（6）: 12-15（in Chinese）.

[4] 韓志勇，李志強(qiáng)，李瑩，等. 汽輪機(jī)調(diào)速系統(tǒng)引起電力電力系統(tǒng)共振機(jī)理低頻振蕩擾動(dòng)分析[J]. 陜西電力，2009，37（7）: 1－5.HAN Zhiyong，LI Zhiqiang，LI Ying，et al. Power system low frequency oscillation of resource mechanism induced by disturbance of turbine-governing system[J]. Shaanxi Electric Power，2009，37（1）: 1－5（in Chinese）.

[5] 楊培宏，馮義，張遇杰. 基于廣域測量系統(tǒng)的電力系統(tǒng)低頻振蕩抑制方法[J]. 陜西電力，2008，36（2）: 11－15.YANG Peihong，F(xiàn)ENG Yi，ZHANG Yujie. Low-frequency oscillation damping method for power system based on wide-area measurement system[J]. Shaanxi Electric Power，2008，36（2）: 11-15（in Chinese）.

[6] 張秀坤，蔣明東. 國產(chǎn)600 MW汽輪機(jī)配汽方式優(yōu)化改造[J]. 節(jié)能技術(shù)，2005，23（1）: 87－90.ZHANG Xiukun，JIANG Mingdong，et al. Optimum governing of 600 MW turbine[J]. Energy Conservation Technology，2005，23（1）: 87－90（in Chinese）.

[7] 楊培宏，劉達(dá). 基于廣域測量系統(tǒng)的電力系統(tǒng)穩(wěn)定器適配器方案[J]. 陜西電力，2009，38（4）: 11－14.YANG Peihong，LIU Da. Study of power system atabilizer configuration based on wide-area measurement system[J].Shaanxi Electric Power，2009，38（4）: 11－14（in Chinese）.

[8] 劉東暉，錢乙衛(wèi)，劉東英. WAMS在陜西電網(wǎng)的應(yīng)用分析[J]. 陜西電力，2012，40（6）: 63－65.LIU Donghui，QIAN Yiwei，LIU Dongying. Application of WAMS in shanxi power network[J]. Shaanxi Electric Power，2008，36（2）: 11－15（in Chinese）.

[9] NWANKPA C O，SHAHIDEHPOUR S M，SCHUSSZ. A stochastic approach to small disturbance stability analysis[J].IEEE Trans on Power Systems，1992，7（4）：1519-15287.

[10] BACHELIER L. Theorie de la speculation[J]. Annales Scitifiques De L.E.N.S. 1990（17）：21-86.

[11] YUE D，HAN Q L. Delay-dependent exponential stability of stochastic systems with time-varying delay，nonlinearity，and Markovian switching[J]. IEEE Trans on Automation Control，2005，52（2）：217-222.

[12] YUE D，WON S. Delay-dependent robust stability of stochastic systems with time delay and nonlinearity uncertainties[J]. IEE Electron Lett，2001，37（15）：992-993.

[13] CHEN W H，GUAN Z H，LU X M. Delay-dependent exponential stability of uncertain stochastic systems with multiple delays：an LMI approach[J]. System & Control，Lett，2005（54）：547-555.

[14] HUA M，DENG F，LIU X，et al. Robust delay-dependent exponential stability of uncertain stochastic system with time-varying delay[J]. Circuits Syst Signal Process，2010（29）：519-526.

[15] XU S Y，CHEN T. Robust H∞control for uncertain stochastic systems with state delay[J]. IEEE Trans on Automatic Control，2002，47（12）：2089-2094.

[16] 倪以信，陳壽孫，張寶霖. 動(dòng)態(tài)電力系統(tǒng)的理論和分析[M]. 北京：清華大學(xué)出版社，2002.

[17] 彭旭東，邱澤晶，肖楚鵬. 電力用戶能效監(jiān)測與需求響應(yīng)研究綜述[J]. 節(jié)能技術(shù)，2013，31（3）: 243-246.PENG Xudong，QIU Zejing，XIAO Chupeng. Overview on power efficiency monitoring and user demand response[J].Energy Conservation Technology，2013，31（3）: 243-246（in Chinese）.

[18] 李明，郭運(yùn)泉，韓英義. 63 kV變電所及其10 kV配電線路最佳無功補(bǔ)償[J]. 節(jié)能技術(shù)，2001，19（2）: 42-43..LI Ming，GUO Yunquan，HAN Yingyi. The optimal reactive power compensation for the 63 kV substation and 10 kV distribution line[J].Energy Conservation Technology，2001，19（2）: 42-43（in Chinese）.

[19] 郭培源. 經(jīng)濟(jì)功率因數(shù)與最佳無功補(bǔ)償[J]. 節(jié)能技術(shù)，1992（3）: 10-13，8 GUO Peiyuan. Economic power factor and optimal reactive power compensation[J]. Energy Conservation Technology，1992（3）: 10-13，8（in Chinese）.

[20] 高寧，趙有東，翟紗，等. 電力行業(yè)中電壓概念淺析[J].西北水電，2013（2）: 57-61.GAO Ning，ZHAO Youdong，ZHAI Sha，et al. The concept of voltage in the electric power industry[J]. Northwest Hydropower，2013（2）: 57-61（in Chinese）.

[21] 肖炏，郭永基. 負(fù)荷隨機(jī)擾動(dòng)對電力系統(tǒng)電壓穩(wěn)定性的影響[J]. 清華大學(xué)學(xué)報(bào)：自然科學(xué)版，2000，40（9）：9-12.XIAO Kai，GUO Yongji. Effects of stochastic load model on power system voltage stability[J]. Journal of Tsinghua University：Science and Technology，2000，40（9）：9-12（in Chinese）.

[22] 邱妍，趙晉泉，朱永忠. 負(fù)荷隨機(jī)擾動(dòng)對電力系統(tǒng)電壓穩(wěn)定性的影響[J]. 電力自動(dòng)化設(shè)備，2009，29（2）：77-81.QIU Yan，ZHAO Jinquan，ZHU Yongzhong. Effect of stochastic load disturbance on power system voltage stability based on bifurcation theory[J]. Electric Power Automation Equipment，2009，29（2）：77-81（in Chinese）.

[23] 劉詠飛，鞠平，薛禹勝，等. 隨機(jī)激勵(lì)下電力系統(tǒng)特性的計(jì)算分析[J]. 電力系統(tǒng)自動(dòng)化，2014，38（9）：137-142.LIU Yongfei，JU Ping，XUE Yusheng，et al. The calculation and analysis of power system under random excitation[J].Automation of Electric Power System，2014，38（9）：137-142（in Chinese）.

[24] HAMDY M A. Euler-maruyama numerical solution of some stochastic functional differential equations[J]. Internation Journal of Applied Mathematics and Computation，2007，11（1）：13-30.