999精品在线视频,手机成人午夜在线视频,久久不卡国产精品无码,中日无码在线观看,成人av手机在线观看,日韩精品亚洲一区中文字幕,亚洲av无码人妻,四虎国产在线观看 ?

計及DG出力相關性的孤島微電網蒙特卡洛法概率潮流

2018-04-08 02:36:51蘇凱森楊家豪鄭澤蔚弓新月
電力工程技術 2018年2期

蘇凱森, 楊家豪, 鄭澤蔚, 弓新月

(1. 廈門大學嘉庚學院,福建 漳州 363105;2. 江蘇省電力有限公司電力科學研究院, 江蘇 南京 211103)

0 引言

近年來,微電網(microgrid,MG)技術得到快速發(fā)展[1-2],其離網孤島運行的能力大大提高了供電可靠性,且對于海島、偏遠地區(qū)等,孤島型微電網也成為了解決供電問題的有效方案。由于無平衡節(jié)點,缺乏主網支撐,通常孤島型微電網采用對等控制方式[3],即通過采取下垂控制策略的分布式電源(distributed generator,DG)、儲能裝置(energy s ̄t ̄o ̄r ̄age,ES)等共同參與微電網的頻率與電壓的調節(jié)。

諸如光伏、風電這一類DG的出力具有較強的隨機性,同時負荷也存在波動,對微電網的頻率及電壓質量造成影響,而孤島微電網的潮流計算是分析基礎。文獻[4—8]建立了孤島微電網的潮流模型,其中文獻[6—7]基于信賴域算法提高潮流計算收斂性,文獻[8]提出類奔德斯分解方法較好地提升潮流收斂速度。但對于單一運行點的潮流結果,并未能對微電網頻率及電壓在不確定性條件下可能的狀態(tài)提供充分的信息。相比之下,概率潮流(probabilistic load flow,PLF)能夠獲得各狀態(tài)變量的完整概率分布信息,更為全面地反映微電網可能的運行狀態(tài)[9-13]。目前關于孤島微電網PLF的研究較少,在現有研究中蒙特卡洛模擬法[11]是求解PLF的常用方法。文獻[12]采用蒙特卡洛法模擬間歇性微電源的隨機特性,從而求解微電網PLF。文獻[13]基于半不變量法求解孤島微電網PLF,計算時將柴油發(fā)電機視為平衡節(jié)點,儲能視為恒功率負荷,相當于采取主從控制模式,由平衡節(jié)點承擔全部功率波動,不符合實際孤島微電網的運行特點。

上述文獻存在2方面不足:(1) 均未考慮對等控制模式下的孤島微電網的情況,即仍是默認頻率為額定值,只考察系統(tǒng)中電壓的概率分布情況。實際中頻率是孤島微電網PLF中重要的待求隨機變量,頻率的概率分布信息對微電網運行決策提供重要的參考依據;(2) 鑒于微電網屬于獨立小型系統(tǒng),微電網中多個光伏電源或風電機組鄰近,其出力往往具有較強的相關性,需要在孤島微電網PLF中計及相關性的影響。

本文的主要貢獻是提出一種計及DG出力相關性的孤島微電網蒙特卡洛概率潮流計算方法。首先簡化孤島微電網節(jié)點分類,建立孤島微電網功率方程;其次建立各類DG及負荷的概率模型,并提出DG出力相關性的處理方法,通過蒙特卡洛模擬計算孤島微電網PLF;最終通過Benchmark 0.4 kV低壓微電網作為算例對文中方法進行驗證。

1 孤島微電網潮流計算

目前孤島微電網潮流模型中的等效節(jié)點類型通常為:PQ節(jié)點、PV節(jié)點及下垂節(jié)點。當出現越限或設備故障脫網時節(jié)點類型可能變化,為簡化節(jié)點分類、減少節(jié)點類型轉換,同時提高計算效率,文中僅設定PQ節(jié)點與PV節(jié)點2種節(jié)點類型,并建立統(tǒng)一的節(jié)點功率方程。

1.1 潮流方程

PQ節(jié)點需列寫有功與無功平衡方程,PV節(jié)點僅需列寫有功平衡方程,方程形式為:

(1)

式中:PCi,QCi分別為節(jié)點i恒功率電源注入的有功和無功,例如光伏、風電等采取最大功率點跟蹤(maximum power point tracking,MPPT)控制方式的DG屬于此類電源;PDi,QDi分別為節(jié)點i具有下垂控制特性的設備注入的有功和無功,例如柴油發(fā)電機、燃氣輪機、儲能等均屬于此類電源;PLi,QLi分別為節(jié)點i的有功與無功負荷;Pi,Qi分別為節(jié)點i注入的總有功功率和無功功率。

節(jié)點注入的有功與無功功率為:

(2)

式中:n為節(jié)點數目;Ui及Uj分別為節(jié)點i與節(jié)點j的電壓;Gij,Bij分別為節(jié)點導納矩陣的實部與虛部;δij為節(jié)點i和節(jié)點j的相角差。

1.2 下垂控制設備建模

具有下垂控制特性的設備注入的功率可統(tǒng)一表示為:

(3)

(4)

(5)

(6)

式中:PDimax,PDimin,QDimax,QDimin分別為節(jié)點i具有下垂控制特性的設備注入的有功和無功上下限,若無此類設備則均取為0;fmax,fmin,Umax,Umin分別為系統(tǒng)頻率、電壓允許上下限;KDfi,KDUi分別為對應的P-f、Q-U下垂系數;f0,f分別為系統(tǒng)頻率的空載值與實際值;U0,Ui分別為電壓的空載值與實際值。其中儲能可以工作在充放電2種工況,在不考慮荷電狀態(tài)的情況下,其PDimin,PDimax取決于儲能的最大充放電功率。

1.3 負荷建模

對負荷進行建模時計及負荷的電壓和頻率靜特性,負荷使用恒阻抗、恒電流和恒功率的組合模型來描述,可表示為:

(7)

式中:PLNi,QLNi分別為節(jié)點i在額定工況下的有功與無功負荷;UNi和fN分別為額定電壓與頻率,取UNi=1 p.u.,fN=1 p.u.;Pi,Qi分別為節(jié)點i注入的總有功功率和無功功率;Api,Bpi,Cpi,Aqi,Bqi,Cqi分別為負荷有功與無功功率中恒阻抗型、恒電流型、恒功率型的百分比系數,分別滿足Api+Bpi+Cpi=1及Aqi+Bqi+Cqi=1;kLpi,kLqi分別為負荷的有功和無功功率的靜態(tài)頻率調節(jié)系數。

1.4 潮流求解

假定孤島微電網中PQ節(jié)點與PV節(jié)點的數目分別為nPQ與nPV,則總計列寫2nPQ+nPV個方程,本文采用牛頓拉夫遜法求解潮流方程組,修正方程簡寫為:

(8)

式中:ΔP,ΔQ為節(jié)點有功與無功不平衡量;Δf,Δδ,ΔU為分別為頻率、相角、電壓的修正量;J為雅克比矩陣,其分塊矩陣分別為:

(9)

由于求解過程中不涉及節(jié)點類型轉換,因此也無需在迭代過程中改變雅克比矩陣的結構,只需相應更新數值,因此能夠為蒙特卡洛模擬提高運算效率。

2 孤島微電網概率潮流

2.1 蒙特卡洛法概率模型

光伏電源、風電機組為了能在每一時刻盡可能輸出最大功率,通常采取MPPT的控制方式,因此DG出力隨著光照與風速的隨機性也呈現不確定性,而負荷同樣存在波動,故而孤島微電網概率潮流計算中具有隨機性的輸入變量主要是間歇性DG以及負荷。文中分析涉及的間歇性DG裝置為光伏電源(photovoltaic,PV)及風力發(fā)電機(wind turbine,WT)裝置,文中假設PV與WT均不發(fā)出無功。

研究表明,一段時間內的光伏出力滿足Beta分布,因此光伏發(fā)電有功出力PPV的概率模型以Beta分布來近似表達:

(10)

式中:Γ為Gamma函數;α,β為Beta分布的形狀參數;Pmax為該時段內光伏出力最大值。

風速概率分布通常滿足Weibull分布,因此風力發(fā)電有功出力PWT的概率模型表示為:

(11)

式中:k為形狀參數;c為尺度參數;a,b為與WT參數有關的常數。

負荷的概率模型用正態(tài)分布描述,其有功和無功功率的概率模型為:

(12)

(13)

式中:μP,μQ分別為有功和無功的均值;σP,σQ分別為有功和無功的標準差。

2.2 DG出力相關性處理

已知m個輸入隨機變量W=[w1,w2,……,wm]的邊際概率分布F(W),并使用相關系數矩陣RW描述出力相關性,相關系數矩陣表達式為:

(14)

矩陣中各元素按下式計算:

(15)

式中:ρwij為隨機變量wi與wj之間的相關系數;Cov(wi,wj)為wi與wj的協(xié)方差;σwi,σwj為wi與wj的標準差。

根據三階多項式正態(tài)變換(third-order p ̄o ̄l ̄y ̄n ̄o ̄m ̄ial normal transformation,TPNT)理論[14],具有相關性的隨機變量空間W可由獨立的服從標準正態(tài)分布隨機變量空間Z的三階多項式表示,即:

(16)

式中:ai(k)(k=0,1,2,3)為各項系數;zi為獨立的服從標準正態(tài)分布的隨機變量。

利用wi_st=(wi-μwi)/σwi將隨機變量wi進行標準化,則有:

(17)

式中:bi(k)(k=0,1,2,3)為標準化后對應的各項系數。顯然各項系數滿足以下關系:

(18)

根據矩法原理,式(17)中等式兩邊的各階原點矩相等,取前四階原點矩的相等關系可展開得到4個方程:

(19)

式中:χ為wi_st的偏度;κ為wi_st的峰度。

通過求解式(19)非線性方程組可得bi(k),再根據式(18)求得ai(k),最后根據式(20)可計算得到具有相關性的服從標準正態(tài)分布的隨機變量空間Y的相關系數矩陣中的各元素ρyij。

(20)

相關系數應為區(qū)間[-1,1]內一實數,同時在上述進行空間映射時,相關系數也應為同號,即ρyij須滿足:

(21)

文中建立滿足給定相關性水平的輸入變量概率模型的步驟為:

(1) 產生滿足獨立正態(tài)分布的樣本矩陣S,其維數為m×N,N為樣本規(guī)模;

(2) 運用TPNT理論以及求解代數方程的方法求取具有相關性的標準正態(tài)分布隨機變量Y的等效相關系數矩陣RY,即通過式(16—21)求取RY中的各元素ρyij;

(3) 根據式(22)對RY進行Cholesky分解[15]

RY=LYLYT

(22)

式中:LY即為Cholesky因式分解所得的下三角矩陣。

隨后可由Y=LYS獲得具有相關性的服從標準正態(tài)分布的樣本矩陣Y;

(4) 通過等概率轉換原則[15-16]得到滿足相關性水平的樣本矩陣W,即對每個樣本利用式(23)由yi求取wi。

F(wi)=Φ(yi)

(23)

式中:F(wi)為wi的累積分布函數,Φ(yi)為yi的累積分布函數。

2.3 算法流程

文所提PLF計算的流程如下:

(1) 輸入原始數據,包含微電網網架數據、電源和負荷的基本參數以及所對應的概率模型;

(2) 設定樣本規(guī)模,對于相互獨立的隨機變量基于概率模型進行抽樣;

(3) 對于相關的隨機變量,由給定的F(W)及RW由式(14—23)建立滿足給定相關性水平的樣本矩陣;

(4) 將全部樣本逐個進行孤島微電網潮流計算,并保存頻率及電壓的計算結果;

(5) 對計算結果進行統(tǒng)計,求取頻率及電壓的的概率密度函數(probability density function,PDF)及累積分布函數(cumulative distribution function,CDF),兩者滿足以下關系:

(24)

式中:X為狀態(tài)變量,即頻率及節(jié)點電壓;FX(x)為狀態(tài)變量X的CDF;fX(x)為狀態(tài)變量X的PDF。算法流程如圖1所示。

圖1 算法流程Fig.1 Algorithm flow chart

3 算例分析

以Benchmark 0.4 kV低壓微電網[17]作為算例系統(tǒng),如圖2所示。其中S1打開,S2閉合,構成孤島微電網系統(tǒng)。系統(tǒng)基準容量取100 kV·A,假定孤島微電網的安全運行范圍為fmax=1.004 p.u.,fmin=0.996 p.u.,Umax=1.05 p.u.,Umin=0.95 p.u.。

圖2 Benchmark 0.4 kV低壓微電網系統(tǒng)Fig.2 Benchmark 0.4 kv low-voltage microgrid system

接入設備參數見表1,其中節(jié)點13—17的電源中,PV及WT采取MPPT的控制方式,在單一樣本潮流計算中視為恒功率電源,而其他的蓄電池組、微型燃氣輪機等采取對等控制,共同參與微電網頻率及電壓的調節(jié)。考慮到負荷情況,在節(jié)點14及17進行必要的電容器固定補償。

表1 算例中接入的設備及參數Tab. 1 Example of access to theequipment and parameters

負荷參數見表2。各節(jié)點恒阻抗、恒電流、恒功率負荷占比統(tǒng)一取為0.3,0.3,0.4,靜態(tài)頻率調節(jié)系數取kLpi=2,kLqi=-2[5]。負荷服從正態(tài)分布,負荷波動的標準差取為期望值的10%。

表2 負荷參數Tab. 2 Load parameters

各PV出力滿足Beta分布,其中形狀參數為α= 1.693,β=5.162。WT出力滿足Weibull分布,切入風速2.7 m/s,額定風速6.7 m/s,可計算得a=-40.5,b=15,另外有形狀參數k=2.94,尺度參數c=3.03。

場景1:假定節(jié)點15的WT退出運行,PV之間的相關系數ρPP取為0.8。

是否計及相關性的2種情況下所求得的頻率及節(jié)點電壓標準差對比如表3所示。

表3 場景1的狀態(tài)變量的標準差對比Tab. 3 Standard deviation comparisonof state variables in scene 1

由表中數據可知,計及相關性后系統(tǒng)頻率可能的波動范圍明顯增大,主要是節(jié)點15與節(jié)點16的PV出力呈現將強的相關性,將引起系統(tǒng)功率的較大幅度波動。若不計相關性的影響將會帶來較大誤差,頻率的相對誤差將近20%。

由于PV只輸出有功功率,有功功率的平衡主要影響頻率。對比節(jié)點15與節(jié)點16的PV出力相互獨立以及計及相關性下的頻率分布情況,頻率的PDF及CDF如圖3所示。

圖3 場景1的頻率PDF及CDFFig.3 Frequency PDF and CDF in scene 1

場景2:節(jié)點15的WT投入運行,PV之間的相關系數ρPP取為0.8,由于風光通常存在互補性,PV與WT之間的相關系數ρPW取為-0.6。

是否計及相關性的2種情況下所求得的頻率及節(jié)點電壓標準差對比如表4所示。

表4 場景2的狀態(tài)變量的標準差對比Tab. 4 Standard deviation comparisonof state variables in scene 2

由表中數據可知,由于PV及WT出力成負相關性,即風光互補,提高了可再生能源輸出的穩(wěn)定性,與不考慮相關性的情況相比頻率波動范圍大大減小。若將風光出力視為相互獨立則會產生極大的誤差,由此可見本文計及DG出力相關性的計算結果更符合實際情況,能夠更準確地評估頻率及電壓質量。

同樣對比各DG出力相互獨立以及計及相關性下的頻率分布情況,頻率的PDF及CDF如圖4所示。

使用盒須圖繪制出場景2下的節(jié)點電壓分布情況,如圖5所示。場景二下各節(jié)點電壓都在一定范圍內波動,均處于合格范圍內,其中節(jié)點16波動區(qū)間較大,主要由于節(jié)點16的PV出力存在隨機性,且容量較大。而節(jié)點15同樣安裝了PV但由于還有WT的存在,風光出力互補因此使得總體出力更為穩(wěn)定,電壓波動也較小些。

圖4 場景2的頻率PDF及CDFFig.4 Frequency PDF and CDF in scene 2

圖5 節(jié)點電壓盒須圖Fig.5 Node voltage box plot

4 結論

文中提出一種計及DG出力相關性的孤島微電網PLF計算方法,可用于求解對等控制模式下的孤島微電網PLF。通過算例分析,獲得以下結論:

(1) 孤島微電網中的可再生能源主要是光伏及風力發(fā)電,由本文場景1的情況可以看出同類DG出力具有強相關性的情況下將加劇孤島微電網的頻率波動;

(2) 一般來說,風光出力存在天然的互補性,由本文場景2對比結果可知,通過風光互補能夠獲得較為穩(wěn)定的出力,使得頻率波動范圍減小,可提高孤島微電網運行的安全性;

(3) 文中計及DG出力相關性的PLF相比于相互獨立的抽樣模擬更符合實際情況,能夠為孤島微電網的運行提供可信的參考數據,存在應用價值。

參考文獻:

[1] HATZIARGYRIOU N,ASAND H,IRAVANI R,et al. Microgrids [J]. IEEE Power and Energy Magazine,2007,5(4):78-94.

[2] KATIRAEI F,IRAVANI R,HATZIARGYRIOU N,et al. Microgrids management[J]. IEEE Power and Energy Magazine,2008,6(3):54-65.

[3] 彭寒梅,曹一家,黃小慶. 對等控制孤島微電網的靜態(tài)安全風險評估[J]. 中國電機工程學報,2016, 36(18):4837-4846.

PENG Hanmei,CAO Yijia,HUANG Xiaoqing. Static security risk assessment of islanded microgrids under peer-peer control[J]. Proceedings of the CSEE,2016,36(18):4837-4846.

[4] MUMTAZ F,SYED M H,HOSANI M A,et al. A novel approach to solve power flow for islanded microgrids using modified Newton Raphson with droop control of DG[J]. IEEE Transactions on Sustainable Energy,2015,6(12):1-11.

[5] 潘忠美,劉健,石夢,等. 計及電壓/頻率靜特性的孤島微電網電壓穩(wěn)定性與薄弱節(jié)點分析[J]. 電網技術,2017,41(7):2214-2221.

PAN Zhongmei,LIU Jian,SHI Meng,et al. Static voltage stability and weak bus analysis of islanded microgrids considering static voltage/frequency characteristics[J]. Power System Technolgy,2017,41(7):2214-2221.

[6] ABDELAZIZ M M,FARAG H E,EL-SAADANY E F,et al. A novel and generalized three-phase power flow algorithm for islanded microgrids using a newton trust region method[J]. IEEE Trans on Power Systems,2013,28(1):190-201.

[7] 彭寒梅,曹一家,黃小慶. 基于 BFGS 信賴域算法的孤島微電網潮流計算[J]. 中國電機工程學報,2014,34(16):2629-2638.

PENG Hanmei,CAO Yijia,HUANG Xiaoqing. Power flow calculation of islanded microgrids based on BFGS trust region method[J]. Proceedings of the CSEE,2014,34(16):2629-2638.

[8] 李培帥,施燁,吳在軍,等. 孤島微電網潮流的類奔德斯分解算法[J]. 電力系統(tǒng)自動化,2017,41(14):119-125.

LI Peishuai,SHI Ye,WU Zaijun,et al. Power flow calculation method similar to benders decomposition for islanded microgrid[J]. Automation of Electric Power Systems,2017,41(14):119-125.

[9] 金楚,黎嘉明,徐沈智,等. 大規(guī)模光伏發(fā)電并網概率潮流計算及對電網的影響[J]. 電力工程技術,2017,36(1):1-8.

JIN Chu,LI Jiaming,XU Shenzhi,et al. Probabilistic load flow calculation and influence analysis for power grid connected with large scale photovoltaic generation system[J]. Electric Power Engineering Technology,2017,36(1):1-8.

[10] 衛(wèi)鵬,劉建坤,周前,等. 基于半不變量和Gram-Charlier級數展開法的隨機潮流算法[J]. 電力工程技術,2017,36(1):34-38.

WEI Peng,LIU Jiankun,ZHOU Qian,et al. A probabilistic power flow algorithm based on semi-variable and gram-charlier series expansion[J]. Electric Power Engineering Technology,2017,36(1):34-38.

[11] 丁明,李生虎,黃凱. 基于蒙特卡羅模擬的概率潮流計算[J]. 電網技術,2001,25(11):10-14, 22.

DING Ming,LI Shenghu,HUANG Kai. Probabilistic load flow analysis based on monte-carlo simulation[J]. Power System Technology,2001,25(11):10-14, 22.

[12] 段玉兵,龔宇雷,譚興國,等. 基于蒙特卡羅模擬的微電網隨機潮流計算方法[J]. 電工技術學報,2011,26(S1):274-278.

DUAN Yubing,GONG Yulei,TAN Xinguog,et al. Probabilistic power flow calculation in microgrid based on monte-carlo simulation. Transactions of China Electrotechnical Society[J].2011,26(S1):274-278.

[13] 劉明波,簡淦楊,董萍. 孤島微電網的隨機潮流計算[J]. 華南理工大學學報(自然科學版),2014,42(4):13-18.

LIU Mingbo,JIAN Ganyang,DONG Ping. Stochastic power flow calculation of islanded microgrids[J]. Journal of South China University of Technology(Natural Science Edition),2014,42(4):13-18.

[14] 劉小團,趙晉泉,羅衛(wèi)華,等. 基于TPNT和半不變量法的考慮輸入量相關性概率潮流算法[J]. 電力系統(tǒng)保護與控制,2013,41(22):13-18.

LIU Xiaotuan,ZHAO Jinquan,LUO Weihua,et al. A TPNT and cumulants based probabilistic load flow approach considering the correlation variables[J]. Power System Protection and Control,2013,41(22):13-18.

[15] 石東源,蔡德福,陳金富,等. 計及輸入變量相關性的半不變量法概率潮流計算[J]. 中國電機工程學報,2012,32(28):104-113.

SHI Dongyuan,CAI Defu,CHEN Jinfu,et al. Probabilistic load flow calculation based on cumulant method considering correlation between input variables[J]. Proceedings of the Chinese Society for Electrical Engineering,2012,32(28):104-113.

[16] 吳巍,汪可友,李國杰. 計及光伏發(fā)電相關性的多重積分法概率潮流計算[J]. 中國電機工程學報,2015, 35(3):568-575.

WU Wei,WANG Keyou,LI Guojie. Probabilistic load flow calculation method based on multiple integral method considering correlation of photovoltaic generation[J]. Proceedings of the CSEE,2015,35(3):568-575.

[17] 賈宏杰,戚艷,穆云飛. 基于家居型溫控負荷的孤立微電網頻率控制方法[J]. 中國科學:技術科學,2013,56(3):693-702.

JIA HongJie,QI Yan,MU YunFei. Frequency response of autonomous microgrid based on family-friendly controllable loads[J]. Science China:Technological Sciences,2013,56(3):693-702.

主站蜘蛛池模板: 欧美在线伊人| 免费一级无码在线网站| 亚洲精品色AV无码看| 青青操视频在线| 69av在线| 91久久青青草原精品国产| 亚洲精品手机在线| 露脸一二三区国语对白| 色欲色欲久久综合网| 特级aaaaaaaaa毛片免费视频| 一级毛片基地| 中国一级特黄大片在线观看| 久久国产黑丝袜视频| 在线视频一区二区三区不卡| 动漫精品中文字幕无码| a级高清毛片| 国产成人综合亚洲欧美在| 国产探花在线视频| 欧美午夜视频| 丰满的熟女一区二区三区l| 日本国产精品一区久久久| 亚洲成a人片7777| 不卡的在线视频免费观看| 色哟哟国产成人精品| 国产粉嫩粉嫩的18在线播放91 | 波多野结衣久久高清免费| 国产成年女人特黄特色大片免费| 免费jjzz在在线播放国产| 国产在线一区二区视频| 日本少妇又色又爽又高潮| www.av男人.com| 久久精品嫩草研究院| 免费av一区二区三区在线| 91 九色视频丝袜| 欧美亚洲激情| 四虎影视8848永久精品| 欧美黄网在线| 国产成人精品男人的天堂下载| 国产精品黄色片| 欧美啪啪视频免码| 天天综合天天综合| 精品欧美日韩国产日漫一区不卡| 亚洲欧美日本国产综合在线| 久草国产在线观看| 无遮挡国产高潮视频免费观看 | 亚洲欧美日本国产专区一区| 国产成人8x视频一区二区| 幺女国产一级毛片| 老司机精品一区在线视频| 日韩精品成人网页视频在线| 成人午夜亚洲影视在线观看| 国产三区二区| 亚洲精品国产乱码不卡| 日韩人妻精品一区| 国产成人欧美| 欧美激情网址| 国产伦精品一区二区三区视频优播| 日韩视频精品在线| 婷婷五月在线| 99热这里只有精品2| 国产精品无码AⅤ在线观看播放| 久久精品午夜视频| 狠狠操夜夜爽| 114级毛片免费观看| 国产网站黄| 青青网在线国产| 国产Av无码精品色午夜| 久久激情影院| 色网站在线视频| 国产精品永久在线| swag国产精品| 真人高潮娇喘嗯啊在线观看| 亚洲欧美精品在线| 精品视频一区二区三区在线播| 国产精品第一区| 亚洲第一成年网| 亚洲AV永久无码精品古装片| 成年人国产视频| 亚洲美女AV免费一区| 亚洲九九视频| 国产女人在线视频| 真人免费一级毛片一区二区|