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

超空泡航行體閉環控制動力學特性研究

2016-01-15 03:13:04熊天紅,包伯成
振動與沖擊 2015年17期

第一作者熊天紅女,講師,1977年生

通信作者包伯成男,博士,教授,1966年生

超空泡航行體閉環控制動力學特性研究

熊天紅1,包伯成2

(1.南京理工大學瞬態物理國家重點實驗室, 南京210094; 2.常州大學信息科學與工程學院, 江蘇常州213164)

摘要:通過對超空泡航行體的動力學描述,采用分段線性滑行力函數擬合復雜非線性滑行力函數,構建了超空泡航行體閉環控制動力學模型,獲得以反饋控制增益為可變參數的四維混沌系統.利用相軌圖、龐加萊映射、分岔圖和Lyapunov指數等動力學分析工具,分析不同反饋控制增益變化時系統復雜的動力學行為。結果表明,超空泡航行體閉環控制動力學行為依賴于各個閉環控制增益,隨著這些參數的變化,系統存在分岔、混沌、周期窗、共存吸引子和不完全費根鮑姆樹等奇異的非線性物理現象;合理選擇反饋增益,能夠實現超空泡航行體的穩定航行。研究結果將對超空泡航行體反饋控制器的設計具有重要的指導意義。

關鍵詞:分岔;混沌;閉環控制;超空泡航行體

基金項目:國家自然科學基金(51277017);江蘇省自然科學基金(BK2012583);國家自然科學基金青年基金(11402116)

收稿日期:2014-05-27修改稿收到日期:2015-01-13

中圖分類號:O322

文獻標志碼:A

DOI:10.13465/j.cnki.jvs.2015.17.028

Abstract:By describing dynamics of supercavitating vehicles and utilizing a piecewise-linear sliding force function to fit a complex nonlinear sliding force function, a closed-loop control dynamic model for a supercavitating vehicle was constructed, a four-dimensional chaotic system with feedback control gains as variable parameters was obtained. By using dynamic analysis tools, such as, phase portrait, Poincaré map, bifurcation diagram and Lyapunov exponent, the complex dynamic behaviors of the system with variation of different feedback control gains were analyzed. The results indicated that the closed-loop control dynamic behaviors of supercavitating vehicles depend on each closed-loop control gain; with these parameters’ varying, the novel nonlinear phenomena, such as, bifurcation, chaos, periodic window, co-existing attractor, imperfect Feigenbaum-tree and so on appear; the stable motion of supercavitating vehicles can be realized by choosing appropriate feedback gains. The study results provided a guidance for feedback controller design of supercavitating vehicles.

Closed-loop control dynamic characteristics for supercavitating vehicles

XIONGTian-hong1,BAOBo-cheng2(1. State Key Laboratory of Transient Physics, Nanjing University of Science and Technology, Nanjing 210094, China;2. College of Information Science and Engineering, Changzhou University, Changzhou 213164, China)

Key words:bifurcation; chaos; closed-loop control; supercavitating vehicles

當航行體在水下高速航行時, 由于伯努利效應,使得航行體周圍的液體汽化,從而產生覆蓋航行體大部分表面的超空泡, 進而降低航行體在水中的阻力,大大提高了航行體的運動速度與航行距離[1-4]。超空泡航行體在水下高速航行時,航行體的尾部與空泡壁接觸時會產生復雜的非線性滑行力,非線性滑行力的出現不僅會增加航行體的摩擦阻力,還會給航行體造成振動與沖擊,進而產生混沌這一復雜的非線性現象[5-6]。

一個非線性動力學系統,當系統參數在一定范圍內變化時,便會出現混沌和分岔等物理現象。混沌和分岔作為一種復雜的非線性物理現象,過去的幾十年里在科學、數學以及工程應用等各個領域得到了研究者極大的關注,關于具體物理系統的動力學建模、非線性物理現象揭示、穩定性和分岔分析等多個方面取得了大量的研究成果[7-9]。目前國內外有關超空泡航行體的非線性動力學研究,主要是針對超空泡航行體開環參數引起的非線性現象和航行體反饋控制的研究[9-16],有關超空泡航行體閉環控制的動力學特性研究,未見有文獻報道。超空泡航行體的閉環控制特性是進行超空泡航行體反饋控制器設計的重要依據,因此完全有必要對其動力學特性進行分析和討論。

以超空泡航行體的動力學描述為出發點,采用分段線性擬合的滑行力函數,構建超空泡航行體閉環控制動力學模型。基于此,利用常規的動力學分析工具,以反饋控制增益為可變量,揭示了超空泡航行體非線性物理現象。通過對航行體閉環控制動力學特性的分析, 合理選擇反饋控制律中各個變量增益,能夠實現超空泡航行體的穩定航行。

1超空泡航行體的動力學描述

航行體在水下高速航行時,周圍液體汽化,從而產生覆蓋航行體大部分表層的超空泡。為研究在一般情況下的空泡形態,引用反映空泡的無量綱空化數,其表達式為[1]

(1)

式中:P∞為外壓,Pc為空泡內壓,ρ為水的密度,V為航行體總速度。對于自然超空泡,空泡內壓力可以近似認為等于水的飽和蒸汽壓力Pc= 2 350 Pa。

超空泡航行體的結構和外形尺寸見圖1(a)。圖中長度L的航行體有兩段構成:后端長度2/3L和半徑R的圓柱段,以及前端長度1/3L的圓錐段。對于該航行體模型,重心(CG)距離頭部為17/28L。航行體的頭部有一個半徑Rn的圓盤空化器。鰭舵采用十字形布局,置于航行體尾端。

航行體與空泡之間的相互作用關系可以通過作用力方程進行建模,在航行體的體坐標系中,作用在航行體上的力見圖1(b), 主要有空化器上的升力Fcavitator、尾翼上的升力Ffins、尾部和空泡壁面之間相互作用產生的滑行力Fplaning,以及航行體質心位置的重力Fgravity。

圖1 超空泡航行體 Fig.1 Supercavtating vehicle

超空泡航行體動力學建模的體坐標系原點位于航行體頭部的圓盤空化器頂端端面的圓心, 超空泡航行體動力學建模采用四個狀態變量來描述超空泡航行體的動力學,分別為垂直位置z、橫向速度w、俯仰角θ和俯仰角速度q。定義橫向速度w與航行體軸線垂直,定義前行速度V與航行體軸線平行。

超空泡航行體的浸沒深度見圖1(c)。超空泡航行體尾端的浸沒深度用參數h′來表示,可表達為[13]

(2)

式中:

f(w)=2w+(w+wt0)tanh[-k(w+wt0)]+

(w-wt0)tanh[k(w-wt0)]

式中:wt0=(Rc-R)V/L為位于過渡點的正w值,k為一個用于選擇控制近似誤差的常數,一般有k=300,Rc為離空化器距離L處的空泡半徑[13]。

按照文獻[14]給出的定義式,航行體中心線與空泡中心線之間的幾何角即超空泡航行體浸沒角α可表達為

(3)

假設在航行過程中各種作用力保持平衡,航行體總速度V保持不變,由上述各參數的定義可得超空泡航行體的動力學模型如下[14]

(4)

其中:m為密度比(ρm/ρ),g為重力加速度,n為尾翼效率;各系數矩陣M0、A0、B0和重力Fgravity分別可表達為

式(4)中超空泡航行體的滑行力Fplaning可表達為

(5)

超空泡航行體設有反饋控制器,其控制輸入分別為δe和δc,一般選擇δe= 0、δc=kzz-kθθ-kqq,kz、kθ和kq均為實常數,分別為控制變量z、θ和q的反饋增益[10-11]。

2滑行力分段線性擬合的動力學模型

2.1滑行力分段線性擬合

超空泡航行體的系統參數值列在表1中所示[14]。基于表1的系統參數值,采用式(5)表達的非線性滑行力Fplaning,則滑行力Fplaning與橫向速度w之間關系曲線如圖2中的虛線所示。

表1 超空泡航行體的參數

為了便于對超空泡航行體閉環控制增益的動力學特性進行定性和定量分析,可采用一個分段線性滑行力Fp函數來擬合式(5)表示的復雜非線性滑行力Fplaning函數。這里采用了一種折衷的五段形式的分段線性滑行力函數。盡管可以采用更多段形式的分段線性滑行力函數來無限逼近非線性滑行力函數,但對系統的動力學特性影響不大,反而又增加了分析難度。

(6)

式(6)表示的分段線性擬合的滑行力Fp與橫向速度w之間關系曲線如圖2中的實線所示。從圖2中可見,非線性滑行力Fplaning和分段線性滑行力Fp之間有著較好的擬合度。

圖2 Fp-w和Fplaning-w的關系曲線Fig.2CurvesofFp-wandFplaning-w圖3 采用不同滑行力函數時系統(4)隨σ變化的分岔圖Fig.3Bifurcationdiagramofsystem(4)withdifferentplaningforcefunctions

當超空泡航行體的反饋控制增益分別為kz= 15、kθ= 30和kq= 0.3,即δe= 0、δc= 15z-30θ-0.3q時,采用非線性滑行力Fplaning函數的系統(4)隨空化數σ變化的關于狀態變量w的分岔圖見圖3(a),而采用分段線性滑行力Fp函數的系統(4)隨空化數σ變化的關于狀態變量w的分岔圖見圖3(b)。

觀察圖3不難發現,兩種分岔圖的變化趨勢是基本一致的,由此表明,采用分段線性滑行力Fp函數來擬合非線性滑行力Fplaning函數后,系統的動力學特性基本保持不變。滑行力函數的分段線性化簡化了超空泡航行體的動力學模型,這對超空泡航行體的動力學特性分析具有重要的意義。

2.2超空泡航行體閉環控制的動力學模型

當系統(4)選用表1參數值,固定σ= 0.0313和δe= 0,且采用分段線性滑行力Fp函數時,可得到超空泡航行體閉環控制動力學模型為

(7)

分段線性滑行力Fp函數可簡化為

(8)

系統(7)的可變參數僅為三個控制增益,因此系統(7)為專用于刻畫超空泡航行體閉環控制動力學行為的簡單模型。

當kz= 15、kθ= 30和kq= 0.3時,系統(7)的運行軌跡即混沌吸引子在各平面上的投影見圖4,相應的龐加萊映射見圖5。利用Jacobi方法計算其Lyapunov指數為L1= 13.705 4,L2=-3.681 6,L3=-29.060 9和L4=-46.340 3。從系統(7)的相軌圖、龐加萊映射以及Lyapunov指數可見,系統(7)為一個四維混沌系統,可生成混沌吸引子。

圖4 混沌吸引子在相平面上的投影 Fig.4 Chaotic attractors

圖5 龐加萊映射在w-θ平面上的投影 Fig.5 Poincaré map in w-θ plane

3閉環控制的動力學特性分析

3.1平衡點和穩定性

代入式(8),可得到系統(7)隨控制增益kz和kθ變化的平衡點。對于典型參數kz= 15、kθ= 30和kq= 0.3,系統(7)只在1.38

S = [1.398 6, 0, 0.019 0, 0.041 6]

把系統(7)在平衡點S處線性化,得到Jacobi矩陣為

平衡點處的特征方程為

det(1λ-JS)=0

可解得四個特征值為

λ1,2=273.59±j344.05,λ3,4=-21.13±j30.90

因此,平衡點S為不穩定鞍焦點,即指數2平衡點,滿足形成混沌吸引子的必要條件。

3.2閉環控制的動力學行為

當反饋控制增益發生變化時,系統(7)有著不同的動力學行為。可利用分岔圖和Lyapunov指數譜等動力學分析工具,分析系統(7)在各個反饋控制增益變化時的動力學行為。

(1)固定kθ= 30和kq= 0.3,kz在[0, 30]區間變化

系統(7)隨反饋增益kz變化的Lyapunov指數譜見圖6(a),相應的狀態變量w的分岔圖見圖6(b)。由圖可見,Lyapunov指數譜與分岔圖是基本一致的,kz∈[0, 18.8]時,系統(7)經倍周期分岔進入混沌狀態,然后由切分岔突變到周期狀態,再歷經倍周期分岔后進入混沌軌道。在此區間從圖5(b)可以觀測到在比較寬的參數范圍內,系統具有一個正值Lyapunov指數,一根零值Lyapunov指數和一根負值Lyapunov指數。在混沌區間內,系統存在多個最大Lyapunov指數小于0的周期窗,主要出現在kz= 1.95、3.45、11.7等附近。在隨著反饋增益kz逐步增大,在kz∈[18.8, 30]時,系統由于混沌危機引發運行軌道狀態突變,形成周期軌道,在此參數區間內最大Lyapunov指數小于或等于0。但在[21.2, 24.5]和[20.2, 21.0]的兩個區間內,其運行軌跡存在不完全費根鮑姆樹現象[17],周期軌道經倍周期正分岔后,又經周期減半逆分岔回到原先的周期軌道。

另外,在一些參數值附近,系統(7)存在吸引子的共存現象。譬如,當kz= 18.6時,系統的運行軌跡在不同的初始值時會發生周期極限環與混沌吸引子的共存現象,分別見圖7(a)和7(b)。

圖7(a)為初始值為(0.5377, 1.833 9,-2.258 8, 0.862 2)時的周期極限環;圖7(b)為初始值為(0.022 9,-0.262 0,-1.750 2,-0.285 7)時的混沌吸引子;而當kz= 23時,系統的運行軌跡在不同的初始值時有兩個共存的不同周期數的極限環,分別見圖7(c)和7(d),圖7(c)為初始值為(-0.001 6, 0.898 5, 0.001 8, 0.145 4)時的周期4極限環; 圖7(d)為初始值為(-0.025 7, -2.329 0, 0.009 0, -0.917 8)時的周期2極限環。

(2)固定kz= 15和kq= 0.3,kθ在[0, 50]區間變化

圖6 隨反饋增益kz變化的動力學行為Fig.6Dynamicalbehaviorswiththevariationofkz圖7 不同kz且不同初始值時w-θ平面上的相軌圖Fig.7Phaseportraitswithdifferentkzanddifferentinitialvalues

系統(7)隨反饋增益kθ變化的Lyapunov指數譜及相應的狀態變量w的分岔圖分別見圖8(a)和8(b)。當kθ∈[0, 10.5]時,系統(7)以周期3軌道運行,相對應圖8(a)的最大Lyapunov指數等于0。kθ逐步增大,在kθ= 10.5附近系統(7)的周期3軌道經過混沌危機直接進入混沌,最大Lyapunov指數突變為正值。并在kθ= 35.2附近系統(7)的運行軌跡經切分岔又回到周期3軌道,最大Lyapunov指數再次突變為0。kθ在[10.5, 35.2]較寬的范圍內,系統(7)處于混沌狀態,并存在豐富的周期窗。在[16, 16.5]和[31, 31.5]區間內,周期窗內系統(7)的運行軌跡存在不完全費根鮑姆樹現象。在[16, 16.5]區間內,出現的不完全費根鮑姆樹現象是周期軌道經倍周期分岔后又經減半逆分岔回到周期軌道[17];而在[31, 31.5]區間內,出現的另一種不完全費根鮑姆樹現象是周期軌道經倍周期分岔后進入混沌狀態,然后又經減半逆分岔回到周期軌道。

(3)固定kz= 15和kθ= 30,kq在[0.2, 0.5]區間變化

系統(7)隨反饋增益kq變化的Lyapunov指數譜及相應的狀態變量w的分岔圖分別見圖9(a)和9(b)。當kq< 0.284時,系統(7)處于周期狀態,最大Lyapunov指數等于0,并在[0.25,0.28]區間內存在不完全費根鮑姆樹現象。當kq∈[0.284,0.308]時,系統(7)基本處于混沌狀態,最大Lyapunov指數基本大于0,且在kq= 0.289和kq= 0.303附近存在兩個較窄的周期窗。當kq>0.308時,系統(7)又回復到最大Lyapunov指數等于0的周期狀態,但其運行軌跡在0.38附近出現了倍周期分岔行為,這時對應的圖9(a)中的第2個Lyapunov指數從負值回到零值然后重新進入負值。

通過超空泡航行體閉環控制動力學特性分析可知,調整反饋控制律中控制增益,能夠實現超空泡航行體的穩定航行。當σ= 0.031 3時,固定控制增益kz、kθ的值,增加反饋增益kq的取值范圍,其狀態變量w隨反饋增益kq變化的分岔圖見圖10(a)。不難觀察到,kq在[1.028 6,1.5]區間內系統平衡點的位置保持不變。如當kq=1.2時平衡點S= (1.398 6, 0, 0.019 0, 0.041 6),系統在平衡點處線性化Jacobi矩陣的特征根為λ1,2=-63.89±j621.86,λ3,4=-9.39±j24.64,這里λ1,2、λ3,4為實部為負的共軛復根,表明平衡點S為穩定的焦點。圖10(b)為系統隨時間演化的Lyapunov指數譜,為清楚顯示,圖10(b)中僅給出三根Lyapunov指數曲線,相應的Lyapunov指數分別為L1= -8.368 2、L2=-8.676 1、L3=-65.208 9和L4=-155.486 9,其最大Lyapunov指數曲線在有限時間尺度內為負值。從系統分岔圖、平衡點處特征根及Lyapunov指數可見,當控制律δc= 15z-30θ-1.2q時,系統處于穩定狀態。

圖11給出了控制律δc=15z-30θ-1.2q時系統響應示意圖,觀察圖11不難發現,各個反饋控制變量初始時刻的值較大, 在反饋控制律的作用下四個狀態變量逐漸穩定到其平衡點上。由此也可表明,航行體在空泡內的位置和姿態固定,處于穩定航行狀態。

圖8 隨反饋增益kθ變化的動力學行為Fig.8Dynamicalbehaviorswiththevariationofkθ圖9 隨反饋增益kq變化的動力學行為Fig.9Dynamicalbehaviorswiththevariationofkq圖10 δc=15z-30θ-1.2q系統動力學行為Fig.10Dynamicalbehaviorsofsysteminδc=15z-30θ-1.2q

圖11 δ c = 15z-30θ-1.2q 時系統響應示意圖 Fig.11 System response in δ c = 15z-30θ-1.2q

4結論

針對超空泡航行體尾部與空泡壁面接觸產生的非線性滑行力,本文采用分段線性擬合的滑行力函數,構建了超空泡航行體閉環控制動力學模型,通過超空泡航行體閉環控制的非線性動力學特性的分析,得到以下結論:

(1)在一定參數范圍內,利用相軌圖、龐加萊映射和Lyapunov指數等數值仿真手段,驗證了超空泡航行體動力學模型是以反饋控制增益為可變參數的四維混沌系統。

(2)閉環控制的超空泡航行體的運行軌跡有著復雜的動力學行為,且依賴于各個閉環控制增益;隨著各個反饋增益的變化,系統出現混沌、分岔、共存吸引子和不完全費根鮑姆樹等非線性現象,從而揭示了超空泡航行體運動規律的不穩定性。

(3)通過合理設置反饋控制律中各個變量增益,可以有效抑制超空泡航行體分岔、混沌等非線性現象,實現超空泡航行體的穩定航行。

參考文獻

[1]Wang G, Ostoja-Starzewski M. Large eddy simulation of a sheet/cloud cavitation on a NACA0015 hydrofoil [J]. Applied Mathematical Modeling, 2007, 31: 417-447.

[2]Singhal A K, LI H Y, Athavale M M, et al. Mathematical basis and validation of the full cavitating high-speed torpedo[J]. Journal of Fluids Engineering, 2003, 125: 459-468.

[3]Li Q T, He Y S, XUE L P. A numerical simulation of pitching motion of the ventilated supercaviting vehicle around its nose[J]. Chinese Journal of Hydrodrnamics, 2011,26(6): 589-685.

[4]李代金,羅凱,黨建軍,等. 超空泡水下航行器空間運動建模與彈道仿真[J]. 兵工學報, 2012, 33(8): 956-961.

LI Dai-jin,LUO Kai,DANG Jian-jun,et al. Kinematic modeling and trajectory simulation for underwater supercavitating vehicles[J]. Acta Armamentarii, 2012, 33(8): 956-961.

[5]魏英杰, 王京華, 張嘉鐘,等. 水下超空泡航行體非線性動力學與控制[J]. 振動與沖擊, 2009, 28(6): 179-204.

WEI Ying-jie, WANG Jing-hua, ZHANG Jia-zhong,et al. Nonlinear dynam ics and control of underwater supercavitating vehicle[J]. Journal of Vibration And Shock, 2009, 28(6): 179-204.

[6]Lin G J, Balachandran B, Abed E H.Nonlinear dynamics and bifurcations of a supercavitating vehicle[J]. IEEE Journal of Oceanic Engineering, 2007, 32(4): 753-761.

[7]包伯成, 胡文, 許建平,等. 憶阻混沌電路的分析與實現[J]. 物理學報, 2011, 60(12): 120502.

BAO Bo-cheng, HU Wen, XU Jian-ping, et al. Analysis and implementation of memristor chaotic circuit[J]. Acta Physica Sinica, 2011, 60 (12): 120502.

[8]李群宏, 閆玉龍, 韋麗梅,等. 非線性傳送帶系統的復雜分岔[J]. 物理學報, 2013, 62(12): 120505.

LI Qun-hong, YAN Yu-long, WEI Li-mei, et al. Complex bifurcations in a nonlinear system of moving belt[J]. Acta Physica Sinica, 2013,62(12):120505.

[9]Bao B C, Xu J P, Zhou G H, et al. Chaotic memristive circuit: equivalent circuit realization and dynamical analysis[J]. Chinese Physics B, 2011, 20(12): 120502.

[10]Lin G J, Balachandran B, Abed E H. Dynamics and control of supercavitating vehicles[J]. Journal of Dynamic Control Systems, Measurement, and Control, 2008,130:021003.

[11]Nguyen V, Balachandran B. Supercavitationvechicles with noncylindrical, nonsymmetric cavities: dynamics and instabilities[J]. Journal of Computational and Nonlinear Dynamics, 2011, 6(4): 041001.

[12]Lin G J, Balachandran B, Abed E. Supercavitating body dynamics, bifurcation and control[C]. 2005 American Control Conference. Portland, OR, USA. 2005, 6: 691-696.

[13]Hassouneh M A, Vincent N, Balakumar B.Stability analysis and control of supercavitating vehicles with advection delay[J]. Journal of Computational and Nonlinear Dynamics, 2013, 8: 021003.

[14]白濤, 孫堯, 莫宏偉. 分叉分析在水下高速運動體穩定控制中的應用[J]. 哈爾濱工程大學學報, 2008, 29(10): 1067-1075.

BAI Tao, SUN Yao, MO Hong-wei. Application of bifurcation analysis to the stability control of underwater high-speed vehicles[J]. Journal of Harbin Engineering University, 2008, 29(10): 1067-1075.

[15]王京華, 魏英杰, 于開平, 等. 基于空泡記憶效應的水下超空泡航行體建模與控制[J]. 振動與沖擊,2010,29(8):160-163.

WANG Jing-hua, WEI Ying-jie, YU Kai-ping, et al. Modeling and control of underwater supercavitating vehicle based on memory effect of cavity[J]. Journal of Vibration And Shock, 2010, 29(8): 160-163.

[16]李魁彬,王安穩,鄧磊.變速超空泡航行體殼結構的動力響應分析[J].振動與沖擊,2013,32(4):138-141.

LI Kui-bin, WANG An-wen, DENG Lei. The influence of cavity shape on maneuvering rotational movement of supercavitating vehicle[J]. Journal of Vibration and Shock, 2013, 32(4):138-141.

[17]包伯成, 胡文, 劉中,等. DOG 小波映射伸縮和平移的動力學分析[J]. 物理學報, 2009, 29(8): 2239-2247.

BAO BO-cheng, HU Wen, LIU Zhong, et al. Dynamical analysis of DOG wavelet mapping with dilation and translation[J]. Acta Physica Sinica, 2009, 58(4): 2239-2247.

主站蜘蛛池模板: 亚洲香蕉久久| 伊人久久综在合线亚洲2019| 久久久黄色片| 久久精品视频一| 性欧美在线| 久久精品嫩草研究院| 美女潮喷出白浆在线观看视频| 亚洲精品麻豆| 女人18毛片一级毛片在线| 国产精品吹潮在线观看中文| 永久免费无码日韩视频| 无码区日韩专区免费系列| 国产精品成人第一区| 伊人久综合| 亚洲一区二区三区国产精品 | 色综合天天综合| 欧美日本中文| 国产精品女熟高潮视频| 亚洲精品成人福利在线电影| 国产乱人伦精品一区二区| 亚洲永久色| 久久久久无码精品国产免费| 国产呦精品一区二区三区下载| a毛片在线| 色婷婷视频在线| 久996视频精品免费观看| 国产精品99久久久久久董美香| 精品久久777| 日韩精品高清自在线| 色欲综合久久中文字幕网| 97青青青国产在线播放| 色综合激情网| 91丝袜乱伦| 美女内射视频WWW网站午夜 | 久久精品无码国产一区二区三区| 国产欧美专区在线观看| 国内丰满少妇猛烈精品播 | 熟妇丰满人妻| 久久精品人人做人人综合试看| 无码专区在线观看| 99激情网| 国产91透明丝袜美腿在线| 九色国产在线| 色视频久久| 久久亚洲中文字幕精品一区| 97影院午夜在线观看视频| 欧美日韩久久综合| 四虎在线高清无码| 91在线免费公开视频| 在线观看网站国产| 日韩欧美国产成人| 日本a级免费| 国产精品女主播| 欧美成人精品在线| 亚洲国产高清精品线久久| 波多野结衣的av一区二区三区| 99热这里只有精品免费| 亚洲欧美精品在线| 免费A∨中文乱码专区| 久久国产精品波多野结衣| 色国产视频| 精品视频免费在线| 中文字幕无码av专区久久| 999精品视频在线| 99久久精品免费看国产电影| AV网站中文| 好吊妞欧美视频免费| 精品剧情v国产在线观看| 亚洲成在人线av品善网好看| 国产自产视频一区二区三区| 婷婷色狠狠干| 久久综合一个色综合网| 亚洲精品自产拍在线观看APP| AV片亚洲国产男人的天堂| 国内精自线i品一区202| 色偷偷综合网| 国产在线自揄拍揄视频网站| 亚洲日本韩在线观看| 色综合a怡红院怡红院首页| 国产一级毛片yw| 国产成人免费观看在线视频| 日韩精品一区二区三区视频免费看|