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

聯合激勵下分數階非線性系統非平穩響應的半解析方法

2024-12-31 00:00:00孔凡廖海君韓仁杰張義洪旭
振動工程學報 2024年8期

摘要: 確定性和隨機激勵聯合作用下的非線性動力系統具有特殊的動力響應特征。本文提出一種用于計算聯合激勵下含分數階阻尼的非線性系統非平穩響應的半解析方法。將系統響應表示為確定性響應和零均值隨機響應之和,則原分數階非線性運動微分方程可等效地化為分數階確定性微分方程和隨機子微分方程的組合。利用時變諧波平衡法處理非線性確定性微分方程,利用統計線性化處理非線性隨機子微分方程。對于后者,結合Prony?SS算法和Laplace變換得到其分數階等效線性方程的半解析解。聯立得到的相關耦合方程,通過數值算法迭代求解響應未知量。蒙特卡羅模擬驗證了此方法的適用性和精度。

關鍵詞: 統計線性化; 時變諧波平衡法; 分數階導數; 非線性系統; Prony?SS算法

中圖分類號: O324; O322" " 文獻標志碼: A" " 文章編號: 1004-4523(2024)08-1339-10

DOI:10.16385/j.cnki.issn.1004-4523.2024.08.008

引" 言

分數階微積分是微積分學的重要分支,距今已有300多年的歷史[1]。近幾十年來,分數階導數模型在工程領域得到了廣泛的應用[2]。土木或機械工程領域的許多情形下,都使用黏彈性阻尼器降低結構振動,而分數階導數模型描述黏彈性材料本構關系具有較大優勢。在這方面,Gemant[3]提出分數階導數模型用以精確描述黏彈性材料的頻率依賴性; Slonimsky[4], Smit等[5]發現黏彈性介質中的應力與應變之間存在分數階微積分關系;Lewandowski等[6]提出黏彈性阻尼器的分數Kelvin?Voigt模型和分數Maxwell模型的參數辨識方法;Bagley等[7?8]通過大量實驗和數據分析得出,分數階導數模型可用于描述黏彈性材料的應力松弛和蠕變現象,且模型簡單、參數少。可見,具有分數階導數阻尼的運動微分方程能很好地描述裝配有黏彈性減/隔振(震)阻尼器結構的動力特性[9]。因此,亟需發展求解分數階運動微分方程的解析或數值方法。

分數階系統的隨機動力響應分析較其確定性響應分析更具挑戰。到目前為止,人們通過若干方法得到了線性系統的隨機動力響應,例如:Pinnola[10]基于復譜矩的概念通過分數階狀態方程的特征向量展開和Melin變換得到了線性分數階系統的二階矩解析解;文獻[11?12]通過Laplace變換和Prony?SS算法得到了分數階線性單自由度系統響應二階矩的半解析解;Di Paola等[13]基于分數階線性系統狀態方程和特征向量展開得到了多自由度線性系統響應的功率譜密度數值解。然而,以上通過特征向量展開或Laplace變換的方法均只適用于分數階線性隨機動力系統,而無法適用于分數階非線性系統。另外,工程隨機激勵(如地震激勵)具有明顯的非平穩特性,而只有基于Laplace變換的方法[11?12]可得到非平穩(半)解析解。

實際工程中的結構會同時受到確定性和隨機激勵聯合作用。例如,在隨機風浪作用下運行的風力發電機[14]。在力學等相關領域,確定性和隨機激勵作用下非線性系統的隨機動力響應得到廣泛關注[15?20]。然而,工程結構中隨機非平穩響應[21]、滯回特性[22]及多自由度系統[23]等問題仍有待進一步研究。

文獻[24?25]發展了非線性系統在聯合激勵作用下非平穩響應的統計線性化方法。該方法的關鍵在于將系統響應分解為確定性和隨機響應分量之和,從而得到兩組耦合的、分別以確定性和隨機響應為未知量的子微分方程;再分別以確定性和隨機動力方法求解從而得到總響應。該方法已被推廣到分數階非線性系統的隨機動力響應分析[26],其中,利用了Yuan?Agrawal(YA)的無記憶方法[27]將分數階導數化為了整數階導數。然而,研究發現:Yuan?Agrawal對分數階導數整數化的處理存在局限,分數階導數接近0或1時方法的精度下降。為此,本文提出使用文獻[11?12]提出的Laplace變換方法計算隨機等效線性方程的近似半解析解,再結合統計線性化方法依照聯合激勵下整數階非線性系統非平穩響應的思路求得分數階系統的響應。

本文結構如下:對分數階非線性運動方程分解,得到等效確定性和隨機子微分方程;利用時變諧波平衡法求解確定性非線性子微分方程;利用統計線性化方法處理隨機非線性子微分方程,并通過Laplace變換和極點?留數方法得到其等效線性方程的半解析解;聯立相關方程并通過數值算法迭代求解響應未知量。最后,通過大量的數值算例驗證本文建議方法的精度和適用性。

1 理論方法

1.1 動力學方程

具有分數階導數阻尼的單自由度非線性系統在隨機與確定性調制諧波聯合激勵作用下的運動方程為:

式中" m,c和k分別為系統的質量、分數階阻尼系數和剛度系數;q為分數階導數的階數;,和分別為結構的位移、速度和加速度;為關于結構位移、速度和加速度的非線性函數;為確定性諧波激勵的頻率,慢變時間調制函數為[28]:

式中" 為諧波激勵的幅值;,為控制調制函數形狀的參數。

為零均值非平穩隨機過程,本文采用調制非平穩過程模型,即

式中" 表示功率譜密度為的白噪聲;為時間調制函數,采用指數函數疊加的形式表示[11]:

式中" 為指數分量的項數;和為實數。可見,式(2)所示調制函數為式(4)的特殊形式;通過調整函數的幅值和指數大小能方便地調整指數函數上升和下降段速率。

式(1)中,表示階分數階導數,采用Caputo定義:

為使分數導數項量綱與經典阻尼項一致,取,其中為結構的自振頻率,為系統阻尼比。

將系統響應分解為確定性諧波與零均值隨機分量之和,即[24?26]

注意到,文獻[28]在處理零均值隨機激勵作用下非零點對稱(如平方)非線性振子時也采用了類似形式,表示響應為非零均值隨機過程;有關該分解的合理性,見文獻[24?26]的評述。對式(6)兩邊求期望可得:

由分數階導數的線性性質[1]可知:

因此,聯合激勵下分數階非線性運動方程式(1)可分解為確定性運動方程(11)和隨機運動方程(12)。由于式(11)包含隨機分量,式(12)包含確定性分量。因此,這兩組方程式相互耦合,須同時考慮才能求得系統響應。

1.2 確定性響應的時變諧波平衡法

1.3 隨機響應的統計線性化方法

可見,其與隨機響應標準差和確定性響應有關。的非平穩性主要源于慢變的隨機激勵調制函數;的時變性來源于慢變的諧波調制函數和快變的被調制諧波/。注意到,以上非平穩性或時變性導致了等效線性剛度為時間的函數。將中的快變分量在一個周期內平均,即

上式中的依賴于等效線性剛度。綜上:式(17)和(18)給出了之間關系的兩個方程;式(43)給出了和之間關系的第三個方程;式(24)建立了和之間關系的第四個方程。可通過數值方法求解以上四個非線性代數方程。

1.4 求解流程

(1) 確定諧波系數和的初值和隨機響應分量的標準差的初值。令非線性強度系數,則,通過式(17)和(18)求得諧波系數的初值和,通過式(43)求得隨機響應分量標準差初值,設迭代步。

(2) 將代入式(48)利用Newton迭代法更新諧波系數和;利用,和通過式(24)求得等效線性剛度系數;根據式(43)求得更新后的標準差。

(3) 判定收斂:如不收斂,,并回到第(2)步;如收斂,結束循環。

值得注意的是,本節的方法同樣也適用于隨機激勵分量為調制色噪聲或/和確定性激勵具有多諧波分量的情形。相應地,需要采用前置濾波器方法對式(12)擴階或/和使用多諧波平衡法處理式(11)。

2 數值算例

2.1 典型響應

本文建議的方法計算得到的響應標準差和均值與MCS的結果對比如圖1所示。其中,MCS中樣本激勵由譜表現方法生成,樣本數為10000個。由圖可見,除確定性響應起始部分外,建議方法得到的響應標準差和均值與MCS估計的結果幾乎吻合,而且可以很好地體現響應的非平穩性。圖1(a)中,由兩種方法得到的響應均值穩態平均功率相差;圖1(b)中,兩種方法得到的響應標準差時間平均則相差,驗證了所建議方法具有較好的精度。

為進一步驗證所建議方法在其他參數設置情況下的適用性,進行以下參數分析。

2.2 確定性諧波激勵頻率

確定性諧波激勵頻率無疑會影響確定性響應分量。定性而言,激勵頻率接近線性系統自振頻率時會增大的幅值。確定性激勵頻率如何影響隨機響應分量,仍有待考察;同時,需研究所建議方法在確定性激勵頻率變化時的適用性。為此,除使確定性諧波激勵頻率變化外,其他系統和激勵參數同2.1節。圖2和3所示為不同頻率引起的系統響應變化曲線。由圖2,3可見,諧波激勵頻率改變時,對確定性響應有較大影響。使確定性響應包絡平均達到峰值的諧波激勵頻率大于自振頻率1.0 rad/s幅?頻響應曲線向高頻傾斜,符合硬化系統的特征。此外,在確定性響應的峰值頻率處,隨機響應標準差的時間平均達到最低。總體而言,所建議方法準確地捕捉了MCS結果顯示的上述響應特征,且和后者預測結果的吻合性較好。值得注意的是,當確定性激勵分量頻率處于1.1~1.2 rad/s之間(主外共振頻率)時,方法對預測響應標準差時間平均的適用性降低。

2.3 確定性諧波激勵幅值

確定性諧波激勵幅值同樣會影響系統確定性響應分量。因此,需考察確定性諧波激勵幅值對隨機響應的影響,同時研究在這種情況下所建議方法的適用性。同樣地,除使諧波激勵幅值變化外,其他系統和激勵參數同2.1節。考察與系統響應之間的關系。如圖4和5所示分別為兩種方法得到的確定性響應包絡平均和隨機響應分量標準差的時間平均與諧波激勵幅值之間的關系曲線。可以看出,隨著諧波激勵幅值增大而增大,兩種方法之間的最大誤差僅為3.37%;反之,隨著諧波激勵幅值增大僅稍有下降,兩種方法之間的最大誤差約為-3.93%。

2.4 隨機激勵譜強度

考察隨機激勵譜強度對系統響應的影響。同樣地,除使變化外,其他系統和激勵參數同2.1節。考察隨機激勵譜強度與系統響應之間的關系。圖6和7所示分別為所建議方法和MCS得到的確定性響應包絡平均和響應標準差時間平均與之間的關系曲線。可見,隨譜強度的增大而增大,兩種方法獲得結果之間的最大誤差約為;隨的增大而減小,兩種方法所獲結果之間的最大誤差約為。

2.5 系統非線性強度

非線性強度對系統響應的影響至關重要。考察系統非線性強度變化時所建議方法的適用性及非線性強度對響應分量的影響。如圖8和9所示為建議方法和MCS得到的確定性響應包絡平均和響應標準差的時間平均隨非線性強度系數變化的曲線。由圖可見,和都隨的增加而降低,其中以更為顯著。圖8中,所建議方法與MCS結果的最大誤差約為3.75%,而圖9中的最大誤差約為,均對應。可見,所建議方法得到的確定性和隨機響應分量的計算精度均隨的增加而降低,但仍在一般統計線性化方法誤差的合理范圍內[28]。

2.6 系統分數階大小

根據Singh等[30]的研究,分數階導數模型以較少的模型參數就能逼近標準線性固體模型(Standard Linear Solid model, SLS, 即多個整數階導數阻尼和剛度單元串并聯)才能模擬的儲能或耗能模量的頻率依賴行為。定性而言,分數階導數增加時,模型阻尼成分增加、剛度成分減小;反之,分數階導數減小時,模型剛度成分增加、阻尼成分減小。本文提出采用Laplace變換方法計算隨機等效線性系統響應,相比Yuan?Agrawal的無記憶方法,適用于分數階導數()較大或較小的情況。因此,有必要考察分數階導數不同時所建議方法的適用性,以及分數階導數大小對響應標準差的時間平均和確定性響應包絡平均的影響。如圖10和11所示分別為和隨分數階變化的曲線,其中給出了所建議方法和MCS得到的結果對比。由圖可見,和均隨增加而逐漸下降,在處達到最大值和最大誤差,均值響應誤差約為,隨機響應誤差約為。

3 結" 論

本文提出了一種求解分數階非線性振子在確定性調制諧波激勵和調制白噪聲作用下非平穩響應的統計線性化方法。首先,系統響應表示為確定性調制諧波響應和零均值響應之和,將原方程等效地分解為耦合的確定性和隨機子微分方程。采用時變諧波平衡和統計線性化方法分別處理了確定性和隨機子微分方程,得到了響應未知量之間的非線性代數關系;其中,基于Laplace變換和極點?留數方法建立了隨機等效線性系統響應方差和等效線性參數之間的關系。最后,利用迭代方法求解了未知量之間耦合的非線性方程。通過蒙特卡羅模擬對該方法進行了驗證和參數分析。

結果表明:該方法繼承了統計線性化方法在處理非線性隨機振動問題方面的普遍適用性,基于Laplace變換的分數階等效線性系統響應的半解析解精度高、適用于分數階的所有情況。注意到,線性分數階系統的脈沖響應函數具有顯式解析表達,采用它而非由頻率響應函數逆變換而來的數值表達(式(29)),可進一步得到線性分數階系統非平穩響應的解析解。本文所建議的方法也可拓展至滯回非線性系統和多自由度非線性系統等情況。

參考文獻:

[1] 吳強, 黃建華. 分數階微積分[M]. 北京: 清華大學出版社, 2016.

[2] Cai Min, Li Changpin. Numerical approaches to fractional integrals and derivatives: a review[J]. Mathematics, 2020, 8(1): 43.

[3] Gemant A. A method of analyzing experimental results obtained from elasto?viscous bodies[J]. Journal of Applied Physics, 1936, 7(8): 311?317.

[4] Slonimsky G L. Laws of mechanical relaxation processes in polymers[J]. Journal of Polymer Science Part C: Polymer Symposia, 1967, 16(3): 1667?1672.

[5] Smit W, de Vries H. Rheological models containing fractional derivatives[J]. Rheologica Acta, 1970, 9(4): 525?534.

[6] Lewandowski R, Chor??yczewski B. Identification of the parameters of the Kelvin?Voigt and the Maxwell fractional models, used to modeling of viscoelastic dampers[J]. Computers amp; Structures, 2010, 88(1?2): 1?17.

[7] Bagley R L, Torvik P J. Fractional calculus—a different approach to the analysis of viscoelastically damped structures[J]. AIAA Journal, 1983, 21(5): 741?748.

[8] Bagley R L, Torvik P J. A theoretical basis for the application of fractional calculus to viscoelasticity[J]. Journal of Rheology, 1983, 27(3): 201?210.

[9] Di Paola M, Pirrotta A, Valenza A. Visco?elastic behavior through fractional calculus: an easier method for best fitting experimental results[J]. Mechanics of Materials, 2011, 43(12): 799?806.

[10] Pinnola F P. Statistical correlation of fractional oscillator response by complex spectral moments and state variable expansion[J]. Communications in Nonlinear Science and Numerical Simulation, 2016, 39: 343?359.

[11] Cao Q Y, Hu S L J, Li H J. Nonstationary response statistics of fractional oscillators to evolutionary stochastic excitation[J]. Communications in Nonlinear Science and Numerical Simulation, 2021, 103: 105962.

[12] Hu S L J, Yang W L, Li H J. Signal decomposition and reconstruction using complex exponential models[J]. Mechanical Systems and Signal Processing, 2013, 40(2): 421?438.

[13] Di Paola M, Pinnola F P, Spanos P D. Analysis of multi?degree?of?freedom systems with fractional derivative elements of rational order[C]∥ ICFDA’14 International Conference on Fractional Differentiation and Its Applications 2014. Catania: IEEE, 2014: 1?6.

[14] Larsen J W, Iwankiewicz R, Nielsen S R K. Nonlinear stochastic stability analysis of wind turbine wings by Monte Carlo simulations[J]. Probabilistic Engineering Mechanics, 2007, 22(2): 181?193.

[15] 陳林聰, 朱位秋. 諧和與寬帶噪聲聯合激勵下含分數導數型阻尼的Duffing振子的平穩響應[J]. 應用力學學報, 2010, 27(3): 517?521.

Chen Lincong, Zhu Weiqiu. Stationary response of Duffing oscillator with fractional derivative damping under combined harmonic and wide band noise excitations[J]. Chinese Journal of Applied Mechanics, 2010, 27(3): 517?521.

[16] Chen Lincong, Zhu Weiqiu. Stochastic jump and bifurcation of Duffing oscillator with fractional derivative damping under combined harmonic and white noise excitations[J]. International Journal of Non?Linear Mechanics, 2011, 46(10): 1324?1329.

[17] Chen Lincong, Zhu Weiqiu. Stochastic averaging of strongly nonlinear oscillators with small fractional derivative damping under combined harmonic and white noise excitations[J]. Nonlinear Dynamics, 2009, 56(3): 231?241.

[18] Chen Lincong, Zhu Weiqiu. Stochastic stability of Duffing oscillator with fractional derivative damping under combined harmonic and white noise parametric excitations[J]. Acta Mechanica, 2009, 207(1-2): 109?120.

[19] Ellermann K. On the determination of nonlinear response distributions for oscillators with combined harmonic and random excitation[J]. Nonlinear Dynamics, 2005, 42(3): 305?318.

[20] Huang Zhilong, Zhu Weiqiu, Suzuki Y. Stochastic averaging of strongly non?linear oscillators under combined harmonic and white?noise excitations[J]. Journal of Sound and Vibration, 2000, 238(2): 233?256.

[21] 孔凡, 沈子恒, 何衛, 等. 周期和隨機聯合激勵作用下滯回系統非平穩隨機響應的一種近似方法[J]. 振動與沖擊, 2022, 41(16): 108?116.

Kong Fan, Shen Ziheng. He Wei,et al. An approximate approach for non?stationary stochastic response of a hysteretic system subjected to combined periodic and stochastic excitation[J]. Journal of Vibration and Shock, 2022, 41(16): 108?116.

[22] Kong Fan, Spanos P D. Stochastic response of hysteresis system under combined periodic and stochastic excitation via the statistical linearization method[J]. Journal of Applied Mechanics, 2021, 88(5): 051008.

[23] Kong Fan, Zhang Huimin, Zhang Yixin, et al. Stationary response determination of MDOF fractional nonlinear systems subjected to combined colored noise and periodic excitation[J]. Communications in Nonlinear Science and Numerical Simulation, 2022, 110: 106392.

[24] 孔凡, 韓仁杰, 張遠進. 確定性周期與隨機激勵聯合作用下非線性系統非平穩響應的統計線性化方法[J]. 振動工程學報, 2022, 35(3): 625?634.

Kong Fan, Han Renjie, Zhang Yuanjin. Non?stationary response of non?linear systems subjected to combined periodic and non?stationary stochastic excitation via the statistical linearization method[J]. Journal of Vibration Engineering, 2022, 35(3): 625?634.

[25] Kong Fan, Han Renjie, Li Shujin, et al. Non-stationary approximate response of non-linear multi-degree-of-freedom systems subjected to combined periodic and stochastic excitation[J]. Mechanical Systems and Signal Processing, 2022, 166: 108420.

[26] 李書進, 張志聰, 孔凡, 等. 周期與色噪聲聯合作用下分數階Duffing振子非平穩響應的無記憶方法[J]. 振動工程學報, 2023, 36(4): 923?933.

Li Shujin, Zhang Zhicong, Kong Fan, et al. A memory?free method for fractional?order Duffing systems subjected to combined periodic and colored excitation[J]. Journal of Vibration Engineering, 2023, 36(4): 923?933.

[27] Yuan Lixia, Agrawal O P. A numerical scheme for dynamic systems containing fractional derivatives[J]. Journal of Vibration and Acoustics, 2002, 124(2): 321?324.

[28] Roberts J B, Spanos P D. Random Vibration and Statistical Linearization[M]. New York, USA: Courier Corporation, 2003.

[29] 陳曉江. 數值分析[M]. 武漢: 武漢理工大學出版社, 2013.

[30] Singh M P, Chang T S. Seismic analysis of structures with viscoelastic dampers[J]. Journal of Engineering Mechanics, 2009, 135(6): 571?580.

A semi?analytical method for non?stationary response determination of nonlinear systems subjected to combined excitation

KONG Fan LIAO Hai?jun HAN Ren?jie ZHANG Yi HONG Xu

(1.School of Civil Engineering amp; Architecture, Wuhan University of Technology, Wuhan 430070, China; 2.College of Civil Engineering, Hefei University of Technology, Hefei 230009, China; 3.College of Civil Engineering, Tongji University, Shanghai 200092, China; 4.China Construction Third Bureau First Engineering Co., Ltd., Wuhan 430040, China)

Abstract: The nonlinear dynamic systems exhibit particular behaviors when subjected to combined deterministic and stochastic excitation. A semi-analytical method for calculating the nonstationary response of a fractional nonlinear oscillator subjected to combined excitation is proposed. Representing the system response as a sum of a deterministic component and zero-mean stochastic component leads to two equivalent sub-equations for the differential equation of motion. The time-varying harmonic balance method is used for the nonstationary solution of the deterministic differential sub-equation, while the statistical linearization method is utilized for obtaining an equivalent linear substitution for the stochastic sub-equation. A semi-analytical solution of the equivalent linear equation is obtained by the Prony-SS and Laplace transform technique. The unknown deterministic/stochastic response components are obtained by solving the derived nonlinear algebraic equations simultaneously. Monte Carlo simulations demonstrate the applicability and accuracy of this method.

Key words: statistical linearization; time-varying harmonic balance; fractional derivative; nonlinear system; Prony-SS algorithm

作者簡介: 孔" 凡(1984―),男,博士,教授。 E?mail: kongfan@hfut.edu.cn。

通訊作者: 洪" 旭(1993—),男,博士,講師。 E?mail: xhong@hfut.edu.cn。

主站蜘蛛池模板: 亚洲成人www| 国产91色在线| 五月婷婷丁香色| 成人中文字幕在线| 国产欧美成人不卡视频| 五月天福利视频| 精品无码视频在线观看| 日本黄色不卡视频| 色婷婷色丁香| 毛片一区二区在线看| 国产一区二区福利| 中国毛片网| 日韩欧美国产精品| 亚洲综合极品香蕉久久网| 国产又粗又猛又爽| 九九九精品视频| 激情乱人伦| 国产高清精品在线91| 日韩123欧美字幕| 精品国产免费观看| 尤物精品国产福利网站| 亚洲综合精品香蕉久久网| 国产精品九九视频| 久久久久九九精品影院| 福利国产微拍广场一区视频在线| 在线观看免费AV网| 国产欧美精品一区二区| 青青青伊人色综合久久| 就去吻亚洲精品国产欧美| 四虎亚洲精品| 色婷婷啪啪| 99青青青精品视频在线| 亚洲高清无在码在线无弹窗| 亚洲系列无码专区偷窥无码| 伊人婷婷色香五月综合缴缴情| Jizz国产色系免费| 亚洲码一区二区三区| 视频一本大道香蕉久在线播放 | 亚洲A∨无码精品午夜在线观看| 久久这里只有精品2| 国产高清不卡| 日韩欧美国产中文| 永久免费无码日韩视频| 精品少妇人妻无码久久| 色妞www精品视频一级下载| 五月婷婷伊人网| 欧美三级视频在线播放| 免费在线a视频| 日本精品αv中文字幕| 日韩无码白| 91精品专区| 国产xxxxx免费视频| 成人毛片免费在线观看| 国产第一页屁屁影院| 亚洲天堂2014| 全免费a级毛片免费看不卡| 制服丝袜国产精品| jizz国产在线| 98超碰在线观看| 中文字幕在线一区二区在线| 欧美三级自拍| 国产aⅴ无码专区亚洲av综合网| 成人国内精品久久久久影院| 国产高清免费午夜在线视频| 热99re99首页精品亚洲五月天| 久久77777| 狠狠综合久久| 国产超碰在线观看| 午夜国产大片免费观看| 欧美成人h精品网站| 中文字幕人妻av一区二区| 91人人妻人人做人人爽男同| 午夜欧美在线| 凹凸国产分类在线观看| 中字无码av在线电影| 色欲国产一区二区日韩欧美| 性欧美在线| 午夜福利无码一区二区| 亚洲日韩Av中文字幕无码| 日韩欧美中文字幕在线韩免费| 试看120秒男女啪啪免费| 无码区日韩专区免费系列|