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

基于改進螢火蟲算法的無人機路徑規劃

2023-06-02 06:47:40龍艦涵許湘揚
計算機測量與控制 2023年5期
關鍵詞:規劃環境

龍艦涵,許湘揚

(瀘州職業技術學院 電氣與電子工程學院,四川 瀘州 646000)

0 引言

無人機(UAV)是一種常見的便攜式遙控設備。UAV由于成本低、適用場景多等優點,廣泛應用于大氣探測行業、電力行業等行業[1-5]。UAV路徑規劃是UAV飛行領域內的關鍵技術之一[6-8]。路徑規劃一般是指確定UAV的起飛位置和目標任務后的一系列優化問題[9-10]。根據UAV的機動性、地理環境威脅和約束條件,規劃出最優或次優飛行路徑。規劃路徑的優劣直接影響無人機能否順利完成任務。因此,無人機路徑規劃技術已成為無人機發展的關鍵。

傳統路徑規劃算法包括人工勢場法[11]、A*算法[14]等,傳統算法需要環境精確建模且運算量較大,存在可行解時可以找到解,但當函數不連續、不可微時,則很難適用。因此在處理解決復雜問題,尤其是離散變量、非連續性問題,智能算法具有更好的適應性。常見的包括蟻群算法[12]、粒子群算法[13],基本蟻群算法全局搜索能力強,穩定可靠,但信息濃度低,精度不高,搜索速度慢[15];差分進化算法,這種算法參數較少,易控易實現,但尋優效率不高[16];魚群算法,基本魚群算法面對復雜環境有較高適應性,但經驗行為學習能力不夠,收斂準確度較低[17];遺傳算法可同時進行多尋優,但后期收斂速度慢,易得局部最優解[18];深度學習法具有較強的搜索能力,但卷積神經網絡較為復雜,三維規劃數據量較大[19]。

螢火蟲算法(FA)是YANG開發的群體啟發算法之一。其受到螢火蟲閃爍的燈光照亮啟發。螢火蟲閃爍是為了尋找配偶或抵抗捕食者的攻擊。當螢火蟲和光源之間的距離增大時,光吸收導致光變得越來越弱。螢火蟲算法作為群體智能算法之一,也容易陷入局部最優。在過去的幾年里,許多學者開發了新的FA變體,以避免過早收斂問題[20-23]。盡管這些改進算法已被證明在解決優化問題上是有效的,但仍存在以下缺點:改進方法僅限于算法的部分優化,沒有從整體上考慮螢火蟲種群的特征。盡管FA算法在個別階段的優化效果較好,但隨著優化問題維數的增加,其性能也會像其他群體智能算法一樣下降。然而,在復雜飛行環境下,無人機路徑規劃問題往往是高維的。另外,通常在路徑規劃階段,如果不考慮運動學約束,算法規劃的路徑在實際情況下無法飛行。

針對以上問題,本文基于多角度改進結合萊維飛行對基本螢火蟲算法做出優化。引入Chebyshev映射優化算法初期粒子分布,使得全局搜索能力得以提升;引入logistic混沌函數改進吸引度系數γ,增強算法局部搜索能力,提高后期收斂速度;利用Levy分布改進螢火蟲個體的搜索步長,得到新的位置更新公式和新的步長更新公式,提高了種群的搜索范圍和有效性。

1 基本螢火蟲算法

螢火蟲算法是一種智能仿生算法,模擬螢火蟲飛向更亮個體的特征。它是基于螢火蟲的閃爍模式和行為,這基本上是基于生物發光現象,這是一個過程,負責閃爍的光和交配伙伴被吸引使用這些光。這些燈的潛在優勢是,更亮的螢火蟲會吸引其他不那么亮的螢火蟲,這種現象使得螢火蟲可以探索搜索空間,即更新它們的位置。利用目標函數計算每只螢火蟲的亮度(適應度)。此外,FA與其他優化算法相比有兩個主要優勢:自動劃分和處理多模態目標函數的能力。這意味著整個螢火蟲種群能夠自主的分裂成若干子群體,各個群體都能夠在局部最優周圍進行聚集,從而在所有的這些解之中得到最優的全局解。此外,如果搜索代理的數量適當大于最優的數量,則分割的性質有利于搜索代理同時找到所有的最優。FA基于兩個概念,即亮度和吸引力,隨著距離的增加而降低。

問題的可行解即螢火蟲的位置,螢火蟲的亮度即適應度函數。可行解不斷接近最優解的過程就是螢火蟲飛向更亮個體的過程,直到滿足預設條件,即完成任務。將這個過程以數學的形式表述如下:

假設有N只螢火蟲,求解問題的維度是D,xi= (xi1,xi2,…,xiD),(i=1,2,…,N)表示第i只螢火蟲的位置,相應的xj= (xj1,xj2,…,xjD),(j=1,2,…,N)表示第j只螢火蟲的位置。

個體i和j之間的距離rij如式(1)計算:

(1)

式中,xid和xjd分別表示d維中,第i只和第j只螢火蟲所在的位置。

螢火蟲的亮度更新公式如式(2),吸引度更新公式如式(3)所示:

(2)

(3)

其中:I0表示螢火蟲初始亮度,β0代表初始吸引度,γ表示吸引系數。

若螢火蟲j更亮,則螢火蟲i向j移動時的位置公式如下:

xid(t+1)=xid(t)+β(xjd(t)-xid(t))+αi(t)ε

(4)

其中:xid(t)和xjd(t)表示在d維時螢火蟲i和j迭代至第t代時的位置,αi(t)表示螢火蟲i第t代時的步長因子,ε服從均勻分布,取值范圍為[-0.5,0.5]。

2 改進的螢火蟲算法

2.1 Chebyshev混沌映射初始化

混沌在自然界的非線性系統中普遍存在,看似無邏輯的混沌過程卻有其規律性。混沌過程的規律性在于其能不重復地遍歷所有狀態,因此遍歷性的優點引入群智能算法,能極大地避免陷入局部極值[24-25]。

通過Chebyshev混沌映射使得螢火蟲隨機混沌的形成初始種群,與其他諸如logistic、sine等混沌映射相比,以Chebyshev映射方式產生的粒子全局分布更廣,信息更加全面,可以更好的避免信息片面帶來的種群易得局部最優的問題。本文采用Chebyshev混沌映射生成混沌變量,如式(5)所示:

xn+1=cos(karccosxn),xn∈[-1,1]

(5)

式中,k為混沌的階次,當k≥2時,初始值任意設置,經過混沌迭代處理后的序列均不相關,本文取k=4,由公式計算生成彼此互不相關,但所有變量具遍歷性質的混沌序列。

利用Chebyshev映射產生初始化種群的步驟如下。

Step 1:若D維空間存在N只螢火蟲,隨機生成D維向量Y= (y1,y2,…,yd),yi∈[-1,1],i∈[1,d]作為第一只螢火蟲的初始信息;

Step 2:用式(5)對Y向量的各維度進行迭代N-1次,得到其他N-1只螢火蟲信息;

Step 3:將N個經過混沌的個體映射到解空間,具體如式(6)所示;

xid=ld+(1+yid)×(ud-ld)/2

(6)

第d維解空間的上界和下界分別用ud和ld表示,式(5)得出的第i個螢火蟲的第d維即yid,xid即經過映射后,第i個螢火蟲的第d維坐標值;

Step 4:將N個螢火蟲的適應度函數分別進行計算并排序,選M個優質個體作初始種群。

2.2 logistic混沌改進吸引度系數γ

經典的logistic混沌映射是一個動態隨機的過程,相比其他方式其具有更好的局部精確性,因此在局部尋優階段對吸引度衰減系數γ采用混沌邏輯處理,可跳出當前狀態逃離局部最優解。混沌序列處理方式如式(7)所示:

xi+1=μxi(1-xi),μ∈[0,4],xi∈[0,1],i=1,2,···,n

(7)

其中:μ是控制參數,本文取μ=4,式(7)此時已進入完全混沌狀態。為使吸引度衰減系數γ和算法運行過程相匹配,引入另一變量hs,其值與適應度函數聯系緊密,這樣處理可以在算法執行初期,當γ的值較小時,I值會很大,此時便可以增強全局搜索能力;算法執行后期,當γ的值較大時,I值會變小,可提高局部搜索能力,hs計算如式(8)所示:

(8)

(9)

其中:Tmax是最大迭代次數,iter是當前迭代次數,通過此公式可調節γ取值,fpi(t)是最佳適應度函數值。

2.3 Levy飛行策略改進步長

Levy flight是法國數學家Levy在20世紀30年代提出的一種隨機行走機制,其行走步長滿足穩定的重尾分布,在局部位置有較大的跳躍概率。Levy飛行的概率密度分布具有尖峰、非對稱和拖尾特征。它的運動模式在頻繁的短距離跳躍和偶爾的長距離跳躍之間交替,可以跳出局部最優,擴大種群搜索區域。許多昆蟲和動物,如蒼蠅和馴鹿,遵循類似自然界Levy飛行的軌跡。

引入Levy飛行的隨機搜索策略代替基本算法的步長因子[25,27]。改進的基本思路是利用Levy分布的概率分布特性,改良個體搜索步長的隨機性,提高了種群的搜索范圍和有效性和全局綜合搜索能力。Levy飛行化簡后的數學描述為:

Levy(s,γ,μ)=

(10)

引入Levy飛行后位置更新公式為:

(11)

其中:⊕是點對點乘法,Levy(λ)是符合Levy分布的步長隨機函數:

Levy(λ)~u=t-λ,1<λ≤3

(12)

改進后步長的計算公式如下所示:

(13)

其中:μ,ν服從正態分布。

(14)

其中:θ取1.5,為方便計算,對該正態分布進行簡化:

(15)

(16)

圖1 不同標度參數下的Levy分布

2.4 算法實現流程

本文算法的實現具體流程如下:

1)建立3維高空作業環境,設置起終點,生成等高度自適應巡航柵格圖。

2)初始化算法中初始亮度I0,初始吸引度β0,最大迭代次數Tmax,螢火蟲數目N,混沌映射后選擇的優質個體M。

3)引入Chebyshev映射,用式(5)~(6)生成螢火蟲初始種群。

4)根據式(7)~(9)計算吸引度系數γ。

5)代入吸引度系數至式(2)~(3),再通過式(1)計算各螢火蟲的相對熒光亮度I,距離rij,吸引度β。

6)根據公式(4)的萊維飛行策略更新位置。

7)對位置更新后的個體重新計算相對熒光亮度I。

8)當前迭代次數iter與最大迭代次數Tmax作比較,若iter>Tmax則退出迭代,若iter

9)算法結束,輸出最優結果。

具體流程圖如圖2所示。

圖2 算法流程圖

2.5 算法時間復雜度分析

為了更清晰地考察和評價改進螢火蟲算法的復雜性,將其與傳統螢火蟲的復雜性進行了比較。

首先,需要分析了FA的時間復雜度。在FA算法中,假設終止條件為Max,其總體大小為N,問題的維數為D,算法的執行時間是t1,求一個解的時間為t2。FA的時間復雜度為O(Max*N*N*D*(t1+t2))。所以它可以表示為O(N2)。

然后,分析改進FA的時間復雜度。在改進算法中,假設終止條件、總體大小和問題的維數分別與FA相同。生成的組的大小為K。計算時間di、生成組M、對xj的評價、對M的評價、對xgbest的評價的時間分別為t3、t4、t5、t6、t7。改進FA的時間復雜度為O(Max*N*D*(N*t3+N*t4+ 4*N*t5+K*K*t6+t7))。因此,它可以表示O(N2+N*k2)。很明顯,改進FA的時間復雜度取決于生成的組的大小。由于自適應多群機制,K的值可以是(0,N)之間的任意值,因此改進FA的時間復雜度在O(N2)和O(N3)之間。

雖然改進FA的時間復雜度比FA要低,但為所提出的自適應多群機制增加計算成本是值得的。原因之一是利用所提出的混沌初始化機制將搜索空間劃分為多個子空間,可以在搜索空間中探索出更滿意的解。另一種是增加對一些較好的螢火蟲的評價,可能有助于FA增加逃脫局部陷阱的概率。

3 無人機飛行環境建模及優化

3.1 飛行空間

UAV在三維空間中執行飛行任務,空間跨度大,環境復雜。此外,無人機受到許多條件的約束。因此,路徑規劃需要平衡搜索時間、路徑代價和計算復雜度,以滿足實際工程的要求。路徑規劃的目標是在復雜環境和各種約束條件下找到到達目標點的可行路徑。在環境模型中,軌跡點的位置用三維坐標(xi,yi,zi)(i= 1,2,…,Ns)表示,UAV飛行空間表示為:

(17)

其中:xmin,xmax,ymin,ymax,zmin,zmax分別表示飛行空間邊界。

3.2 約束條件

1)距離成本約束:

在UAV執行飛行任務的過程中,距離較短的飛行軌跡(Mfirst,M2,M3,…,Mlast)比距離較遠的軌跡具有更好的性價比。通過縮短飛行距離、降低無人機燃油消耗和縮短飛行時間,降低了飛行成本。考慮到無人機自身攜帶的燃料有限,需要將飛行距離控制在一定限度內,可用以下表達式來表征距離成本約束。

(18)

(19)

式中,G(G≤Gmax)為距離代價,N為構成一條完整飛行軌跡的軌跡段數,fi為軌跡i的距離代價。

2)飛行高度限制:

氣壓、風速等因素對無人機的飛行有著巨大的影響。考慮到無人機飛行的安全,飛行高度應控制在一定范圍內,飛行高度不能過高或過低。設環境模型在航跡點(xMi,yMi,zMi)處的地形高度Tmap為,fh表示最小飛行高度,F定義為最小飛行高度可調信息。

(20)

3)最大偏角約束:

在飛行過程中,無人機需要及時調整飛行姿態,這受到其機動性的限制。只有當偏角小于最大偏角時,無人機才能保持穩定的飛行姿態。設某航跡段與其相鄰航跡段的向量為pi=(xMi+1-xMi,yMi+1-yMi,zMi+1-zMi),最大偏角滿足約束條件:

(21)

4)最大俯仰角約束:

影響無人機俯仰角的因素包括風向、氣候和自身機動性。此約束限制 UAV在特定角度范圍內向上或向下飛行。最大俯仰角可以用數學表達式描述:

(22)

3.3 平滑操作

為滿足無人機特殊的飛行機動特性,需要對算法生成的飛行路徑進行平滑處理。普遍采用的方法是對生成的航跡曲線進行擬合,增加了工程的工作量。三次樣條插值是一種分段的插值方法,它在原軌跡點中利用若干個插值點形成一條平滑的軌跡。三次樣條插值法處理過的無人機飛行軌跡更加平滑,無需對算法規劃的路徑進行曲線擬合。既可以使得無人機在緊急制動或者緊急轉向時具有良好的動態性能,又可以減少路徑生成過程的冗余度。因此,本文采用三次樣條插值曲線來構建和表示無人機飛行路徑。

3.3.1 三次樣條插值

設Δ為區間[a,b]的一部分,若函數S(x)滿足:

1)S(x)∈C2[a,b];

2)S(xi)=f(xi)i= 0,1,2,…,n;

3)S(x)是每個子區間[xi,xi+1](i= 0,1,2,…,n-1)上的三次多項式,那么S(x)被稱為分區Δ的三次樣條函數。

如果在節點xi上給定一個函數值,且S(xi)=yi,則S(x)稱為三次樣條插值函數。為了確定唯一,需要添加3個常見的邊界條件:

1)S(a)=y0,S(b)=yn;

2)S(a)=y0,S(b)=yn;

3)S(a+0)=S(b-0),S(a+0)=S(b-0),S(a+0)=S(b-0)。

3.3.2 代碼設計

三次樣條插值是一種分段插值方法。每段之間的樣條是不一樣的。全部軌跡曲線都是一階連續,而在軌跡節點是二階連續的。因此,軌跡節點的數量表示整條飛行軌跡的轉彎次數。雖然環境條件很復雜,但一般只需要三至五次就可以躲避全部的障礙。根據環境變化選擇飛行軌跡節點的位置。

本文將一條路徑上的所有節點視為螢火蟲的編碼,即一條路徑上n個路徑節點的三維坐標構成螢火蟲個體的n-3維坐標。假設已經確定了n個路徑節點的坐標,分別為r1= (xr1,yr1,zr1),...,rn= (xrn,yrn,zrn),路徑的起點和終點坐標分別為(xs,ys,zs)和(xt,yt,zt)。分別在區間[xs,xr1,xr2,...,xrn,xt],[ys,yr1,yr2,...,yrn,yt]和[zs,zr1,zr2,…,zrn,zt],插值點坐標為(xS1,xS2,…,xSm),(yS1,yS2,…,ySm),(zS1,zS2,…,zSm),并通過三次樣條插值。那么將軌跡點、插值點以及起終點的連接起來就是UAV規劃的飛行路線。

3.3.3 適應度函數與碰撞檢測

無人機的飛行軌跡不能與障礙物產生沖突,路徑長度要求盡可能短。碰撞檢測的目的是檢查規劃的飛行路線是否穿過障礙物。對三次樣條插值后的曲線按路徑節點進行分段。每個節點產生{(x,y)| (xri,yri),(xS1,yS1),…,(xSm,ySm)}和{z|zri,zS1,zS2,…,zSm}集合,兩個集合構成路徑曲線上的一點。

本文將障礙定義為三元函數F=aX+bY+cZ。(a,b,c)表示權重插值。對函數F數據進行插值,得到集合{(x,y)}值對應的插值函數值。將函數值與集合{(z)}進行比較,以確定路徑是否通過障礙物。本文構造的適應度函數為:

z=L??

(23)

?=?1f1+?2f2+?3f3

(24)

(25)

式中,f1、f2、f3分別為飛行高度、最大偏航角、最大仰角的代價。滿足?1+ ?2+ ?3= 1。

4 仿真及分析

4.1 算法性能實驗

為了研究改進FA 的性能,對來自 CEC2013的 28 個測試函數分別利用傳統FA和本文算法獨立執行30 次實驗。在Windows環境下使用MATLAB2018b進行了測試。參數設為隨機參數α=0.2,固定光吸收系數-315,初始亮度I0= 1,系數γ=1,r=0,最大的迭代次數為800,問題維數為10。結果如表1所示。

表1 算法性能對比結果

將本文算法的精度與傳統FA算法進行比較,其中結果是使用30次運行的平均值得到的。從表中可以看出,結果表明,在CEC2013的28個測試函數中,本文算法在測試函數中精度均優于傳統FA算法。從表中的結果可以看出,基于多角度改進策略的螢火蟲算法在單峰函數上表現良好,而對于多峰函數以及組合函數,搜索步長的設置會影響本文提出的算法。另外引入Wilcoxon秩檢驗,在5%的顯著水平下進行試驗,當P≤0.05時,則存在顯著性差異。反之當P>0.05時,則表示差異不顯著,用N/A表示,具體結果如表2所示。由其可以看出,改進算法與基本FA算法存在顯著差異,改進算法更具優勢。

表2 Wilconxon秩檢驗P值

4.2 環境建立

本文以數學模型建立模擬高空環境,主要研究山地地形的UAV的自適應巡航路徑規劃問題。以山峰的形式模擬無人機高空工作面臨的地形狀況,將無人機實際工作環境用三維的山地環境表示,如圖3所示。整個模型空間為150×150×250,x∈(0,150),y∈(0,150),z∈(0,250)。

圖3 高空三維環境建模

4.3 參數設置

為檢驗改進后螢火蟲算法的規劃性能,基于Windows11系統、Matlab R2017a平臺、內存16 GB等硬件環境,分別利用改進前后FA算法進行規劃對比。基于上節中建立的模擬山地環境,并且設置無人機的飛行起始點為(0,0,160),終點為(150,150,160),將無人機視為質點,尺寸忽略不計。算法中主要參數初始值需通過多次仿真實驗獲取,具體如表3所示。

4.4 實例仿真分析

基于模擬山地環境,并按上節參數取值初始化算法中的相關參數,分別利用多角度改進的FA算法和基本FA算法對無人機的航線進行規劃。圖3為基本FA算法規劃的飛行路徑,圖4為本文算法規劃的飛行路徑。兩種算法各項指標對比如表3所示。

圖4 算法改進前規劃路徑

從圖4、圖5及表4可以看出,相比于基本FA算法,多角度改進的FA算法規劃路線的飛行距離更短,轉向次數更少且可以更快尋得最佳飛行路線。具體表現在路徑長度減少7.47%,節點減少31.57%,收斂時間減少18.54%。另外飛行軌跡中沒有尖銳的急轉向角度,整個航線平滑度得以提升,滿足無人機在真實環境下的飛行約束。因此在對無人機飛行路線進行規劃時,本文算法對無人機真實作業場景具有更強的適應力,能夠幫助無人機安全高效的完成飛行作業。而基本螢火蟲算法規劃的飛行路線中存在一些無用路徑,且存在尖銳的轉角,這可能導致無人機硬件損耗增加,無法很好適應復雜環境。

表4 指標對比

圖5 算法改進后規劃路徑

為進一步驗證改進算法的性能,將改進算法與其他規劃算法進行對比,仿真環境與圖2保持一致,各算法結果指標對比如表5所示。

表5 指標對比

由表5可知,本文改進算法相較其他對比算法,在路徑長度、節點數、收斂時間等指標上,都有不同程度的提升,具體表現為路徑長度縮短約為1.14%~7.89%,節點數減少約7.14%~38.10%,收斂時間減少約1.09%~20.53%。因此,驗證了本文算法的優良性能。

4.5 算法消融實驗

為檢驗本文提出的改進策略之間的獨立性及必要性,對本文算法進行消融實驗。基于相同環境,對僅使用Chebyshev映射改進的算法(FA-C)、僅使用logistic映射改進的算法(FA-l)、僅使用Levy飛行策略改進的FA算法(FA-L)與本文改進FA進行比較,結果如表6所示。

表6 消融實驗結果

從表6可看出,相較于基本FA算法,單一的改進策略在路徑長度及收斂時間上均有不同提升,引入Chebyshev映射可以是擴展尋優空間,logistic映射則可以提升算法局部搜索性能及尋優效率,Levy飛行策略能夠擴大搜索范圍,使得結果更加精確。結果表明,同時引入3種策略可有效提高算法性能及效率。

5 結束語

本文針對無人機的最優路徑規劃問題,提出一種多角度改進的螢火蟲算法。首先利用Chebyshev混沌特性初始化種群,改善了初始種群不易產生的問題。針對步長因子過于固定的問題,引入Levy飛行算法,利用萊維分布得到新的位置更新公式和新的步長更新公式,提高了種群的搜索范圍和有效性。最后針對搜索過程中避免陷入局部極值的問題,引入logistic混沌變異改進吸引度系數,使算法更容易跳出當前狀態而逃出局部最優區域,增強了算法局部搜索的能力,保證了算法后期的收斂速度。在模擬高空山地環境仿真后對比數據得出結論,相較改進前螢火蟲算法,本文算法路徑長度減少7.47%,節點減少31.57%,平順度優于改進前,收斂速度更快,收斂時間減少18.54%,相較其他對比算法路徑長度縮短約為1.14%~7.89%,節點數減少約7.14%~38.10%,收斂時間減少約1.09%~20.53%,因此本文具有良好的規劃性能。

猜你喜歡
規劃環境
長期鍛煉創造體內抑癌環境
一種用于自主學習的虛擬仿真環境
發揮人大在五年規劃編制中的積極作用
孕期遠離容易致畸的環境
不能改變環境,那就改變心境
環境
孕期遠離容易致畸的環境
規劃引領把握未來
快遞業十三五規劃發布
商周刊(2017年5期)2017-08-22 03:35:26
多管齊下落實規劃
中國衛生(2016年2期)2016-11-12 13:22:16
主站蜘蛛池模板: 亚洲无码37.| 国产一级α片| 天堂av综合网| 18黑白丝水手服自慰喷水网站| 欧美一区二区福利视频| 最新加勒比隔壁人妻| 亚洲视频免费在线| 久久免费视频6| 久久女人网| 欧美第九页| 国产精品无码久久久久久| 亚洲精品无码抽插日韩| 国产人人干| 亚洲午夜片| 久久久久亚洲Av片无码观看| 在线国产91| 一级成人a毛片免费播放| 好吊色妇女免费视频免费| 亚洲国产精品日韩av专区| 天天综合网在线| 无码'专区第一页| 四虎影视无码永久免费观看| 亚洲精品无码AV电影在线播放| 亚洲国产91人成在线| 国产精品香蕉| 国产成年无码AⅤ片在线| 亚洲天堂精品在线观看| 久久香蕉国产线看精品| 性欧美久久| 久久国产精品77777| 天堂亚洲网| 国产精品一线天| 在线播放国产99re| 午夜欧美理论2019理论| 日韩a级毛片| 亚洲av无码人妻| 成人va亚洲va欧美天堂| 国产大片喷水在线在线视频| 伊人久久大香线蕉综合影视| 原味小视频在线www国产| 国产精欧美一区二区三区| 免费在线a视频| 亚洲区第一页| 曰AV在线无码| 99久久国产综合精品2023| 久久久久久国产精品mv| 中日韩一区二区三区中文免费视频 | 久久亚洲综合伊人| 日本精品视频一区二区| 男人的天堂久久精品激情| 成人年鲁鲁在线观看视频| 亚洲午夜福利在线| 99久久成人国产精品免费| 在线精品亚洲一区二区古装| 亚洲性视频网站| 99re精彩视频| 天天综合网亚洲网站| 69国产精品视频免费| 免费人成又黄又爽的视频网站| 亚洲日韩高清在线亚洲专区| 中文字幕不卡免费高清视频| 在线不卡免费视频| 亚洲综合亚洲国产尤物| 色呦呦手机在线精品| 波多野结衣爽到高潮漏水大喷| 欧美成人在线免费| 中文字幕亚洲电影| 91国内在线观看| 精品国产成人av免费| 日本久久久久久免费网络| 综合色88| 亚洲天堂久久久| 好吊妞欧美视频免费| 欧美精品色视频| 亚洲综合一区国产精品| 免费又爽又刺激高潮网址| 亚洲欧美日本国产综合在线 | 动漫精品中文字幕无码| 国产亚洲精品自在线| 97视频精品全国免费观看| 国产理论一区| 亚洲国产清纯|