亚洲免费av电影一区二区三区,日韩爱爱视频,51精品视频一区二区三区,91视频爱爱,日韩欧美在线播放视频,中文字幕少妇AV,亚洲电影中文字幕,久久久久亚洲av成人网址,久久综合视频网站,国产在线不卡免费播放

        ?

        具有不確定性連接剛度的太陽電池陣基頻分析

        2018-01-08 07:25:35龐夢非朱春艷張美艷唐國安
        上海航天 2017年6期
        關(guān)鍵詞:基頻太陽電池鉸鏈

        龐夢非,朱春艷,張美艷,唐國安

        (1.復(fù)旦大學(xué) 航空航天系,上海 200433; 2.上海宇航系統(tǒng)工程研究所,上海 201109)

        具有不確定性連接剛度的太陽電池陣基頻分析

        龐夢非1,朱春艷2,張美艷1,唐國安1

        (1.復(fù)旦大學(xué) 航空航天系,上海 200433; 2.上海宇航系統(tǒng)工程研究所,上海 201109)

        以鉸鏈副抗彎剛度具不確定性的太陽電池陣為研究對象,考慮加工和實驗誤差等隨機因素,用數(shù)值計算方法仿真電池陣的基頻實驗,對安裝在同一電池陣上的不同組、抗彎剛度滿足一定分布規(guī)律的鉸鏈副進行多次仿真實驗,建立了鉸鏈副抗彎剛度和電池陣基頻的仿真樣本,用高斯過程回歸方法訓(xùn)練出兩者間的映射關(guān)系,獲得了電池陣基頻期望值及置信區(qū)間關(guān)于鉸鏈副抗彎剛度的變化關(guān)系。給出了電池陣基頻關(guān)于鉸鏈副剛度參數(shù)不確定性的分析流程。分析結(jié)果表明:在鉸鏈副剛度不確定條件下,高斯過程回歸方法可確定電池陣基頻期望值和協(xié)方差表達式中的超參數(shù),與其他回歸或插值方法等相比,不僅能預(yù)示電池陣基頻的期望值,而且可估計其置信區(qū)間。在加工和實驗誤差不可避免的實際情況下,有必要生產(chǎn)多組鉸鏈副進行模態(tài)實驗并對同一組鉸鏈副適當增加重復(fù)實驗次數(shù),使其能合理預(yù)示電池陣基頻的不確定性,減小實驗誤差的影響,預(yù)示合理的基頻。

        太陽電池陣; 鉸鏈副; 仿真實驗; 剛度; 固有頻率; 不確定性; 高斯過程回歸; 置信區(qū)間

        0 引言

        不確定性問題廣泛存在于工程結(jié)構(gòu)和裝備的裝配、制造中,如受加工工藝限制,結(jié)構(gòu)的剛度、質(zhì)量和阻尼等力學(xué)參數(shù)總表現(xiàn)出一定的隨機性。正確估計參數(shù)不確定性對結(jié)構(gòu)力學(xué)特性和行為的影響有重要的工程應(yīng)用價值。模糊方法、區(qū)間方法和概率方法是目前不確定性建模的三種主要方法。用模糊集合理論描述不確定性,由此發(fā)展出模糊有限元方法。文獻[1]用模糊方法對薄壁組合梁進行了不確定性分析。文獻[2]求解了具不確定參數(shù)的結(jié)構(gòu)動力學(xué)分析的模糊特征值求解問題。用基于模糊集理論的方法描述不確定性時,需獲知不確定性結(jié)構(gòu)參數(shù)的隸屬度函數(shù)。多數(shù)情況下,確定隸屬度函數(shù)非常困難,使用者常不得不主觀性地選取相應(yīng)的隸屬度函數(shù)。這樣的分析結(jié)果可靠性并無保障。區(qū)間分析方法是在信息不充分條件下,描述工程問題不確定性的一種方法,將這些不確定性結(jié)構(gòu)參數(shù)視為未知變量在具已知邊界的區(qū)間內(nèi)取值。文獻[3]較早將區(qū)間分析法引入工程結(jié)構(gòu)中不確定性問題的研究。文獻[4-5]用區(qū)間有限元方法計算了不確定系統(tǒng)的頻響函數(shù)包絡(luò),提出了不確定結(jié)構(gòu)頻響函數(shù)包絡(luò)分析的區(qū)間靈敏度理論。文獻[6]綜合了區(qū)間分析法與有限元方法,求解傳遞函數(shù)的包絡(luò)。區(qū)間分析方法無需事先獲知不確定參數(shù)的很多信息,只需知道不確定參數(shù)的上下限。對直接采用區(qū)間數(shù)學(xué)運算法則求解線性區(qū)間方程組的方法進行了研究,但該法難以用于求解大規(guī)模問題,更重要的是其預(yù)測的范圍大于實際值,結(jié)果過于保守。文獻[7-9]最先引入貝葉斯概率框架估計不確定性結(jié)構(gòu)模型參數(shù),貝葉斯概率框架在解決模型修正過程中不確定性的效果顯著,但由于計算量過大,在工程實際中的應(yīng)用受限。文獻[10]綜合了攝動和有限元法,研究了一般概率攝動有限元方法,并用隨機攝動法分析了有不確定參數(shù)和不確定載荷的振動機械系統(tǒng)響應(yīng)[11]。文獻[12]用攝動有限元法對不確定結(jié)構(gòu)靜動態(tài)特性的穩(wěn)健優(yōu)化設(shè)計進行了研究。高斯過程回歸方法建立在貝葉斯概率理論的基礎(chǔ)上,在機器學(xué)習(xí)等領(lǐng)域中發(fā)揮了重要作用,但目前用其解決工程中不確定性問題的研究較少。本文將采用高斯過程回歸方法,探討其在解決工程不確定性問題中的適用性。

        本文以鉸鏈副抗彎剛度具不確定性的太陽電池陣為研究對象,用數(shù)值計算方法仿真電池陣的基頻實驗,考慮加工和實驗誤差等隨機因素,對安裝在同一電池陣上的不同組、抗彎剛度滿足一定分布規(guī)律的鉸鏈副進行多次仿真實驗,建立鉸鏈副抗彎剛度和電池陣基頻的仿真樣本,通過高斯過程回歸方法訓(xùn)練兩者間的映射關(guān)系。在此基礎(chǔ)上,希望得到電池陣基頻期望值及置信區(qū)間關(guān)于鉸鏈副抗彎剛度的變化關(guān)系,使給定一個抗彎剛度實驗值就能預(yù)測出電池陣基頻的期望值,并對方差作出估計。

        1 單輸出高斯過程回歸基本原理

        高斯過程回歸問題可描述為:記x1,…,xn∈Rd為n個輸入向量;X=[x1…xn]∈Rd×n為輸入矩陣;y1,…,yn∈R1為n個輸出標量;y=[y1…yn]T∈Rn為輸出向量。此處:Rd,Rd×n分別為d維和d×n維線性空間。再令D={X,y}為已存在的訓(xùn)練樣本集,高斯過程回歸的基本任務(wù)是:根據(jù)樣本集D,學(xué)習(xí)輸入X與輸出y間的映射關(guān)系f(x):Rd→R1,使對新的輸入向量x*∈Rd能預(yù)測出對應(yīng)的映射值f(x*)∈R1。

        對單輸出高斯回歸問題,考慮以下模型

        yi=f(xi)+εi

        式中:f(xi)為當輸入為xi時對應(yīng)的隨機變量;εi為滿足N(0,(σn)2)的獨立同分布的噪聲;yi為f(xi)受到噪聲污染后的觀察值;i=1,2,…,n[13]。此處:N(0,(σn)2)表示均值為0、方差為(σn)2的高斯分布。

        因高斯過程的任意有限維分布均滿足高斯分布,故可得觀察值y的先驗分布為

        y~N(0,C)

        (1)

        即y滿足均值為0、方差為C的高斯分布[14]。此處:C=K(X,X)+(σn)2In;K(X,X)∈Rn×n為對稱正定協(xié)方差矩陣。

        觀察值y與預(yù)測值f*=f(x*)的聯(lián)合先驗分布為

        (2)

        式中:K(X,x*)∈Rn×1為預(yù)測點x*與訓(xùn)練樣本集的輸入X間的n×1階協(xié)方差矩陣;K(x*,X)=(K(X,x*))T;k(x*,x*)為預(yù)測點x*自身的協(xié)方差。

        由式(2)可得預(yù)測值f*的后驗分布為

        (3)

        式中:

        (4)

        cov(f*)=k(x*,x*)-K(x*,X)·

        C-1K(X,x*)

        (5)

        則可得以概率形式表示的預(yù)測結(jié)果。

        高斯過程的性質(zhì)取決于其均值函數(shù)和協(xié)方差函數(shù)。常用的協(xié)方差函數(shù)有平方指數(shù)協(xié)方差函數(shù)、Matern類協(xié)方差函數(shù)和神經(jīng)網(wǎng)絡(luò)協(xié)方差函數(shù)。平方指數(shù)協(xié)方差的形式為

        k(xp,xq)=

        (6)

        式中:xp,xq為任意兩個d維向量;(σf)2=covf為隨機變量f的方差,表示隨機變量偏離期望值的程度。若令M=λ-2I,則式(6)可簡化為

        (7)

        式中:λ為輸入量的特征距離。協(xié)方差中的未知數(shù)σf,λ被稱為超參數(shù)。單輸出高斯過程回歸方法中的超參數(shù)還包括噪聲的標準差σn。令θ為預(yù)測模型中所有超參數(shù)組成的向量,對式(7)表示的協(xié)方差函數(shù),θ=[λ(σf)2(σn)2],一般可用極大似然法估計超參數(shù)的初始值[14]。同時,根據(jù)概率論中的貝葉斯原理,可得超參數(shù)的概率表達形式并將問題轉(zhuǎn)換為對其極大值的求解,有

        (8)

        設(shè)θ接近于均勻分布,即概率密度p(θ)可近似為常量,則尋找函數(shù)p(θ|y,X)的最大值等價于尋找p(y|X,θ)的最大值問題

        (9)

        優(yōu)化設(shè)計理論證明了梯度下降法速度較快的優(yōu)點,但其缺點是易陷入局部極值點。為此,文獻[15]提出了一種超參數(shù)初始值的優(yōu)化方法。本文采取了該優(yōu)化方法,可得到較好的結(jié)果。

        2 電池陣基頻關(guān)于鉸鏈副剛度參數(shù)不確定性分析

        航天器上使用的電池陣很多是由若干太陽電池板通過壓緊-鎖定的鉸鏈副連接的,航天器姿態(tài)控制律設(shè)計等需要電池陣的基頻數(shù)據(jù)。鉸鏈副是一種機構(gòu),因加工、裝配等各種因素的差異,對電池陣整體結(jié)構(gòu)提供的連接剛度存在相當大的不確定性?;l的實驗測試結(jié)果也必定存在隨機性,最終得到的基頻數(shù)據(jù)必具有不確定性?;谝延械你q鏈副剛度參數(shù)和電池陣基頻實測結(jié)果,預(yù)示一個給定鉸鏈副剛度實驗值下的電池陣基頻并估計其不確定性有重要的工程應(yīng)用價值。

        2.1 太陽電池陣模型及其不確定性分析過程

        為說明基頻及其不確定性預(yù)示的過程和結(jié)果,本文用數(shù)值計算方法仿真電池陣的基頻實驗。電池陣的基頻實驗如圖1(a)所示。實驗時用有機玻璃板替代實際的電池板,局部放大的是連接兩塊電池板的鉸鏈副[16]。仿真所用電池板的有限元模型如圖1(b)所示,模型中鉸鏈副等效為梁單元。在此狀態(tài)下,基頻對應(yīng)的模態(tài)是垂直于XY平面的彎曲。有機玻璃板及等效梁模型的屬性參考了文獻[17]。

        電池陣基頻關(guān)于鉸鏈副剛度參數(shù)的不確定性分析流程為:

        a)設(shè)置多組鉸鏈副剛度系數(shù),構(gòu)造輸入矩陣X;

        b)電池陣基頻測試實驗仿真,獲得輸出向量y;

        e)預(yù)示給定鉸鏈副剛度實驗值下的電池陣基頻并估計其不確定性。

        2.2 高斯過程回歸方法分析結(jié)果

        僅討論單輸入情形,連接有機玻璃板的一對鉸鏈副具有相同的剛度。基頻對應(yīng)的模態(tài)主要受到繞y軸的抗彎剛度的影響,將鉸鏈副等效為梁單元后即為繞y軸的抗彎剛度EI,文獻[17]用實驗方法辨識得到EI=316 N·m2。由于加工和裝配工差等因素,即使是同一生產(chǎn)批次的不同鉸鏈副,其抗彎剛度也存在差異。這種差異具隨機性,根據(jù)加工誤差統(tǒng)計特性可認為其服從正態(tài)分布。為模擬此情形,將文獻[17]的實測結(jié)果EI≈316 N·m2作為期望值,取標準差為60 N·m2,正態(tài)分布的三個抗彎剛度作為對三個鉸鏈副進行仿真實驗,由得到的抗彎剛度值可構(gòu)造出輸入矩陣

        對上述三個鉸鏈副的抗彎剛度樣本分別進行基頻實驗。仿真時用有限元計算程序依次計算出鉸鏈副剛度分別為X中各分量時對應(yīng)的電池陣基頻,即基頻分別為

        f= [f1f2f3]=

        [9.139 5 9.677 5 9.878 9] Hz

        為模擬實驗過程存在的各種誤差,將f中的每個元素疊加服從高斯分布的隨機向量ε~N(0,0.12),作為由仿真實驗得到的輸出樣本

        y= [y1y2y3]=

        [9.249 5 9.832 0 9.887 4] Hz

        圖2中:藍線即為高斯過程回歸模型的訓(xùn)練結(jié)果;綠線為理論曲線;紅線包圍的部分為99%置信區(qū)間,說明將一個已測得剛度參數(shù)的鉸鏈副安裝在電池陣上進行模態(tài)實驗,電池陣的基頻有99%的可能性處于紅色線范圍內(nèi)。由圖2可知:所得結(jié)果曲線受采樣點影響很大,與理論曲線趨勢相差較大,卻給出了很小的置信區(qū)間,這與實際情況不符。

        為避免實驗誤差導(dǎo)致的不合理預(yù)測結(jié)果,用相同方法構(gòu)造出輸入矩陣

        X=[EI1EI2EI3]=

        [346.648 5 276.623 6 308.499 0] N·m2

        獲得鉸鏈副抗彎剛度值后進行基頻仿真實驗,每個抗彎剛度樣本分別進行基頻實驗5次。仿真時,用有限元計算程序依次計算鉸鏈副剛度分別為X中各分量時對應(yīng)的電池陣基頻,所得基頻分別為

        f= [f1f2f3]=

        [9.863 4 9.366 6 9.611 3] Hz

        同樣,為模擬實驗過程存在的各種誤差,對f中的每個元素疊加服從高斯分布的隨機向量ε~N(0,0.12),每個點實驗5次即疊加5次。記由仿真實驗得到的輸出樣本為y=[y1…y15]。此處:y1=9.810 4 Hz;y2=9.874 0 Hz;y3=9.976 3 Hz;y4=9.937 7 Hz;y5=9.977 8 Hz;y6=9.275 1 Hz;y7=9.384 6 Hz;y8=9.268 3 Hz;y9=9.405 1 Hz;y10=9.399 2 Hz;y11=9.740 9 Hz;y12=9.721 2 Hz;y13=9.676 6 Hz;y14=9.560 8 Hz;y15=9.563 7 Hz。其中y1~y5,y6~y10,y11~y15分別對應(yīng)同一個抗彎剛度值。

        不同抗彎剛度下兩塊有機玻璃板和一對鉸鏈副組成的電池陣的基頻如圖3所示。

        由圖3可知:計算曲線與理論曲線非常接近,預(yù)測效果較好;模型還能給出樣本點輸出的置信區(qū)間。圖3中疊加了計算曲線的99%置信區(qū)間(紅色線包圍部分)。本文研究結(jié)果表明:對此單輸出高斯過程回歸問題,若每個抗彎剛度值對應(yīng)只進行1次基頻實驗,則實驗誤差可能會導(dǎo)致結(jié)果曲線與理論曲線偏差較大。經(jīng)過多次的仿真和計算,發(fā)現(xiàn)在同一個抗彎剛度值下進行5次實驗就可獲得較好的預(yù)測效果。

        3 結(jié)束語

        本文對具不確定性連接剛度的太陽電池陣基頻仿真分析方法進行了研究。高斯過程回歸方法能通過少量的樣本點訓(xùn)練出大量的預(yù)測值,將其用于工程實際問題的不確定性分析,可減少試驗工況?;谔栯姵仃嚹B(tài)實驗的仿真結(jié)果表明:在鉸鏈副剛度不確定條件下,高斯過程回歸方法可確定電池陣基頻期望值和協(xié)方差表達式中的超參數(shù),與其他回歸或插值方法等相比,不僅能預(yù)示電池陣基頻的期望值,并且能估計其置信區(qū)間。分析結(jié)果表明:太陽電池陣的同一鉸鏈副必須進行多次模態(tài)實驗,才能減少實驗誤差造成的影響、預(yù)示出合理的基頻結(jié)果。對真實的模態(tài)實驗,應(yīng)用高斯過程回歸方法預(yù)示太陽電池陣基頻有待進一步研究。如實際情況中,同一電池陣中不同鉸鏈副的剛度存在差異,要求預(yù)示的頻率不僅限于基頻,這些多輸入多輸出問題也需拓展。

        [1] PAWAR P M, NAM JUNG S, RONGE B P. Fuzzy approach for uncertainty analysis of thin walled composite beams[J]. Aircraft Engineering and Aerospace Technology, 2012, 84(1): 13-22.

        [2] XIA Y, FRISWELL M. Efficient solution of the fuzzy eigenvalue problem in structural dynamics[J]. Engineering Computations, 2014, 31(5): 864-878.

        [3] RAO S S, BERKE L. Analysis of uncertain structural systems using interval analysis[J]. AIAA Journal, 1997, 35(4): 727-735.

        [4] MOENS D, VANDEPITTE D. An interval finite element approach for the calculation of envelope frequency response functions[J]. International Journal for Numerical Methods in Engineering, 2004, 61(14): 2480-2507.

        [5] MOENS D, VANDEPITTE D. Interval sensitivity theory and its application to frequency response envelope analysis of uncertain structures[J]. Computer Methods in Applied Mechanics and Engineering, 2007, 196(21): 2486-2496.

        [6] DESSOMBZ O, THOUVEREZ F, LANé J P, et al. Analysis of mechanical systems using interval computations applied to finite element methods[J]. Journal of Sound and Vibration, 2001, 239(5): 949-968.

        [7] BECK J L, KATAFYGIOTIS L S. Updating models and their uncertainties I: Bayesian statistical framework[J]. Journal of Engineering Mechanics, 1998, 124(4): 455-461.

        [8] KATAFYGIOTIS L S, BECK J L. Updating models and their uncertainties. II: model identifiability[J]. Journal of Engineering Mechanics, 1998, 124(4): 463-467.

        [9] KATAFYGIOTIS L S, PAPADIMITRIOU C, LAM H F. A probabilistic approach to structural model updating[J]. Soil Dynamics and Earthquake Engineering, 1998, 17(7): 495-507.

        [10] 張義民, 劉巧伶, 聞邦椿. 多自由度非線性隨機參數(shù)振動系統(tǒng)響應(yīng)分析的概率攝動有限元法[J]. 計算力學(xué)學(xué)報, 2003, 20(1): 8-11.

        [11] 張義民, 聞邦椿. 非線性隨機振動機的動力學(xué)研究[J]. 工程力學(xué), 1999, 16(3): 100-104.

        [12] 亢戰(zhàn), 程耿東. 非確定性結(jié)構(gòu)靜動態(tài)特性穩(wěn)健優(yōu)化設(shè)計[J]. 力學(xué)學(xué)報, 2006, 23(1): 57-65.

        [13] RASMUSSEN C E, WILLIAMS C K I. Gaussian process for machine learning[M]. Cambridge: MIT Press, 2006: 8-24.

        [14] 盛驟, 謝式千, 潘承毅. 概率論與數(shù)理統(tǒng)計[M]. 北京: 高等教育出版社, 2008: 14-20+46+152-156.

        [15] 單志輝. 基于高斯過程回歸的翼型快速設(shè)計研究[D]. 南京: 南京航空航天大學(xué), 2011.

        [16] 馬越峰, 江禮俐, 蘇里, 等. 太陽電池陣的鉸鏈副剛度參數(shù)識別[J]. 上海航天, 2006, 23(3): 23-27.

        [17] 馬越峰. 振動結(jié)構(gòu)的參數(shù)識別與模型修正及其應(yīng)用[D]. 上海: 復(fù)旦大學(xué), 2005.

        [18] MAGRAB E B. An engineer's guide to MATLAB: with applications from mechanical, aerospace, electrical, civil, and biological systems engineering[M]. Englewood Cliffs: Prentice Hall, 2011: 710-712.

        FundamentalFrequencyAnalysisofSolarCellArraywithUncertainConnectionStiffness

        PANG Meng-fei1, ZHU Chun-yan2, ZHANG Mei-yan1, TANG Guo-an1

        (1. Department of Aeronautics and Astronautics, Fudan University, Shanghai 200433, China;2. Aerospace System Engineering Shanghai, Shanghai 201109, China)

        With the consideration of the random factors such as machining and experimental error, solar cell array with a hinge pair which had uncertain bending stiffness was taken as an example and the fundamental frequency test was simulated by numerical method in this paper. The different hinge pairs and their bending stiffness satisfying certain distribution law which were amounted on the same solar cell were simulated many times. The simulation model of the hinge pair bending stiffness and the cell array fundamental frequency were established. The mapping between them was trained by the Gaussian process regression method. On this basis, the relationship between the expected value and confidence interval of the cell array fundamental frequency and the hinge pair bending stiffness was obtained. The uncertainty analysis process of cell array fundamental frequency with the respect to the stiffness parameters of hinge pair was given. The results show that the fundamental frequency expected value of the cell and super parameters in the covariance expression can be determined by the Gaussian process regression methods under the condition of uncertainty of hinge pair stiffness. It can not only indicate the expected value of the fundamental frequency but also estimate its confidence interval, which is different from other regression method or interpolating method. In order to reasonably predict the uncertainty of the cell array fundamental frequency, it is necessary to produce multiple sets of hinge pairs to do modal experiments and to increase the number of repetitive experiments on the same hinge pair in the actual case where the processing and experimental errors are inevitable.

        solar cell array; hinge pair; simulation; stiffness; natural frequency; uncertainty; Gaussian process regression; confidence interval

        1006-1630(2017)06-0103-06

        V414.33

        A

        10.19328/j.cnki.1006-1630.2017.06.016

        2017-05-20;

        2017-07-12

        國家自然科學(xué)基金項目資助(11572089);上海市科學(xué)技術(shù)委員會揚帆計劃項目資助(15YF1411900)

        龐夢非(1993—),女,碩士生,主要研究方向為結(jié)構(gòu)動力學(xué)。

        唐國安(1962—),男,教授,主要研究方向為結(jié)構(gòu)動力學(xué)。

        猜你喜歡
        基頻太陽電池鉸鏈
        語音同一認定中音段長度對基頻分析的影響
        基于時域的基頻感知語音分離方法?
        橋面鋪裝層對中小跨徑橋梁基頻影響分析
        基于虛擬鉸鏈打開機構(gòu)的艙門提升機構(gòu)研究
        球鉸鏈防塵罩抱緊力優(yōu)化
        北京汽車(2019年6期)2019-02-12 05:19:50
        汽車連接器帶鉸鏈護殼產(chǎn)品的塑料模具設(shè)計改進
        汽車電器(2018年10期)2018-11-01 01:17:26
        45000kHz基頻晶體濾波器
        電子制作(2017年20期)2017-04-26 06:57:35
        幾種新型鈣鈦礦太陽電池的概述
        光對聚合物太陽電池的影響
        柔性砷化鎵太陽電池
        亚洲日本中文字幕乱码在线| 日韩狼人精品在线观看| 亚洲国产一区二区视频| 亚洲国产黄色在线观看| 一区二区三区激情免费视频| 欧美色欧美亚洲另类二区| 大地资源中文第三页| 国产成人精品蜜芽视频| 天天色天天操天天日天天射| 亚洲精品久久久av无码专区| 在线播放国产一区二区三区| 人妻少妇中文字幕久久69堂| 国产片在线一区二区三区| 亚洲av日韩av天堂久久| 国产精品美女久久久久久2018| 最新福利姬在线视频国产观看| 亚洲女同人妻在线播放| 亚洲国产av自拍一区| 精品无码久久久久成人漫画 | 国模私拍福利一区二区| 蜜桃网站在线免费观看视频| 日韩av毛片在线观看| 日本japanese丰满多毛| 国产一区二区三区啪| 日韩精品一区二区三区av| 亚洲精品一区二区三区麻豆| 18禁无遮挡羞羞污污污污网站| 欧美日韩国产专区| 国产av一区二区内射| 无码精品一区二区三区在线| 三级4级全黄60分钟| 国产亚洲精品国看不卡| 日本精品一区二区三区在线观看| 亚洲熟妇久久国产精品| 欧美a级在线现免费观看| 亚洲va精品va国产va| 国产一区二区av免费在线观看| 狠狠人妻久久久久久综合蜜桃| 久久精品国产99国产精2020丨 | 中文字幕日本女优在线观看| 亚洲综合日韩精品一区二区|