亚洲免费av电影一区二区三区,日韩爱爱视频,51精品视频一区二区三区,91视频爱爱,日韩欧美在线播放视频,中文字幕少妇AV,亚洲电影中文字幕,久久久久亚洲av成人网址,久久综合视频网站,国产在线不卡免费播放

        ?

        考慮道路約束的應(yīng)急物資調(diào)度優(yōu)化模型與算法

        2022-05-24 01:01:20王付宇
        關(guān)鍵詞:算子交叉變異

        王付宇,張 康

        (安徽工業(yè)大學(xué) a.管理科學(xué)與工程學(xué)院;b.復(fù)雜系統(tǒng)多學(xué)科管理與控制安徽普通高校重點(diǎn)實(shí)驗(yàn)室,安徽 馬鞍山 243002)

        0 引言

        近年來(lái),災(zāi)害時(shí)常發(fā)生,如2008年汶川大地震、2020年新冠疫情和2020年南方特大洪澇災(zāi)害等,這些災(zāi)害對(duì)人們的生命財(cái)產(chǎn)和安全帶來(lái)了巨大的威脅。因此,在重大災(zāi)害發(fā)生后開(kāi)展及時(shí)、有效的救援工作十分重要,而應(yīng)急物資和車輛調(diào)度對(duì)災(zāi)后救援和恢復(fù)有著重要意義。

        針對(duì)災(zāi)后應(yīng)急物資和車輛調(diào)度問(wèn)題,國(guó)內(nèi)外學(xué)者進(jìn)行了大量的研究。1998年,List等[1]首次在放射性危險(xiǎn)物品運(yùn)輸優(yōu)化模型中引入了應(yīng)急問(wèn)題,為應(yīng)急資源調(diào)度問(wèn)題的研究奠定了基礎(chǔ);王海軍等[2]考慮到不同應(yīng)急響應(yīng)階段和目標(biāo)重要程度不同,通過(guò)決策者對(duì)總運(yùn)輸時(shí)間和應(yīng)急成本的動(dòng)態(tài)賦權(quán),實(shí)現(xiàn)了應(yīng)急調(diào)度的動(dòng)態(tài)決策;段曉紅等[3]針對(duì)城市路網(wǎng)下多事故應(yīng)急救援工作的特點(diǎn),構(gòu)建雙層規(guī)劃模型,對(duì)應(yīng)急車輛調(diào)度和交通疏散進(jìn)行協(xié)同決策;王付宇等[4]等針對(duì)災(zāi)害初期存在道路受損和救援物資不足的情況,以救援的公平性作為調(diào)度目標(biāo),使用步長(zhǎng)遞減的天牛須算法進(jìn)行求解;Yi W等[5]針對(duì)人員疏散問(wèn)題和應(yīng)急物資調(diào)度的運(yùn)輸問(wèn)題,設(shè)計(jì)了集成分布模型;Moreno A等[6]針對(duì)災(zāi)害救援中車輛運(yùn)力的可重復(fù)利用問(wèn)題,設(shè)計(jì)了混合整數(shù)規(guī)劃啟發(fā)式方法。

        上述研究大部分是基于單種配送工具的應(yīng)急資源調(diào)度問(wèn)題,但是在道路發(fā)生損毀的情況下,單種配送工具很難滿足現(xiàn)實(shí)需求。2010年,Hu Z H[7]使用改進(jìn)免疫算法,求解基于救援網(wǎng)絡(luò)的集裝箱多式聯(lián)運(yùn)問(wèn)題;Fikar C等[8]針對(duì)選取陸運(yùn)和空運(yùn)協(xié)調(diào)運(yùn)輸中轉(zhuǎn)點(diǎn)選擇問(wèn)題,設(shè)計(jì)了一種決策支持系統(tǒng),打通了應(yīng)急救援工作的最后一公里;陳雷雷[9]建立了大規(guī)模突發(fā)事件下以受災(zāi)點(diǎn)滿意度為目標(biāo)的多車輛和多物資的物資優(yōu)化調(diào)度模型;阮俊虎等[10]對(duì)應(yīng)急響應(yīng)中的“直升飛機(jī)+車輛”醫(yī)療物資聯(lián)合運(yùn)送問(wèn)題進(jìn)行研究;Ruan J H等[11]為解決大型自然災(zāi)害中的車輛和直升機(jī)聯(lián)運(yùn)問(wèn)題,構(gòu)建了平衡的車輛-直升機(jī)聯(lián)運(yùn)網(wǎng)絡(luò);Erdemir E T等[12]為解決突發(fā)事件下陸-空聯(lián)合救援問(wèn)題,提出了一種位置覆蓋模型和貪婪啟發(fā)式算法。

        以上文獻(xiàn)大多未考慮災(zāi)害初期運(yùn)力不足和道路通行受約束等情況。而在現(xiàn)實(shí)情景下,由于災(zāi)害的突發(fā)性,導(dǎo)致災(zāi)害發(fā)生初期救援工具很難在短時(shí)間內(nèi)集結(jié)。2013年,王旭平等[13]針對(duì)運(yùn)力受限情況下的救援車輛路徑選擇和應(yīng)急物資分配問(wèn)題,構(gòu)建了混合整數(shù)模型,使用遺傳算法求解模型。2020年,薛星群等[14]針對(duì)震后應(yīng)急資源調(diào)度的特點(diǎn),考慮道路通行受約束和運(yùn)力受限等條件,構(gòu)建多目標(biāo)規(guī)劃模型,并使用改進(jìn)的NSGA-II算法求解模型。但是該研究并未考慮物資的裝卸時(shí)間或者準(zhǔn)備時(shí)間,在特殊情況下,救援物資準(zhǔn)備時(shí)間以及裝載時(shí)間對(duì)于整個(gè)救援過(guò)程來(lái)講不可忽視,該問(wèn)題仍有改進(jìn)的空間。因此本文在考慮道路通行受約束和運(yùn)輸能力不足的條件下,將救援物資的裝卸或者救援工具的準(zhǔn)備時(shí)間嵌入模型中,建立了考慮道路約束和多式聯(lián)運(yùn)的應(yīng)急物資調(diào)度模型。

        對(duì)于多目標(biāo)問(wèn)題,求解方法大多使用NSGA-II算法,但NSGA-II算法存在收斂速度慢、種群多樣性差等[15-16]缺點(diǎn)。多年來(lái),學(xué)者對(duì)NSGA-II算法進(jìn)行了大量的研究。張國(guó)富等[17]提出一種NSGA-II與蟻群優(yōu)化的混合智能搜索算法,并應(yīng)用到多受災(zāi)點(diǎn)、多需求點(diǎn)的應(yīng)急物資調(diào)度問(wèn)題中;李燕等[18]利用改進(jìn)的NSGA-II對(duì)交叉口車輛延誤及機(jī)動(dòng)車碳排放兩方面進(jìn)行優(yōu)化;王付宇等[19]提出了多目標(biāo)蜂群算法,并將其應(yīng)用到重大疫情事件下的應(yīng)急資源調(diào)度問(wèn)題中。從現(xiàn)有文獻(xiàn)分析可知,對(duì)NSGA-II算法仍有改進(jìn)的空間,例如,如何根據(jù)研究問(wèn)題的特殊性,設(shè)計(jì)相應(yīng)算子以增加算法搜索空間等。

        從檢索到的文獻(xiàn)可知,很少有使用自適應(yīng)的NSGA-II算法求解多式聯(lián)運(yùn)問(wèn)題。為求解多目標(biāo)應(yīng)急調(diào)度模型,本文使用帶有自適應(yīng)機(jī)制的NSGA-II算法對(duì)模型進(jìn)行求解,所提算法對(duì)NSGA-II算法做出如下改進(jìn):1)設(shè)計(jì)自適應(yīng)機(jī)制,綜合種群進(jìn)化的橫向信息和縱向信息引導(dǎo)種群進(jìn)化;2)針對(duì)聯(lián)合調(diào)度中多車輛配送的特點(diǎn),設(shè)計(jì)隨機(jī)變鄰域搜索算子,對(duì)解空間充分搜索;3)將貪婪思想嵌入鄰域探索過(guò)程。

        1 模型構(gòu)建

        1.1 假設(shè)條件

        本文針對(duì)所建立的應(yīng)急資源調(diào)度模型,作出如下假設(shè):

        1)受災(zāi)點(diǎn)的位置和兩兩受災(zāi)點(diǎn)之間的距離已知;

        2)物資集散中心的物資量充足而且種類單一;

        3)受災(zāi)點(diǎn)的種類以及需求已知,并且各受災(zāi)點(diǎn)僅能被訪問(wèn)一次;

        4)配送工具平均速度不變,運(yùn)行時(shí)長(zhǎng)不受限制。

        1.2 模型構(gòu)建

        物資救援體系使用完全有向圖構(gòu)建G=(V,A),V表示所有節(jié)點(diǎn)的集合,V={0,1,…,n},其中v=0表示集散中心。V′=V/{0}表示所有受災(zāi)點(diǎn)的集合。A={(i,j)|i,j∈V,且i≠j}表示救援物資網(wǎng)絡(luò)的弧集;K=K1∪K2表示所有運(yùn)輸工具的集合,K1表示直升機(jī)集合,K2表示汽車運(yùn)輸工具集合;Vc表示受交通約束的受災(zāi)點(diǎn)的集合。

        因此,救援物資調(diào)度模型為

        (1)

        F2=min(W+U+V)

        (2)

        (3)

        (4)

        (5)

        (6)

        (7)

        (8)

        (9)

        (10)

        2 改進(jìn)的NSGA-II算法原理

        NSGA-II算法是解決多目標(biāo)優(yōu)化問(wèn)題常用的算法之一,但經(jīng)典算法搜索能力較弱并且收斂性較差,為了提高算法的收斂性和搜索范圍,提出了基于種群熵的自適應(yīng)機(jī)制,充分利用進(jìn)化的橫向和縱向信息。種群進(jìn)化的橫向信息指種群層級(jí)的個(gè)體水平,某一階段種群進(jìn)化的橫向信息包含所有個(gè)體目標(biāo)函數(shù)值、Pareto層級(jí)數(shù)和不同Pareto層級(jí)中個(gè)體的數(shù)量。種群進(jìn)化的縱向信息指以種群進(jìn)化的代數(shù)為縱向時(shí)間軸,種群中非支配解集的分布狀況。

        2.1 自適應(yīng)機(jī)制

        2.1.1 自適應(yīng)交叉概率

        遺傳算法中交叉概率和變異概率對(duì)算法收斂性和路經(jīng)求解質(zhì)量起著至關(guān)重要的作用[20],而固定的交叉概率對(duì)種群進(jìn)化很不利,眾多學(xué)者對(duì)交叉概率做了大量研究[21-22]。本文借鑒文獻(xiàn)[23-25]的設(shè)計(jì)思想,設(shè)計(jì)了自適應(yīng)交叉概率公式。

        pc=mean(pci)

        (11)

        (12)

        其中,pc1為最大交叉概率;pc2為最小交叉概率;f′i為待交叉染色體中第i個(gè)目標(biāo)函數(shù)值較大的個(gè)體;pc為個(gè)體的交叉概率,具體表示為所有目標(biāo)函數(shù)適應(yīng)度值的平均值;favgi,fmini,fmaxi為當(dāng)前種群中第i個(gè)目標(biāo)的平均值、最小值和最大值。

        2.1.2 基于種群熵的變異概率

        變異概率的大小一定程度上影響了種群進(jìn)化的速度,因此本文借鑒文獻(xiàn)[26]中種群熵的概念,設(shè)計(jì)了自適應(yīng)的變異概率,變異概率如式(16)所示。

        (13)

        (14)

        (15)

        pm=P0G(i)

        (16)

        當(dāng)種群中所有個(gè)體都分布在某個(gè)子空間中時(shí),由式(15)可知:Gimax最大值為ε。一般來(lái)說(shuō)pmmax最大不超過(guò)0.2,因此本文取pmmax=0.2。因此ε=pmmax/p0max=0.33。不同方差下的變異概率如圖1所示。

        圖1 不同方差下變異概率對(duì)比圖Fig.1 Comparison of variance probabilities under different variances

        從圖1中任一曲線可以看出,變異概率在算法的前期處于較小的值,在算法的后期,變異概率的值穩(wěn)定維持在0.08左右。通過(guò)變異概率對(duì)比圖可以看出,當(dāng)δ=0.1時(shí),變異概率在進(jìn)化中期有著較大的探索空間,因此本文δ=0.1。

        2.2 改進(jìn)的交叉和變異算子

        2.2.1 變異算子

        設(shè)計(jì)合理的變異算子同樣是遺傳算法中的難點(diǎn),本文結(jié)合配送工具資源可重復(fù)利用和多式聯(lián)運(yùn)等特點(diǎn),設(shè)計(jì)了隨機(jī)變異算子。隨機(jī)變異算子的描述如下:

        第1種變異(random_v):隨機(jī)從染色體中選擇一個(gè)未受到通行約束的受災(zāi)點(diǎn),隨機(jī)插入到該染色體中,判斷插入之后染色體路徑長(zhǎng)度是否優(yōu)化,若未優(yōu)化,則重復(fù)選擇插入點(diǎn),直到產(chǎn)生一個(gè)較優(yōu)的解。

        第2種變異(diff_v):從兩個(gè)不同配送工具中選擇兩個(gè)受災(zāi)點(diǎn)交換位置,若受災(zāi)點(diǎn)是受到通行約束的,則在直升機(jī)配送路徑中選擇一個(gè)不受約束的受災(zāi)點(diǎn),然后進(jìn)行交換位置。

        第3種變異(road_v):隨機(jī)選擇一個(gè)配送工具,對(duì)該配送工具的不同配送路徑中的兩個(gè)受災(zāi)點(diǎn)進(jìn)行交換,若交換之后路徑未得到優(yōu)化,則重新選擇,直到路徑得到優(yōu)化為止。

        2.2.2 交叉算子

        本文在順序交叉算子基礎(chǔ)上進(jìn)行了改進(jìn),改進(jìn)的順序交叉算子如下所示。

        2.3 混合選擇策略

        遺傳算法中如何從交配池中選擇合適的父代,以及如何選擇算子是值得關(guān)注的問(wèn)題,本文使用式(17)和式(18)選擇父代個(gè)體。

        (17)

        (18)

        其中,poprank為種群經(jīng)過(guò)排序后第rank層級(jí)群體;T為最大迭代次數(shù),即T=tmax;T1=β*T;T2=(1-β)T,通常β=0.382或者β=0.258;random_v、diff_v和road_v為設(shè)計(jì)的3種變異算子。

        2.4 改進(jìn)算法流程

        改進(jìn)算法的流程如下:

        步驟1:初始化參數(shù),最大迭代次數(shù)T_max、最大最小交叉概率pc1和pc2以及最大貪婪搜索次數(shù)t_max,產(chǎn)生初始種群;

        步驟2:計(jì)算個(gè)體適應(yīng)度值;

        步驟3:非支配排序,并選擇精英個(gè)體,形成父代種群;

        步驟4:根據(jù)式(17)選擇交配個(gè)體,并進(jìn)行交配,形成子代種群;

        步驟5:根據(jù)式(16)計(jì)算變異概率;

        步驟6:通過(guò)式(18)從子代種群中選擇變異個(gè)體,并使用隨機(jī)變異算子形成新個(gè)體;

        步驟7:選出較優(yōu)的個(gè)體,判斷貪婪搜索次數(shù)是否大于t_max,若滿足,則保存較優(yōu)個(gè)體到新子代種群中,并轉(zhuǎn)到步驟8;否則,返回步驟6;

        步驟8:合并父代和子代種群,并判斷進(jìn)化次數(shù)是否大于T_max,若不滿足,返回步驟3;

        步驟9:進(jìn)行非支配排序,輸出非支配解集。

        具體算法流程圖如圖2所示。

        圖2 改進(jìn)非支配快速排序算法流程圖Fig.2 flow chart of the improved non dominated quick sort algorithm

        3 算例分析

        3.1 算例數(shù)據(jù)和算法參數(shù)

        本文采用的數(shù)據(jù)[14]如下,應(yīng)急物資集散中心坐標(biāo)(50,50),20個(gè)受災(zāi)點(diǎn)的坐標(biāo)如表1所示,加黑字體為受約束的受災(zāi)點(diǎn)。隨機(jī)指定3個(gè)受約束的受災(zāi)點(diǎn),兩架直升飛機(jī),3輛貨車,配送工具的相關(guān)參數(shù)如表2所示。由于文獻(xiàn)[14]未考慮物資的裝卸時(shí)間和成本,因此本文對(duì)相關(guān)數(shù)據(jù)進(jìn)行修改,在文獻(xiàn)[14]數(shù)據(jù)的基礎(chǔ)上增加了裝卸時(shí)間和成本。算法的測(cè)試環(huán)境為windows 10操作系統(tǒng),4 G運(yùn)行內(nèi)存,intel四核2.5 Ghz,測(cè)試軟件選擇python 3.7版本。

        表1 受災(zāi)點(diǎn)位置數(shù)據(jù)Tab.1 Disaster location data

        表2 配送工具參數(shù)Tab.2 Distribution tool parameters

        本文種群中染色體的個(gè)數(shù)為100,演化的代數(shù)定為500,最大貪婪搜索次數(shù)設(shè)為100,pc1=0.9,pc2=0.7,pc3=0.5。獨(dú)立運(yùn)行程序15次,任取其中一次的實(shí)驗(yàn)結(jié)果,Pareto前沿如圖3所示,平均等待時(shí)間最小為50.37 min,總成本最小為2 142.67萬(wàn)元。從圖3中的變化趨勢(shì)可以看出,當(dāng)平均等待時(shí)間減少時(shí),總成本將會(huì)增加。并且當(dāng)平均時(shí)間變化很小時(shí),總成本會(huì)大幅度攀升。

        圖3 改進(jìn)NSGA-II算法Pareto前沿Fig.3 Pareto frontier of improved NSGA-II

        3.2 實(shí)驗(yàn)結(jié)果與分析

        針對(duì)決策人的不同偏好和模型的兩個(gè)目標(biāo),時(shí)間偏好型下的時(shí)間最小調(diào)度結(jié)果如表3,經(jīng)濟(jì)成本型下的調(diào)度結(jié)果如表4所示,從表3和表4可以看出,當(dāng)平均等待時(shí)間最小時(shí),總的調(diào)度成本高達(dá)4 472.96萬(wàn)元;當(dāng)總成本最小時(shí),平均等待時(shí)間為145.38 min。

        表3 平均等待時(shí)間最小下的調(diào)度計(jì)劃Tab.3 Schedule with minimum average total waiting time

        表4 總成本最小下的調(diào)度計(jì)劃Tab.4 Schedule with minimum cost

        為檢驗(yàn)所提算法的有效性,通過(guò)與文獻(xiàn)算法和NSGA-II算法進(jìn)行對(duì)比,將3個(gè)算法的種群數(shù)量和迭代次數(shù)均設(shè)為500,NSGA-II和文獻(xiàn)算法的變異和交叉概率分別設(shè)為0.1和0.8。對(duì)比結(jié)果如圖4所示。

        圖4 算法的對(duì)比圖Fig.4 Comparison of different algorithms

        為避免求解的偶然性,進(jìn)行15組重復(fù)測(cè)試,表5中F1的單位為min,F(xiàn)2的單位為萬(wàn)元,由表5的數(shù)據(jù)可知:所提算法得到的平均等待時(shí)間的均值為44 min,比文獻(xiàn)[14]和NSGA-II分別優(yōu)化了6.3%和8.4%;總成本方面,比文獻(xiàn)[14]和NSGA分別優(yōu)化了1.6%和18.21%;穩(wěn)定性方面,本文所提算法都要優(yōu)于另外兩種算法。

        表5 不同算法的對(duì)比結(jié)果Tab.5 Comparison results of different algorithms

        3.3 算法的性能分析

        為測(cè)試改進(jìn)算法相關(guān)性能[27],本文采用分布性指標(biāo)(SP)和覆蓋范圍(MS)來(lái)衡量算法的性能。

        多樣性SP指標(biāo)是衡量種群多樣性的一個(gè)指標(biāo)。種群分布越均勻,則SP指標(biāo)越小。SP指標(biāo)如式(19)所示。

        (19)

        非支配解集覆蓋范圍MS指標(biāo),衡量Pareto前沿的覆蓋范圍,如式(20)所示。

        (20)

        為了說(shuō)明本文算法在多樣性和分布性能的效果,重復(fù)運(yùn)行每個(gè)算法30次,得到MS和SP指標(biāo)的均值和方差,具體數(shù)據(jù)如表6所示。

        表6 SP和MS指標(biāo)比較Tab.6 Comparison of SP and MS indexes

        表6中SP和MS指標(biāo)數(shù)據(jù)表明,所提算法結(jié)果比傳統(tǒng)算法更優(yōu),與文獻(xiàn)[14]算法相比,結(jié)果不比其差。

        圖5a為文獻(xiàn)[14]和所提算法中非被占優(yōu)解集的對(duì)比圖,圖5b為所提算法與文獻(xiàn)[14]中非被占優(yōu)解集真實(shí)占比對(duì)比圖。由圖5a和圖5b分析可知,文獻(xiàn)[14]的種群含有大量的重復(fù)解。從圖5中非被占優(yōu)解集占比收斂點(diǎn)可以看出,所提算法在進(jìn)化后期有著較強(qiáng)的擾動(dòng)能力。證明了所設(shè)計(jì)的交叉算子和消除重復(fù)解策略取得了效果。

        圖5 非被占優(yōu)解變化趨勢(shì)Fig.5 trend of non-dominant solutions

        圖6顯示所提算法得到的MS指標(biāo)值大于文獻(xiàn)[14]數(shù)據(jù),表明在種群進(jìn)化過(guò)程中,改進(jìn)NSGA-II算法所得種群覆蓋范圍更大,搜索范圍更大。種群覆蓋范圍更大是因?yàn)樵O(shè)計(jì)的變異和交叉算子大大增強(qiáng)了算法的搜索能力。圖7表明,與文獻(xiàn)[14]相比,所提算法所得SP值更小、一致性更高。

        圖6 不同算法MS值對(duì)比圖Fig.6 Comparison of MS values of different algorithms

        圖7 SP值對(duì)比圖Fig.7 Comparison of SP values

        4 總結(jié)

        本文主要研究在重大突發(fā)事件時(shí),運(yùn)輸能力受限制和道路通行受到約束的情境下的應(yīng)急調(diào)度問(wèn)題。在求解模型時(shí),本文結(jié)合多式聯(lián)運(yùn)的特征,通過(guò)自適應(yīng)機(jī)制引導(dǎo)種群進(jìn)化。根據(jù)種群進(jìn)化的橫向信息和縱向信息設(shè)計(jì)了自適應(yīng)交叉和變異公式,以加速種群進(jìn)化,提高算法效率。本文結(jié)合配送工具可重復(fù)利用以及多車輛調(diào)度的特點(diǎn),設(shè)計(jì)了多種變異算子以增加算法的鄰域搜索能力。最后通過(guò)算例和相關(guān)指標(biāo)表明,與其他算法相比較,所提算法在改進(jìn)種群多樣性和分布性方面有著較好的結(jié)果。

        在重大突發(fā)事件下,需求的不確定性、受災(zāi)點(diǎn)群眾的服務(wù)滿意度和公平性等問(wèn)題同樣是決策者考慮的重要因素,因此,在突發(fā)事件下,考慮受災(zāi)點(diǎn)的社會(huì)公平性、滿意度以及綠色調(diào)度等可能是未來(lái)的研究重點(diǎn)。算法方面,未來(lái)調(diào)度問(wèn)題往往是大規(guī)模、高維的,因此,使用多目標(biāo)算法與AI的混合算法解決大規(guī)模高維問(wèn)題,將是未來(lái)研究的熱點(diǎn)。

        猜你喜歡
        算子交叉變異
        擬微分算子在Hp(ω)上的有界性
        變異危機(jī)
        變異
        各向異性次Laplace算子和擬p-次Laplace算子的Picone恒等式及其應(yīng)用
        “六法”巧解分式方程
        一類Markov模算子半群與相應(yīng)的算子值Dirichlet型刻畫
        Roper-Suffridge延拓算子與Loewner鏈
        連一連
        變異的蚊子
        基于Fast-ICA的Wigner-Ville分布交叉項(xiàng)消除方法
        国产乱码精品一区二区三区久久| 亚洲狼人社区av在线观看| 精品国产又大又黄又粗av| 男女调情视频在线观看| 岛国av无码免费无禁网站 | 成人无码区免费a片www| 亚洲无AV码一区二区三区| 亚洲一二三四五中文字幕| 亚洲乱码无人区卡1卡2卡3| 亚洲午夜精品久久久久久人妖| 一级一级毛片无码免费视频 | 亚洲熟妇色xxxxx欧美老妇| 极品少妇一区二区三区四区| 中文字幕一区二区中出后入| 久久久国产精品免费a片3d| 国产特级全黄一级毛片不卡| 男女男在线精品免费观看| 亚洲国产美女高潮久久久| 99在线精品免费视频九九视| 人伦片无码中文字幕| 麻豆成人久久精品二区三区91| 久久久久99精品成人片欧美| 真实单亲乱l仑对白视频| 亚洲精品成人国产av| 亚洲一区二区自偷自拍另类| 综合色就爱涩涩涩综合婷婷| 国产午夜无码视频免费网站| 久久伊人网久久伊人网| 国产日韩精品中文字幕| 另类老妇奶性生bbwbbw| 国产亚洲sss在线观看| 精品不卡视频在线网址| 97精品国产一区二区三区| 久久人妻公开中文字幕| 日韩精品极品免费观看| 蜜桃一区二区三区视频网址| 国产精品对白刺激久久久| 亚洲专区一区二区在线观看| 精品人妻一区二区三区狼人| 好吊妞无缓冲视频观看| 波多野结衣一区|