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

基于最優(yōu)插入子集的動態(tài)規(guī)劃算法求解旅行商問題

2023-01-31 08:56:00段隆振
計算機(jī)應(yīng)用與軟件 2022年12期

但 開 段隆振

(南昌大學(xué)信息工程學(xué)院 江西 南昌 330031)

0 引 言

經(jīng)典的旅行商問題(Travelling Salesman Problem,TSP)可以描述為:給定一組城市,這些城市之間的距離已知,尋找一條經(jīng)過所有城市一次且僅一次并返回起始城市的最短路線。旅行商問題有著廣泛的實(shí)際應(yīng)用,生活中經(jīng)常出現(xiàn)的各類組合優(yōu)化問題都可以歸為旅行商這類問題。如物資運(yùn)輸路線中,汽車應(yīng)走怎樣的路線使路程最短[1];大型工業(yè)制件在增材制造中的路徑規(guī)劃[2];CNC在線檢測中,如何盡可能縮短檢測路徑長度[3];在航海相關(guān)問題中,無人駕駛船舶路徑規(guī)劃或船舶航線自動生成[4];多無人機(jī)協(xié)同作業(yè)中,規(guī)劃出符合無人機(jī)機(jī)動性能約束和安全要求的最優(yōu)航跡的問題[5]等。

三角不等式是TSP的一個很自然的限制:對任意三座城市A、B和C,從A到C的路程加上從C到B的路程應(yīng)該大于等于從A到B的路程。旅行成本是對稱的,而且滿足三角不等式的題目稱為TSP的度量(metric)題目。當(dāng)把城市視為平面上的點(diǎn),許多距離計算方式都滿足這個自然條件,例如城市距離為歐幾里得距離(Euclidean distance)、曼哈頓距離(Manhattan distance)、切比雪夫距離(Chebyshev distance)的TSP都屬此類。如果兩個城市之間的距離是相應(yīng)的歐幾里得距離,且兩個城市之間的距離在兩個方向上都是相同的,這樣的TSP稱為對稱的Euclidean TSP。許多自然實(shí)例都是對稱的Euclidean TSP。

旅行商問題的最優(yōu)化求解非常困難,其最根本的原因是求解這些問題的現(xiàn)有算法需要極長的運(yùn)行時間以及極大的存儲空間,以至于根本不可能在現(xiàn)有的計算機(jī)上實(shí)現(xiàn),即存在所謂的“組合爆炸”。例如31個城市的旅行商問題,使用窮舉法需要遍歷30!/2種排列。如果使用目前我國排名第一的超級計算機(jī)——神威太湖之光(Sunway TaihuLight,含有10 649 600個處理器,運(yùn)算速度可達(dá)每秒125 436萬億次),假設(shè)檢驗(yàn)每條路線只需要一次算術(shù)運(yùn)算,那么解決這道31個城市的TSP問題預(yù)計需要超過3 352萬年。

通過Cook[6]和Karp[7]的理論,可以將旅行商問題和許多其他問題聯(lián)系起來。很多的組合優(yōu)化問題,比如背包問題、車間調(diào)度問題、分配問題都和旅行商問題同屬于NP完全問題。如果其中任意一個能用多項(xiàng)式確定性算法解決,那么其他所有NP完全問題也能用多項(xiàng)式算法解決。事實(shí)上,有很多方法本來是研究旅行商問題時提出來的,經(jīng)過不斷的發(fā)展,后來又推廣到了其他NP完全問題上去。

動態(tài)規(guī)劃(Dynamic Programming)是解決多階段決策過程最優(yōu)化的一種數(shù)學(xué)方法。在多階段決策過程中,動態(tài)規(guī)劃方法將問題的過程分成幾個相互聯(lián)系的階段,恰當(dāng)?shù)剡x取狀態(tài)變量和決策變量及定義最優(yōu)值函數(shù),從而把一個大問題化成一組同類型的子問題,然后逐個求解。即從邊界條件開始,逐段遞推尋優(yōu),在每一個子問題的求解中,均利用了它前面的子問題的最優(yōu)化結(jié)果,依次進(jìn)行,最后一個子問題所得的最優(yōu)解,就是整個問題的最優(yōu)解。動態(tài)規(guī)劃的方法在電網(wǎng)調(diào)峰[8]、航天航空[9]和保險金融[10]等領(lǐng)域中都有廣泛的應(yīng)用。

已知求解旅行商問題的精確算法的運(yùn)行時間界限是Held-Karp[11]動態(tài)規(guī)劃法的O(n22n),這個紀(jì)錄已經(jīng)保持了58年。如果能突破這個界限,哪怕只是一點(diǎn)點(diǎn)改進(jìn),也有可能獲得更快的時間界限,就有可能適用于實(shí)際應(yīng)用,從而推進(jìn)旅行商問題的實(shí)戰(zhàn)前線。

1 簡單插入與兩個猜想

1.1 簡單插入

插入算法(Insertion Algorithm)的思路是從一條周游幾個城市的子路線出發(fā),逐個增加新的城市并插入到合適位置,直至得到一條包含所有城市的新路線。算法核心部分的偽代碼為:

Begin

設(shè)算法運(yùn)行過程中M為當(dāng)前子路線;

根據(jù)新城市選取函數(shù)選取待插入城市i;

根據(jù)插入位點(diǎn)選取函數(shù)選取插入位點(diǎn)p;

M.insert(p,i);

if所有城市均在M中:

返回M;

else:

進(jìn)行下一步插入操作;

End

其中插入位點(diǎn)選取的規(guī)則是遍歷所有可插入位點(diǎn),計算插入后的路線長度,其中最短路線長度所對應(yīng)的位點(diǎn)作為插入位點(diǎn)。將這條規(guī)則稱作簡單選取規(guī)則,通過簡單選取規(guī)則進(jìn)行的插入稱為簡單插入。

1.2 一個插入法的實(shí)例

例1十三個城市的TSP,城市從0到12依次編號,城市橫坐標(biāo)依次為{466,826,555,79,781,251,329,336,978,793,489,496,735},城市縱坐標(biāo)依次為{407,389,485,264,618,766,640,670,195,15,714,218,438}。

其最優(yōu)路線為[7,6,2,0,3,11,9,8,1,12,4,10,5],最優(yōu)路線長為3 082個單位長度,如圖1所示。

圖1 例1的最優(yōu)路線

對于旅行商問題,有兩個最基本定理,分別是交叉避免和凸包規(guī)則。

交叉避免:最優(yōu)路線不會自交,因此在解題過程中應(yīng)避免出現(xiàn)交叉的情況。

用一步2-opt可以證明這個事實(shí),若巡回路徑中存在交叉路徑,則打開交叉路徑得到的路徑長度一定小于原路徑長度。

凸包規(guī)則:最優(yōu)路線中經(jīng)過邊界各點(diǎn)的順序和凸包邊界經(jīng)過邊界各點(diǎn)的順序相同。它們在邊界上的順序與在最優(yōu)路線上的順序相同,其余城市的順序則是由邊界向內(nèi)側(cè)拐入短的子路徑而得。

最優(yōu)路線需要避免自交,所以經(jīng)過邊界各點(diǎn)的順序必須和凸包邊界經(jīng)過邊界各點(diǎn)的順序相同。

對于例1,使用凸核插入法,即初始子路線從城市點(diǎn)集的凸包開始,插入位點(diǎn)選取規(guī)則不變,遍歷所有可能的插入順序,可能得到的結(jié)果只有4種,如圖2所示。其中路線長度最小為3 106個單位長度,大于最優(yōu)路線的3 082個單位長度。

圖2 凸核插入法求解例1的所有可能結(jié)果

作為近似算法而言,在插入法中加入能夠把握整體形狀的初始子路線可以在插入城市的過程不斷完善細(xì)節(jié),往往能得到比不加入初始子路線更優(yōu)質(zhì)的路線。但在尋找最優(yōu)解方面,由于插入法在插入過程中不會改變已插入的城市排列,所以加入的初始子路線首先必須保證其排列與其在最優(yōu)路線上的排列一致。即使找到凸包這種特殊的初始子路線滿足排列一致的條件,例1也說明了其有可能無法達(dá)到最優(yōu)解。

1.3 簡單插入最優(yōu)性猜想

例1已經(jīng)指出增加初始子路線的做法會使得最優(yōu)解被排除在算法的解空間之外,那么,不指定初始子路線,用簡單插入的方式,遍歷所有的城市插入順序,是否必定能得到最優(yōu)解呢?對于例1,按照插入順序5→10→4→12→11→1→8→6→2→0→9→3→7,可以得到最優(yōu)解,如圖3所示。

圖3 不指定初始子路線基于簡單插入的例1最優(yōu)解插入過程

針對對稱的Euclidean TSP,考慮如下猜想。

簡單插入最優(yōu)性猜想:對于任意TSP,依據(jù)簡單選取規(guī)則的插入法的所有可能中,必定包括TSP的最優(yōu)解。

1.4 插入過程與插座

定義1對于任意的TSP,從其最優(yōu)路線中選擇任一城市a剔除,剩下的城市排列為城市a的插座。

定理1通過將a簡單插入a的插座中得到的必定是最優(yōu)路線。

證明:用反證法。

若將a簡單插入a的插座中得到的不是最優(yōu)路線,設(shè)a插入的位點(diǎn)為i,最優(yōu)路線中a的插入位點(diǎn)為j,根據(jù)簡單插入的定義,有插入i后的路線長度tour(i)小于插入j后的路線長度tour(j),這與最優(yōu)路徑長度最短的定義矛盾。證畢。

插座本身可以看作原問題規(guī)模減一的TSP,如果插座中的城市排列恰好是這個子問題的最優(yōu)解,這樣的插座稱為同型插座,其他的不是子問題最優(yōu)解的則稱為異型插座。城市數(shù)量為n的TSP,插座總數(shù)為n。

1.5 同型插座猜想

易知,存在異型插座數(shù)量為0情形(n個城市的凸包),那么存在同型插座數(shù)量為0的情形嗎?

針對對稱的Euclidean TSP,考慮如下猜想。

同型插座猜想:對于任意TSP,必然存在同型插座。

根據(jù)定理1,對于任意TSP,其最優(yōu)解可以由其中一個城市通過簡單插入操作插入到與其對應(yīng)的插座中得到。而根據(jù)同型插座猜想,對于任意TSP,必然存在同型插座,同型插座即為子集的最優(yōu)解。結(jié)合定理1以及同型插座猜想,由數(shù)學(xué)歸納法即可推出下面的同型插座猜想的推論。

同型插座猜想的推論:對于任意TSP,依據(jù)簡單選取規(guī)則,必然存在這樣插入過程:插入中每一步的子路線都是當(dāng)前的TSP最優(yōu)解。

這樣的插入過程中可以根據(jù)子路線中不同城市的數(shù)量分為n個階段,每個階段所有子路線的集合稱為該階段的最優(yōu)插入子集。

圖3中的插入過程就是這一推論的一個例子,在每一步的插入過程中,其插入后形成的路線都是當(dāng)前已插入城市組成的子問題的最優(yōu)解。

1.6 同型插座數(shù)量分布規(guī)律的統(tǒng)計分析

在1 000×1 000正方形區(qū)域內(nèi)隨機(jī)取m個點(diǎn)的坐標(biāo)數(shù)據(jù),這些數(shù)據(jù)構(gòu)成了城市規(guī)模為m的旅行商問題,用X表示這個旅行商問題的同型插座數(shù)量,X是一個隨機(jī)變量,它的概率密度記為f(x),其均值為μ,標(biāo)準(zhǔn)差為σ。設(shè)X1,X2,…,Xn表示來自f(x)的一列獨(dú)立的隨機(jī)變量。根據(jù)中心極限定理,當(dāng)n很大時(在實(shí)際的抽樣中,一般當(dāng)樣本容量大于等于30時就可以應(yīng)用中心極限定理),隨機(jī)變量Yn近似地服從正態(tài)分布N(nμ,nσ2)。

已知Yn的概率分布模型,可以利用樣本的密度構(gòu)造似然函數(shù)求出μ和σ2的極大似然估計值。Yn的似然函數(shù)為:

(1)

對式(1)兩邊取對數(shù):

(2)

對式(2)求偏導(dǎo):

解得:

考慮一般情況,假設(shè)城市個數(shù)為n的所有TSP組成一個集合Un,用Cm表示Un中異型插座的數(shù)量為m個的所有元素集合,m=0,1,…,n。如果同型插座猜想成立,那么對于任意Un,可知Cn=?。考慮三個城市的TSP,其所有排列方式都是等價的,換言之,三個城市的TSP只有一條路線,這條路線即為最佳路線。根據(jù)插座理論可知,四個城市的TSP有四個插座,這些插座都是三個城市,又因?yàn)槿齻€城市只有一種可能路線(即最佳路線),所以四個插座皆為同型插座。如果同型插座猜想不成立,則存在城市個數(shù)為k的TSP不存在同型插座。綜上所述,k應(yīng)該大于等于5。

表1 插座型號數(shù)量分布實(shí)驗(yàn)的數(shù)據(jù)

從表1同型插座數(shù)量頻數(shù)分布的實(shí)驗(yàn)數(shù)據(jù)中可以看出,同型插座為1的情形,在城市數(shù)量大于7時開始出現(xiàn),且不同城市數(shù)量的3萬次獨(dú)立實(shí)驗(yàn)中出現(xiàn)的頻數(shù)均為個位數(shù),計算其平均頻率為0.000 103;同型插座為0的情形,在實(shí)驗(yàn)中沒有出現(xiàn)過一次。這樣的實(shí)驗(yàn)結(jié)果有兩種可能:(1) 同型插座猜想成立,任意TSP至少存在一個同型插座;(2) 存在同型插座猜想的反例,即存在同型插座數(shù)量為0的TSP,但反例的出現(xiàn)是極小概率的事件,以至于在本次實(shí)驗(yàn)中未出現(xiàn)。

假設(shè)同型插座猜想存在反例。根據(jù)伯努利大數(shù)定律,當(dāng)n足夠大時,事件A出現(xiàn)的頻率將接近于其發(fā)生的概率。比如同型插座數(shù)量為1的情形,用頻率估計其概率,概率約為萬分之一。同型插座數(shù)量為0的情形,頻率為0,根據(jù)假設(shè)其概率不為0,那么,其概率必然遠(yuǎn)小于萬分之一。

從圖4中可以看出,同型插座的分布形態(tài)隨著城市個數(shù)的增加,越來越接近兩邊少、中間多的鐘形分布概率模型。在典型的鐘形分布即正態(tài)分布中,如果a小于標(biāo)準(zhǔn)差μ,那么P{x

圖4 同型插座的頻數(shù)分布

如果考慮同型插座數(shù)量小于4的情形,從圖5中可以看到隨著城市個數(shù)的增加,其出現(xiàn)的次數(shù)是遞減的,這意味著隨著問題規(guī)模增加,同型插座數(shù)量小于4出現(xiàn)的概率在下降。

圖5 同型插座數(shù)量小于4的頻數(shù)隨城市規(guī)模的變化

綜上所述,同型插座猜想的反例即使存在,也是小概率事件,而且隨著城市規(guī)模越大,這個概率呈減小的趨勢。

2 算法設(shè)計與實(shí)現(xiàn)

2.1 基于簡單插入子集的動態(tài)規(guī)劃法

考慮n座城市的TSP題目,城市從1到n依次編號命名。依據(jù)簡單插入子集中的城市數(shù)量可以劃分為不同的階段,用k表示,k分別等于1,2,…,n。

以城市集{1,2,3,4,5}為例,城市數(shù)量為3的子集有{1,2,3}、{1,2,4}、{1,2,5}、{1,3,4}、{1,3,5}、{1,4,5}、{2,3,4}、{2,3,5}、{2,4,5}、{3,4,5}共10個,因?yàn)槿齻€城市的路線實(shí)際只有一條,這里每一個城市子集對應(yīng)一個子路線,所以第三階段有10個狀態(tài)。

其城市數(shù)量為4的子集有{1,2,3,4}、{1,2,3,5}、{1,2,4,5}、{1,3,4,5}、{2,3,4,5}共5個。其中{1,2,3,4}可由{1,2,3}插入4、{1,2,4}插入3、{1,3,4}插入2、{2,3,4}插入1這四種情況得到相應(yīng)的簡單插入子路線,這四條路線有可能重復(fù),記錄下不同的成為第四階段的一個狀態(tài)。一般第四階段以后,每個城市子集都有可能對應(yīng)多個子路線,至少為一個,且不同題目對應(yīng)不同的數(shù)量。

第五階段的城市子集為{1,2,3,4,5},可由{1,2,3,4}插入5、{1,2,3,5}插入4、{1,2,4,5}插入3、{1,3,4,5}插入2、{2,3,4,5}插入1得到。以{1,2,3,4}插入5的情況為例,在第四階段找到城市子集{1,2,3,4}對應(yīng)的簡單插入子路線,假如有a條子路線,那么就將城市5逐個插入其中,記錄下不重復(fù)的新的子路線。

在最后一個階段的狀態(tài)轉(zhuǎn)移完成后,會得到{1,2,3,4,5}所對應(yīng)的多條簡單插入子路線,根據(jù)簡單插入最優(yōu)性猜想,最優(yōu)解必在其中,此時只要計算這些路線的長度,最短的那條即為該旅行商問題的最短旅行路線。

2.2 基于最優(yōu)插入子集的新動態(tài)規(guī)劃法

依據(jù)同型插座猜想的推論,TSP的最優(yōu)解可以通過同型插座進(jìn)行簡單插入得到,城市個數(shù)為n的TSP有n個插座,這些插座可能為同型,也可能為異型。同型插座猜想指出這n個插座中必然存在一個同型插座。n個插座對應(yīng)n個TSP的子問題,這些城市個數(shù)為n的TSP都有最優(yōu)解,它們與插座的關(guān)系是,同型插座等價于子問題的最優(yōu)解,異型插座不是子問題的最優(yōu)解。

假設(shè)我們知道所有城市個數(shù)為n的子問題的最優(yōu)解,那么這些最優(yōu)解中必然至少有一個會是原問題的同型插座,對所有這些最優(yōu)解進(jìn)行簡單插入,假設(shè)a為插入城市,那么可能的結(jié)果只有兩種,(1) (城市a)+(a的同型插座)=最優(yōu)解;(2) (城市a)+(非插座)=非最優(yōu)解,由前文得知第一種結(jié)果必定存在,因此可以遍歷所有結(jié)果,其中路線長度最小的即為最優(yōu)解。同時,城市個數(shù)為n-1的子問題的最優(yōu)解可以如法炮制,從城市個數(shù)為n-2的子問題中得到,即滿足無后效性(馬爾可夫性)。因此,可以從城市個數(shù)為3(3個城市無須排列即為最優(yōu)路線)開始,逐步求得城市個數(shù)為4,5,…,n的最優(yōu)解。

根據(jù)同型插座猜想的推論,必然存在這樣插入過程:插入中每一步的子集都是當(dāng)前的TSP最優(yōu)解,這樣的子集稱為最優(yōu)插入子集。這是簡單插入最優(yōu)性猜想的一個特例。與基于簡單插入子集的動態(tài)規(guī)劃法相比,新動態(tài)規(guī)劃法只需關(guān)注簡單插入子集中最優(yōu)的部分,有效減少了每個階段的狀態(tài)數(shù)量,因此,算法的效率得到了提升。

城市個數(shù)為n的TSP,共有2n個子集(包括它自己和空集,事實(shí)上算法是從城市數(shù)目為3的子集開始的,為方便計算,此處增加到最大值2n),在每一個子集中,至多有n個插座,對每一個插座,主要操作為復(fù)雜度為n的簡單插入。2n乘以n再乘以n,所以新動態(tài)規(guī)劃法的時間復(fù)雜度為O(n22n)。

由于動態(tài)規(guī)劃法本身存在的“維數(shù)障礙”,在電子計算機(jī)上計算時,每遞推一段,都必須把前一段算出的最優(yōu)值函數(shù)在相應(yīng)的狀態(tài)集合上的全部值存入內(nèi)存中,因而當(dāng)維數(shù)增大時,所需的內(nèi)存量會呈指數(shù)倍增長。在新動態(tài)規(guī)劃法中,每一個城市子集都對應(yīng)一個狀態(tài),一共有2n個城市子集,每個城市子集最多有n個城市,因此空間復(fù)雜度為O(n2n)。

2.3 算法實(shí)現(xiàn)

算法運(yùn)行的硬件環(huán)境為Intel Core i5- 3320M CPU,2.6 GHz,4 GB內(nèi)存,開發(fā)環(huán)境為Python 3.7。實(shí)驗(yàn)采用TSPLIB上的不同數(shù)據(jù)集,驗(yàn)證算法在各種潛在樣本上的有效性。實(shí)驗(yàn)采用了eg7146、kz9976、mo14185、xqg237和xql662的實(shí)驗(yàn)數(shù)據(jù),每次實(shí)驗(yàn)隨機(jī)選擇數(shù)據(jù)集中的12個城市,用Held-Karp動態(tài)規(guī)劃法(HK),基于簡單插入子集的動態(tài)規(guī)劃法(簡插法),基于最優(yōu)插入子集的動態(tài)規(guī)劃法(優(yōu)插法)計算得到結(jié)果,表2是一次實(shí)驗(yàn)的數(shù)據(jù)。

表2 在不同數(shù)據(jù)集上的一次實(shí)驗(yàn)結(jié)果

這樣的實(shí)驗(yàn)重復(fù)進(jìn)行一萬次,實(shí)驗(yàn)結(jié)果解異同均為TRUE,即三種算法的解都是最優(yōu)解。這初步表明在不同數(shù)據(jù)集中,基于最優(yōu)插入子集的動態(tài)規(guī)劃法能有效求得最優(yōu)解。

3 結(jié) 語

基于簡單插入子集的動態(tài)規(guī)劃法的理論依據(jù)是簡單插入最優(yōu)性猜想。這一猜想理論上的證明或許能帶來對TSP結(jié)構(gòu)上的深入認(rèn)識。

基于最優(yōu)插入子集的動態(tài)規(guī)劃法的理論依據(jù)是同型插座猜想。對于同型插座猜想,實(shí)驗(yàn)測試的統(tǒng)計分析表明它有可能存在反例,但是即使存在反例,其在實(shí)際情況中出現(xiàn)的概率極小。而且城市規(guī)模越大,反例出現(xiàn)的概率越接近0。

如果無數(shù)多的猴子在無數(shù)多的打字機(jī)上隨機(jī)地打字,并持續(xù)無限久的時間,那么在某個時候,它們必然會打出莎士比亞的全部著作。這個設(shè)想本身在現(xiàn)實(shí)生活中是不可能重現(xiàn)的,就像一個房間內(nèi)的空氣總是均勻分布的。同型插座猜想可能也有類似的性質(zhì),理論上存在反例,但實(shí)際遇到的TSP都滿足同型插座猜想。算法實(shí)現(xiàn)中采用了TSPLIB上的不同數(shù)據(jù)集,基于最優(yōu)插入子集的動態(tài)規(guī)劃法均能找到最優(yōu)解,實(shí)驗(yàn)結(jié)果在一定程度上支持了這種觀點(diǎn)。

新動態(tài)規(guī)劃算法為求解旅行商問題提供了新的思路,對于認(rèn)識、分析眾多復(fù)雜問題具有一定的啟示,對促進(jìn)組合優(yōu)化問題的研究有著積極的作用。如果將其與其他啟發(fā)式算法融合或許會催生出性能優(yōu)異的新算法,這將給實(shí)際解決大規(guī)模TSP問題提供新的方案。

主站蜘蛛池模板: 欧美一级视频免费| 欧美有码在线观看| 尤物精品视频一区二区三区| 98精品全国免费观看视频| www亚洲精品| 露脸国产精品自产在线播| 日韩精品一区二区三区swag| 欧美一级一级做性视频| 中文字幕在线一区二区在线| 国产91久久久久久| 在线欧美一区| 99r在线精品视频在线播放 | 国产美女叼嘿视频免费看| 日韩美女福利视频| 精品国产美女福到在线不卡f| 国产无码高清视频不卡| 这里只有精品在线| 欲色天天综合网| 国产精品短篇二区| 国产自在线拍| 欧美成人精品一级在线观看| 好紧好深好大乳无码中文字幕| 综合色区亚洲熟妇在线| 国产嫖妓91东北老熟女久久一| 无码人妻热线精品视频| 亚洲无码精彩视频在线观看| 久久综合伊人 六十路| 亚洲无线一二三四区男男| 91精品国产自产在线观看| 久久性妇女精品免费| 午夜无码一区二区三区在线app| 特级欧美视频aaaaaa| 成人一区专区在线观看| 福利在线不卡| 国产一区二区福利| 国产极品美女在线| 91美女在线| 久久精品亚洲中文字幕乱码| 午夜精品久久久久久久99热下载| 欧美综合激情| 国产91高跟丝袜| 国产九九精品视频| 国产打屁股免费区网站| 999国内精品久久免费视频| 在线观看亚洲人成网站| 精品亚洲欧美中文字幕在线看 | 日韩最新中文字幕| 亚洲视频一区| 久久综合九九亚洲一区| 伊人久久综在合线亚洲2019| 久久精品国产免费观看频道| 四虎亚洲精品| 久久久久亚洲精品成人网| 无码中字出轨中文人妻中文中| 五月天香蕉视频国产亚| 亚洲国产成人久久77| 青青久久91| 久久国产高潮流白浆免费观看| 成人午夜久久| 制服丝袜无码每日更新| 欧美精品不卡| 国产又色又爽又黄| 久久熟女AV| 久久久久久午夜精品| 国产精品99在线观看| 亚洲国产精品久久久久秋霞影院| 亚洲人成在线免费观看| 无码AV日韩一二三区| 2022国产91精品久久久久久| 欧美一区中文字幕| 婷婷色婷婷| 高清色本在线www| 99热最新网址| 全裸无码专区| 国产白浆一区二区三区视频在线| 国产喷水视频| 男人天堂伊人网| 免费一极毛片| 一本一道波多野结衣av黑人在线| 97超级碰碰碰碰精品| 午夜激情婷婷| 波多野结衣久久高清免费|