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

正交各向異性層合板殼耦合結構的能量傳遞特性

2022-03-18 00:46:04陳海波
振動與沖擊 2022年5期
關鍵詞:模態結構

李 翱, 陳海波, 鐘 強, 王 幸

(中國科學技術大學 近代力學系 中國科學院材料力學行為與設計重點實驗室,合肥 230026)

板殼結構廣泛應用于汽車、船舶、航空航天等工程領域,其中很多采用的是復合材料層合結構,且在許多服役狀態下都需要分析其高頻振動響應。板殼高頻振動的一個顯著特征是結構變形的波長遠小于結構的整體尺寸,這意味著通過標準有限元建模計算成本高,非確定性因數較為敏感,因此計算魯棒性不強[1]。為了避免這兩個問題,多種能量方法已被提出,其中最流行的是統計能量分析(statistical energy analysis,SEA)以及它的擴展方法能量有限元法(energy finite element method,EFEM)、能量輻射傳遞法(radiative energy transfer method,RETM)等[2],它們大多基于彈性波在結構中的運動進行研究[3-4]。彈性波在傳播過程中,如遇結構不連續處,即多結構耦合處以及邊界處,會產生折反射且伴隨著入射功率在不同波場形式間的重新分配,這也意味著能量的重新分配。研究耦合處波的傳遞是一個傳統問題,它體現了能量在耦合處的分配,而能量傳遞系數則從數值上量化了這種分配[5]。因此,討論彈性波在結構不連續處的透反射情況,研究能量在板殼結構中的耦合傳遞系數是極為必要的,它是計算系統耦合損耗因子(coupling loss factor,CLF)的關鍵[6]。實際上采用功率流有限元等能量方法研究耦合結構的振動傳遞特性時,能量傳遞系數的計算是必須首先要解決的關鍵問題之一[7-8]。

目前,對耦合結構能量傳遞特性的研究主要針對板結構,如Le Bot推導了多板沿一條邊界耦合的能量傳遞系數,Langley等[9]針對多板耦合于梁的系統,研究了彈性波的傳遞特性。江民圣等[10]對任意夾角的耦合板能量傳遞系數進行研究,討論了板厚及激勵頻率對L型耦合板能量傳遞系數的影響。葛月[11]計算了耦合板的能量傳遞系數,并將其應用于對應耦合結構的能量邊界元分析中。Yan等[12-13]對復合材料板的耦合能量傳遞系數進行推導,然后基于此對系統進行能量有限元分析。Xie等[14]計算了只考慮彎曲波時正交各向異性加筋板的能量傳遞系數,對耦合系統的彎曲振動能量有限元進行了驗證分析。

曲板與平板耦合、曲板與曲板耦合常見于飛機機身結構中,例如,機翼-機身連接、地板-機身連接,但對于彈性波在此類耦合系統中的傳遞特性研究較少。Langley計算了彈性波在曲板耦合系統中的能量傳遞系數和耦合損耗因子,Tso等[15]研究了波在封閉圓柱殼與肋板周向耦合系統中的傳播,但均未涉及復合材料層合結構。對層合曲板的研究目前主要針對其單個子系統的波數、模態密度和傳聲損失[16],只有Lee討論了波在復合材料圓柱殼軸向和環向上的傳播。復合材料層合結構在實際工程中的應用日益廣泛,占比也越來越重,復合材料耦合結構中的能量傳遞系數計算已成為耦合系統分析必須解決的基本問題之一[17]。簡言之,研究波在復合材料層合板殼耦合系統中的能量傳遞特性,極具工程意義。

本文基于復合材料的力學性能對層合板殼結構的平衡微分方程進行推導,進而根據波法的位移假設求解波數曲線,然后借鑒Langley提出的波動剛度矩陣方法對正交各向異性層合板殼耦合結構的能量傳遞系數進行求解,并在此基礎上計算了層合板殼耦合系統的耦合損耗因子。本文的工作完善了層合板殼耦合結構高頻振動相關理論,為后續進一步進行能量預報奠定了基礎。

1 正交各向異性層合板殼的平衡微分方程

平衡微分方程是研究層合結構耦合系統能量傳遞特性的關鍵,本文基于薄殼理論[18]和經典層合板理論[19]推導正交各向異性層合板殼的平衡微分方程。

1.1 層合板殼的應力應變關系

圖1為正交各向異性單層板單元體。圖1中:x-y為偏軸坐標系;1-2為材料主方向即正軸坐標系;x軸與材料主方向1之間的夾角α為復合材料單層方向角。

圖1 一般正交各向異性板的兩個坐標系

其正軸模量矩陣為

(1)

其偏軸模量矩陣為

(2)

(3)

對于一般復合材料層合結構如圖2所示。

圖2 一般層合板殼的幾何結構

其內力-應變關系為

(4)

(5)

式中:N、M為內力矩陣;Nx、Ny為軸力;Nxy為切向力;Mx、My為彎矩;Mxy為扭矩;ε為中面應變矩陣;κ為中面曲率矩陣。

1.2 層合板殼的平衡微分方程

對于圖3所示的層合殼結構。

圖3 環向開敞層合柱殼單元

其中面應變和曲率如下

(6)

式中:u、v、w為中面位移;R為柱殼的曲率半徑。通過式(5)和(6)可求解內力矩陣式(4),結合式(4)和式(7)

(7)

便可計算出剪力Qx、Qy及加上扭矩等效剪力后的總剪力Qxt、Qyt。

將內力與剪力代入柱殼微元體的力平衡方程

(8)

式中,qx、qy、qz為層合柱殼所受載荷在x、y、z三個方向上的分量。

得層合柱殼的平衡微分方程

L(u)=0

(9)

(10)

(11)

(12)

2 正交各向異性層合板殼耦合能量傳遞

針對各向異性結構無法利用亥姆霍茲方程對平衡微分方程進行解耦,這里直接基于波法給出位移表達式,由Snell折反射定律[20]得波數在耦合邊界上的投影相等,將位移代入平衡微分方程求解特征值和特征向量,最后基于耦合處力的平衡條件和位移連續性條件,對能量傳遞系數進行求解。

2.1 正交各向異性層合板殼波數曲線

在圖3中,考慮曲率半徑的影響,對于x邊界,坐標y為常數,其轉角為v/R-?w/?y。所以四個邊界位移可以寫成

d=S(u)

(13)

S是一個4×3的微分算子,其非零元素如下

(14)

同理,對應的四個邊界力Nyx、Ny、Qyt、My可以寫成如下形式

t=T(u)

(15)

T的元素如下

(16)

(17)

(18)

(19)

基于波法,位移矩陣u的解可以假設為以下形式

u=αe-ikxx-ikyy+iωt

(20)

式中:ω為圓頻率;kx=kcosθ、ky=ksinθ分別為波數k在x和y方向的分量;θ為波的入射角;α為位移的復波幅向量。將式(20)代入平衡微分方程式(11)得

H(kx,ky,ω)α=0

(21)

H是一個3×3的矩陣,給定ω,可以計算出相應的波數曲線。如圖4所示,曲線1,2代表彎曲波,曲線3代表剪切波,曲線4代表縱波(也叫拉壓波)。對層合柱殼只有當ω增大到一定值時才會出現縱波,否則只存在剪切波和彎曲波。

2.2 耦合系統的能量傳遞系數

令λ=-iky,式(21)可以寫為

H(kx,λ,ω)α=0

(22)

由Snell定律可知,kimcosθim=kjncosθjn(m、n為結構編號;i、j=1,2,3,4為圖4中波的編號),即每個子結構中的波沿耦合邊界具有相同的運動,它們在x方向的波數分量kx相等。 當給定ω和kx,沿y軸正向可以解出四個有效特征值λj,及對應的四個特征向量αj。每個λj值對應一個有效波,即位移受多個波場的影響,且每個λj所對應的單位幅值的波的功率為

(23)

其中,Sj和Tj分別由式(14)和式(16)~(19)給出,式中?/?x由-ikx替代,?/?y由λj替代。一旦四個有效λj值確定,位移向量便可以寫成以下形式

(24)

式中:Aj對應第j個有效波的幅值;A0為入射波的幅值;λ0,α0為入射波在式(22)中的解。聯立等式(13)和(14),則邊界位移可以寫成

d=P1A+d0

(25)

t=P2A+t0

(26)

聯立式(25)、(26),邊界位移和力的關系式如下

(27)

式中:K為結構的波動剛度矩陣;f為入射波引起的動力載荷。

在圖5所示的半無限層合板殼耦合系統中,第n個子系統在連接處的切線方向與全局坐標軸yg的夾角為Φn,其局部坐標系中的耦合處邊界位移為dn,邊界力為tn。 設全局坐標系中連接處的位移dg,合力為tg。

圖5 半無限層合板殼耦合系統示意圖

根據耦合處位移連續性,dn與dg的轉換關系如下

(28)

(29)

式中,Rn為坐標轉換矩陣。根據耦合處力的平衡可得

(30)

式中,N為耦合系統中子系統總數。根據式(27)可知tn=Kndn-fn。當入射波作用在第m個子系統上時,式(30)可寫為

(31)

(32)

所以產生的各個波在傳遞時攜帶的能量為

(33)

式中,Pyj由式(23)給出。從而可以計算出入射波傳遞到結構耦合處所產生的各個波的能量傳遞系數為

τijmn=Pijmn/Pim

(34)

式中:Pim為在結構m中入射i型波的功率;Pijmn為結構m中的i型波傳遞到結構n中產生j型波的功率。當m≠n時,τijmn表示結構m中的i型波傳遞到結構n中產生j型波的能量透射系數;當m=n時,τijmm表示結構m中的i型波經過耦合處反射到結構m中產生j型波的能量反射系數。所有的能量傳遞系數總和為1。

2.3 耦合損耗因子CLF

式(23)代表波沿y方向傳播的功率流,同理可以計算出波沿x方向傳播的功率流Pxj。所以能量傳播角為

θe=arctan(Pyj/Pxj)

(35)

又總能量密度

(36)

所以波群速度為

(37)

方向為波數曲線的外法線方向,如圖6所示。

圖6 波與能量的傳播方向

群速度在y方向和θ方向的投影分別為

cgy=cgsinθe,cgθ=cgcos(θ-θe)

(38)

對于面積為Lx×Ly的層合板殼結構,每個模態頻率對應的網格面積為π2/(LxLy),所以圖6所示波數曲線與坐標軸圍成的面積包含的模態數為

(39)

模態密度為

(40)

又波速c=ω/k,式(38)中θ方向的波群速度也可以由cgθ=dω/dk得到,所以式(40)可變換為如下形式

(41)

所以給定ω時,模態密度隨入射角的變化為

(42)

相應的模態能量為

(43)

式中,E為結構振動總能量。所以在整個耦合邊界上,由結構m中入射i型波傳遞到結構n中轉換為j型波的總能量流為

(44)

式中,Ljunction為耦合邊界的長度。

在SEA理論中,qijmn=ηijmnωEim,ηijmn為耦合損耗因子。所以,對于圖5所示的線耦合結構,由式(44)可計算出能量傳遞系數和耦合損耗因子之間存在如下關系

(45)

式中,τijmn由式(34)給出。

當我們將能量傳遞系數應用于耦合系統的能量預報時,基于混響場假設需要對能量傳遞系數進行如下平均

(46)

當材料為各向同性時,式(46)便可簡化為我們常見的標準各向同性混響場平均傳遞系數計算公式

(47)

3 數值算例分析

3.1 頻散曲線

以石墨/環氧樹脂復合材料為例,其材料參數如表1所示。計算圖3所示的層合柱殼的頻散曲線并與Ghinet的結果進行對比驗證。各層材料方向角:[0°/45°/-45°/90°/-45°/45°/0°]。每層材料厚度為2 mm,單曲率環向開敞層合柱殼的曲率半徑為2 m。厚徑比為0.007,小于0.05,可視為薄殼。

表1 石墨/環氧復合材料參數

(a) θ=0°

3.2 模態密度

本節以三層碳纖維/環氧樹脂復合材料為例(每層鋪設角度如圖8所示),其材料參數如表2所示。分別就層合板和環向開敞層合柱殼的模態密度進行計算,并討論兩者的異同之處。每層材料厚度為1 mm,面積為1 m×1 m,柱殼的厚徑比為0.01。

圖8 三層碳纖維/環氧方向示意圖

表2 碳纖維/環氧復合材料參數

圖9中,在頻率較低時,環向開敞層合柱殼的彎曲波與剪切波模態密度要比相同表面積和厚度的層合板小,此時層合柱殼不存在縱波。曲率的引入相當于增大了層合柱殼的剛度,使得其模態數減少,即出現曲率的剛化效應,這與各向同性柱殼的特性類似。

(a) 彎曲波

隨著頻率增加,層合柱殼的彎曲波及剪切波模態密度也隨之急劇增加,并超過相對應的層合板模態密度。當頻率增大至環頻率fr=1 698 Hz時,層合柱殼的縱波開始出現,且遠大于層合板的縱波模態密度。在環頻率fr附近層合柱殼的彎曲波和剪切波模態密度出現峰值,此后層合柱殼的模態密度開始向層合板過渡,隨著頻率的增大兩者趨于一致。對于層合板的彎曲波,其模態密度基本不隨頻率變化,保持在0.041 67 (rad·s-1)-1,而層合板的面內剪切波和縱波模態密度則隨頻率呈線性增加,即式(41)中對應的面內波速度與群速度的乘積關于入射角的積分在頻率增大到一定值后開始收斂。此外,從以上三個對比結果來看,曲率半徑對剪切波的影響相較彎曲波和縱波小。

3.3 能量傳遞系數

以層合板與環向開敞層合柱殼沿x軸耦合為例,計算在層合板上分別入射彎曲波、剪切波、縱波時的能量傳遞系數,如圖10所示。耦合角度φ=90°,φ=|Φm-Φn|,結構都為三層碳纖維/環氧樹脂復合材料,每層材料厚度為1 mm,材料方向角為[0°/90°/0°],結構2的厚徑比為0.01。激勵頻率4 000 Hz,結構1,2的波數曲線分別如圖11和圖12所示。兩個結構的臨界值kx1、kx2基本重合,kx1=3.252 3為縱波的臨界值,kx2=11.930 8為剪切波的臨界值。kx4為結構2彎曲波在x方向的最大波數值,即曲線1,2的分界點。能量傳遞系數隨kx的變化情況如圖13所示。

圖10 層合板殼耦合系統

圖11 結構1的波數曲線

圖12 結構2的波數曲線

(a) 入射彎曲波

圖13(a)中入射波為彎曲波,為方便起見,采用雙y軸坐標進行繪圖,后續部分繪圖亦如此。f1、f2分別為圖12中結構2的彎曲波1、2。當kx接近kx1(θ=88.13°)時,縱波的透射和反射能量系數τfl12、τfl11快速衰減并在kx1處完全終止,而彎曲波能量系數處于相對平穩狀態。與此同時開始發生彎曲波向剪切波的波型轉換,并隨著kx的增大這種波型轉換逐漸增強。當kx增大到kx2(θ=83.09°)時,剪切波終止,僅反射和透射彎曲波,且反射強于透射。當kx為kx3=45.496 8(θ=58.16°)時,對應兩個不同的ky,即開始在兩個不同的方向同時透射彎曲波τff112、τff212(當頻率較大時不發生此現象)。在kx4=52.798 7(θ=42.87°)以后,不發生透射,為全反射狀態。

圖13(b)中入射波為剪切波,當kx

圖13(c)中入射波為縱波,可以發現整體上反射回層合板內的縱波占主導地位。當入射角較大時,透射波以彎曲波為主,在很大范圍內維持平穩。彎曲波的反射系數與透射系數變化規律一致,但從右側坐標刻度可以看出彎曲波反射系數很小。隨著kx的增大,反射的縱波逐漸衰減,當kx>3后,反射回層合板內的剪切波、透射到層合柱殼內的剪切波和縱波先后到達峰值,且峰值點依次增大。另一方面,也反映出當入射角較小時,透射能量占比較大,以面內波能量為主。

從圖14中可以看到能量傳遞系數之間滿足互易性,即τijmn=τjinm。

圖14 能量傳遞系數滿足互易性

綜上所述,對于耦合角度為90°且材料相同的層合板殼耦合系統,在連接處波的能量反射強于透射。另外,可以發現三個圖中的能量傳遞系數總和恒為1,滿足保守耦合系統的能量守恒,這從側面證明了計算的正確性,且能量傳遞系數滿足互易性。

將圖10所示耦合結構的材料退化為鋁合金進行計算,鋁合金參數如表3所示。所得結果如圖15所示,這與Langley的計算結果圖10完全一致。

表3 鋁合金材料參數

圖15 在結構1上入射彎曲波,ω=3 076.8 rad/s

本文還計算了如圖16所示三結構耦合的情況,結構材料及厚度與圖10所示的算例一致。激勵頻率為4 000 Hz,以在結構1上輸入彎曲波為例,對各子系統間的能量傳遞系數進行計算。如圖17所示,能量傳遞系數總和恒為1,同樣滿足了保守耦合系統的能量守恒。

圖16 三結構耦合

圖17 能量傳遞系數

3.4 耦合角度和激勵頻率對CLF的影響

本節在3.3節分析的基礎上計算系統的耦合損耗因子,層合板和柱殼的面積尺寸均為1 m×1 m,分別討論其隨耦合角度和頻率的變化。

圖18為耦合損耗因子隨耦合角度的變化,頻率為4 000 Hz。如圖18(a)所示入射波為彎曲波,以ηff11和ηff112為主。當φ≤15.75°和φ≥150.8°時,ηff11隨耦合角的增大而減小,當15.75°<φ<150.8°時,ηff11整體趨勢單調遞增,透射到層合柱殼中的彎曲波耦合損耗因子ηff112變化情況則與之相反。并且在此頻率下,出現彎曲波多角度透射,ηff212分別在耦合角φ=4.5°、171°出現峰值,耦合角接近90°時,其值較小。此外與面內剪切波、縱波對應的耦合損耗因子占比很小,剪切波的透反射耦合損耗因子ηfs12、ηfs11大于縱波的ηfl12、ηfl11,在0°~180°間都出現三個極值點,且在0°~180°值為0。

(a) 入射彎曲波

如圖18(b)所示,入射剪切波時,以ηss11和ηss12為主。 隨著耦合角度的增加ηss12先急劇下降后緩慢增加,在61.2°到達最小值,在180°值最大,此時除ηss11外其余耦合損耗因子近似為0。耦合角度為0°時,面內波之間存在波型轉換,不存在面內波到離面彎曲波的轉換。隨著耦合角度的增加,ηsl12逐漸減小至0,ηsl11占比最小,ηss11先增后減在28.8°~133.2°之間占比最大,ηsf11、ηsf112存在兩個波峰,且先后在75.6°和79.2°達到最低。

如圖18(c)所示,入射縱波時,可以發現各耦合損耗因子的變換情況與入射剪切波時類似,當φ≤18°和φ≥147.6°時以ηll12為主,其余耦合角度以ηll11為主,且ηlf112較大占比僅次于ηll11。

另外還發現無論耦合角度如何變化耦合損耗因子的總和都恒為定值,且入射彎曲波時耦合損耗因子比入射剪切波及縱波時小,即入射彎曲波時層合板殼之間的耦合作用明顯弱于入射面內波時的情況,這與各向同性材料耦合結構的規律一致。

圖19給出耦合損耗因子隨激勵頻率的變化,耦合角度φ=90°保持不變。從圖19三個結果圖來看,耦合損耗因子變化較為復雜的區域大致集中在5 000 Hz以前,且在這區域多數占比較大的耦合損耗因子衰減較快,在5 000 Hz以后衰減較為緩慢。當頻率低于環向開敞層合柱殼的環頻率fr=3 032.5 Hz時,不存在透射的縱波耦合損耗因子。

如圖19(a)入射彎曲波時,當頻率大于fr時,ηfl12開始出現隨后趨于穩定,ηfs12也開始出現,并在6 760 Hz處出現突變,隨后緩慢遞減。ηfl11基本保持恒定值,ηff212存在較大波動且主要集中在環頻率之前。總體上以彎曲波的透反射耦合損耗因子ηff11和ηff112為主,且ηff11占比最大。

入射剪切波時,如圖19(b)所示,主要以ηss11、ηss12、ηsf112為主,且反射的耦合損耗因子占比同樣大于透射。入射縱波時,如圖19(c)所示,ηll11全程占主導地位。入射面內波時,耦合損耗因子的曲線較為光滑,而入射彎曲波時,曲線存在波折。

整體而言無論入射何種波,耦合損耗因子總和伴隨頻率增大而逐漸衰減,頻率越大衰減越緩慢。簡言之,頻率越大子系統間的耦合作用越弱,且入射彎曲波時層合板殼之間的耦合作用同樣明顯弱于入射面內波時,這同樣與各向同性材料耦合結構的規律一致。

3.5 不同鋪設方式對能量總透射系數的影響

前文已對各種波的能量傳遞系數進行了詳細討論,故本節重點分析鋪設方式對總的能量傳遞的影響。同樣以3.3節中層合板與環向開敞層合柱殼90°耦合為例,設置結構1的材料鋪設角度為[0°/0°/0°]保持不變。分別入射彎曲波、剪切波和縱波,探究結構2在不同鋪設方式下對能量總透射系數的影響。結構1為三層碳纖維/環氧樹脂復合材料,結構2為四層碳纖維/環氧樹脂復合材料,選擇對稱鋪設為[α/-α/-α/α],反對稱鋪設為[α/-α/α/-α],鋪設角度為45°,每層材料厚度都為1 mm。

圖20~圖22計算了兩種鋪設方式下,三種入射波的能量總透射系數。從圖中可以發現,雖然不同鋪設方式對彎曲波能量的透射較面內剪切波和縱波的影響略大,但兩種鋪設方式的計算結果整體來看相差甚小,即能量總透射系數對鋪設方式的敏感度較低。計算還表明,鋪設角度為30°和60°時的計算結果與45°結果的規律相同。

4 結 論

本文以正交各向異性層合板和柱殼為研究對象,基于薄殼和層合板理論推導了位移為基本變量的平衡微分方程,通過計算環向開敞層合柱殼的頻散曲線與已有研究進行對比,驗證本文理論推導的正確性。然后根據波動剛度矩陣方法對耦合系統間的能量傳遞系數進行求解,并結合模態密度推導耦合處的耦合損耗因子(CLF)。數值算例分析表明:

(1) 在頻率較低時,曲率的引入相當于增大了層合柱殼的剛度,使得其模態數減少,即出現曲率的剛化效應。但在環頻率以后,隨著頻率的增大層合板與柱殼的模態密度趨于一致。此外曲率半徑對剪切波的影響相較彎曲波和縱波小。

(2) 能量傳遞系數總和恒為1,滿足保守耦合系統的能量守恒,且能量傳遞系數具有互易性。對于材料相同的層合板殼90°耦合系統,在連接處波的能量反射強于透射。此外發現能量總透射系數對鋪設方式的敏感度較低。

(3) 耦合損耗因子CLF的總和不受耦合角度影響,且隨頻率增大而逐漸衰減,當入射彎曲波時耦合損耗因子比入射剪切波及縱波時小,即入射彎曲波時層合板殼之間的耦合作用明顯弱于入射面內波的時候。

猜你喜歡
模態結構
《形而上學》△卷的結構和位置
哲學評論(2021年2期)2021-08-22 01:53:34
論結構
中華詩詞(2019年7期)2019-11-25 01:43:04
新型平衡塊結構的應用
模具制造(2019年3期)2019-06-06 02:10:54
論《日出》的結構
車輛CAE分析中自由模態和約束模態的應用與對比
國內多模態教學研究回顧與展望
高速顫振模型設計中顫振主要模態的判斷
航空學報(2015年4期)2015-05-07 06:43:35
創新治理結構促進中小企業持續成長
現代企業(2015年9期)2015-02-28 18:56:50
基于HHT和Prony算法的電力系統低頻振蕩模態識別
由單個模態構造對稱簡支梁的抗彎剛度
計算物理(2014年2期)2014-03-11 17:01:39
主站蜘蛛池模板: www精品久久| 狠狠色综合网| 国产欧美高清| 在线观看免费黄色网址| 天天综合色网| 国产香蕉97碰碰视频VA碰碰看| 国产视频a| 亚洲人成色在线观看| 毛片免费在线视频| 专干老肥熟女视频网站| AV天堂资源福利在线观看| 久久婷婷五月综合色一区二区| 欧美日韩第二页| 热热久久狠狠偷偷色男同| 亚洲AⅤ波多系列中文字幕| 97超爽成人免费视频在线播放| 国产丝袜无码一区二区视频| 色爽网免费视频| 亚洲午夜片| 熟妇丰满人妻| 好久久免费视频高清| 九九九国产| 亚洲男人的天堂视频| 久久国产高清视频| 亚洲一区二区三区在线视频| 日本午夜三级| 精品国产成人三级在线观看| 国产成人精品一区二区| 一区二区自拍| 91精品啪在线观看国产91| 亚洲娇小与黑人巨大交| 91人妻日韩人妻无码专区精品| 伊在人亚洲香蕉精品播放| 欧美精品一区在线看| www.狠狠| 国产精品55夜色66夜色| 99久久亚洲综合精品TS| 国产91在线免费视频| 成年人福利视频| 丁香亚洲综合五月天婷婷| 亚洲日本中文字幕乱码中文| 欧洲日本亚洲中文字幕| 国产精品手机在线播放| 国产不卡一级毛片视频| 3D动漫精品啪啪一区二区下载| 国产一级小视频| 亚洲精品大秀视频| 凹凸精品免费精品视频| 999精品视频在线| 久久精品最新免费国产成人| 国产尤物视频在线| 久久永久免费人妻精品| 久久久波多野结衣av一区二区| 欧日韩在线不卡视频| 免费无码又爽又刺激高| 青青青视频免费一区二区| 国产日韩欧美在线播放| 国产欧美日韩另类| 国产一区成人| 国产精品区视频中文字幕| 亚洲精品午夜天堂网页| 中文字幕久久波多野结衣| 中文字幕有乳无码| 熟女日韩精品2区| 国产在线拍偷自揄拍精品| 国产女人水多毛片18| 国产精品视频白浆免费视频| 国产菊爆视频在线观看| 久热精品免费| 亚洲成人在线免费观看| 亚洲人成电影在线播放| 99热这里只有成人精品国产| 国产精品久久精品| 毛片网站免费在线观看| 欧美在线一二区| 国产亚洲精久久久久久久91| 一级全黄毛片| 无码免费的亚洲视频| 国产精品偷伦视频免费观看国产| 亚洲免费人成影院| 免费毛片全部不收费的| 2021国产乱人伦在线播放|