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

面向離散優(yōu)化問題的改進(jìn)二元粒子群算法

2015-06-24 13:29:17印桂生崔曉暉董宇欣楊雪
關(guān)鍵詞:優(yōu)化

印桂生,崔曉暉,董宇欣,楊雪

(哈爾濱工程大學(xué)計算機(jī)科學(xué)與技術(shù)學(xué)院,黑龍江哈爾濱150001)

面向離散優(yōu)化問題的改進(jìn)二元粒子群算法

印桂生,崔曉暉,董宇欣,楊雪

(哈爾濱工程大學(xué)計算機(jī)科學(xué)與技術(shù)學(xué)院,黑龍江哈爾濱150001)

二元粒子群算法被廣泛用于求解離散組合優(yōu)化問題。在求解離散優(yōu)化問題時,二元粒子群算法會出現(xiàn)解空間利用率低,速度和狀態(tài)趨同以及退化和波動等演化問題。針對這些問題,提出一種改進(jìn)的二元粒子群算法。算法使用Gray碼演化基編碼,混沌初始化過程,改進(jìn)速度和狀態(tài)調(diào)整方法以及子代處理方法用于提高種群利用率和種群多樣性。在不同類型的檢驗函數(shù)以及多選擇背包問題上,和現(xiàn)有優(yōu)化算法及其他二元粒子群算法相比,改進(jìn)算法能夠獲得較高的收斂精度以及較快的收斂速度,體現(xiàn)出多離散優(yōu)化問題的實際效用。

二元粒子群;Gray碼;混沌;子代處理;離散優(yōu)化

粒子群優(yōu)化(particle swarm optimization,PSO)用于求解連續(xù)解空間的優(yōu)化問題[1]。在實際中,存在這大量離散優(yōu)化問題[2?4],為解決這類問題,學(xué)者們提出了用于離散問題的二元PSO[5?10]。

隨問題維度的上升,現(xiàn)有二元PSO通常增大種群規(guī)模或增加迭代次數(shù)來提高求解精度,導(dǎo)致算法的時間復(fù)雜度以及空間復(fù)雜度隨之升高,喪失了優(yōu)化算法在處理離散問題上的實際效用。為此,本文提出一種面向離散問題的改進(jìn)二元PSO算法。該算法采用的Gray碼演化基編碼方式,混沌初始化過程,改進(jìn)的速度和狀態(tài)調(diào)整過程以及子代處理過程,旨在問題維度上升的情況下,仍能夠高效且穩(wěn)定的收斂到高精度的演化狀態(tài),提高算法對離散問題的實際效用。

1 經(jīng)典二元粒子群算法及其問題

Kennedy借鑒經(jīng)典連續(xù)PSO的速度調(diào)整和狀態(tài)調(diào)整思路,設(shè)計了用于離散問題的二元粒子群算法BPSO(binary particle swarm optimization)[5],其速度和狀態(tài)調(diào)整過程表示為

式中:ωvi,j(t)為第t代的第i個粒子在第j維的速度,t)∈{0,1}分別表示第i個粒子的個體最優(yōu)位置和種群最優(yōu)位置。U(0,φ1)和U(0,φ2)生成(0,1)間的隨機(jī)數(shù)。vmax表示速度的上限。

式中:sig(x)為sigmoid函數(shù),用于將速度變量轉(zhuǎn)換成(0,1)內(nèi)的概率。

BPSO以及基于BPSO的改進(jìn)算法[5?6,8?10]主要存在以下問題:解空間的利用率較差。由于缺乏有效的演化基初始化策略,導(dǎo)致算法無法充分利用有限的種群規(guī)模。速度趨同現(xiàn)象。迭代后期,BPSO隨著各個維速度的不斷累加,演化基上各維狀態(tài)只能無限趨近于1,從而出現(xiàn)趨同現(xiàn)象[8]。狀態(tài)遺忘特征。BPSO現(xiàn)有狀態(tài)完全依賴于速度轉(zhuǎn)換后的概率,導(dǎo)致狀態(tài)調(diào)整過程呈現(xiàn)遺忘特征。演化波動以及退化現(xiàn)象。BPSO的收斂精度曲線隨迭代次數(shù)增加而不斷波動且出現(xiàn)退化的現(xiàn)象。

針對上述問題,提出一種面向離散優(yōu)化問題的改進(jìn)二元粒子群算法(improved binary particle swarm opti?mization,IBPSO),采用混沌初始化種群,引入二維速度變量并設(shè)計子代處理過程。

2 高效離散粒子群算法

2.1 IBPSO的編碼方式

IBPSO的演化基編碼形式包括:用于適應(yīng)度計算以及滿足約束判斷的整數(shù)演化基形式I;用于中間狀態(tài)的二進(jìn)制BCD碼形式B;用于演化過程的Gray碼形式。各形式的表示以及轉(zhuǎn)換方法如下:

針對種群規(guī)模為s和問題規(guī)模為n的離散問題,種群中某一整數(shù)演化基Ik表示為

假設(shè)Ik中每一維ikl的取值約束為[al,bl],則可通過ikl-al將ikl轉(zhuǎn)換為無符號的正整數(shù)i'kl進(jìn)而轉(zhuǎn)換為無符號的二進(jìn)制BCD碼Bkl,即

為便于IBPSO在約束范圍內(nèi)有效演化,通過左補0的方式使Bkl的寬度等于「lb(bl-al)。

式(4)的BCD編碼雖然是整數(shù)演化基的直觀二元表示,但存在Hamming懸崖問題[9],因此,IBPSO采用式(6)將BCD編碼的演化基轉(zhuǎn)換為對應(yīng)的Gray碼Gkl=(gkl,1,gkl,2,…,gkl,d),并使用Gkl作為種群中速度和狀態(tài)調(diào)整的操作對象。

2.2 算法改進(jìn)思路

2.2.1 混沌初始化方法

可利用混沌過程產(chǎn)生初始種群,充分利用解空間,形成分布均勻的初始種群。IBPSO的混沌初始化的步驟為:

1)產(chǎn)生與問題規(guī)模n一致的(0,1)隨機(jī)種子集合{Z0l}(l=1,2,…,n)。

2)將Z0l作為輸入,通過式(6)所示的混沌發(fā)生器,產(chǎn)生2倍種群規(guī)模的{Zkl}(k=1,2,…,2s)。

3)將具有同樣下標(biāo)k的ikl構(gòu)成候選初始整數(shù)演化基Ik=(ik1,ik2,…,ikn)。根據(jù)適應(yīng)度,對2s個候選初始整數(shù)演化基進(jìn)行排序,選擇適應(yīng)度較優(yōu)的前s個不同個體作為初始種群。

4)根據(jù)l維的約束條件[al,bl],通過式(7),將Zkl映射為相應(yīng)的整數(shù)變量ikl。

5)利用2.1節(jié)中的介紹的方法,將整數(shù)演化基轉(zhuǎn)換為Gray碼演化基,并初始化速度變量Vk。對于Vkl(l=1,2,…,n)的每一維變量vkl,o,初始化為[vkl,o(t,0),vkl,o(t,1)](o∈[1,2,…,d]),其中,vkl,o(t,0)為第t代第k個Gray碼演化基的第l維上所對應(yīng)二進(jìn)制的第o個分量趨近于0的程度,vkl,o(t,1)為趨近1的程度。vkl,o的初值為vkl,o(0,0)=vkl,o(0,1)=vmax/2。

2.2.2 改進(jìn)的速度調(diào)整過程

以二維速度向量為基礎(chǔ),改進(jìn)BPSO的速度調(diào)整過程,采用式(8)調(diào)整演化基中每一分量gkl,o(o∈[1,2,…,d])的速度。

式中:vkl,o(t-1,s)表示在t-1代演化基接近s的程度,pkl,o(t-1)表示在t-1代,個體k的最優(yōu)演化基在該維的取值,gkl,o(t-1)表示在t-1代,種群最優(yōu)演化基在該維的取值。ω為非慣性權(quán)重。當(dāng)s取0或1時,IBPSO依據(jù)vkl,o(t-1,s),pkl,o(t-1)和gkl,o(t-1)的取值情況,同時調(diào)整vkl,o(t,0)和vkl,o(t,1),如式(9)。

由于pkl,o(t-1),gkl,o(t-1)∈{0,1},式(9)還可繼續(xù)轉(zhuǎn)換為4種形式,以其中pkl,o(t-1)=1和gkl,o(t-1)=1為例。

當(dāng)pkl,o(t-1)=1時,說明個體k的最優(yōu)演化基在該維的狀態(tài)為1,此時,為促進(jìn)個體k的快速收斂到最優(yōu)狀態(tài),應(yīng)增強趨近于1的程度,同時,減弱趨近于0的程度。當(dāng)gkl,o(t-1)=1時,說明種群中最優(yōu)演化基在該維的狀態(tài)為1,為促進(jìn)個體k向全局狀態(tài)收斂,應(yīng)增強趨近于1的程度,同時,減弱趨近于0的程度。根據(jù)上述思路,粒子在解空間中可自由移動,避免速度趨同現(xiàn)象。

2.2.3 改進(jìn)的狀態(tài)調(diào)整過程

IBPSO采用式(10)調(diào)整Gray碼演化基中每一分量gkl,o(o∈[1,2,…,d])的狀態(tài)。

式中:U(0,1)為(0,1)的隨機(jī)數(shù),sgn(x)為符號函數(shù)。當(dāng)x>0時,sgn(x)=1;否則,sgn(x)=0,sig(x)為sig?moid函數(shù)。

以gkl,o(t-1)=0為例,考察gkl,o(t-1)的相反取值gkl,o(t-1)所對應(yīng)的速度vkl,o(t-1,1)和隨機(jī)產(chǎn)生數(shù)值U(0,1)的大小關(guān)系,確定符號函數(shù)sgn的取值。當(dāng)sgn=0時,說明gkl,o(t-1)的取值更趨近于1,則利用式(10)將gkl,o(t-1)取反,當(dāng)sgn=1時,說明gkl,o(t-1)的取值更趨近于0,保持gkl,o(t-1)不改變。在gkl,o(t-1)=1,式(10)的狀態(tài)調(diào)整思路同gkl,o(t-1)=0。上述過程中IBPSO的狀態(tài)調(diào)整參考前一代的狀態(tài),緩解了迭代過程中的狀態(tài)趨同問題。

2.2.4 子代處理過程

將演化基上不滿足約束的維度調(diào)整到約束范圍內(nèi),如式(11)所示:

篩選并形成子代的過程:將調(diào)整后的候選演化基I'(t)和父代演化基I(t-1)構(gòu)成數(shù)量為2s的集合IS={I'(t),I(t-1)},按適應(yīng)度排序,將排序中適應(yīng)度較大且不同的前s個演化基形成子代I(t)。

2.3 IBPSO的主要步驟

1)執(zhí)行混沌初始化操作。生成初始Gray碼演化基G、速度變量V、種群和個體最優(yōu)位置。

2)速度以及狀態(tài)調(diào)整。更新速度變量,生成候選子代G'(t)。

3)按約束調(diào)整演化基。將不滿足約束的整數(shù)演化基調(diào)整為有效演化基。

4)篩選并形成子代。依據(jù)子代篩選過程,從IS中生成子代I(t)。

5)更新操作。更新演化基每一維的最優(yōu)位置以及種群最優(yōu)位置。

6)結(jié)束判斷。

2.4 IBPSO的時間復(fù)雜度

IBPSO的時間復(fù)雜度包括初始化、速度調(diào)整、狀態(tài)調(diào)整及子代處理等過程。初始化的時間復(fù)雜度為

O(n(2+2s)+2sz(n)+2slb2s+

速度和狀態(tài)調(diào)整過程時間復(fù)雜度為

子代處理的時間復(fù)雜度為

如果假設(shè)IBPSO設(shè)定的結(jié)束代數(shù)為tmax,則IBPSO時間復(fù)雜度為

2.5 IBPSO的收斂性分析

在任意代t,種群中每個的粒子,通過式(8)和(10)確定Ik轉(zhuǎn)換的Gray碼表示Gk上的每一維轉(zhuǎn)移概率pkl,o(t)。由于pkl,o(t)∈(0,1),因此,粒子的每一步轉(zhuǎn)移概率為(k)>0,轉(zhuǎn)移過程構(gòu)成Markov鏈,即任意解間依概率可達(dá)。同時,IBPSO鄰代種群序列是單調(diào)的。由上述理論分析,IBPSO滿足依概率1收斂到全局最優(yōu)解的充要條件。

3 仿真分析

以檢驗函數(shù)以及背包問題為基礎(chǔ),分析算法在經(jīng)典離散問題和實際離散問題上的尋優(yōu)效果。參數(shù)設(shè)置如表1所示。限定s=10以及tmax=100。

表1 PSO相關(guān)參數(shù)設(shè)定情況Table1 Experimental parameters setting of PSO

3.1 檢驗函數(shù)的實驗及結(jié)果分析

3.1.1 實驗設(shè)計

本文選擇具有單峰、多峰等典型特征的檢驗函數(shù)[11]構(gòu)成測試環(huán)境,如表2所示。

表2 檢驗函數(shù)列表Table2 List of test functions

將BPSO[5],NBPSO[7]和EBPSO?作為第1組對比算法,分析IBPSO的改進(jìn)效果。NBPSO:采用與IBPSO類似的速度及狀態(tài)調(diào)整,但未采用混沌初始化以及子代處理過程。EBPSO?:采用混沌初始化但未采用子代處理過程。

將SPSO[13]、EnBPSO[14]以及MBPSO[15]作為第2組對比算法,分析IBPSO與其他改進(jìn)的二元粒子的求解差異。SPSO:與經(jīng)典BPSO的狀態(tài)調(diào)整方法相同,僅在速度調(diào)整過程中,速度變量偏向全局最優(yōu)的比例隨迭代次數(shù)的增加而增大。EnBPSO:在SPSO基礎(chǔ)上,引入非線性慣性權(quán)重。MBPSO:使用線性函數(shù)而非sig?moid函數(shù)。

3.1.2 結(jié)果分析

在檢驗函數(shù)f1~f5上,分別取維度n={10,30,50},計算IBPSO以及各算法執(zhí)行100次的實驗結(jié)果,如表3所示。AVG為平均演化精度。ENUM為計算代價,z(n)為適應(yīng)度函數(shù)執(zhí)行一次計算的單位時間。

對于AVG指標(biāo),BPSO無論在單峰檢驗函數(shù)還是多峰檢驗函數(shù)下,都無法收斂到全局最優(yōu)。NBPSO僅在高維檢驗函數(shù)上,陷入局部最優(yōu)解。EBPSO?一定程度的緩解NBPSO陷入局部最優(yōu)解的問題,但仍然無法收斂到全局最優(yōu)解。SPSO在低維度問題中能夠收斂到全局最優(yōu)解,但在高維度問題上難以收斂到全局最優(yōu)解。和SPSO相比,MBPSO對單峰問題具有較好求解能力,但對多峰問題較差。IBPSO無論是單峰還是多峰函數(shù),均可有效收斂。

在衡量指標(biāo)ENUM上,SPSO、EnBPSO和MBPSO與BPSO的算法結(jié)構(gòu)類似,因此,它們的ENUM與BP?SO相同。由于IBPSO和EBPSO?采用的混沌初始化過程,因此IBPSO和EBPSO?的ENUM數(shù)量增加了10個單位的z(n)。雖然增加了較少的計算代價,但I(xiàn)BP?SO和EBPSO?的收斂精度AVG普遍優(yōu)于BPSO和NB?PSO。

表3 各二元PSO算法在檢驗函數(shù)上的實驗結(jié)果Table3 Experimental results of the different binary PSO algorithms on the test functions

以多峰函數(shù)f3為例分析各算法在處理50維問題時的收斂情況,如圖1所示。

圖1 各算法求解50維f3時的收斂情況Fig.1 Convergence of each algorithm on f2with 50 di?mensions

圖1 中,BPSO與IBPSO在收斂速度方面的差異較單峰函數(shù)f1更加明顯,在100代演化結(jié)束時,IBPSO能夠獲取f3的全局最優(yōu)值,但BPSO仍未收斂。SPSO和EnBPSO的收斂情況類似,雖然波動較低,但仍無法避免退化問題。

上述分析表明:和其他算法相比,IBPSO能夠在有限迭代次數(shù)內(nèi),收斂到理想狀態(tài)。

3.2 多選擇背包問題實驗及結(jié)果分析

多維多選擇背包問題(multi?dimensional multi?choice knapsack problem,MMKP)是經(jīng)典0?1背包問題的擴(kuò)展,其約束包括:物品分為M類,指定第m個類中最多有Nm個物品,第i類中第j個物品的重量和價值分別為wij和vij,背包總?cè)萘繛镃,背包中對于一類物品至多有一個,其目標(biāo)是最大化背包內(nèi)組合物品數(shù)量的收益。

實驗數(shù)據(jù)來源自Brunel University提供的測試集[12],選取該數(shù)據(jù)集中13個測試用例,每個測試用例中物品分類數(shù)量從5個到400個不等,且各類數(shù)量中物品數(shù)量相同。記錄運行結(jié)果的平均值A(chǔ)VG以及適應(yīng)度最大值MAX,如表4所示,其中,M為物品類別數(shù),N為物品總個數(shù),R為每類物品數(shù)。

根據(jù)表4,在低維度多選擇背包問題上,貪心算法能夠獲得較優(yōu)的收益平均值,但隨著物品種類數(shù)量以及每類物品數(shù)量的上升(M>15和R>10),GA,PSO,IB?PSO和EnBPSO和MPSO表現(xiàn)出較高的平均收益和最大收益。在智能優(yōu)化算法中,二元粒子群系列改進(jìn)算法(IBPSO,EnBPSO和MBPSO)所獲得的收益高于經(jīng)典優(yōu)化算法(GA和PSO)。

表4 各二元PSO算法在背包問題上的實驗結(jié)果Table4Experimental results of the different binary PSO algorithms on knapsack problem

綜合背包問題的實驗分析結(jié)果,在多選擇背包問題維度上升時,適合采用IBPSO算法進(jìn)行求解。

4 結(jié)束語

針對現(xiàn)有二元PSO求解離散優(yōu)化問題時存在的問題,提出一種改進(jìn)的二元PSO算法(IBPSO)。實驗結(jié)果表明,在經(jīng)典檢驗函數(shù)以及實際背包問題中,IBPSO較相同參數(shù)設(shè)定下其他優(yōu)化算法具有更高的收斂精度以及更有效的收斂過程,體現(xiàn)了對離散問題的求解優(yōu)勢。

在今后的研究中,將IBPSO與其他離散問題相結(jié)合,研究與問題相關(guān)的特定啟發(fā)式過程,指導(dǎo)IBPSO求解相關(guān)領(lǐng)域的離散問題。同時,基于策略融合是提高優(yōu)化算法的有效手段,可將IBPSO與其他離散演化策略相互融合,借鑒其他演化策略的收斂優(yōu)勢,設(shè)計混合策略的離散優(yōu)化算法。

[1]ZHU Hanhong,WANG Yi,WANG Kesheng,et al.Particle swarm optimization(PSO)for the constrained portfolio opti?mization problem[J].Expert Systems with Applications,2011,38(8):10161?10169.

[2]ZHANG Zike,ZHOU Tao,ZHANG Yicheng.Personalized recommendation via integrated diffusion on user?item?tag tri?partite graphs[J].Physica A:Statistical Mechanics and its Applications,2010,389(1):179?186.

[3]WU Changjun,KALYANARAMAN A.An efficient parallel approach for identifying protein families in large?scale met?agenomic data sets[C]//Proceedings of the 2008 ACM/IEEE Conference on Supercomputing[S.l.],2008:1?10.

[4]CAVUSLU M,KARAKUZU C,KARAKAYA F.Neural i?dentification of dynamic systems on FPGA with improved PSO learning[J].Applied Soft Computing,2012,12(9):2707?2718.

[5]KENNEDY J,EBERHART R C.A discrete binary version of the particle swarm algorithm[C]//Proceedings of the Con?ference on Systems,Man,and Cybernetics.Piscataway,USA,1997:4104?4109.

[6]QIN Jin,LI Xin,YIN Yixin.An algorithmic framework of discrete particle swarm optimization[J].Applied Soft Com?puting,2012,12:1125?1130.

[7]MOJTABA A K,TOOSI K N.A novel binary particle swarm op?timization[C]//Proceedings of the 15th Mediterranean Confer?ence on Control and Automation.[S.l.],2007:33?41.

[8]RATHI A,RATHI P,VIJAY R.Optimization of MSA with swift particle swarm optimization[J].International Journal of Computer Allication,2010,12(8):28?33.

[9]GANESH M,KRISHNA R,MANIKANTAN K,et al.Entro?py based binary particle swarm optimization and classification for ear detection original research article[J].Engineering Applications of Artificial Intelligence,2014,27:115?128.

[10]BANSAL J,DEEP K.A modified binary praticle swarm op?timization for knapsack problems[J].Applied Mathematics and Computation,2012,218(22):11042?11069.[11]CHEN Enxiu,LI Jianqin,LIU Xiyu.In search of the essential binary discrete particle swarm[J].Applied Soft Computing,2011,11:3260?3269.

[12]YANG Xiaohua,YANG Zhifeng,YIN Xinan,et al.Chaos gray?coded genetic algorithm and its application for pollu?tion source identifications in convection?diffusion equation[J].Communications in Nonlinear Science and Numerical Simulation,2008,13(8):1678?1688.

[13]AIHARA K.Chaos and its applications[J].Procedia IU?TAM,2012,5:199?203.

[14]De JONG K A.An analysis of the behavior of a class of ge?netic adaptive systems[D].Ann Arbor:University of Michi?gan,1975:196?249.

[15]KHAN S,LI K F,MANNING E G,et al,Solving the knapsack problem for adaptive miltimedia system[J].Stu?dia Informatica:Special Issue on Cutting,Packing and Knapsacking Problems,2001,9:157?178.

An Improved binary particle swarm optimization for discrete optimization problems

YIN Guisheng,CUI Xiaohui,DONG Yuxin,YANG Xue
(College of Computer Science and Technology,Harbin Engineering University,Harbin 150001,China)

Binary particle swarm optimization(BPSO)is wildly used to solve discrete combinational optimization problems.With the low amount of population size and limit iterations,BPSO would have the evolutional problems such as low utilization of solution space,convergence of the speed and status of particles,as well as the degradation and volatility during the iterations.To solve these problems,an improved binary PSO(IBPSO)is designed,which uses the Gray code evolution based coding,chaos initialization process of population,improved modification of the speed and status of particles and the off?spring processing to increase the diversity and utilization of the population.According to the experimental results on the test functions with different types and multiple choice knapsack prob?lems,IBPSO outperforms the existing optimization algorithms and other binary algorithms with higher precision solu?tion and faster convergence speed,which shows the practicality of multiple discrete optimization problems.

binary particle swarm optimization;Gray code;chaos;off?spring processing;discrete optimization

10.3969/j.issn.1006?7043.201306024

http://www.cnki.net/kcms/doi/10.3969/j.issn.1006?7043.201306024.html

TP181

A

1006?7043(2015)02?0191?05

2013?06?05.網(wǎng)絡(luò)出版時間:2014?11?27.基金項目:國家自然科學(xué)基金資助項目(61272186,61100007);黑龍

江省自然科學(xué)基金資助項目(F200937,F(xiàn)201110);黑龍江省博士后基金資助項目(LBH?Z12068);中央高校基本科研業(yè)務(wù)費專項資金資助項目(HEUCF100608).

印桂生(1964?),男,教授,博士生導(dǎo)師;

崔曉暉(1984?),男,博士研究生.

印桂生,E?mail:yinguisheng@hrbeu.edu.cn.

猜你喜歡
優(yōu)化
超限高層建筑結(jié)構(gòu)設(shè)計與優(yōu)化思考
PEMFC流道的多目標(biāo)優(yōu)化
能源工程(2022年1期)2022-03-29 01:06:28
民用建筑防煙排煙設(shè)計優(yōu)化探討
關(guān)于優(yōu)化消防安全告知承諾的一些思考
一道優(yōu)化題的幾何解法
由“形”啟“數(shù)”優(yōu)化運算——以2021年解析幾何高考題為例
圍繞“地、業(yè)、人”優(yōu)化產(chǎn)業(yè)扶貧
事業(yè)單位中固定資產(chǎn)會計處理的優(yōu)化
4K HDR性能大幅度優(yōu)化 JVC DLA-X8 18 BC
幾種常見的負(fù)載均衡算法的優(yōu)化
電子制作(2017年20期)2017-04-26 06:57:45
主站蜘蛛池模板: 色屁屁一区二区三区视频国产| 人与鲁专区| 国产成人啪视频一区二区三区| 在线一级毛片| 国产欧美中文字幕| 丁香婷婷激情网| 97国产在线视频| 2022精品国偷自产免费观看| 色综合久久综合网| 最近最新中文字幕在线第一页| 国产综合日韩另类一区二区| 中美日韩在线网免费毛片视频| 69av免费视频| 国产经典在线观看一区| 丁香六月激情婷婷| 一级爆乳无码av| 谁有在线观看日韩亚洲最新视频| WWW丫丫国产成人精品| 新SSS无码手机在线观看| 中国特黄美女一级视频| 国产福利免费视频| 色视频久久| 色爽网免费视频| 国产永久在线观看| 亚洲热线99精品视频| 青青青草国产| 免费AV在线播放观看18禁强制| 538国产视频| 国产精品久久久久久久久| 国产精品高清国产三级囯产AV| 欧美三級片黃色三級片黃色1| 国产成人精品免费av| 亚洲,国产,日韩,综合一区| 中文字幕人妻无码系列第三区| 欧美在线伊人| 91视频青青草| 色婷婷狠狠干| 啪啪免费视频一区二区| 国产精品xxx| 99re经典视频在线| 97青草最新免费精品视频| 一区二区三区成人| 国产精品亚洲一区二区在线观看| 欧美激情二区三区| 国产在线拍偷自揄观看视频网站| 国产色伊人| 成人毛片免费在线观看| 国产香蕉在线| 国产性猛交XXXX免费看| 亚洲另类第一页| 99在线视频免费| 日本欧美一二三区色视频| 亚洲国产欧美中日韩成人综合视频| 九色在线观看视频| 亚洲第一成年人网站| 精品自窥自偷在线看| 国产成人高清精品免费| 婷婷色丁香综合激情| 在线国产欧美| 婷婷综合色| 成人欧美日韩| 日韩AV无码免费一二三区| 波多野结衣国产精品| 思思热在线视频精品| 在线观看91精品国产剧情免费| 国产呦视频免费视频在线观看| 成人午夜视频免费看欧美| 日韩在线欧美在线| 欧美中文字幕在线二区| 国产97视频在线| 欧美不卡视频在线观看| 91区国产福利在线观看午夜 | 精品亚洲国产成人AV| 亚洲精品视频免费| 91精品网站| 四虎综合网| 久久这里只有精品66| 二级毛片免费观看全程| 天天综合网色| 一级香蕉视频在线观看| 伊人色婷婷| 国产男人的天堂|