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

基于適應性動態步長的變異果蠅優化算法

2016-07-19 19:39:43王行甫陳靜王琳
計算機應用 2016年7期
關鍵詞:優化

王行甫 陳靜 王琳

摘要:針對基本果蠅優化算法(FOA)容易陷入局部最優值、后期收斂速度變慢和收斂精度較低的缺點,提出了一種基于適應性動態步長的變異果蠅優化算法(MFOAADS)。首先,利用佳點集法選取種群初始位置,降低算法初始點選取的隨機性和陷入局部最優值的概率;然后,采用適應性動態步長優化策略,提高收斂速度和求解精度;最后,若算法陷入了早熟,則對種群最優個體按一定概率執行柯西變異擾動,賦予其跳出局部最優的能力。經5個經典函數測試表明,固定迭代次數時MFOAADS的收斂精度與收斂速度明顯優于FOA;固定目標精度時,MFOAADS相對于FOA平均迭代次數有著大幅下降且成功率達97%以上。實驗結果表明,所提算法求解精度、運行效率以及可靠性相對于基本FOA算法都有著顯著提高。

關鍵詞:

果蠅優化算法;早熟收斂;佳點集;適應性動態步長;柯西變異

中圖分類號: TP181 文獻標志碼:A

0引言

群智能優化算法常用于求解過于復雜的優化問題,已廣泛應用于工程和科學領域。果蠅優化算法(Fruit Fly Optimization Algorithm, FOA)是2011年由潘文超[1]提出的一種尋求全局最優解的群智能算法。該算法主要模擬果蠅通過嗅覺和視覺能力發現食物的過程來實現尋優,與遺傳算法、蟻群算法、粒子群算法等典型群智能算法相比,該算法簡單易懂、所需調參數量少、尋優效率較高,因此較容易應用到實際問題的解決當中[2]。目前FOA已成功將其應用于求解數學函數極值[3]、微調ZSCORE模型系數優化[3]、廣義回歸神經網絡參數優化[4]與支持向量機參數優化[5]等,具有潛在的研究價值。

但與其他群智能算法一樣,FOA也存在著陷入早熟的風險。在近幾年的研究中,許多關于FOA的改進和應用逐步涌現。文獻[6]中將混合蛙跳算法的局部深度搜索策略融入到果蠅優化算法中來平衡種群的探索和開發能力。文獻[7]通過向算法中引入一個新的控制參數自適應調整搜索范圍內的群位置來解決連續函數優化問題。文獻[8]則提出一種基于云模型的自適應調整參數的果蠅優化算法。這些算法都是從一個側面進行了改進,在一定程度上優化了算法性能,但都未能完全解決算法陷入早熟的問題。為此,本文提出一種基于適應性動態步長的變異果蠅優化算法(Mutation Fruit Fly Optimization Algorithm based on Adaptive Dynamic Step Size, MFOAADS),其中采用弱化種群位置初始化的隨機性、自適應動態變化步長、陷入早熟防控機制3個措施對算法進行了更加全面的改進。實驗表明該算法能在一定程度上有效擺脫局部極值的干擾,且所達到的收斂速度和求解精度都顯著優于基本FOA。

1基本FOA

果蠅是一種嗅覺和視覺都優于大多數動物的一種生物,它利用靈敏的嗅覺器官搜集空氣中的各種氣味并向食物飛近,再利用視覺能力發現食物或同伴聚集的位置,并向該方向飛去。FOA就是基于果蠅的這種自然行為設計而成的一種新型全局尋優方法。

FOA基本步驟歸納如下。

1)初始化參數。參數包括種群規模N、種群迭代次數Maxgen以及果蠅群體初始位置X_axis、Y_axis。

2)根據果蠅群體位置,初始化群體中每個果蠅個體的位置坐標,使該個體利用嗅覺尋找食物,h為果蠅的搜索步長,rand取0到1之間的隨機值:

Xi=X_axis+2h*rand-h

Yi=Y_axis+2h*rand-h (1)

3)由于無法得知食物位置,因此先計算果蠅個體與原點的距離Disti,在計算味道濃度判定值Si,此值為距離的倒數。

Disti=X2i+Y2i

Si=1/Disti (2)

4)將味道濃度判定值代入味道濃度判定函數Function(或稱為適應度函數),用來求出果蠅個體位置的味道濃度Smelli:

Smelli=Function(Si)(3)

找出該果蠅群體中味道濃度最高的果蠅個體:

[bestSmellbestIndex]=max(Smell)(4)

其中:bestSmell為最高的味道濃度;bestIndex為果蠅物體中具有最高味道濃度的果蠅個體序號;Smell為果蠅群體的味道濃度集合。

5)判斷味道濃度是否優于前一代的最優濃度:若是,則保留最佳味道濃度值與最佳果蠅個體的位置坐標,使果蠅群體利用視覺向該位置飛去;否則,重復執行第2)步至第6)步。

smellBest=bestSmell

Xbest=X(bestIndex)

Ybest=Y(bestIndex) (5)

6)判斷當前迭代次數是否達到最大迭代數Maxgen:若是則輸出果蠅最優個體位置;否則,執行第2)步。

FOA有著結構簡單、控制參數少、易調節、計算量較小、全局尋優能力強、收斂速度快的優點;但由于種群初始位置由隨機方式決定,算法搜索步長固定且沒有充分利用目標函數所提供的信息,使得迭代后期收斂速度變慢,容易陷入局部最優值,在復雜高維多極值優化問題更是如此。

2適應性動態步長的變異FOA

針對基本FOA的缺點,MFOAADS分別從3個方面進行改進:通過佳點集法選取種群初始位置,降低初始位置選取的隨機性;通過種群適應值信息動態調整算法的步長大小來平衡全局搜索能力和局部搜索能力,進而提高算法的執行效果和效率;在算法中加入早熟判斷和柯西變異機制,進一步避免算法陷入局部最優值。

2.1佳點集法選取種群初始位置

在基本FOA中大多數都是隨機選擇果蠅群體初始點來初始化果蠅群體,這可能導致群體初始點距離最優位置較遠或容易陷入局部最優點,從而影響算法的全局尋優效率和搜索精度。佳點集方法最初由華羅庚等[9]提出的一種有效的、可以減少實驗次數的實驗方法。佳點集的定義與構造[10]:設Gs是S維歐氏空間中的單位立方體,如果r∈Gs,形為pn(m)={({r1(n)*m},…,{r(n)i*m},…,{r(n)s*m})|1≤m≤n},其偏差φ(n)滿足φ(n)=C(r,ε)n-1+ε,其中C(r,ε)是只與r,ε(ε>0)有關的常數,則稱pn(m)為佳點集,r稱為佳點。取rm={2 cos (2πk/p,1≤m≤s}(p是滿足(p-s)/2≤s的最小素數),或rm={exp(m)},1≤m≤s,{a}表示a的小數部分。

佳點集產生的點無重疊點且分布均勻具有較好的多樣性,如圖1所示。圖1中展示的是橫縱坐標區間均為[0,1]的二維空間中佳點集法所產生的點。在相同取點個數的條件下,佳點序列要比其他方法選取的點序列更均勻[11]。 因此,本文采用佳點集方法在事先規定的區間內產生點,并從中選擇味道濃度最高的點作為果蠅群體初始點,然后初始化種群。

2.2適應性動態步長調整策略

基本FOA中果蠅個體以群體位置為中心,以固定步長進行隨機搜索,當步長過大時全局尋優能力較強但是容易跳過最優值,穩定性較差,當步長過小時局部尋優能力較強但是容易陷入局部最優解,收斂速度較慢且求解精度較低,因此步長的選擇直接影響算法的執行效率和效果。FOA具有很大的隨機性,并且搜索過程因問題而異,所以不論固定步長還是事先規定好的步長變化都不能很好地適應搜索空間中的不確定性。在基本FOA中,種群迭代進化中只是單純利用了適應度函數對應的適應度值,沒有充分利用其中包含的種群信息,使得最優值搜索方向的啟發性不強。因此,本文提出一種適應性動態步長調整策略,基于種群反饋的適應值變化信息來調整算法中步長h的大小。

第t代種群的適應度平均值Mt計算公式如下:

Mt=(∑Ni=1fi(t))/N(6)

其中:N為種群規模; fi(t)表示第i個果蠅在第t代時的適應度值。假設此時是最小化情況,種群反饋的適應值變化信息由下面的定義表示。

定義1種群平均適應值的相對變化率為:

k=(Mt-1-Mt)/Mt-1(7)

其中t≤2。

當變化率k較大時,表明搜索空間相對步長較為平滑和單調,表明了當前步長對搜索空間的開發率較高,所以此時適當地增加步長可以提高種群的探索性能,有效減少算法的執行次數,提高算法的搜索效率。當k較小時,表明搜索空間相對步長較為復雜,找到更好的解的概率降低,所以此時適當地減小步長可增加種群對搜索空間的開發性,提高算法的搜索精度,因此,本文提出的步長調整公式如下:

ht+1=min [ht(1+α ln(1+k)),hmax],k>1

ht,1≤k≤0.01

max [ht/(1+βe-k),hmin],k<0.01(8)

其中:ht為種群第t代搜索步長;α和β均為步長變化幅度調節參數。由于在不同的迭代進程中k的值變化較大,本文在增加步長時使用自然對數,在減小步長時使用e作為指數。它們在工程中較常使用,可以明顯降低其變化幅度,改善ht+1的平滑性。同時保證步長值在最大步長和最小步長之間。從式(8)可看出,當k較大時,k越大算法搜索步長越大,趨向極值點的速度加快,全局搜索能力增強;當k較小時,k越小越接近極值點,對應步長減小幅度越大,局部搜索能力也越強。

2.3早熟判斷和柯西變異策略

本文通過使用佳點集方法篩選得到較好的果蠅群體初始點,且對算法中的搜索步長進行適應性動態調整優化,雖然能提高算法的執行效率和降低算法發生早熟收斂的可能性,卻無法保證避免出現早熟現象,算法仍然有可能陷入局部最優值。

本文采用種群適應度方差σ2[12]對算法是否陷入早熟進行判斷,種群適應度方差反映了種群的收斂程度,當σ2越小,說明種群的聚集程度越高,種群多樣性越低,尋優進程發生停滯。在本文中,當σ2≤δ(δ為適應度方差閾值),判斷此時算法趨于收斂;同時為排除全局收斂的情況設置當前全局最優適應值大于理論最優適應度或目標精度γ,則判定算法陷入了早熟。種群適應度方差σ2計算公式如下:

σ2=∑Ni=1(fi-Mf)2(9)

其中: fi表示當前第i個果蠅的適應度值;M為當前種群的平均適應度值; f為用于限制σ2大小的歸一化定標因子。 f的取值需符合兩個條件:1)歸一化后,要保證種群中所有個體的|fi-M|的最大值小于等于1;2)f的值隨算法變化而變化。 f的具體取值公式如下:

f=max(fi-M),max(fi-M)>1

1,其他(10)

當算法滿足陷入早熟的條件時,算法陷入僵局,種群的果蠅個體無法找到更優解,此時本文給算法添加一種擾動機制,賦予果蠅種群飛出局部最優值的能力。本文中以一定概P執行柯西變異對種群中最優個體進行V次擾動,增加種群的多樣性,然后對變異后的果蠅群體進行二次尋優,避免算法陷入局部最優。柯西變異操作如下:

其中:Xbest,、Ybest為當前最優果蠅個體位置, ρ為擾動幅度調節參數,Cauchy(0,1)為標準柯西分布的隨機數生成器。Cauchy(0,1)的隨機變量生成函數為ω=tan[π(η-0.5)](其中η是[0,1]上服從均勻分布的隨機變量)。

2.4MFOAADS流程

基于以上基本原理和改進方法,本文提出的適應性步長的變異果蠅優化算法的基本流程如下。

1)初始化種群規模N,種群迭代次數Maxgen,適應度方差閾值δ,變異概率P,采用佳點集方法產生果蠅群體初始位置X_axis、Y_axis。

2)執行基本FOA的第2)步到第6)步。

3)根據式(6)~(8)計算出種群平均適應值的相對變化率,從而確定算法步長h的值。

4)由式(9)、(10)計算果蠅群體適應度方差σ2,如果σ2≤δ且當前種群適應度最優值smellBest大于理論最優適應度或目標精度γ,則判斷算法陷入局部最優,重新初始化步長,然后轉第5)步;否則,轉第6)步。

5)若在[0,1]上服從均勻分布的隨機數C小于P,使用式(11)對果蠅種群進行柯西變異操作,然后執行基本FOA的第3)步到第6)步。

6)若當前迭代次數等于最大迭代數Maxgen或當前種群適應度最優值達到理論最優適應度或目標精度要求,則輸出全局最優適應值和最優果蠅位置;否則,轉向第2)步,進入下一代尋優。

在MFOAADS流程中將能反映種群狀況的適應值信息用于下一次迭代過程的搜索步長的確定,提高算法尋優效果和效率。同時為監督算法是否陷入局部最優值,在算法流程加入了早熟判斷機制,且一旦發現早熟情況,就采用柯西變異策略使算法跳出局部最優,重新進入全局最優值的搜索過程,從而大大減少算法陷入早熟的幾率。

3實驗與結果分析

3.1實驗設計

實驗時采用Matlab進行算法實現,使用聯想筆記本電腦,雙核2.4GHz,4GB內存,Windows 8操作系統。為了驗證本文提出的MFOAADS的性能,本文設計了FOA優化實驗和MFOAADS優化實驗這兩類測試實驗。為測試算法的尋優性能,選取5個經典函數進行算法測試(本實驗為求解最小值的情況)。這些函數具有函數優化時所碰到的普遍特性:多峰性、多維性和凹凸性等,所以它們也常常用來測試智能算法的尋優性能。函數形式、維數、搜索區間、理論極值和函數峰值特征如表1所示。

待優化函數的維數越高,自變量范圍越大,目標精度越高,優化過程的難度就越大[13]。為了便于比較和突出MFOAADS的性能,本文選用較為嚴苛的實驗參數,具體參數設置為:群體規模N=30,最大迭代數Maxgen=2000,初始步長h=1, f1、 f2、 f3適應度方差閾值δ=0, f4、 f5適應度方差閾值為δ=1E-8,種群搜索區間如表1所示,步長調節參數α=1、 β=0.1,擾動幅度調節參數為ρ=0.5,變異概率為P=0.3。

同時在本次實驗中,將實驗結果與參考文獻[14-16]中的算法進行對比。文獻[14]中混合布谷鳥算法的差分進化算法(hybrid optimization algorithm of Cuckoo Search and Differential Evolution, CSDE)是在向差分進化算法中引入布谷鳥搜索算法來增加粒子的搜索活力,有效地克服了差分進化算法的缺陷,使尋優精度有較大提高。參考文獻[15]中一種新的自適應粒子群優化(new Adaptive Particle Swarm Optimization, APSO)算法中的全局最優位置與個體最優位置分別替換為相關個體最優位置的加權平均,更好地平均了算法的全局與局部搜索能力,極大地提高了算法的多樣性和搜索效率。參考文獻[16]中全局版人工魚群算法(Global edition Artificial Fish Swarm Algorithm, GAFSA)用整個魚群的中心位置和全局極值位置代替人工魚的鄰域中心位置和鄰域極值位置;同時采用動態調整視野和步長的方法,提高搜索精度和速度。這些算法表現出的較高的尋優性能,通過與它們進行對比,從而在一定程度上證明本文所提出的算法相對于其他智能優化算法的優越性。

3.2實驗結果分析與比較

將5個測試函數進化次數固定為2000次,分別采用FOA和MFOAADS進行50次獨立運行實驗,并將實驗結果與CSDE算法和APSO算法相比較,如表2所示(注:表2中的“—”表示參考文獻中沒有介紹)。從表2中可看出,MFOAADS在所有函數測試中的性能表現均優于基本FOA,其中f1、 f2的平均最優值都達到了全局最優解0,對于f3、 f4、 f5的平均最優值相對于基本FOA分別降低了9、2和8個數量級,且所有測試函數的標準差均優于基本FOA,說明MFOAADS顯著地提高函數優化過程中的尋優效果和穩定性。MFOAADS在5個測試函數中的平均最優值和標準差都優于CSDE算法。同時除去擁有連續高峰的函數Ackley(f5)性能表現相對于APSO算法較差,MFOAADS其他測試函數的平均最優值和標準差均優于APSO算法。由此可見,與文獻[10-11]中的算法相比,MFOAADS雖然在Ackley函數上的尋優性能稍差;但在其他函數上它的尋優效果和穩定性還是具有很大優勢的。

圖2是本文兩類算法的3個函數適應度對數值進化曲線,其實線MFOAADS,虛線表示基本FOA(注:由于篇幅問題僅給出3個函數圖像,為了方便觀察進化曲線,本文對函數適應度取以10為底的對數,其中由于Schaffers極值為-1,將其加1再取其以10為底的對數來觀察)。由圖2可看出,在f2中MFOAADS的最優值收斂曲線幾乎成直線下降,成功收斂到全局最優值, f1也取得了類似的效果。在其他幾個函數的最優值尋找過程中的后期,MFOAADS也比FOA更容易跳出局部最優值。

對于5個測試函數,分別采用FOA和MFOAADS經過50次獨立運行后,記錄達到各函數目標精度的平均迭代次數和成功率 (成功率=達到精度的運行次數/實驗總次數),且與APSO算法和GAFSA算法比較,如表3所示(注:表3中的“—”表示參考文獻中沒有介紹)。MFOAADS在所有測試函數中的達到目標精度的平均迭代次數均小于基本FOA、APSO算法以及GAFSA,且相應成功率均達到了97%以上??梢婋m然MFOAADS在Schaffers函數上的成功率比GAFSA稍差,它的尋優效率和成功率的總體表現比其他算法更加突出。

4結語

針對果蠅優化算法中的缺點,本文提出了一種適應性動態步長的變異果蠅優化算法,分別從種群初始化、算法搜索步長、擺脫局部最優機制3個方面作出改進。通過對5個經典測試函數進行實驗對比,可看出本文所提出的算法在收斂速度、收斂進度以及收斂可靠性上比基本FOA有較大提高,其中部分函數成功擺脫局部最優值的干擾,找到對應的全局最優值;但同時本文提出的算法仍無法完全解決陷入早熟的問題,對如何更有效地平衡算法的全局搜索能力和局部搜索能力還有待于進一步的解決。

參考文獻:

[1]

潘文超.果蠅最佳化演算法[M].臺北:滄海書局,2011:10-12.(PAN W T. Fruit Fly Optimization Algorithm [M]. Taipei: Tsang Hai Book Publishing, 2011: 10-12.)

[2]

吳小文,李擎.果蠅算法和5種群智能算法的尋優性能研究[J].火力與指揮控制,2013,38(4):17-20.(WU X W, LI Q. Research of optimizing performance of fruit fly optimization algorithm and five kinds of intelligent algorithm [J]. Fire Control & Command Control, 2013, 38(4): 17-20.)

[3]

PAN W T. A new fruit fly optimization algorithm: taking the financial distress model as an example [J]. KnowledgeBased Systems, 2012, 26: 69-74.無期

[4]

潘文超.應用果蠅優化算法優化廣義回歸神經網絡進行企業經營績效評估[J].太原理工大學學報:社會科學版,2011,29(4):1-5.(PAN W T. Using fruit fly optimization algorithm optimized general regression neural network to construct the operating performance of enterprises model [J]. Journal of Taiyuan University of Technology (Social Sciences Edition), 2011, 29(4): 1-5.)

[5]

王欣,杜康,秦斌,等.基于果蠅優化算法的LSSVR干燥速率建模[J].控制工程,2012,19(4):630-633.(WANG X, DU K, QIN B, et al. Drying rate modeling based on FOALSSVR [J]. Control Engineering of China, 2012, 19(4): 630-633.)

[6]

劉成忠,黃高寶,張仁陟,等.局部深度搜索的混合果蠅優化算法[J].計算機應用,2014,34(4):1060-1064.(LIU C Z, HUANG G B, ZHANG R Z, et al. Shuffled fruit fly optimization algorithm with local deep search [J]. Journal of Computer Applications, 2014, 34(4): 1060-1064.)

[7]

PAN Q K, SANG H Y, DUAN J H, et al. An improved fruit fly optimization algorithm for continuous function optimization problems [J]. KnowledgeBased Systems, 2014, 62: 69-83.無期

[8]

韓俊英,劉成忠.自適應調整參數的果蠅優化算法[J].計算機工程與應用,2014,50(7):50-55.(HAN J Y, LIU C Z, Fruit fly optimization algorithm with adaptive parameter [J]. Computer Engineering and Applications, 2014, 50(7): 50-55.)

[9]

華羅庚,王元.數論在近似分析中的應用[M].北京:科學出版,1978:83-84.(HUA L G, WANG Y. Number Theory in the Application of Approximate Analysis [M]. Beijing: Science Press, 1978: 83-84.)

[10]

張鈴,張鈸.佳點集遺傳算法[J].計算機學報,2001,24(9):917-922.(ZHANG L, ZHANG B. Good point set based genetic algorithm [J]. Chinese Journal of Computers, 2001, 24(9): 917-922.)

[11]

龍文,梁昔明,徐松金,等.聚類佳點集交叉的約束優化混合進化算法[J].計算機研究與發展,2015,49(8):1753-1761.(LONG W, LIANG X M, XU S J, et al. A hybrid evolutionary algorithm based on clustering goodpoint set crossover for constrained optimization [J]. Journal of Computer Research and Development, 2015, 49(8): 1753-1761.)

[12]

呂振肅,侯志榮.自適應變異的粒子群優化算法[J].電子學報,2004,32(3):416-420.(LYU Z S, HOU Z R. Particle swarm optimization with adaptive mutation [J]. Acta Electronica Sinica, 2004, 32(3): 416-420.)

[13]

王凌.智能優化算法及其應用[M].北京:淸華大學出版社,2001:148-149.(WANG L. Intelligent Optimization Algorithm and Its Application [M]. Beijing: Tsinghua University Press, 2001: 148-149.)

[14]

李明,曹德欣.混合CS算法的DE算法[J].計算機工程與應用,2013,49(9):57-60.(LI M, CAO D X. Hybrid optimization algorithm of cuckoo search and DE [J]. Computer Engineering and Applications, 2013, 49(9): 57-60.)

[15]

林川,馮全源.一種新的自適應粒子群優化算法[J].計算機工程,2008,34(7):181-183.(LIN C, FENG Q Y. New adaptive particle swarm optimization algorithm [J]. Computer Engineering, 2008, 34(7): 181-183.

[16]

王聯國,洪毅,施秋紅.全局版人工魚群算法[J].系統仿真學報,2009,21(23):7483-7502.(WANG L G, HONG Y, SHI Q H. Global edition artificial fish swarm algorithm [J]. Journal of System Simulation, 2009, 21(23): 7483-7502.)

猜你喜歡
優化
超限高層建筑結構設計與優化思考
房地產導刊(2022年5期)2022-06-01 06:20:14
PEMFC流道的多目標優化
能源工程(2022年1期)2022-03-29 01:06:28
民用建筑防煙排煙設計優化探討
關于優化消防安全告知承諾的一些思考
一道優化題的幾何解法
由“形”啟“數”優化運算——以2021年解析幾何高考題為例
圍繞“地、業、人”優化產業扶貧
今日農業(2020年16期)2020-12-14 15:04:59
事業單位中固定資產會計處理的優化
消費導刊(2018年8期)2018-05-25 13:20:08
4K HDR性能大幅度優化 JVC DLA-X8 18 BC
幾種常見的負載均衡算法的優化
電子制作(2017年20期)2017-04-26 06:57:45
主站蜘蛛池模板: 亚洲国产中文在线二区三区免| 老色鬼欧美精品| 婷婷六月天激情| 福利一区在线| 日韩欧美高清视频| 国产麻豆永久视频| 国产不卡在线看| 久久一色本道亚洲| 免费va国产在线观看| 亚洲美女操| 国产性生交xxxxx免费| 伊人色在线视频| 亚洲国产成人精品无码区性色| 国产一级毛片yw| 亚洲Va中文字幕久久一区| 亚洲欧洲自拍拍偷午夜色| 免费看a级毛片| 亚洲午夜久久久精品电影院| 欧美日韩国产在线观看一区二区三区 | 免费 国产 无码久久久| 国产精品成人啪精品视频| 国产成人一区二区| 国产自在线拍| 国产成人精品视频一区二区电影| 免费A∨中文乱码专区| 26uuu国产精品视频| 欧美啪啪视频免码| 欧美爱爱网| 91无码人妻精品一区二区蜜桃| 国产成人高清精品免费软件| 欧美午夜性视频| 免费毛片网站在线观看| 欧美激情成人网| 国产一区二区精品福利| 亚洲三级成人| 亚洲精品无码不卡在线播放| 国产中文一区二区苍井空| 中字无码精油按摩中出视频| 欧美第九页| 一级毛片网| 国产无吗一区二区三区在线欢| 国产网站免费看| 国产喷水视频| 91热爆在线| 中文字幕丝袜一区二区| 欧美在线天堂| 91国内外精品自在线播放| 国产丝袜第一页| 久青草免费在线视频| 久久精品亚洲热综合一区二区| 中文字幕欧美成人免费| 国产精品久久久久无码网站| 91欧美在线| 天天色天天综合| 久久这里只有精品8| 国产呦视频免费视频在线观看| 东京热一区二区三区无码视频| 久久香蕉国产线看精品| 91精品国产一区自在线拍| 亚洲高清中文字幕| 人妻精品全国免费视频| 国产视频入口| 欧美特黄一级大黄录像| 午夜免费视频网站| 欧美一级黄片一区2区| 国产精品亚欧美一区二区三区 | 亚洲首页在线观看| 国产精品久线在线观看| 国产精品亚洲综合久久小说| 亚洲福利视频一区二区| 国产高清在线观看| 无码免费试看| 久久91精品牛牛| 欧美在线一级片| 亚洲精品福利视频| 九色视频在线免费观看| 亚洲国产精品久久久久秋霞影院| 成人毛片免费在线观看| 美美女高清毛片视频免费观看| 欧美日韩成人| 国产综合无码一区二区色蜜蜜| 久久久久九九精品影院|