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

基于改進遺傳算法的磁流變阻尼器多目標空間 優(yōu)化布置

2023-06-03 07:58:44張香成徐宏輝趙軍楊洋
湖南大學學報(自然科學版) 2023年5期
關鍵詞:優(yōu)化結構

張香成,徐宏輝,趙軍,楊洋

(1.鄭州大學 力學與安全工程學院,河南 鄭州 450001;2.鄭州大學 水利與土木工程學院,河南 鄭州 450001)

磁流變 阻尼器(Magnetorheological Damper,MRD)作為結構振動控制中最具有前景的半主動控制裝置之一,因結構簡單、響應迅速、阻尼力大且連續(xù)順逆可調(diào)等特點而廣泛應用于土木結構減震控制中[1-5].阻尼器對結構的減震控制效果不僅與阻尼器的數(shù)量、阻尼力大小有關,而且與阻尼器在結構中的位置密切相關,為了降低控制成本并且提高控制效率,需要對結構中的阻尼器進行優(yōu)化配置.

1997 年Takewaki[6]和Wu 等[7]在結構振動控制領域引入遺傳算法,并進行了大量的研究工作,至此,遺傳算法在結構阻尼器優(yōu)化布置中得到普遍應用.針對具體的工程優(yōu)化問題,選用適當?shù)木幋a和相應的遺傳算子進行優(yōu)化,是遺傳算法的一個重要研究方向.

很多學者采用遺傳算法“0-1”編碼方式對阻尼器進行優(yōu)化布置研究,0 表示該處不設置阻尼器,1表示設置阻尼器.貝偉明等[8]基于遺傳算法采用等效二次型性能指標對MRD 進行優(yōu)化布置,證明了遺傳算法的有效性;閤東東等[9]采用H2范數(shù)控制理論和改進遺傳算法對相鄰結構間阻尼器布置位置進行優(yōu)化研究,結果表明該方法可有效減小結構的二次型性能指標;孫彤等[10]基于遺傳算法對軌道式負剛度裝置提出一種優(yōu)化布置數(shù)學模型,并以10 層鋼筋混凝土(Reinforced Concrete,RC)結構為例研究了該裝置優(yōu)化布置基本原則;燕樂緯等[11]提出一種相對適應度函數(shù),提高了遺傳算法的進化效率;金波等[12]基于遺傳算法將黏滯阻尼器替換桿件的模態(tài)應變能作為優(yōu)化函數(shù),實現(xiàn)了對網(wǎng)架結構的優(yōu)化控制.“0-1”編碼多將結構簡化為層間剪切模型,忽略了阻尼器在結構平面和空間上分布方式的影響;對于大型空間結構,若實現(xiàn)精準定位則會造成“維度災難”,計算量和計算時間難以承受.

實數(shù)編碼是用具體實數(shù)表示每層阻尼器數(shù)量,郭勇等[13]采用實數(shù)編碼方式將遺傳算法與劣出優(yōu)入算法相結合,對輸電塔阻尼器的布置方式進行了多目標優(yōu)化研究.數(shù)字序列編碼中每個基因位數(shù)字表示樓層號,基因位數(shù)表示阻尼器個數(shù).燕樂緯等[14]基于遺傳算法提出數(shù)字序列編碼,并驗證了該方法的有效性;馬宏偉等[15]提出粗粒度并行遺傳算法,提高了優(yōu)化結果收斂速度及精度.實數(shù)編碼和數(shù)字序列編碼大大縮短了基因維度,解決了二進制編碼不能完備表達求解空間的問題,但無法實現(xiàn)阻尼器的精準定位.

為實現(xiàn)空間結構中阻尼器優(yōu)化布置的精準定位,基于遺傳算法提出一種新型編碼方式——“AB”型數(shù)字編碼,推導了MRD 附加剛度和附加阻尼矩陣,采取H2范數(shù)優(yōu)化控制理論、相對適應度函數(shù)、多目標優(yōu)化、代間比較權重等理論提出了改進遺傳算法.采用MATLAB 軟件開發(fā)了改進遺傳算法MRD 多目標空間優(yōu)化布置分析程序.基于文獻[16-17]的建模理論,以10 層框剪偏心結構為算例驗證了本文改進遺傳算法的正確性和程序的有效性.

1 改進遺傳算法

1.1 編碼方式

為實現(xiàn)高層空間結構阻尼器優(yōu)化布置的精確定位,對遺傳算法進行合理改進,提出一種新型編碼方式——“A-B”型序列編碼.其中,A 對應建筑結構的樓層數(shù),B 對應該樓層的具體位置,例如基因型為[2-1,2-4,8-7,5-3,6-10,7-5]的個體對應表現(xiàn)型為該結構共布置6 個阻尼器,分別布置在第二層第一個位置、第四個位置,第八層第七個位置等.假設10層結構擬布置18個阻尼器,共100個可選位置,采用遺傳算法對其進行優(yōu)化布置研究,4種編碼方式對比如表1所示.

“A-B”型序列編碼方式僅需明確定義每層可安設阻尼器位置的編號,相比于實數(shù)編碼和數(shù)字序列編碼可精確定位空間結構阻尼器具體安設位置,相比于“0-1”編碼大大縮短基因長度,從而使得遺傳算法更容易尋找最優(yōu)解,加快尋優(yōu)速度.

1.2 遺傳操作

遺傳操作包括選擇、交叉、變異操作,是遺傳算法的主體模塊,其運算方法的設計決定了種群進化的方向.

采用輪盤賭選擇方式對各個體進行優(yōu)勝劣汰操作,采用單點交叉方式,對A、B 基因分開進行交叉操作,分別給予特定的交叉概率,互不影響.其中交叉位置隨機生成.例如:

a1 A:3 5 4 6-8 5 6;B:6-1 9 6 3 8 7.

a2 A:3 6 6 5-4 5 8;B:6-5 9 6 3 8 1.

a1 和a2 表示兩個不同基因型的個體,其中“-”表示交叉位置,完成交叉操作后得到兩個全新個體為:

a1 A:3 5 4 6-4 5 8;B:6-5 9 6 3 8 1.

a2 A:3 6 6 5-8 6 5;B:6-1 9 6 3 8 7.

變異操作和交叉一致,A、B 基因互不影響,其中變異基因位、變異基因編碼在允許范圍內(nèi)隨機生成.

1.3 約束性條件及精英保留策略

相比于“0-1”編碼及數(shù)字序列編碼方式,“A-B”型編碼方式無須考慮交叉、變異操所造成的阻尼器數(shù)量變化約束處理問題,但會出現(xiàn)基因重復現(xiàn)象,即在交叉、變異過程中有可能出現(xiàn)兩個基因位基因代碼相同的情況.基于“A-B”型編碼方式的特點,提出一種基因代碼轉換方式,將個體各基因位代碼轉換為十進制數(shù)字,每個數(shù)字在交叉、變異過程中進行甄別,若出現(xiàn)相同數(shù)字則對其進行“變異”操作,直至無重復數(shù)字為止.

在每代的選擇、交叉、變異操作之后,將當代種群中最小適應度值與上述操作前最大適應度值進行對比,保留適應度值高的個體進入下一代,即“末尾淘汰、精英保留”.隨著種群的不斷迭代,種群整體適應度值不斷提高,即基因優(yōu)良性不斷提高.

2 優(yōu)化控制理論

2.1 優(yōu)化目標函數(shù)

采用H2范數(shù)優(yōu)化控制理論對MRD 進行優(yōu)化布置研究,MRD受控結構減震運動微分方程如下:

式中:M、C、K分別是結構的總質量矩陣、阻尼矩陣、總剛度矩陣;分別為受控框剪結構的加速度、速度、位移響應;為施加在結構上的地震波加速度;R是地震作用列向量;H是MRD 的位置矩陣;F(t)是MRD 所提供的控制力向量.阻尼矩陣C采用Rayleigh正交阻尼模型,詳見文獻[16].

對阻尼力進行變形重組,可得等效剛度和等效阻尼,如下:

結合公式(1)、(2)可得:

式中:=C+HG2為結構的總阻尼矩陣;=K+HG1為結構的總剛度矩陣.

式(3)的狀態(tài)空間的狀態(tài)函數(shù)為:

對受控系統(tǒng)狀態(tài)方程進行拉普拉斯變換,并求取其傳遞函數(shù):

式中:s為拉普拉斯算子[18].

傳遞函數(shù)的H2范數(shù)可用式(6)計算:

式中:Lc為Lyapunov 方程ALc+LcAΤ=-BBΤ的非負定解.

2.2 多目標優(yōu)化與代間比較權重

多目標優(yōu)化最優(yōu)解稱為Pareto 最優(yōu)解,其特點是隨著優(yōu)化進程各目標函數(shù)齊頭并進,且在無法提升任何目標函數(shù)的同時保證不削弱任何目標函數(shù).多目標優(yōu)化算法以MRD 受控結構層間位移的H2范數(shù)指標和加速度的H2范數(shù)指標聯(lián)合控制對MRD 布置位置進行優(yōu)化計算,優(yōu)化目標函數(shù)如式(7)所示.

代間比較權重定義如下:

棗棉間作條件下,不同灌水量對棉花產(chǎn)量存在顯著影響(圖3),適量供水和輕度水分脅迫下產(chǎn)量表現(xiàn)較好,充分供水次之,中度水分脅迫產(chǎn)量最低。M3W3處理產(chǎn)量較M3W1、M3W4分別高62. 5%和13. 9%,M3W2處理產(chǎn)量較M3W1、M3W4分別高55. 6%和9. 0%,M3W4處理產(chǎn)量較M3W1高42. 7%,M3W3處理產(chǎn)量與M3W2處理無顯著差異。

代間比較權重使得進步慢的目標函數(shù)將獲得較大的權值,以促進各目標函數(shù)間的均衡發(fā)展.相比于人為定義權重值以及適應性權重的方法,其不僅排除了人為因素定義兩個權重值不同造成的影響,而且提供了種群向最優(yōu)目標進化的方向和搜索壓力,同時又使得這種壓力在各指標間均勻地進行分配.

2.3 適應度函數(shù)

適應度函數(shù)值的大小影響著個體被遺傳到下一代的概率,個體適應度函數(shù)值越大,遺傳到下一代的概率越高.為實現(xiàn)目標函數(shù)到適應度函數(shù)的良性映射,提高優(yōu)良個體遺傳到下一代的概率,采用相對適應度函數(shù)進行遺傳操作.

式中:為相對適應度函數(shù);f(xi)為個體目標函數(shù);{f(x)}為當代種群目標函數(shù)集合.

2.4 程序編制

基于MRD 受控結構三維計算模型、經(jīng)典遺傳算法、“A-B”型編碼方式、H2范數(shù)優(yōu)化控制理論、多目標優(yōu)化與代間比較權重等理論,采用MATLAB 軟件開發(fā)了改進遺傳算法阻尼器空間位置尋優(yōu)程序.程序計算流程如圖1所示.

圖1 改進遺傳算法阻尼器空間位置尋優(yōu)程序流程圖Fig.1 The flow chart of the optimization procedure for the spatial position of the damper in the improved genetic algorithm

3 優(yōu)化模型及地震波選取

3.1 優(yōu)化模型

以10 層RC 框剪結構為例,結構的混凝土強度等級為C35,彈性模量為3.25×104MPa,泊松比為0.2,鋼筋混凝土密度為2 500 kg/m3,樓板厚度為0.12 m,剪力墻厚度為0.20 m.結構的三維模型如圖2 所示,一至三層框架柱的截面尺寸為0.7 m×0.7 m,四層及以上框架柱的截面尺寸為0.6 m×0.6 m;跨度為9 m梁的截面尺寸為0.3 m×0.7 m,跨度為4.5 m 梁的截面尺寸為0.25 m×0.5 m,跨度為6.3 m 和6.0 m 梁的截面尺寸為0.25 m×0.6 m.RC 框剪結構的前兩階振型阻尼比假定為0.05,并將底層柱下端視為固結.圖2 中數(shù)字為節(jié)點編號值及阻尼器位置編號.建模理論、程序中所采用的MRD 及其在結構中的設置方式和文獻[17]中一致,此處不贅述.

圖2 三維框剪偏心模型(單位:m)Fig.2 3D frame shear eccentric model(unit:m)

3.2 地震波選取

根據(jù)《建筑抗震設計規(guī)范》(GB 50011—2010)[19],擬建場地建筑抗震設防烈度為7(0.1g)度、Ⅲ類場地、設計地震分組為第一組,地震動加速度反應譜特征周期為0.45 s,水平地震影響系數(shù)為0.5.根據(jù)場地類別和設計地震分組,選取ATC-63 FEMA P-695 規(guī)范[20]推薦的22條遠場地震記錄中的5條水平分量作為時程分析的實震輸入,另外選取2 條人工波作為補充,將7 條地震波加速度時程最大值均設置為220 gal,分別繪制加速度響應譜曲線、平均譜曲線與設計譜曲線,如圖3所示.與設計反應譜相比,7組地震波平均反應譜在結構主周期點上誤差為6.48%(設計反應譜:0.209 3,平均反應譜:0.223 8),滿足誤差不超過20%的規(guī)范要求.

圖3 輸入地震波各時程譜及平均譜與設計譜對比Fig.3 Comparison of each time-history spectrum and average spectrum of the input seismic wave with the design spectrum

4 優(yōu)化方案驗證

結構共100 個可選位置(每層10 個),擬設18 個MRD.經(jīng)過程序試算,初始種群大小設為300,交叉概率設置為0.6,A、B 變異概率均設置為0.1,程序可迅速得出最優(yōu)解.程序設置為連續(xù)進化20 代適應度函數(shù)值不發(fā)生改變或者進化到50 代時程序自動停止搜索.對于超大型建筑結構可參照文獻[15],能提高尋優(yōu)速度和最終優(yōu)化結果的準確度.

4.1 位置優(yōu)化結果

程序運行235 min 計算完畢,多目標優(yōu)化進程如圖4 所示.圖4(a)~(c)分別為以層間位移為輸出的H2范數(shù)指標、以加速度為輸出的H2范數(shù)指標和以多目標評價函數(shù)隨進化代數(shù)的變化.從圖4(c)中可以看出多目標評價函數(shù)指標隨著進化進程前期迅速增大,中期緩慢上升,后期趨于穩(wěn)定.為驗證優(yōu)化結果的正確性,對4 種工況進行對比,即:工況一,不設置MRD;工況二,隨機布置;工況三,底層均勻布置;工況四,優(yōu)化布置方案.優(yōu)化方案中X向MRD共12個,Y向6 個,為排除各工況X、Y向阻尼器數(shù)量不同的影響,工況二、工況三中X、Y向阻尼器數(shù)量調(diào)整為12個、6 個,并分別計算各工況下以層間位移為輸出結構的H2范數(shù)指標值H2d和以加速度為輸出結構的H2范數(shù)指標值H2a,如表2所示.

圖4 改進遺傳算法多目標優(yōu)化進程Fig.4 Improved genetic algorithm multi-objective optimization process

表2 各工況MRD布置位置及H2范數(shù)值Tab.2 MRD layout position and H2 norm value of each scheme

4.2 時程響應對比

為了更為明顯地對比4種工況下MRD對結構的減震控制效果,以頂層88 號節(jié)點為例,7 條地震波激勵下的水平雙向位移時程響應分別取前30 s 進行繪制(RSN-848 地震波總時長28 s,取全程),并分別在位移峰值處局部放大,如圖5 所示.從圖5 中可以明顯地看出,無論是在哪一條地震波激勵下,工況二、工況三、工況四對應結構頂層88 號節(jié)點的X、Y向位移均明顯小于工況一條件下對應時間點下的位移響應,說明在地震波激勵下MRD 可明顯減小結構的位移響應.工況二、工況三、工況四在7 條地震波激勵下結構頂層88號節(jié)點X向最大位移的平均減震率分別為31.97%、34.81%、39.45%,Y向最大位移的平均減震率分別為22.93%、29.95%、32.68%.該結果與表2各工況范數(shù)指標值相對應,說明了基于范數(shù)優(yōu)化控制理論的改進遺傳算法空間位置尋優(yōu)程序的有效性,實現(xiàn)了對阻尼器最佳布置位置的精確定位.

圖5 7條地震波激勵下結構頂層88號節(jié)點水平雙向位移時程響應Fig.5 Time history response of horizontal bidirectional displacement of No.88 node on the top floor of the structure under the excitation of seven seismic waves

4.3 各層最大位移、加速度對比

圖6 為7 條地震波作用下框剪結構各層水平雙向位移、加速度平均值包絡圖.由圖6 可知,MRD 受控結構各層最大位移、最大加速度響應均小于未控結構.其中,工況四X向各層位移明顯小于工況二、工況三,Y向略小于工況二、工況三.以中間樓層 第五層為例,7條地震波激勵下工況一至工況四X向平均位移分別為54.52 mm、36.84 mm、36.12 mm、34.88 mm,相比于工況一,工況二、工況三、工況四分別減少了32.43%、33.75、36.05%;Y向平均位移分別為23.18 mm、18.34 mm、15.18 mm、14.96 mm,相比于工況一,工況二、工況三、工況四分別減少了20.88%、34.51%、35.46%.Y向一至三層位移略小于其他各層,原因是該結構一至三層在Y向布置了剪力墻,該位置自身剛度足夠大,抗剪能力足夠強.

圖6 框剪結構各層水平雙向位移、加速度包絡圖Fig.6 Horizontal bidirectional displacement and acceleration envelopes of each layer of frame-shear structure

在地震波激勵下,結構剛度和阻尼的增加均會減小位移響應,對于加速度而言,雖然阻尼的增加會減小結構的加速度響應,但剛度的增加會加大結構加速度響應.因此,MRD 受控結構的加速度減震效果不如位移減震效果明顯,甚至在某些時刻加速度會出現(xiàn)增大的現(xiàn)象.因結構加速度響應結果的復雜性,單層加速度對比有大有小,故采用式(10)對工況二、工況三、工況四的加速度響應進行綜合評價.

式 中:J1為加速 度指 標;n為 樓層 數(shù);aixmax、aiymax、aix0max、aiy0max分別為X、Y向受控和未控狀態(tài)下各層各節(jié)點最大加速度均值.

經(jīng)過計算,工況二、三、四加速度綜合指標值J1分別為1.794、1.778、1.771,說明3 種工況相比,工況四的加速度綜合減震率最大,工況三次之,工況二相對最小.

4.4 結構扭轉控制

7 條地震波激勵下結構各層水平雙向最大位移對比如表3所示.從表3可看出工況二、工況三、工況四與工況一相比,X向扭轉略有增加但不明顯,其中工況四各層位移比值更接近于1;Y向工況二扭轉加劇,工況三、工況四明顯減小,其中工況四各層位移比值更接近于1,最大位移比值頂層1.064,遠遠小于工況一頂層最大位移比值1.314.對比分析結果表明:MRD 布置位置不當會加劇結構的扭轉響應;空間位置尋優(yōu)算法計算結果對結構的扭轉控制作用顯著.

表3 7條地震波作用下各層水平雙向最大位移比平均值Tab.3 The average value of the horizontal bidirectional maximum displacement ratio of each layer under the action of seven seismic waves

4.5 數(shù)量優(yōu)化研究

基于結構安全性和設置MRD 帶來的經(jīng)濟支出,本節(jié)對安設MRD 數(shù)量與結構減震性能的關系進行分析.提出3種評價指標,分別為加速度指標J1、層間位移角指標J2、綜合評價指標J3,從多角度綜合分析不同數(shù)量的阻尼器對結構整體的減震控制效果.

式中:θixmax、θiymax、θix0max、θiy0max分別為X、Y向受控和未控狀態(tài)下各層各節(jié)點最大層間位移角均值.

為充分說明不同數(shù)量的阻尼器對結構綜合減震率的影響,基于MRD 受控框剪結構三維框剪計算模型和改進遺傳算法,對該結構不同數(shù)量的阻尼器布置位置進行優(yōu)化.優(yōu)化程序參數(shù)設置和上文一致,分別計算設置12、16、20、24、28、32、36 個阻尼器的最優(yōu)分布位置,并計算在7 條地震波激勵下各工況以及滿置下結構的綜合評價指標,并取其平均值繪制阻尼器數(shù)量對應綜合評價指標影響關系曲線,如圖7所示.從圖7 中可看出:J1與J2相比,J2即層間位移角指標值較小,說明MRD 對結構的位移減震效果強于加速度的減震效果;3 種評價指標均隨著MRD 數(shù)量的增多總體上逐漸減小,減小幅值趨于平緩,說明隨著結構中MRD 的不斷增加,結構的層間位移角、加速度減震率雖逐漸增大,但減震性價比逐漸減小;通過對比3 種評價指標走勢發(fā)現(xiàn),當MRD 數(shù)量超過18時,3 種指標下降率均趨于平緩,在不考慮其他因素的情況下,說明該結構設置18個MRD減震控制效果最佳.

圖7 各減震評價指標隨MRD數(shù)量變化趨勢圖Fig.7 Trend chart of each shock absorption evaluation index with the number of MRDs

4.6 阻尼器耗能分析

圖8 為RSN-1116 地震波作用 下X向2~4 號和Y向6~7號MRD 的阻尼力-位移滯回曲線.從圖8中可以看出,在地震波的作用下,MRD 的滯回環(huán)呈橢圓形,且阻尼力幅值隨位移幅值的增大而增大,說明MRD 能夠穩(wěn)定地耗能;此外,位移和阻尼力不會在同一時刻達到最大值.當滯回環(huán)中位移幅值小于 2 mm 時,MRD 的阻尼力為10 kN;當滯回環(huán)的位移幅值大于 14 mm 時,滯回環(huán)中MRD 的阻尼力最大可以達到200 kN.

圖8 RSN-1116波作用下MRD的阻尼力-位移滯回曲線Fig.8 Damping force-displacement hysteresis curve of MRD under the action of RSN-1116 wave

5 結論

1)基于遺傳算法提出了“A-B”型數(shù)字編碼,推導了MRD 附加剛度和附加阻尼矩陣,采用H2范數(shù)優(yōu)化控制、代間比較權重等理論提出了改進遺傳算法,并開發(fā)了MATLAB 程序,實現(xiàn)了對空間結構中MRD的精確定位和多目標優(yōu)化布置.

2)基于MRD 受控RC 結構三維計算模型程序,采用該算法對10 層框剪結構MRD 的布置進行優(yōu)化研究,結果表明:改進遺傳算法優(yōu)化結果迅速收斂;4種工況中優(yōu)化方案控制效果最佳,頂層88號節(jié)點X、Y向減震率分別達39.45%、32.68%,結構扭轉得到有效控制.

3)對不同數(shù)量MRD 優(yōu)化布置進行分析,結果表明其綜合減震指標雖然隨著數(shù)量的增加而減小,但降低趨勢也逐漸減小,該指標可為阻尼器的數(shù)量優(yōu)化提供依據(jù).

4)MRD 的滯回環(huán)呈橢圓形,且阻尼力幅值隨位移幅值的增大而增大,說明MRD 能夠穩(wěn)定地發(fā)揮耗能作用.

猜你喜歡
優(yōu)化結構
超限高層建筑結構設計與優(yōu)化思考
《形而上學》△卷的結構和位置
哲學評論(2021年2期)2021-08-22 01:53:34
民用建筑防煙排煙設計優(yōu)化探討
關于優(yōu)化消防安全告知承諾的一些思考
一道優(yōu)化題的幾何解法
由“形”啟“數(shù)”優(yōu)化運算——以2021年解析幾何高考題為例
論結構
中華詩詞(2019年7期)2019-11-25 01:43:04
新型平衡塊結構的應用
模具制造(2019年3期)2019-06-06 02:10:54
論《日出》的結構
創(chuàng)新治理結構促進中小企業(yè)持續(xù)成長
主站蜘蛛池模板: 亚洲欧洲日产国码无码av喷潮| 超薄丝袜足j国产在线视频| 不卡网亚洲无码| 久久综合丝袜日本网| 99在线视频免费| 欧美精品另类| 欧美成人综合在线| 国产一级二级三级毛片| 波多野结衣在线se| 国产精品女人呻吟在线观看| 九色视频一区| 色偷偷男人的天堂亚洲av| 免费精品一区二区h| 黄色a一级视频| 丁香六月综合网| 欧美成人国产| 国产成人三级| 国产永久在线视频| 国产成人高清精品免费5388| 伊人成人在线| 欧美中文字幕在线二区| 国产成人高清精品免费软件| www.youjizz.com久久| 九九这里只有精品视频| 精品超清无码视频在线观看| 国产真实二区一区在线亚洲| 亚洲熟女中文字幕男人总站| 2020最新国产精品视频| 日韩乱码免费一区二区三区| 色妺妺在线视频喷水| 亚洲精品视频免费| 亚洲AⅤ波多系列中文字幕| 99国产精品免费观看视频| 精品三级在线| 曰AV在线无码| 精品三级在线| 亚洲人网站| 日本一区二区三区精品国产| 国产aⅴ无码专区亚洲av综合网| 日本亚洲成高清一区二区三区| 国产大全韩国亚洲一区二区三区| a免费毛片在线播放| 中文字幕在线日韩91| 中文字幕在线播放不卡| 免费高清a毛片| 亚洲国产91人成在线| 国产综合在线观看视频| 婷婷亚洲天堂| 日韩av高清无码一区二区三区| 国内精品一区二区在线观看| 91网址在线播放| jijzzizz老师出水喷水喷出| 九九线精品视频在线观看| 色综合中文| 国产99久久亚洲综合精品西瓜tv| 伊人久久久久久久久久| 国产精品自在在线午夜区app| 久久 午夜福利 张柏芝| 曰韩人妻一区二区三区| 91免费国产在线观看尤物| 又爽又大又光又色的午夜视频| 日本三级精品| 黄色网页在线播放| 亚洲日韩欧美在线观看| 91九色视频网| 欧美日韩另类在线| 亚洲一级毛片在线观播放| 在线精品视频成人网| 久久永久精品免费视频| 亚洲人成日本在线观看| 日本在线欧美在线| 99这里只有精品在线| 亚洲一区二区在线无码| 播五月综合| 欧美一区福利| 国产一级毛片yw| 成人国产精品视频频| 亚洲最新在线| 国产69囗曝护士吞精在线视频| 91欧美在线| 欧美成人h精品网站| 日韩在线成年视频人网站观看|