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

壓電纖維致動柔性結構的仿魚體波振動特性及流場分布研究

2024-06-17 00:00:00溫志偉婁軍強陳特歡崔玉國魏燕定李國平
中國機械工程 2024年3期

摘要:提出利用壓電纖維(MFC)致動柔性結構的二階振動模態來模仿鲹科魚類的魚體波運動。基于假設模態法求解智能柔性結構的前三階振型,實驗測得水下柔性結構的二階振型,驗證了智能柔性結構的二階振型逼近鲹科鱸魚的魚體波。采用計算流體動力學(CFD)分析了柔性結構周圍流線及壓力分布情況,結果表明:二階振動行為下,柔性結構節點前后流線方向相反;兩側流場始終存在兩組高、低壓集中區域,且節點前后的壓力集中區域分布相反;柔性結構兩側始終存在正壓梯度,其中向前的壓力分量為柔性結構提供持續的推進力,而節點前后壓力分量產生的側向力方向相反,部分抵消,增強了柔性結構的側向穩定性。研究結果為水下仿生推進器的設計和性能分析提供了重要參考。

關鍵詞:水下柔性結構;仿魚體波振動;壓電纖維;流場分布特性;計算流體動力學

中圖分類號:TP242;O351.2

DOI:10.3969/j.issn.1004132X.2024.03.003

0引言

由于與海洋工程相關的需求逐漸增加,如海底光纜檢修、海底資源勘探、水下隱蔽軍事偵察以及海洋環境保護等[1],水下推進器或航行器已經引起了廣泛關注。在長期的自然進化過程中,魚類進化出了高效的游動性能、良好的靈活性和機動性。與傳統螺旋槳推進的水下航行器相比,模仿魚類游動方式的仿生推進器具有高靈活性、低噪聲、高能效等優勢,成為水下航行器領域備受關注的研究熱點[2]。

近85%的魚類采用身體/尾鰭(bodyorcaudalfin,BCF)游動模式[3],因此,大部分仿生推進器都圍繞BCF游動模式展開。電機驅動的BCF仿生推進器配備高功率的電機驅動系統,具有速度快、機動性能高及負載能力強等優點。然而復雜的結構為此類仿生推進器帶來能耗高、噪聲大、柔度低、控制復雜以及難以小型化等缺點。近年來,具有低噪聲、體積小、柔順度高等優點的智能材料如離子基聚合物材料(ionicpolymermetalcomposite,IPMC)、形狀記憶合金(shapememoryalloys,SMA)、壓電材料(piezoelectricceramictransducer,PZT)以及介電高彈體(dielectricelastomer,DE)等為仿生推進器的設計提供了新思路[4-5]。目前,越來越多的研究者開始關注使用智能材料驅動的微型水下推進系統[6]。為了模擬BCF魚類的靈活游動,許多研究人員將智能材料驅動的懸臂梁作為水下推進系統的核心組件。雖然懸臂梁模型無法完全復制BCF魚類的復雜運動,但相關研究工作對微型水下推進器的設計、實現和評估具有重要意義。薄片式壓電纖維復合材料(macrofibercomposite,MFC)作為一種新型智能材料,具有柔度高、響應快和驅動位移大等優勢[7],越來越多地應用于BCF仿生推進器的設計。CEN等[8]開發了雙MFC驅動的BCF推進器,并展現了其優越的游泳性能。官源林等[9]提出了一種MFC作為驅動材料的復合型魚尾,該魚尾具有較大的擺動位移。HU等[10]提出了一種使用MFC驅動的柔性尾鰭微型機器魚,并通過實驗測量了其推進力,證實了MFC致動柔性結構在機器魚和仿生推進器領域的優勢和前景。但是,目前大部分微小型仿生推進器仍局限于利用智能材料驅動柔性結構往復擺動來模擬魚類的BCF運動模式[11],而對BCF魚類游動過程中的魚體波行為和運動姿態缺乏充分關注,導致目前的仿生推進器與被模仿魚類的游動性能還存在較大差距,需要深入研究。

本文根據鲹科鱸魚的游動方式和特點,提出利用壓電纖維致動柔性結構的二階振動模態來模仿鲹科魚類的魚體波運動。基于假設模態法分析了MFC致動柔性結構的多階振動特性,然后通過多點實驗測得水下MFC致動柔性結構的實際二階振型。采用計算流體動力學(CFD)分析了仿魚體波振動柔性結構周圍流場的速度流線和壓力集中區域的分布特性和演化過程,最后計算得到MFC致動柔性結構仿魚體波振動下產生的推進力和側向力變化情況。

1BCF推進模式及魚體波

生物學家通過大量的魚類游動行為觀察實驗,根據軀干參與波動的長度比例將魚類BCF游動模式劃分為鰻鱺科模式、亞鲹科模式、鲹科模式和鮪科模式四類。進一步發現BCF魚類在游動過程中主要依靠其肌肉組織使脊椎發生反復彎曲變形,從而形成一個從魚體頭部傳遞到尾部的行波[12]。魚類游動過程中的脊椎彎曲曲線即為“魚體波”。BCF鲹科模式具有游動速度高、穩定性好且柔度適中的優點。鲹科推進模式魚類的游動過程如圖1所示,鲹科魚類在游動推進過程中頭部基本保持不動,可視為剛體,主要靠身體后半部分的柔性軀干彎曲變形獲得較大推進力和較快推進速度,是近年來水下仿生推進器重點研究和模仿的推進方式[13]。

魚體波曲線近似可由魚體波幅值包絡線和傅里葉級數合成得到,魚類的游動性能與魚體波波數、擺動幅值以及擺動頻率密切相關。鲹科魚類的魚體波函數可表示為[14]

其中,x軸坐標原點可視為魚類剛體頭部和柔性軀體的連接點,a(x)為魚體波包絡方程,a0、a1、a2為包絡線系數,k為魚體波數,ω=2πf為擺動角頻率,f為擺動自然頻率,相關魚體波參數主要取決于魚類身體形狀特征、游速游姿以及其游動狀態,主要通過生物學測量得到。圖2所示為鲹科魚類中的典型代表——鱸魚在一個周期T內不同時刻的魚體波曲線,圖2中BL為柔性軀干體長(bodylength,BL)。

2MFC致動柔性結構的振動特性

由于鲹科魚類在推進過程中頭部基本保持不動,主要靠尾部帶動尾鰭擺動產生推進作用,故許多研究者都采用經典的懸臂式Euler-Bernoulli梁模型來模擬鲹科魚類的BCF推進模式[15]。本文中MFC致動水下柔性結構示意圖見圖3,柔性鋁梁基體一端固定,一端自由,且梁基體根部上下兩面對稱粘貼一對MFC致動器。本文采用型號為M5628-P1的壓電纖維復合材料作為MFC致動柔性結構的致動器,其中P1型對應d33壓電效應,即電場方向與壓電纖維變形方向一致。基于d33壓電效應,MFC致動器在外部驅動電壓作用下進行伸長或收縮從而帶動鋁梁基體彎曲振動。在致動器與柔性梁基體完美粘貼且小變形的假設前提下,可將MFC致動柔性結構簡化為Euler-Bernoulli梁模型,其動力學方程可表示為

式中,w(x,t)為柔性結構的橫向振動位移;Cv為柔性結構的總阻尼比,包括結構阻尼和流體黏性阻尼;Fdrive(x,t)為MFC致動器施加在柔性結構上的驅動力;Fhydro(x,t)為柔性結構運動過程中和周圍流體的相互作用力;Y(x)I(x)、ρ(x)A(x)分別為MFC致動柔性結構在x處橫截面的等效抗彎剛度和單位質量密度。

由于MFC致動柔性結構可分為粘貼MFC致動器的主動部分和未粘貼MFC的從動部分,簡化后的MFC致動柔性結構的梁模型如圖4所示,兩部分的相應參數可表示為

式中,l為結構長度;Y、I、ρ、A分別為彈性模量、截面極慣性矩、體密度和橫截面面積;下標ac和pa分別表示主動部分和從動部分;下標m和b分別表示MFC致動器和鋁基體。

復合梁在水下和空氣中的模態振型基本一致[16],本文將重點研究自由振動狀態下MFC致動柔性結構的振動特性。采用變量分離法,MFC致動柔性結構主被動部分的振型函數可表示為[17-18]

式中,Φi(x)為柔性結構的第i階模態振型函數;Qi(t)為第i階模態下的廣義模態坐標。

在無阻尼自由振動條件下,MFC致動柔性結構的動力學方程(式(2))的解可表示為

式中,ωi為柔性結構的第i階固有頻率。

MFC致動柔性結構一端固定、一端自由,其邊界條件可表示為

相應地,MFC致動柔性結構在主動、被動部分連接處的連續條件可表示為

施加質量歸一化條件,并令特征系數矩陣Mi的行列式為零,即可得到MFC致動柔性結構的各階固有頻率以及對應的主被動部分模態振型函數。

MFC致動器和柔性鋁梁的基本參數見表1,圖5所示為MFC致動柔性結構的前三階歸一化模態振型,圖5中還給出了勻質同一截面梁的前三階振型作為對比。由圖5可以看出,MFC致動柔性結構的第二階振型與圖2中鲹科鱸魚的魚體波曲線較為相似,故采用MFC致動柔性結構在二階固有頻率下的振動行為來模仿和研究鲹科鱸魚的魚體波特性。

3測控平臺的搭建及實驗結果

為了測試MFC致動柔性結構的水下振動特性,搭建圖6所示的測控平臺。MFC致動柔性結構一端夾緊,置于760mm×460mm×500mm的水箱中。實驗過程中控制中心PC機發出驅動信號,經功率放大器(Trek-PZD700A)放大傳輸至MFC致動器上,驅動柔性結構產生振動。同時,激光位移傳感器(KeyenceLK-G80)實時測量MFC致動柔性結構的水下振動位移,并將獲得的柔性結構振動信息傳回至控制中心。整個實驗平臺基于LabVIEW軟件構建,并使用NIcDAQ-9178數據采集機箱和D/A模塊(NI-AO9263)與A/D(NI-AI9205)模塊完成數據采集和通信。

為獲得MFC致動柔性結構的水下固有頻率,對MFC致動器施加峰峰值800V、頻率范圍0.1~35Hz的正弦掃頻信號。整個掃頻時間持續20s,采樣頻率設置為2000Hz。實驗測得MFC致動柔性結構末端振動位移的時域、頻域響應如圖7所示,初步確定MFC致動柔性結構的前兩階水下固有頻率分別為0.78Hz和5.4Hz。

為了精確獲得MFC致動柔性結構在二階固有頻率下的振型曲線,實驗中以柔性結構水下二階固有頻率為激勵頻率,測得了MFC致動柔性結構上不同監測點的水下振動位移響應,其中9個監測點P1~P9的水下振動位移峰峰值依次為0.720,1.822,3.194,4.262,3.394,1.578,0.650,1.778,5.390mm。圖8所示為9個監測點水下振動位移的時域響應曲線,可以看出監測點P1~P9的水下振幅呈現出先增大后減小、最后又逐漸增大的特點。

以t=20.089s為柔性結構水下振動的起始時刻,圖9進一步給出了各監測點在半個周期內、五個時序時刻振動位移的變化情況,五個時序時刻分別為t0=20.089s、t1=t0+(1/8)T、t2=t0+(2/8)T、t3=t0+(3/8)T、t4=t0+(4/8)T。由圖9可以看出,在t3即1/4周期時刻,MFC致動柔性結構上各監測點的振動位移基本為零,且在約160mm處的監測點,其振動位移在任意時刻都基本接近零,即為MFC致動柔性結構二階振型節點。采用三次樣條曲線對采集得到的監測點在五個時序時刻的振動位移分別進行擬合,圖9還給出了擬合得到的MFC致動柔性結構在不同時刻的連續變形曲線。可以看出實驗測得的MFC致動柔性結構的變形曲線與鱸魚不同時刻的魚體波曲線較為接近,證實了利用MFC致動柔性結構在二階固有頻率下的振動變形來模擬鲹科魚類魚體波變形的可行性。

4柔性結構二階振動下流場CFD分析

為分析MFC致動柔性結構在二階固有頻率下周期振蕩行為引起的周圍流場變化情況,進一步采用ANSYS/Fluent對水下二階振蕩柔性結構進行CFD流固耦合仿真分析。仿真過程中采用低雷諾數的k-ε湍流模型進行Navier-Stokes方程的數值求解[19]。在建立MFC致動柔性結構初始幾何模型的基礎上,為減少壁面效應并提高計算精度,以MFC致動柔性結構末端的實測振幅為基準參考量,設置流體計算域為650mm×400mm×300mm。網格劃分采用適應性強的非結構化網格和動網格技術,模型運動通過加載UDF實現,計算網格劃分結果如圖10所示,最終計算域劃分為441537個網格單元和83836個網格節點。

仿真過程中采用MFC致動柔性結構在二階水下固有頻率5.4Hz下的實測結果開展CFD流固耦合分析。為了描述柔性結構周圍流場的分布演化情況,定義柔性結構末端擺動到下側極限位置的時刻相位為0°;而當柔性結構向上擺動到中間平衡位置時,向上速度最大時刻對應的相位為90°,上極限位置時刻對應的相位為180°。另外,設定柔性結構末端從下極限位置擺動到上極限位置的過程為外擺階段(0°→180°),反之,則為回擺階段(180°→360°)。由于MFC致動柔性結構在啟動和加速階段的流場處于不穩定狀態,流場特征復雜紊亂,故選取第20個周期的穩定流線矢量圖來描述MFC致動柔性結構周圍流場結構的分布特征和變化過程,并截取流域中間平面(y=0)進行流場特性分析。

圖11所示為柔性結構周圍流場瞬態速度流線在一個擺動周期內的穩定變化情況,從0°相位開始,以45°為間隔,共8個典型相位。可以看出:由于MFC致動柔性結構做往復連續振動,其表面為運動邊界,故觀察到既有流線從其表面流出,也存在流線運動終止在其表面上的情況。MFC致動柔性結構往復振動過程中,周圍流體施加在振蕩柔性結構上的反作用力可分解為作用在柔性結構表面法線方向上的流體壓力和沿著其表面切線方向的黏滯曳力[20],柔性結構擺動導致了周圍流體的速度流線發生改變,從圖11的8幅分圖中觀察發現運動柔性結構表面向外推水,表現為其表面流體流線流出;而其另外一側表面吸水,表現為流體流線流入。需要指出的是,由于水下柔性結構二階共振狀態下存在著節點,而節點前后柔性結構上各點振動方向相反,導致柔性結構周圍流線以節點為界,前后流動方向相反。

MFC致動柔性結構的往復振蕩行為對周圍流體產生了交替變化的吸力和推力[21],誘導了流體的非定向運動,從而在柔性結構周圍出現了分布交替變化的壓力集中區域。圖11中流場速度流線分布時刻對應的流場壓力分布云圖見圖12。可以看出:在圖中任意時刻,柔性結構兩側存在著兩組高、低壓集中區域,且節點前后的高/低壓力集中區域分布剛好相反,這是因為柔性結構節點前后各點振動方向相反。顯然,在圖12中的任意擺動時刻,MFC致動柔性結構兩側始終存在著正壓梯度。具體而言,在前半個擺動周期內(270°→90°),擺動柔性周圍壓力分布仍以節點為界,前段壓力梯度向前向上,而后段存在著向前向下的壓力分量。向前的壓力分量為MFC致動柔性結構提供持續的推進力。而后半個擺動周期內(90°→270°),節點前段出現了向前向下的壓力分量,而后端的壓力分量向前向上。顯然,方向呈周期變化的壓力分量提供了側向力促使柔性結構向兩側平動,而且同一時刻節點前后壓力分量產生的側向力方向相反,部分抵消,增強了柔性結構擺動過程中的側向穩定性。

圖13進一步給出了在二階固有頻率下MFC致動柔性結構振動過程中所受的推進力Fx和側向力Fz隨時間的變化情況。由圖13可見,推進力和側向力呈典型周期性變化,且都在t=0和t=T/2兩個時刻附近達到最大值。結合壓力云圖12發現:這兩個相位時刻附近(0°和180°)流速最低,形成的壓力集中區域最大,故而誘導產生的水動力最大,得到的瞬時最大推進力和側向力分別為3.2mN和182.2mN,見圖13。對比圖13a和圖13b后發現,由于MFC致動柔性結構在外擺和回擺過程中交替產生旋向相反的渦結構,產生的推進力在一個穩定擺動周期內出現兩次波峰起伏,故推進力Fx的變化頻率是側向力Fz變化頻率的兩倍。

由圖13可以看出:Fx在大部分擺動周期內取值均小于零(方向為-x方向),故可以持續產生推進效果,如圖13a中藍色區域所示。而圖13b中側向力Fz關于x軸對稱,在一個擺動周期內積分為零,保證了柔性結構擺動過程中的側向穩定性[22]。為了定量表征MFC致動柔性結構產生的推進力大小,利用梯形數值積分計算得到其在一個穩定擺動周期T內產生的平均推進力:

顯然,數值積分的結果就是圖13a中一個周期內藍色區域與紅色區域面積之差,最終得到MFC致動柔性結構在二階固有頻率下的周期平均推進力為1.5mN,單位時間內產生的平均推進力為8.1mN,證實了MFC致動柔性結構的推進效果。

5結語

本文借鑒鲹科鱸魚的推進方式和游動特點,提出利用壓電纖維致動柔性結構的二階振動模態來模仿鲹科魚類的魚體波運動。通過測量水下MFC致動柔性結構上8個監測點在二階水下固有頻率下的振動響應,得到MFC致動柔性結構的水下二階振型。最后采用CFD數值計算分析了柔性結構二階振動下的振動特性和周圍流場分布情況。

CFD仿真結果表明:在仿魚體波振動下,柔性結構表面既有流線流出也有流線終止,且節點前后流線方向相反。在任一時刻柔性結構兩側流場始終存在兩組高壓、低壓集中區域,且節點前、后的壓力集中區域分布相反;同時柔性結構兩側始終存在正壓梯度,其中向前的壓力分量為柔性結構提供持續的推進力,而節點前后壓力分量產生的側向力方向相反,部分抵消,增強了柔性結構的側向穩定性。研究結果為仿生推進器的仿魚體波設計和性能分析提供了參考。

參考文獻:

[1]HES,WANGN,HOM,etal.DesignofaNewStressWaveCommunicationMethodforUnderwaterCommunication[J].IEEETransactionsonIndustrialElectronics,2021,68(8):7370-7379.

[2]ZHONGQ,ZHUJ,FISHFE,etal.TunableStiffnessEnablesFastandEfficientSwimminginFish-likeRobots[J].ScienceRobotics,2021,6(57):eabe4088.

[3]XIEF,LIZ,DINGY,etal.AnExperimentalStudyontheFishBodyFlappingPatternsbyUsingaBiomimeticRobotFish[J].IEEERoboticsandAutomationLetters,2020,5(1):64-71.

[4]ZHAOQ,LIUS,CHENJ,etal.Fast-movingPiezoelectricMicro-roboticFishwithDoubleCaudalFins[J].RoboticsandAutonomousSystems,2021,140:103733.

[5]邱志成,陳思文.移動雙柔性梁系統的振動主動控制[J].振動測試與診斷,2022,42(1):62-67

QIUZhicheng,CHENSiwen.ActiveVibrationControlofaTranslationalDoubleFlexibleBeamSystem[J].JournalofVibration,Measurementamp;Diagnosis,2022,42(1):62-67.

[6]郭健,潘彬彬,崔維成,等.基于智能材料的深海執行器及海洋仿生機器人研究綜述[J].船舶力學,2022,26(2):301-313.

GUOJian,PANBinbin,CUIWeicheng,etal.ReviewofDeep-seaActuatorsandMarineBionicRobotsBasedonIntelligentMaterials[J].JournalofShipMechanics,2022,26(2):301-313.

[7]KUANGY,ZHUM.EvaluationandValidationofEquivalentPropertiesofMacroFibreCompositesforPiezoelectricTransducerModelling[J].CompositesPartB:Engineering,2019,158:189-197.

[8]CENL,ERTURKA.Bio-inspiredAquaticRoboticsbyUntetheredPiezohydroelasticActuation[J].Bioinspirationamp;Biomimetics,2013,8(1):016006.

[9]官源林,李華峰,楊熙鑫.基于MFC的復合型仿生魚尾的振動性能[J].振動測試與診斷,2013(增刊2):165-168.

GUANYuanlin,LIHuafeng,YANGXixin.VibrationPerformanceofCompoundBionicFishtailBasedonMFC[J].JournalofVibration,Measurementamp;Diagnosis,2013(S2):165-168.

[10]HUD,LOUJ,CHENT,etal.MicroThrustMeasurementExperimentandPressureFieldEvolutionofBionicRoboticFishwithHarmonicActuationofMacroFiberComposites[J].MechanicalSystemsandSignalProcessing,2021,153:107538.

[11]DURAISAMYP,KUMARSR,NAGARAJANSM.Design,Modeling,andControlofBiomimeticFishRobot:aReview[J].JournalofBionicEngineering,2019,16(6):967-993.

[12]CUIZ,YANGZ,SHENL,etal.ComplexModalAnalysisoftheMovementsofSwimmingFishPropelledbyBodyand/orCaudalFin[J].WaveMotion,2018,78:83-97.

[13]ZHOUC,NAHAVANDIS,GUN,etal.AnalysisandImplementationofSwimmingBackwardforBiomimeticCarangiformRobotFish[J].IFACProceedingsVolumes,2008,41(2):3082-3086.

[14]MACIASMM,SOUZAIF,BRASILACPJr,etal.Three-dimensionalViscousWakeFlowinFishSwimming-aCFDStudy[J].MechanicsResearchCommunications,2020,107:103547.

[15]YANB,YUN,MAH,etal.ATheoryforBistableVibrationIsolators[J].MechanicalSystemsandSignalProcessing,2022,167:108507.

[16]KWONYW,PRIESTEM,GORDISJH.InvestigationofVibrationalCharacteristicsofCompositeBeamswithFluid-structureInteraction[J].CompositeStructures,2013,105:269-278.

[17]楊斌強,徐文潭,陸國麗,等.寬頻壓電振動能量采集器的分布參數模型與實驗[J].中國機械工程,2017,28(2):127-134.

YANGBinqiang,XUWentan,LUGuoli,etal.DistributedParameterModelandExperimentsofaBroadbandPiezoelectricVibrationEnergyHarvester[J].ChinaMechanicalEngineering,2017,28(2):127-134.

[18]胡凱明,李華.軸向受壓梁非線性隨機最優電壓有界控制[J].浙江大學學報(工學版),2020,54(5):940-946.

HUKaiming,LIHua.NonlinearStochasticOptimalVoltageBoundedControlforAxialCompressedBeam[J].JournalofZhejiangUniversity(EngineeringScience),2020,54(5):940-946.

[19]袁巧玲,計時鳴,文東輝,等.基于改進的低雷諾數湍流模型的軟性磨粒流加工仿真與實驗[J].中國機械工程,2014,25(6):800-807.

YUANQiaoling,JIShiming,WENDonghui,etal.SimulationandExperimentofSoftAbrasiveFlowMachiningBasedonImprovedLow-Reynolds-numberTurbulenceModel[J].ChinaMechanicalEngineering,2014,25(6):800-807.

[20]GUPTAS,AGRAWALA,HOURIGANK,etal.AnguilliformandCarangiformFish-inspiredHydrodynamicStudyforanUndulatingHydrofoil:EffectofShapeandAdaptiveKinematics[J].PhysicalReviewFluids,2022,7(9):094102.

[21]THEKKETHILN,SHARMAA,AGRAWALA.Self-propulsionofFishes-likeUndulatingHydrofoil:aUnifiedKinematicsBasedUnsteadyHydrodynamicsStudy[J].JournalofFluidsandStructures,2020,93:102875.

[22]ZHOUH,HUT,LOWKH,etal.Bio-inspiredFlowSensingandPredictionforFish-likeUndulatingLocomotion:aCFD-AidedApproach[J].JournalofBionicEngineering,2015,12(3):406-417.

主站蜘蛛池模板: 国产麻豆精品久久一二三| 久久免费精品琪琪| 婷婷久久综合九色综合88| 国内精品手机在线观看视频| 欧美在线免费| 久久99精品久久久久纯品| 免费黄色国产视频| 福利国产在线| 国产精品va免费视频| 亚洲欧洲天堂色AV| a欧美在线| 曰AV在线无码| 91在线精品麻豆欧美在线| 亚洲欧美日韩久久精品| 麻豆国产原创视频在线播放| AV片亚洲国产男人的天堂| 在线亚洲精品自拍| 少妇精品网站| 国产黑丝视频在线观看| 亚洲伦理一区二区| 亚洲国产中文精品va在线播放| 国产原创自拍不卡第一页| 成AV人片一区二区三区久久| 中文字幕伦视频| 日韩成人午夜| 久久精品中文字幕少妇| 国产精品高清国产三级囯产AV| jizz国产视频| 国产特级毛片aaaaaa| 国产精品浪潮Av| 黄色网址手机国内免费在线观看| 国产精品美人久久久久久AV| 朝桐光一区二区| 亚洲成a人片在线观看88| 国产精品网址在线观看你懂的| 91在线中文| 第一页亚洲| swag国产精品| 精品国产成人av免费| 精品少妇人妻av无码久久| 五月激激激综合网色播免费| 欧洲亚洲欧美国产日本高清| 日韩a级毛片| 精品视频一区在线观看| 亚洲综合激情另类专区| 婷婷色中文| 爱做久久久久久| 精品国产Av电影无码久久久| 国产男人的天堂| 中国成人在线视频| 午夜国产精品视频| 色婷婷亚洲综合五月| 欧美成人一级| 国产最新无码专区在线| 国产69精品久久久久孕妇大杂乱 | 狠狠色丁香婷婷综合| 亚洲一区二区无码视频| 国产一区亚洲一区| 国产a网站| 久久综合婷婷| 亚洲二区视频| 国产午夜无码专区喷水| 中文字幕首页系列人妻| 国产精品制服| 亚洲一区二区视频在线观看| 亚洲久悠悠色悠在线播放| 国产福利观看| 91香蕉视频下载网站| 亚洲综合久久成人AV| 精品99在线观看| 国产精品无码一区二区桃花视频| 国产精品嫩草影院视频| 亚洲黄色成人| 88国产经典欧美一区二区三区| 国产精品女主播| 国产精品网址在线观看你懂的| 午夜福利在线观看入口| 国产精品任我爽爆在线播放6080| 女人18毛片一级毛片在线 | 国产日本欧美亚洲精品视| 国产精品香蕉在线| 国产综合另类小说色区色噜噜|