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

基于拓展傅里葉水電機組軸系參數敏感性研究

2024-05-12 06:42:26
水電站機電技術 2024年4期
關鍵詞:方向系統

劉 姮

(國網白銀供電公司,甘肅 白銀 730900)

隨著我國社會與經濟的快速發展,水力發電行業也取得了飛速的進步。水輪發電機組的單機容量日益增大,且隨著能源結構的改變,水電機組與風電、光伏等新能源的多能互補運行成為重要發展方向。但隨之而來的就是由水力、機械、電氣和其它各種因素共同耦合引起的機組運行穩定性問題,機組運行失穩會直接影響機組使用壽命與發電質量水平,嚴重時甚至還會影響水電站的安全運行[1-3]。

水電機組運行工況包括啟動、停機、工況轉換等過程,特別是在與風電、光伏等新能源互補運行時,工況調節更加頻繁,機組在水力、機械和電磁等因素共同作用下,呈現出復雜的動力響應特征,由于這些外部作用因素并不是相互獨立的,且不同因素引起的機組部位振動往往是相互耦合的,故將機組作為一個耦合的動態系統進行研究已經成為一種趨勢,但其難度和工作量也隨之大大增加[4-6]。機組軸系的振動往往是導致水電機組振動事故的主要原因,如導軸承間隙不均勻,大軸弓狀回旋導致水輪機密封間隙偏斜,引發水力不平衡;機組軸系的動態響應導致定轉子間氣隙不均,造成發電機定轉子間不平衡磁拉力加大,加劇機組軸系的振動幅值等[7-9]。因此,對水電機組軸系的各個系統參數影響規律進行針對性研究,對掌握機組振動誘發機理和減振避振措施具有重要意義。

本文以水輪發電機組軸系系統為研究對象,考慮不平衡磁拉力、非線性油膜力和碰摩力對轉子軸承系統的耦合作用,并利用擴展傅里葉敏感性分析方法進行機組軸系系統參數敏感性分析,得出軸系系統參數對其運行穩定性的影響規律。

1 軸系外激勵模型

1.1 不平衡磁拉力

在機組正常運行時,發電機轉子在均勻的磁場中旋轉,轉子徑向受均勻的磁拉力影響。若出現水力不平衡、轉子質量不平衡或發電機的定轉子不圓等因素則會造成定、轉子間氣隙不均,從而產生不平衡磁拉力。不平衡磁拉力主要作用于發電機轉子上,且與轉子偏心同向,這會進一步加劇機組的震動。

目前,不平衡磁拉力的求解方法主要分為3 種情況,分別是求解線性解析解、求解非線性解析解和數值計算。其中求解線性解析解由于使用的是線性模型,沒有考慮到實際運行中的非線性關系,往往會產生較大的誤差。數值計算方法又可分為有限差分法和有限元法,在求解電機內部非線性電磁場時常用這種方法,這種方法的精度非常高,而且考慮的影響因素也非常全面,但計算過程會非常復雜繁瑣。求解非線性解析解的方法則剛好介于二者之間,有較好的精度,求解也相對容易。由于在水輪發電機組中不平衡磁拉力引起的機組振動有很強的非線性特性,并且結合實際對于計算精度的要求,通常都選用求解非線性解析解的方法[10]。

如圖1 所示,O'表示轉子外圓中心,在O-xyz坐標系中,設O'坐標為(X1,Y1),偏心轉子氣隙可表示為:

圖1 發電機偏心轉子氣隙示意圖

式(1)中,δ0表示發電機轉子不偏心時的平均氣隙寬度,表示發電機轉子徑向位移,α表示氣隙處于x軸夾角,γ表示發電機轉子轉角。

氣隙磁導的級數展開形式可以表示為:

式中,Λn表示氣隙磁導的Fourier 系數,具體表示為:

由電機學原理可知,三相同步發電機在空載工況下的氣隙基波磁動勢為:

式(4)中,Fj表示轉子勵磁電流的基波磁動勢,Fj=kjIj,Ij表示發電機轉子的勵磁電流,kj表示氣隙基波磁動勢系數,p表示發電機的磁極對數,ω表示轉動頻率。

氣隙磁密的分布表達式為:

在一般情況下,可將磁密的切向分量忽略不計。假設鐵心的磁導無限大,則可將垂直于鐵心或空氣邊界位置的Mɑxwell應力表示為:

將轉子表面積分可以得到不平衡磁拉力的解析表達式為:

從式(7)可知,不平衡磁拉力可以分為兩個部分:與時間無關部分和與時間相關部分。與時間無關部分的幅值大小為f1,且方向指向氣隙最小的位置;與時間相關部分的波動頻率是機組轉頻的二倍,幅值大小分別為f2、f3、f4。當磁極對數p>3 時,只存在與時間無關部分的不平衡磁拉力。對水輪發電機組而言,負載時會改變負載電流和功率因數,但對于不平衡磁拉力的負載與空載時差別并不大。因此,負載時的不平衡磁拉力可用空載磁場時得到的不平衡磁拉力來近似表示。水輪發電機組磁極對數一般較多,當p>3 時,機組不平衡磁拉力的解析表達式如下:

其中Rr表示發電機轉子半徑,Lr表示轉子長度。

1.2 碰摩力模型

工業上為了提高旋轉機械的效率,定轉子間的間隙設計往往很小,但這也導致有較大幾率會產生定轉子碰摩現象。碰摩力主要發生在定轉子間發生碰撞的情況下,而碰摩可能導致燒瓦現象,嚴重時更可能造成機組毀壞等破壞性事故,并且碰摩現象可能導致系統響應出現擬周期、混沌、分岔等非線性運動形式。碰摩工況主要有在一段時間內發生碰摩力連續但不光滑的碰摩和忽略碰摩細節,認為碰摩整個過程瞬間完成這兩種情況,這兩種碰摩過程有著不同的力學機理和力學意義。

目前,對于定轉子的碰摩非線性動力學問題主要有兩種模型進行處理,包括雙線性剛度模型和基于Hertz 接觸理論的連續碰摩模型。第一種模型主要針對的是剛體或者非常硬的碰撞體;第二種則有相對較大的應用范圍,可以考慮碰撞體的材料和結構等,也可以包含碰撞速度、阻尼和能量損失等信息,但建模過程也會相對更加復雜。由于水輪發電機組的轉速相對較低,且定轉子質量較大,可在忽略變形的條件下認為其保持剛性,所以在本文中綜合考慮,采用雙線性剛度模型[11-13]。

如圖2 所示,忽略摩擦熱效應的影響且認為定子為線性變形,定轉子間的摩擦符合庫倫定律。當發生碰摩過程,碰摩力可表示為:

圖2 定轉子碰摩示意圖

式(9)中,FN表示法向碰摩力,FT表示切向摩擦力,β表示徑向碰摩角,f表示摩擦系數,kr表示定子徑向剛度。

碰摩力在x,y方向的分量可以表示為:

式(10)中,H表示Heɑviside函數,且

1.3 非線性油膜力模型

非線性油膜力作用于軸承上,現今常用的非線性油膜力模型主要有Capone 模型、穩態短軸承油膜力模型、穩態長軸承油膜力模型和非穩態短軸承油膜力模型四種。其中第二種和第三種是不考慮瞬態擾動的,在本文中不適用,第四種和第一種相比,收斂性和精度都略有不足,因此,本文在處理非線性油膜力時采用Capone 模型[14-16]。

非線性油膜力的表達式為:

1.4 機組軸系動力學模型

本文研究的水輪發電機組軸系系統結構示意圖如圖3 所示。水輪發電機組軸系模型由發電機轉子和水輪機轉輪以及3 個導軸承構成。考慮到機組軸系系統在不平衡磁拉力、非線性油膜力和碰摩力作用下的動力特性較為復雜,本文只考慮軸系的橫向振動,并在轉子和聯軸器為剛性連接情況下進行分析。

圖3 機組軸系系統結構示意圖

忽略轉軸重量以及轉子重心回轉體的極慣性距,假設機組轉軸兩端的剛度均為Ke,上、下導軸承距轉子垂直方向距離均為ɑ。通過拉格朗日方程可得水輪發電機組軸系系統的運動微分方程為:

其中,m1、m2分別為轉子及軸頸質量,c1、c2分別為轉子及軸承處阻尼。X1、Y1為轉子的軸心位移,X2、Y2為軸承軸頸位移。

2 敏感性分析方法

在進行機組軸系系統參數的敏感性分析時,首先根據各個系統參數的取值范圍,分析系統參數在選定取值范圍內的取值變化對機組軸系系統輸出響應影響機理。目前,參數敏感性分析方法主要分為兩種,一種是局部敏感性分析,即在保持其他給定參數不變的條件下,分析選定的一個參數的變化對系統輸出響應的影響規律;另一種是全局敏感性分析,全局敏感性分析方法可以分析多個參數同步變化時,其中一個選定參數和其他參數之間相互作用對系統的輸出響應的影響程度[17]。

目前常用的全局敏感性分析方法有3 類,分別為:回歸分析法、Sobol 指數分析法、拓展傅里葉敏感性分析法(EFAST)[18-21]。其中,結合了傅里葉幅度敏感性檢驗法的優勢提出來的拓展傅里葉敏感性分析法不但可以分析系統某一個參數變化對輸出響應的影響,而且也可以用于多個參數同時變化交互作用時對系統輸出響應影響的分析。因此,本章采用拓展傅里葉敏感性分析法(EFAST)對水輪發電機組軸系系統參數進行全局敏感性分析。

EFAST 用的是方差分析的思想,即假定各個輸入參數和參數間的相互作用是由模型輸出的方差引起的[13,14]。因此可以反映出模型輸出響應對于輸入參數的敏感性。故可通過模型的方差的分解來獲取各個參數和參數相互耦合作用對總方差的影響的占比,稱為參數的敏感性指數。參數的敏感性指標可分為單個輸入參數對于模型輸出響應的影響和多個輸入參數間的互相耦合作用對于模型輸出響應的影響兩大類,我們分別稱其為主敏感性指數(Main effect)和總敏感性指數(Total effect)。在采用EFAST進行參數敏感性分析前,需先確定出各個輸入參數的取值范圍和幾何分布。

EFAST 的原理及計算過程如下所示:

首先,在所有待分析參數所組成的多維空間里選擇出合適的搜索函數Gi。然后用搜索函數Gi把模型Y=f(x1,x2,…,xn)轉化成Y=f(s)。其中搜索函數Gi與輸入參數xi的關系如下所示:

上式中,i∈(1,n),n表示不確定輸入參數的個數,ωi為搜索頻率。

搜索函數Gi的表達式為:

其中,Pi表示輸入參數xi的概率密度函數。

對式(14)進行積分和傅里葉變換可得:

其中,輸入參數xi通過xj作用耦合方差用Vij表示,輸入參數xi通過xj和xk作用的耦合方差用Vijk表示,輸入參數xi通過x1,x2,…,xn作用的耦合方差用xn表示。通過歸一化處理可獲得輸入參數xi的一階敏感性指數Si,Si的表達式為Si=Vi/V,輸入參數xi的總敏感性指數STi可以表示為STi=(V-V-i)/V,其中V-i表示不包括參數xi的所有參數方差和。

水輪發電機組軸系系統在運行的過程中會受到水力、機械和電磁等多方面因素互相耦合的影響。本文中挑選出了10 個軸系系統參數來進行敏感性分析,具體參數如表1 所示:機組旋轉角速度(ω),轉子偏心(e0),勵磁電流(Ij),軸承長度(Lb),均勻氣隙大小(δ0),轉子長度(Lr),軸承半徑(Rb),轉子半徑(Rr),軸承間隙(Cz),定子徑向剛度(kr)。

表1 軸系系統參數及其取值范圍

研究對象為式(12)中的水輪發電機組軸系系統模型,利用EFAST 分析上述10 個軸系系統參數的變化對軸系系統模型輸出響應的敏感性。令e=[x1,x2,x3,x4,x5,x6,x7,x8,x9,x10]為方差函數。其中,[x1,x2,x3,x4,x5,x6,x7,x8,x9,x10]=[ω,e0,Ij,Lb,δ0,Lr,Rb,Rr,Cz,kr]。

令Ns∈(-π,π),其中Ns=2Mωmax+1,M=6,ωmax是ωi的最大值。可得,y=f(x1,x2,…,x10)可以轉化為y=f(s)。

對函數f(s)進行傅里葉變換可得到:

首先給輸入參數xi分配一個較大的頻率ωi,然后再給其他輸入參數分配一組較小頻率{ωi},且參數頻率滿足ωi≥2M·max{ωi}。然后將頻域分為[1,M.max({ωi})]和[M.max{ωi}+1,(Ns-1)/2]兩部分。由此可得,輸入參數xi的總敏感性指數的計算公式如下所示:

3 軸系系統參數敏感性分析

3.1 轉子X/Y 軸方向參數敏感性分析

假設水輪發電機組軸系系統模型的輸入參數分布規律均為正態分布的假設,且根據模型機制通過尺度比例設置輸入參數方差。軸系系統模型輸出響應包含了轉子在X方向和Y方向的位移(X1,Y1)以及軸頸在X方向和Y方向的位移(X2,Y2)。輸入參數xi采用Monte Carlo 方法來采樣,采樣次數設為5 000 次。使用上述條件來計算,我們可以得出水輪發電機組軸系系統的10 個輸入參數對轉子X/Y 軸方向輸出響應的敏感性指數如圖4 所示。

圖4 軸系系統參數對轉子X 軸方向位移敏感性指數

由圖4 可得轉子X 軸方向位移的主敏感性和總敏感性分析結果。當以轉子X 軸方向位移為水輪發電機組軸系系統模型的輸出響應時,軸系系統參數主敏感性指數排名前4 分別為均勻氣隙大小(δ0)、機組旋轉角速度(ω)、定子徑向剛度(kr)、軸承長度(Lb),他們的平均值分別為0.31、0.09、0.08 和0.07,其他軸系系統參數的主敏感性指數均小于0.05。說明在這10 個軸系系統參數中,這4 個軸系系統參數對轉子X 軸方向位移的動態輸出響應的直接影響最為顯著,而其他6 個軸系系統參數的主敏感性指數均較小,表明這些系統參數對轉子X 軸方向位移影響較弱。另一方面,可得軸系系統參數的總敏感性指數排名前4 的分別為定子徑向剛度(kr),均勻氣隙大小(δ0),轉子半徑(Rr),機組旋轉角速度(ω),他們的平均值分別為0.57、0.46、0.41 和0.29,其他6 個軸系系統參數的總敏感性指數均小于0.20,說明這4 個軸系系統參數與其他參數相互間作用時對轉子X 軸方向的位移輸出響應有比較顯著的影響。

綜上所述,均勻氣隙大小(δ0)、機組旋轉角速度(ω)和定子徑向剛度(kr)這3 個系統參數單獨變化時會對轉子X 軸方向位移響應產生較大影響,其中均勻氣隙大小(δ0)對轉子X 軸方向位移響應影響最為顯著;當機組軸系系統在多個參數同時變化時,定子徑向剛度(kr)、均勻氣隙大小(δ0)、轉子半徑(Rr)和機組旋轉角速度(ω)這4 個系統參數可以通過參數間相互作用而對轉子X軸方向位移響應產生較大影響,其中定子徑向剛度(kr)產生的參數間相互作用最強。

由圖5 可得轉子Y 軸方向位移的主敏感性和總敏感性分析結果。當以轉子Y 軸方向位移為水輪發電機組軸系系統模型的輸出響應時,軸系系統參數主敏感性指數排名前2 位分別為轉子偏心(e0)和均勻氣隙大小(δ0),他們的平均值分別為0.89和0.06,其他軸系系統參數的主敏感性指數均小于0.01。說明在這10 個軸系系統參數中,這2 個軸系系統參數對轉子Y 軸方向位移的動態輸出響應直接影響最為顯著,而其他8 個軸系系統參數的主敏感性指數均較小,表明這些系統參數對轉子X 軸方向位移影響較弱。另一方面,可得軸系系統參數的總敏感性指數排名前2 位的仍為轉子偏心(e0)和均勻氣隙大小(δ0),平均值分別為0.92 和0.08,其他8個軸系系統參數的總敏感性指數均小于0.03,說明這2 個軸系系統參數與其他參數相互間作用時對轉子Y 軸方向的位移輸出響應有比較顯著的影響。

圖5 軸系系統參數對轉子Y 軸方向位移敏感性指數

綜上所述,轉子偏心(e0)和均勻氣隙大小(δ0)這2 個系統參數單獨變化時會對轉子Y 軸方向位移響應產生較大影響,其中轉子偏心(e0)對轉子Y 軸方向位移響應影響最為顯著;當機組軸系系統在多個參數同時變化時,也還是轉子偏心(e0)和均勻氣隙大小(δ0)這2 個系統參數可以通過參數間相互作用而對轉子Y 軸方向位移響應產生較大影響,其中轉子偏心(e0)產生的參數間相互作用最強。

3.2 軸頸X/Y 軸方向參數敏感性分析

同上所示,帶入計算過程,我們可以得到水輪發電機組軸系系統的10 個輸入參數對軸頸X/Y 軸方向輸出響應的敏感性指數如圖6 所示。

圖6 軸系系統參數對軸頸X 軸方向位移敏感性指數

由圖6 可得軸頸X 軸方向位移的主敏感性和總敏感性分析結果。當以軸頸X 軸方向位移為水輪發電機組軸系系統模型的輸出響應時,軸系系統參數主敏感性指數排名前5 位分別為軸承間距(Cz)、軸承長度(Lb)、軸承半徑(Rb)、均勻氣隙大小(δ0)和轉子偏心(e0),他們的平均值分別為0.35、0.29、0.11、0.06 和0.05,其他軸系系統參數的主敏感性指數均小于0.05。說明在這10 個軸系系統參數中,這5 個軸系系統參數對軸頸X 軸方向位移的動態輸出響應的直接影響最為顯著,而其他5 個軸系系統參數的主敏感性指數均較小,表明這些系統參數對軸頸X 軸方向位移影響較弱。另一方面,可得軸系系統參數的總敏感性指數排名前5 位的分別為軸承間距(Cz)、軸承長度(Lb)、軸承半徑(Rb)、均勻氣隙大小(δ0)和轉子偏心(e0),他們的平均值分別為0.37、0.33、0.13、0.07 和0.06,其他5 個軸系系統參數的總敏感性指數均小于0.01,說明這5 個軸系系統參數與其他參數相互間作用時對軸頸X 軸方向位移輸出響應有比較顯著的影響。

綜上所述,軸承間距(Cz)、軸承長度(Lb)、軸承半徑(Rb)、均勻氣隙大小(δ0)和轉子偏心(e0)這5 個系統參數單獨變化時會對軸頸X 軸方向位移響應產生較大影響,其中軸承間距(Cz)對軸頸X 軸方向位移響應影響最為顯著;當機組軸系系統在多個參數同時變化時,仍是上述5 個系統參數可通過參數間相互作用而對軸頸X 軸方向位移響應產生較大影響,其中軸承間距(Cz)產生的參數間相互作用最強。

由圖7 可得軸頸Y 軸方向位移的主敏感性和總敏感性分析結果。當以軸頸Y 軸方向位移為水輪發電機組軸系系統模型輸出響應時,軸系系統參數主敏感性指數排名前4 位分別為軸承長度(Lb)、軸承間距(Cz)、軸承半徑(Rb)和轉子偏心(e0),他們的平均值分別為0.51、0.23、0.16 和0.12,其他軸系系統參數的主敏感性指數均小于0.02。說明在這10個軸系系統參數中,這4 個軸系系統參數對軸頸Y軸方向位移的動態輸出響應的直接影響最為顯著,而其他6 個軸系系統參數的主敏感性指數均較小,表明這些系統參數對軸頸Y 軸方向位移影響較弱。另一方面,可得軸系系統參數的總敏感性指數排名前4 位仍為軸承長度(Lb)、軸承間距(Cz)、軸承半徑(Rb)和轉子偏心(e0),他們的平均值分別為0.55、0.24、0.18 和0.13,其他6 個軸系系統參數的總敏感性指數均小于0.02,說明這4 個軸系系統參數與其他參數相互間作用時對軸頸Y 軸方向位移輸出響應有比較顯著的影響。

圖7 軸系系統參數對軸頸Y 軸方向位移敏感性指數

綜上所述,軸承長度(Lb)、軸承間距(Cz)、軸承半徑(Rb)和轉子偏心(e0)這4 個系統參數單獨變化時會對軸頸Y 軸方向位移響應產生較大的影響,其中軸承長度(Lb)對軸頸Y 軸方向位移響應影響最為顯著;當機組軸系系統在多個參數同時變化時,仍是上述4 個系統參數可以通過參數間相互作用而對軸頸Y 軸方向位移響應產生較大影響,其中軸承長度(Lb)產生的參數間相互作用最強。

4 結論

為了全面系統地分析機組軸系系統參數對轉子中心和軸頸中心橫向位移的影響程度,在機組軸系系統動力學模型的基礎上,應用拓展傅里葉敏感性分析方法,推導并計算出軸系系統參數的全局敏感性。現將結果總結如下:

(1)基于所建立的水輪發電機組軸系系統動力學模型,通過給出軸系系統10 個模型參數取值范圍和分布特征,利用拓展傅里葉敏感性分析方法,對上述10 個系統參數進行了全局敏感性分析,不僅揭示某一參數單一變化時對系統輸出影響程度,而且給出了多參數同時變化時參數間相互作用效應對系統輸出影響敏感性指數。

(2)均勻氣隙大小、機組旋轉角速度和定子徑向剛度這3 個系統參數單獨變化時會對轉子X 軸方向位移響應產生較大影響,其中均勻氣隙大小對轉子X 軸方向位移響應影響最為顯著;當機組軸系系統在多個參數同時變化時,定子徑向剛度、均勻氣隙大小、轉子半徑和機組旋轉角速度這4 個系統參數可以通過參數間相互作用而對轉子X 軸方向位移響應產生較大影響,其中定子徑向剛度產生的參數間相互作用最強。轉子偏心和均勻氣隙大小這2 個系統參數單獨變化時會對轉子Y 軸方向位移響應產生較大影響,其中轉子偏心對轉子Y 軸方向位移響應影響最為顯著;當機組軸系系統在多個參數同時變化時,也還是轉子偏心和均勻氣隙大小這2 個系統參數可以通過參數間相互作用而對轉子Y 軸方向位移響應產生較大影響,其中轉子偏心產生的參數間相互作用最強。

(3)軸承間距、軸承長度、軸承半徑、均勻氣隙大小和轉子偏心這5 個系統參數單獨變化時會對軸頸X 軸方向位移響應產生較大影響,其中軸承間距對軸頸X 軸方向位移響應影響最為顯著;當機組軸系系統在多個參數同時變化時,仍是上述5 個系統參數可通過參數間相互作用而對軸頸X 軸方向位移響應產生較大影響,其中軸承間距產生的參數間相互作用最強。軸承長度、軸承間距、軸承半徑和轉子偏心這4 個系統參數單獨變化時會對軸頸Y 軸方向位移響應產生較大的影響,其中軸承長度對軸頸Y 軸方向位移響應影響最為顯著;當機組軸系系統在多個參數同時變化時,仍是上述4 個系統參數可以通過參數間相互作用而對軸頸Y 軸方向位移響應產生較大影響,其中軸承長度產生的參數間相互作用最強。

猜你喜歡
方向系統
Smartflower POP 一體式光伏系統
工業設計(2022年8期)2022-09-09 07:43:20
2022年組稿方向
計算機應用(2022年2期)2022-03-01 12:33:42
2022年組稿方向
計算機應用(2022年1期)2022-02-26 06:57:42
2021年組稿方向
計算機應用(2021年4期)2021-04-20 14:06:36
2021年組稿方向
計算機應用(2021年3期)2021-03-18 13:44:48
WJ-700無人機系統
2021年組稿方向
計算機應用(2021年1期)2021-01-21 03:22:38
ZC系列無人機遙感系統
北京測繪(2020年12期)2020-12-29 01:33:58
基于PowerPC+FPGA顯示系統
半沸制皂系統(下)
主站蜘蛛池模板: 欧美激情一区二区三区成人| av午夜福利一片免费看| 亚洲视频四区| 日韩 欧美 小说 综合网 另类| 国产一区免费在线观看| 丝袜久久剧情精品国产| av午夜福利一片免费看| 九九热免费在线视频| 美女被操黄色视频网站| 欧美三级不卡在线观看视频| 美女被躁出白浆视频播放| 无码区日韩专区免费系列| 91亚洲国产视频| 日韩无码视频专区| 五月激激激综合网色播免费| 亚洲熟女中文字幕男人总站| 红杏AV在线无码| 国产一级一级毛片永久| 午夜少妇精品视频小电影| 香蕉国产精品视频| 91亚洲精选| 亚洲精品国产成人7777| 成人福利在线观看| 中文字幕亚洲乱码熟女1区2区| 欧美综合激情| 国产视频自拍一区| 9966国产精品视频| 国产成人艳妇AA视频在线| 亚洲毛片在线看| 欧美成人午夜影院| 国产成人91精品| 日本成人不卡视频| 国产在线拍偷自揄拍精品| 国产女人水多毛片18| 日韩精品无码免费专网站| 亚洲全网成人资源在线观看| 亚洲综合片| 99偷拍视频精品一区二区| 国产爽妇精品| 午夜毛片免费观看视频 | 亚洲综合经典在线一区二区| 亚洲精品777| 久久国产精品娇妻素人| 国产97区一区二区三区无码| 亚洲无线一二三四区男男| 欧美日韩资源| 亚洲日韩精品无码专区| 天堂亚洲网| 黄色一级视频欧美| 国产亚洲精品97AA片在线播放| 国模粉嫩小泬视频在线观看| 国产香蕉在线视频| 手机看片1024久久精品你懂的| 国产精品亚洲va在线观看 | 成年人国产网站| 在线观看亚洲精品福利片 | 亚洲三级成人| 亚洲精品无码抽插日韩| 国产成人做受免费视频| 国产综合无码一区二区色蜜蜜| 日韩123欧美字幕| 国产麻豆aⅴ精品无码| 在线免费无码视频| 为你提供最新久久精品久久综合| 久久综合结合久久狠狠狠97色| 小蝌蚪亚洲精品国产| 婷婷六月综合网| 激情成人综合网| 天堂在线www网亚洲| 国产一区二区影院| 国产精品污视频| 伊人国产无码高清视频| 国产精品成人久久| 亚洲伊人电影| 黄片在线永久| 国产综合色在线视频播放线视| www亚洲天堂| 国产色图在线观看| 57pao国产成视频免费播放| 国产噜噜在线视频观看| 国产午夜人做人免费视频中文| 亚洲国产精品成人久久综合影院|