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

復合壓電水平懸臂梁振動俘能特性試驗及機電耦合模型研究

2023-01-03 04:37:00薛曉敏馬強力王佳佳
振動與沖擊 2022年24期
關鍵詞:振動模型

薛曉敏, 馬強力, 王佳佳

(西安交通大學 人居環境與建筑工程學院,西安 710049)

近年來,由于壓電材料優良的機電轉換性能,其振動能量回收能力受到了極大的關注[1]。為了盡可能提高壓電振動能量回收的能力,眾多學者提出了各種提高壓電俘能效率的方案,包括研發新型壓電材料與結構、能量收集接口電路和存儲電路優化設計等[2-4]。針對壓電材料俘能結構而言,最為常見的即為單、雙晶壓電陶瓷懸臂梁形式[5]。例如,Fang等[6]設計了一種自由端部設置鎳質量塊的壓電陶瓷懸臂梁俘能器,可以從1g振動環境下俘獲2.16 μW的電能,且顯示振動共振頻率為609 Hz。Shen等[7]設計了類似的壓電陶瓷懸臂梁俘能器,為了降低共振頻率,在自由端部設置了硅質量塊,平均收集功率和功率密度分別為0.32 μW和416 μW/cm3,共振頻率降低到183.8 Hz。總體而言,壓電懸臂梁俘能結構振動能量回收能力還非常微弱,僅為若干毫瓦;且結構振動主頻過高,這與環境振動頻率的低頻現象相互矛盾(振動環境一般一百赫茲以內甚至更低)[8]。此外,現有試驗研究主要著眼于壓電陶瓷俘能器的能量回收能力,主要關注采用何種方法或技術使得俘能器俘獲最大電能,然而,有關于力-電物理特性研究甚少涉及,實際上,機電耦合特性是振動能轉換電能的本質內在屬性,其深入研究對了解壓電俘能結構的內在特性并進一步提高能量轉換效率至關重要。

另外,值得注意的是,俘能器的電學能量輸出與振動能量輸入的關系復雜,且隨環境振動頻率、電路阻抗、阻尼比等因素均有密切關系[9]。為清楚描述壓電懸臂梁俘能結構機電關系,很多學者開展了相關模型研究。早期模型研究采用了集中參數下單自由度解法[10]預測壓電能量收集器耦合系統的動力學特性,該方法僅局限單個振動模態,無法精確描述多物理場耦合問題。作為一種改進的方法,Sodano等[11-12]基于Euler-Bernoulli理論提出Rayleigh-Ritz離散形式可以對無限自由度分布系統描述為有限自由度離散系統。最近,Erturk等[13]改進了雙壓電晶片懸臂模型,并找出了雙壓電晶片配置的串聯和并聯連接的閉合形式解析解,并通過與有限元仿真對比,結果顯示與解析解吻合。經試驗數據對比顯示,以上模型驗證可較好預測能量回收能力,即,電功率的峰值或總體均值。實際上,模型驗證還應該對壓電俘能器的本質特性——機電耦合特性進行對比驗證,應當反映實際系統機電耦合特性的重要現象。唯有此,才能充分說明模型的正確性和精確性。然而現有模型還甚少進行機電耦合特性試驗驗證,其精確性和有效性還有待于在試驗或實際工程中獲得進一步驗證。

鑒于以上背景,本文將對壓電單晶懸臂梁俘能片的能量轉換特性進行相關試驗和理論研究。設計并搭建復合壓電懸臂水平梁振動電能回收電路系統,并在懸臂梁自由端配重質量塊,以此改變俘能結構振動特性。通過試驗觀測,探討外接電路負載、激勵頻率和質量塊質量對俘能效率的影響規律,為壓電俘能器優化設計提供基礎依據;通過試驗觀測,探究不同加載條件下的壓電俘能系統的機電耦合特性及其變化,為后續模型驗證提供重要基礎;以試驗研究結果為基礎,提出有限元模型和機電耦合數學模型用以模擬壓電懸臂梁的能量收集重要耦合特性;最后,將仿真結果與試驗結果進行對比,驗證了所提機電耦合數學模型的有效性。

本文的試驗和理論研究對進一步挖掘壓電材料振動-電能收集能力具有重要意義,為壓電俘能結構及優化電路設計提供重要參考。

1 壓電單晶懸臂薄梁俘能性能試驗

PPA(piezo protection advantage)是美國美德公司以壓電陶瓷為原材料制成的商業化系列產品,包括作動器、傳感器、俘能器等。其中1001型號是最為典型的壓電單晶薄片的振動俘能器產品之一,如圖1所示。該壓電片為多種材料復合黏結而成,除壓電陶瓷外還包括聚酯、金屬銅、不銹鋼和聚酰亞胺,如表1所示。以上多種材料的復合有效改善了壓電陶瓷脆性特性并適當提高結構剛度。

圖1 PPA壓電俘能片實物Fig.1 PPA piezoelectric harvester patch

表1 PPA壓電俘能片結構及參數

1.1 試驗方法

為探究壓電俘能片的俘能特性,本文開展其懸臂梁結構的振動能量回收電能的特性試驗研究。試驗裝置實物如圖2(a)所示。激振路徑為:通過信號發生器輸送激勵信號,經功率放大器對信號放大,繼而輸入到激振器用以激振俘能片進行電能俘能。試驗數據測試對象包括振動應變(懸臂梁固定端中心處)、俘獲電壓和電流,均由數據采集器進行采集。本試驗中采用的儀器如表2所示。

表2 壓電懸臂梁俘能試驗儀器

1.1.1 電路設計及測試

根據電功率公式(即P=i·v,其中:v和i分別為電壓和電流;P為功率)設計振動能量俘獲電路:PPA與電阻箱串聯形成回路,數據采集器電壓通道與電阻箱并聯用以測量回路電壓,數據采集器電流通道與電阻箱串聯用以測量回路電流。根據數據采集器采集的電流和電壓即可換算得出PPA俘獲的電能功率。

1.1.2 結構設計及測試

在激振器頂桿通過夾具將PPA俘能片夾持,形成懸臂梁俘能器結構,如圖2(b)所示。此外,為了探究結構振動特性(自振頻率)的變化與電能回收能力的關系,還在懸臂梁自由端配重不同質量的質量塊(0,12.5 g,18.5 g和25 g)。

圖2 PPA懸臂梁振動俘能試驗系統Fig.2 Test system for energy harvest from PPA cantilever

1.2 試驗測試結果

基于1.1節試驗系統,本文考慮不同的加載工況:電路負載為0~90 kΩ(間隔為5 kΩ),簡諧激勵頻率為0~150 Hz(間隔5 Hz),自由端配重0,12.5 g,18.75 g,25 g。典型工況(電路負載30 kΩ,激勵頻率30 Hz和配重12.5g)試驗測試結果,如圖3所示。圖3(a)~圖3(c)為應變、電壓和電流時程響應,其響應頻率均與激勵頻率一致。其中,應變峰值1.905×10-9、電壓峰值17.99 V和電流峰值0.619 4 mA。圖3(d)~圖3(f))為特性關系曲線,其中,圖3(d)和圖3(e)顯示電壓-應變和電流-應變關系曲線呈現近似橢圓形回線,這是因為能量轉換存在時滯現象致使應變與電壓(電流)同頻異相導致的。圖3(f)顯示功率-應變關系曲線呈現近似“蝴蝶形”回線,這是因為歐姆定律導致應變與功率存在倍頻關系。

需要指出的是,圖3(a)應變幅值向上明顯偏置了0.362×10-9,為其半幅的23%。應變的偏移導致電壓和電流也發生偏置,其中電壓偏置7.4%,電流偏置6.4%;以上應變、電學響應的偏置現象普遍存在于所有試驗結果中,僅表現為不同工況的偏置程度有差別而已。比如,配重為0的工況(激勵頻率120 Hz和負載10 kΩ),其應變偏置11%,電壓偏置2.5%和電流偏置2.2%。總體而言,當懸臂梁自由端配重加大時偏置現象更顯著。由此判斷,懸臂梁自身和配重的重力效應是導致俘能器在簡諧振動環境下,應變響應、電學響應仍存在一定程度偏置的主要原因。

圖3 典型工況的試驗結果Fig.3 Typical test results

為進一步探究不同加載工況(電路負載、激勵頻率和配重)對壓電俘能器機電轉換關系的影響規律,本文對各種加載條件下電能回收能力進行了觀測分析。

1.2.1 變激勵頻率分析

預設激振器一定加速度幅值(增益),由0均勻增大激振頻率至150 Hz,觀察俘能器應變響應的變化規律,即掃頻試驗。一般而言,當頻率接近自振頻率時將產生諧頻共振現象,此時的應變、電學響應幅值將達到峰值,其對應的頻率即為結構主頻。當懸臂梁無配重時,掃頻至123.3 Hz,單變應變振幅譜峰值為2.696 2 Hz-1,說明結構配重為零時俘能片諧振頻率為123.3 Hz。

匯總試驗結果,結構配重為0,12.5 g,18.75 g和25 g時,其諧振頻率分別為123.2 Hz,31.2 Hz,25.75 Hz和21.94 Hz,顯然較大的配重可顯著降低結構主頻,如圖4所示。圖4的縱坐標功率譜幅值進行歸一化處理在0~100%內,便于各工況的對比分析。

圖4 激勵頻率對俘能性能影響Fig.4 Influence of excitation frequency on harvesting capacity

1.2.2 變電路負載分析

除了振動環境的影響,俘能器內部電路負載的大小也會對能量采集效率產生影響。于是此次設定特定激勵頻率和末端配重,僅通過調整負載電阻值觀察其對壓電系統能量俘獲能力的影響規律。

試驗結果顯示電阻負載對應變、電學響應幅值具有明顯影響,與激振頻率類似,俘能器俘能效率與內在電阻數值顯著相關,只有將電路負載調整到最佳數值附近,才能為俘能器挖掘潛在最大電能提供重要保證。比如,當激振頻率20 Hz、配重25 g條件下,電阻調整到50 kΩ時獲得功率峰值為2.249 mW。

匯總試驗結果,當配質量為0,12.5 g,18.5 g和25 g,電阻負載分別在10 kΩ,30 kΩ, 35 kΩ和50 kΩ時,俘能器達到功率峰值,如圖5所示。圖5的縱坐標功率進行歸一化處理在0~100%內,便于各工況的對比分析。

圖5 電阻負載對俘能性能影響Fig.5 Influence of excitation frequency on harvesting capacity

1.2.3 最佳配置分析

通過研究,發現壓電單晶懸臂梁俘能器的激勵頻率、電路負載對懸臂梁俘能結構的能量俘獲能量有重要影響。基于以上結果,接下來將探討俘能結構優化(配重)問題。值得一提的是,為了保證試驗條件的統一性,激振器增益設定了相同數值,但由于開環電路特性使得各個工況加速度激勵幅值并非保持同一水平。鑒于此,本文將各工況下的電功率換算為單位加速度幅值(1 m/s2)下的數值,便于比較分析。

單位加速度幅值加載條件下,不同配重俘能結構在各自最佳電阻負載設置下,獲得的電功率隨激勵頻率的變化,如圖6所示。圖6中不同配重俘能結構在各工況下的轉換功率峰值和功率均值,如表3所示。通過對比分析,發現配重為18.75 g的俘能結構,其在近似諧頻25 Hz和最佳電阻負載35 kΩ下,獲得了最大的功率峰值35.42 mW和最大的功率均值5.589 8 mW。說明了壓電懸臂梁結構俘能系統的最佳配重應為18.75 g,而并非是越大的配重應獲得更大的電能。

圖6 俘能性能綜合對比Fig.6 Comprehensive comparison of harvesting capacity

表3 俘能性能綜合對比

總而言之,電路負載、激勵頻率和結構配重需要綜合考慮,尋求系統的最佳匹配是激發其潛在電能俘獲能力的重要途徑。

2 壓電懸臂梁俘能有限元建模

有限元分析軟件ABAQUS具有強大的多物理場分析能力,同時還可做系統的分析和研究,本文采用該軟件對壓電單晶復合懸臂梁進行實體建模。通過仿真分析探討壓電懸臂梁的機電轉換關系和特性。

依據實際情況設定分層材料特性參數、幾何尺寸,設定一端固定一端自由邊界條件,并于固定端施加加速度簡諧激勵。在模型中,采用六面體八節點壓電單元(C3D8E)模擬壓電陶瓷層(PZT-5H);采用六面體八節點實體單元模擬非壓電材料(聚酯、金屬銅、不銹鋼和聚酰亞胺)。層間連接采用“TIE”約束命令以此模擬封裝的PPA制品。本次仿真可實現無質量塊壓電懸臂梁俘能特性計算。ABAQUS的PPA懸臂梁模型,如圖7所示。

圖7 PPA壓電懸臂梁ABAQUS模型Fig.7 ABAQUS model of PPA piezoelectric cantilever

為了驗證仿真方法的正確性,本文對工況一致的仿真與試驗結果進行對比。頻率50 Hz簡諧激勵下的俘能壓電片的電壓時程響應,其中試驗電路負載為10 kΩ,如圖8所示。該工況下的仿真電壓時程響應與試驗測試結果幅值和頻率均吻合較好,幅值誤差為8%(見圖8(a))。該工況下的電壓-應變曲線關系,試驗曲線呈現橢圓形,而仿真結果為線性關系(見圖8(b))。對比所有工況,仿真電壓幅值預測與試驗結果吻合良好,誤差保證在15%以內;仿真電壓與應變關系均為線性,與試驗回線不同,但其斜率與試驗回線的平均斜率基本保持一致。

出現以上結果的原因解釋如下:①仿真計算獲得的電壓為俘能片內部電壓vp,而試驗中很難直接測得內部電學響應,通常需通過連接外部電路測得,因此,試驗結果實際上為外部(電路)電壓v。根據本文試驗電路情況,vp與數據采集獲得v為并聯關系,即v=vp,因此內、外電壓的幅值和頻率應一致。仿真和試驗的電壓響應幅值和頻率一致性吻合良好,驗證了有限元仿真計算的正確性,見圖8(a);②有限元仿真對象是俘能片在振動環境下由于電荷發生偏移集中繼而產生電位差的過程,該vp與結構應變同頻同相,而試驗中由于外接了電路,使得v與結構應變存在同頻異相的現象;內外部電壓相位差為0.005 s(=0.087~0.074),見圖8(a);仿真(內部)電壓-應變的同頻同相,使得特性曲線為線形關系,而試驗的(外部)電壓-應變的同頻異相,使得特性曲線呈現橢圓形,見圖8(b)。

總之,有限元模型法可較好的預測俘能器的內部電壓,但要獲得俘能器外接電路的電學響應和俘能效果,還需借助壓電懸臂梁俘能數學模型進行預測。

圖8 典型工況的ABAQUS模型與試驗結果對比Fig.8 Typical response comparison of ABAQUS model with test data

3 壓電懸臂梁俘能數學建模

鑒于有限元模型模擬壓電單晶懸臂梁俘能特性的缺陷,提出一種機電耦合數學模型。

3.1 不考慮重力效應的機電模型(原模型)

3.1.1 耦合振動方程

基于Euler-Bernoulli假設,在基礎激勵下,有末端配重的均質單晶懸臂梁受迫振動方程如下

(1)

圖9 PPA壓電懸臂梁俘能機電模型示意Fig.9 Schematic diagram of energy harvest system for PPA cantilever

(2)

式中:qk(t)為模態坐標;φk(x)為模態函數;vp為壓電材料產生的內部電壓;ωk為第k階模態無阻尼固有頻率;ξk為模態機械阻尼比;fk為激勵模態力函數和λk為耦合項[14]。

3.1.2 耦合電路方程

壓電本構方程定義了應變(S)或應力(T),電荷密度位移(D)和電場(E)之間的相互作用。于是,壓電本構的張量(應力-電荷形式)可簡化表示如下:

(3)

(4)

式(3)中壓電的軸向應變分量在復合梁中的子結構層與梁在位置x處的曲率成正比,由此可以計算出由于該應變引起的電場分量。由于外部電路導納為1/R,根據Gauss定律積分形式可獲得輸出電壓為

(5)

式中:R為電阻負載兩端的電壓;A3為x-z平面的電極面積。將式(4)代入式(5),并結合壓電層平均彎曲應變曲率關系,獲得式(6)

(6)

然后,根據Kirchhoff定律,內部電流ip表達式如下

(7)

式(6)和式(7)右端相同,于是合并兩式求解內部電流。一般而言,試驗中很難直接測得內部電流ip,通常都是測試獲得外部回路的電流i。于是,根據歐姆定律則可獲得外部電路電流

i(t)=vp(t)/R

(8)

根據電路情況,內部電壓vp與數據采集獲得外部電壓v為并聯關系,即v=vp。

此外,為便于求解,式(6)按模態展開形式如下

(9)

3.2 考慮重力效應的機電模型(修正模型)

機電耦合模型以簡諧振動為假設,依此求解獲得的電壓和電流響應也為簡諧振動。然而,在1.2節特性試驗過程中,發現簡諧激勵振動下的懸臂梁俘能電學響應存在明顯偏置現象,而且該現象會隨著自由端配重的增加而更加顯著。鑒于此,本文將充分考慮配重重力效應對壓電懸臂梁俘能響應的影響,對機電耦合模型進行改進,進一步提高模型精度,使其適應更廣泛的懸臂梁俘能系統。

基于Gauss理論,自由端質量塊m對懸臂梁產生附加應變,見圖9。接著引發的額外電壓如下

(10)

(11)

(12)

(13)

(14)

3.3 機電模型求解

耦合振動方程和耦合電路方程綜合形成機電耦合方程,即式(2)和式(9)的聯立。實際上,它們為一對微分方程組,但要獲得該微分方程組的解析解非常困難。頻率響應函數(frequency response function,FRF)可用來表征測試系統在給定頻率下的穩態輸入個輸出的關系,其具有直觀反映測試系統在各個不同頻率簡諧輸入信號的響應特性。為獲得穩態解,本文將采用FRF對機電耦合微分方程組進行數值求解。

假設懸臂梁固定端激勵簡諧力函數fk(t)=Fkeωtj,其幅值Fk可由式(15)表達

(15)

式中:A0為固定端橫向位移幅值;ω為激勵頻率;j為虛部數值。

基于線性系統假設,模態坐標寫成幅值為Qk的復數qk(t)=Qkeωtj,電壓響應寫成幅值為V的復數v(t)=Veωtj,于是機電耦合方程可以表達如下

(16)

(17)

求解式(16)、式(17),即可獲得外部/內部電壓v/vp和內部電流ip響應

v(t)/vp(t)=eωtj×

(18)

ip(t)=eωtj×

(19)

4 機電模型驗證

為便于誤差分析,本文對所有算例根據配重組合特定電阻分類4組。第一組配重為0和電路負載10 kΩ下的算例;第二組配重為12.5 g和電路負載為45 kΩ下的算例;第三組配重為18.75 g和電路負載45 kΩ的算例;第四組配重為25 g和電路負載50 kΩ的算例。每組則通過變化激勵頻率獲得各種不同工況。

4.1 誤差分析

為了驗證模型對特性曲線模擬的精確性和有效性,定義特性曲線誤差百分比如下

Err(k)=

(20)

式中:Err(k)為特定電路電阻和質量塊配置下,第k工況(激勵頻率為變量)模型誤差;vexp/vsim為測試/模型結果的電壓響應;iexp/isim則表示測試/模型結構的電流響應,由于電流和電壓對功率的貢獻為等效的,因此公式中它們的權重各占50%。

采用FRF法,對同樣工況的試驗算例進行了對應機電耦合方程的求解,并獲得一系列仿真結果。負載50 kΩ,頻率10 Hz,配重25 g工況下的原模型和修正模型仿真效果對比,如圖10所示。發現兩種模型均可對應變-電壓橢圓曲線有效描述,其中,修正模型更為有效,還可對該回線的偏置(向右下)趨勢進行準確預測(見圖10(a)和圖10(b));類似地,修正模型除有效描述應變-功率蝴蝶曲線之外,還可對該回線的偏置(向右)現象有效模擬(見圖10(c)和圖10(d))。直觀而言,修正模型對壓電單晶懸臂梁機電特性描述更靈活、更有效。

圖10 模型與修正模型擬合典型試驗數據的效果對比Fig.10 Simulation comparison between model and improved model with typical test data

通過觀察所有工況下的仿真結果,采用原模型及修正模型模擬試驗數據的特性曲線誤差百分比,如圖11所示。總體而言,誤差基本集中在20%以內,但不可否認,誤差范圍隨著工況不同有所變化,越接近最佳工況(最佳頻率),其誤差呈現遞減趨勢,均小于10%;當工況逐漸遠離最佳工況時,誤差也逐漸增大。以第一組為例(見圖11(a)),在120 Hz(最接近諧振頻率123.2 Hz)下,采用原模型和其修正模型的誤差比較接近,誤差分別為8.09%和7.80%;而當20 Hz(最遠離共振頻率)下,采用原模型和修正模型誤差則差別加大,誤差分別為84.63%和49%。為了觀測模型整體精度,將第一組的10個工況的模型誤差進行平均,則顯示原模型和修正模型平均誤差分別為29.37%和23.2%,采用修正模型平均誤差比原模型明顯提高21%。類似地,第二組、第三組、第四組的誤差分布規律基本類似,(見圖11(b)~圖11(d))。

圖11 模型與修正模型誤差對比Fig.11 Error comparison between model and improved model

4.2 等效誤差分析

對于振動俘能器而言,振動環境為其能量轉換的主要來源,而振動環境非簡單特定頻率和幅值的簡諧激勵,而是具有不可預測的各種幅值和頻率的隨機激勵。當環境振動頻率與結構主頻接近時,俘能能力大大提升,反之,當環境振動頻率遠離結構主頻時,俘能能力大大減弱。因此,模型的有效性不但和特性曲線描述誤差有關,還應和相應工況的能量占比相關。基于此,本文提出“等效誤差百分比EErr”,用以衡量模型預測能量回收能力的評判指標。

以前述不同配重分類組別為例,定義組內第k工況下獲得的功率響應幅值為N(k),組內所有工況獲得的總功率響應幅值為∑N(k),于是第k工況的能量占比NR(k)為

(21)

于是,該g組采用模型獲得的等效誤差百分比為

EErr(g)={∑NR(k)·Err(k)}×100%

(22)

等效誤差百分比EErr(g)綜合體現了模型對g結構在各頻率外激勵下的能量俘獲能力的有效評判。原模型和修正模型在4個分組中的等效誤差,如表4所示。其中原模型平均等效誤差為16.93%,模型經過改進后,其平均等效誤差減小到12.98%,精度整體提高了23.33%。驗證了修正模型的精確性和有效性,說明改進模型可用于精確預測壓電單晶懸臂梁回收能量數值以及機電特性曲線的描述。

表4 原模型與修正模型等效誤差對比

5 結 論

本文對壓電單晶懸臂水平梁的俘能特性進行了試驗和模型研究。

(1)通過壓電俘能懸臂梁俘能性能試驗,發現其力-電轉換呈現復雜特性;發現懸臂梁的自由端配重、外部激勵頻率以及內部電路負載對系統能力回收均具有重要影響,其最佳匹配下才可為俘能器挖掘潛在最大電能提供重要保證。

(2)提出壓電俘能懸臂俘能力電耦合數學模型用以預測其振動能量俘獲的電能。通過與有限元模型對比,驗證了其不但可以精確模擬力、電響應幅值和頻率,還可有效擬合應變-電壓、應變-電流和應變-功率等復雜特性曲線;通過與試驗結果的誤差對比分析,同時驗證了考慮重力效應的修正模型還可較好預測特性曲線的偏置現象,充分驗證所提修正模型描述俘能特性更為全面、更為有效。

以上試驗研究和模型研究成果對探究壓電懸臂梁俘能系統振動能量回收能力及其在工程上的有效應用提供重要的理論基礎。

猜你喜歡
振動模型
一半模型
振動的思考
科學大眾(2023年17期)2023-10-26 07:39:14
噴水推進高速艇尾部振動響應分析
重要模型『一線三等角』
This “Singing Highway”plays music
重尾非線性自回歸模型自加權M-估計的漸近分布
振動攪拌 震動創新
中國公路(2017年18期)2018-01-23 03:00:38
中立型Emden-Fowler微分方程的振動性
3D打印中的模型分割與打包
FLUKA幾何模型到CAD幾何模型轉換方法初步研究
主站蜘蛛池模板: 亚洲一区二区三区香蕉| 中文字幕av一区二区三区欲色| 欧美va亚洲va香蕉在线| 亚洲欧美自拍中文| 露脸一二三区国语对白| 色国产视频| 国产AV无码专区亚洲A∨毛片| 亚洲日韩每日更新| 国产乱肥老妇精品视频| 刘亦菲一区二区在线观看| 亚洲欧美综合精品久久成人网| 国产尤物jk自慰制服喷水| 国产手机在线ΑⅤ片无码观看| 丝袜高跟美脚国产1区| 老色鬼久久亚洲AV综合| 国产成人精品视频一区视频二区| 中日韩一区二区三区中文免费视频 | 亚洲最大福利网站| 一区二区三区精品视频在线观看| a毛片在线| jizz在线免费播放| 韩日免费小视频| 无码中文AⅤ在线观看| 沈阳少妇高潮在线| 青草91视频免费观看| 无码内射中文字幕岛国片| 亚洲一级毛片免费观看| 内射人妻无码色AV天堂| 国产成人成人一区二区| 亚洲五月激情网| 亚欧成人无码AV在线播放| 伊人激情综合网| 人妻无码AⅤ中文字| 九九热免费在线视频| 99热这里只有免费国产精品| 国产真实二区一区在线亚洲| 天天激情综合| 日韩午夜片| 伊人无码视屏| 国产精品黄色片| 亚洲精品国产成人7777| 中文字幕人成人乱码亚洲电影| 亚洲成人一区二区| 久久精品娱乐亚洲领先| 热热久久狠狠偷偷色男同| 亚洲欧美人成人让影院| 香蕉久久国产精品免| 国产欧美日韩在线在线不卡视频| 欧美在线三级| 亚洲视频免费播放| AV不卡无码免费一区二区三区| 日韩在线视频网站| 亚洲日韩AV无码一区二区三区人| 在线播放91| 一区二区午夜| 草逼视频国产| 午夜精品久久久久久久无码软件| 欧美一区二区自偷自拍视频| 亚洲天堂777| 91啦中文字幕| 精品无码专区亚洲| 日韩视频精品在线| 亚洲欧美另类专区| 亚洲欧洲自拍拍偷午夜色| 波多野结衣久久高清免费| 亚洲天堂免费观看| 国产精品欧美亚洲韩国日本不卡| 男人天堂亚洲天堂| 国产黄色视频综合| 欧美特黄一级大黄录像| 婷婷五月在线视频| 久996视频精品免费观看| 国产成人精品优优av| 欧美国产在线看| 欧洲日本亚洲中文字幕| 亚洲福利网址| 美女被躁出白浆视频播放| 天天色综网| 99久久无色码中文字幕| 国产日韩欧美成人| 国产精品任我爽爆在线播放6080| 亚洲视屏在线观看|