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

動態圖上基于2-HOP COVER的TOP-K最短路徑算法

2019-04-15 07:44:38
計算機應用與軟件 2019年4期

施 琴 兒

(復旦大學計算機科學技術學院復旦-眾安區塊鏈與信息安全聯合實驗室 上海 200433)(上海市區塊鏈工程技術中心 上海 200433)

0 引 言

圖論中最短路徑問題是一個經典的問題。在一個有向帶權圖中,最短路徑查詢問題是查找圖中兩個節點之間的最短路徑,最短距離問題是求出兩個節點之間的最短的路徑長度。top-k最短路徑查詢是對前k條最短路徑的查詢。近年來top-k最短路徑問題在實際中有諸多的應用,如多目標跟蹤[1]、基因網絡[2]和多序列分析[3]等等。

很多現實中的網絡是在隨著時間而改變的。對于網絡的實時變化,這些針對靜態的圖的top-k算法需要對每次變化的圖進行重新的計算。在實際生活中的網絡點之間的關系很復雜,網絡的數據也在日益增多。對于大規模的圖,每一次的變化都會引起靜態算法的重新計算。對于微小的變化而言,絕大部分計算是重復的,在實際運用中效率低下。同時,我們也希望有方法分析圖的改變趨勢。例如,在社交網絡圖中,一個用戶在一個時刻關注某件事的幾個焦點,下個時刻關注的焦點又可能發生了改變。通過記錄top-k距離的改變,可以分析用戶行為。再例如,在道路網中,通過追蹤兩個城市之間道路的改變,可以分析道路建造的趨勢。

傳統的最短路徑算法不能滿足大規模圖的top-k最短路徑的計算,需要預處理原始圖的數據關系,然后在圖進行動態更新后,通過更新原始圖的部分數據,動態圖的更新算法能快速地得到新的top-k最短路徑。

1 相關研究

1959年Hoffman和Pavley[4]提出了KSP(k shortest path)問題后,一些研究人員提出了解決該問題的算法。

KSP問題是基于最短路徑問題提出的,一部分解決該問題的算法也是基于最短路徑算法的。經典的最短路徑算法有Dijkstra算法[5]、Floyd-Warshall[6]算法等。其中,Dijkstra算法對于每次查詢給定的節點對都需要重新計算最短路徑;Floyd算法對于節點對的查詢耗時很少,但需要預先計算路徑表。這些算法用于大規模圖時,其時間和空間消耗不能滿足實際需要。為了提高大規模圖最短路徑計算的性能,一些研究人員提出了兩步驟算法,這類算法在求解最短路徑問題時分為了兩個步驟:預處理和查詢。目前該類算法中的主流是基于樹分解或2-hop cover[7]的。Wei在2010年[8]提出了一種高效的樹分解方法:TEDI。在此之后Akiba等[9]基于TEDI做了一些優化。Xu[10]在2016年提出的BBQ優化了預處理時間以及提出了批量查詢算法。Akiba等[11]在2013年提出了一個剪枝的2-hop cover最短路徑算法:PLL。為了彌補PLL局限于靜態圖的不足,Akiba等[12]又提出了基于PLL的動態最短路徑算法。

近年來主要的KSP算法分為兩類:一類是不可重復路徑的top-k最短路徑算法,以Yen[13]提出的Yen′s算法為代表;另一類是可重復路徑的top-k最短路徑算法,以Eppstein[14]提出的算法為代表,該算法允許環的存在。Akiba在2015年[15]根據2-hop cover的思想提出了在大規模圖中高效求解top-k的算法,在預處理得到索引集后,對于每一次查詢能利用索引集快速得到top-k最短路徑。

這些對KSP問題的研究,針對的都是靜態圖上KSP的計算。如果圖進行了少許的更新改變,對這些算法而言就相當于是另一幅圖了。在實際中,關系網等網絡都是在實時更新的。對于靜態圖中的KSP算法,這些路徑計算的算法都需要重新再執行一遍,極大地影響了效率。本文借鑒了兩步驟算法中的2-hop cover的思想,對圖進行預處理,建立索引集。如果圖進行了更新,只需要更改原始圖中預處理得到的部分索引集,就能得到更新后圖的索引集。

2 預處理算法分析

圖論中圖分為兩種:有向圖和無向圖。由于在大部分的社交網絡關系圖中,兩個人不一定是相互關注的,并且在好友中親密值也會不同,對應的是有向帶權圖,因此本文中主要研究的是有向帶權圖上的更新算法。動態top-k算法分為兩個步驟:預處理步驟和查詢步驟。然后根據預處理得到的索引集進行更新算法的操作。

2.1 2-hop cover

本文用的是基于2-hop cover的框架,具體的定義(參照文獻[7])如下:

定義1在一個有向圖G=(V,E)中,其對應距離標記集為LG,LG={L(v)},其中v∈V。距離標記L(v)=(Lin(v),Lout(v)),其中:

(1) 對于所有的v∈V,Lin(v)是所有(u,δ(u,v))的集合,δ(u,v)=d(u,v);

(2) 對于所有的v∈V,Lout(v)是所有(u,δ(v,u))的集合,δ(v,u)=d(v,u)。

任意兩點x和y之間的距離定義為:

δ(x,y)=min{δxv+δvy|(v,δxv)∈Lout(x),
(v,δvy)∈Lin(y)}

(1)

如果對于圖G中任意節點x和y,通過上述的距離標記集LG得到的距離等于它們在圖中的實際距離,那么LG就是圖G的2-hop cover。由上述的定義可知,計算任意兩個節點最短距離的時間為O(|Lout(x)+Lin(y)|)。

2.2 索引算法

動態top-k算法需要對原始圖進行預處理計算,對于給定的一個圖G=(V,E),得到索引集IN=(L,Lr,C)。記Gr是圖G的反轉圖,Gr=(V,Er)。其中V是節點的集合,E是邊的集合,L是圖的距離標記表集合,Lr是反轉圖的距離標記表集合,C是自循環標記表集合。

索引算法是對圖的預處理過程。在算法中,需要分別計算自循環標記表和距離標記表。默認節點編號按照度的大小給出(為了方便后面能剪枝更多的節點)。具體過程見算法1。

算法1索引算法

Input: G=(V,E),k

Output: index of G

1 procedure CreateIndex()

2 Gr=Reverse(G)

3 for i=1 to n do

4 ComputeCircle(G,vi,k)

5 end for

6 L(v),Lr(v)=empty for all v∈V

7 for i=1 to n do

8 PrunedBFS(G,Gr,L,Lr,C,vi,k)

9 end for

10 return (L,Lr,C)

算法1包含兩個子程序:ComputeCircle和PrunedBFS,分別計算自循環標記表和距離標記表。

計算自循環標記表的基本思路為:對于每個節點v,用BFS來計算它到自身的距離;在每個節點的BFS過程中,剪枝掉已經訪問過的節點;每個節點能被訪問最多k次。這樣得到一個自循環標記表。具體過程見算法2。

算法2自循環標記表算法

Input: G,v,k

Output: C(v)

1 procedure ComputeCircle(G,v,k)

2 a priority queue Q←(v,0), C(v)←empty

3 while Q is not empty do

4 pop up (u,δ) from Q

5 if u=v then C(v)=C(v)∪{δ}

6 if |C(v)|≥k and δ≥C(v)kreturn C(v)

7 for all (u,w)∈E and w≥v do

8 push (w,δ+e(u,w)) in Q

9 end for

10 end while

11 return C(v)

計算距離標記表的基本思路為:對于每個節點v,同樣用剪枝的BFS來計算;在訪問節點v時,如果該節點已經被訪問k次,或者δvu≥mink(L,C,v,u),將節點u剪枝,之后的搜索將不會再訪問;否則,即找到一條v到u的新路徑,該路徑的長度小于現有索引集得到第k條路徑的長度,這時將二元組(u,δvu)放入L(u)中。具體過程見算法3。

算法3距離標記表算法

Input: L,Lr,C,v,k

Output: L,Lr

1 procedure PrunedBFS(G,Gr,L,Lr,C,v,k)

2 a priority queue Q ←(v,0)

3 while Q is not empty do

4 pop up (u,δ) from Q

5 if maxk(querydis(IN,v,u))>δthen

6 add (v,δ) to L(u)

7 for all (u,w)∈E and w≥v do

8 push (w,δ+e(u,w)) in Q

9 end for

10 end if

11 end while

12 same method for Lr

13 return (L,Lr)

2.3 查詢算法

通過索引算法得到的索引集IN=(L,Lr,C),計算任意兩點u和v的距離公式為:

dis(u,v)={δuw+δww+δwv}

(2)

式中:(w,δuw)∈Lr(u),δww∈C(w),(w,δwv)∈L(v)。

由式(2)可知,在計算兩個節點的最短top-k距離時,需要查找到Lr(u)和L(v)中的都可到達的候選節點,即wi∈Lr(u)∩L(v)。先通過(Lr(u),C)得到到達(w1,w2,…,wt)的k個最短距離δu,wi={δ1,δ2,…,δk},再通過L(v)中的二元組 (wi,δwi,u)得到節點u到v的top-k最短距離。

引理1對于每個點對(x,y),通過IN=(L,Lr,C),其k個最短路徑dis(x,y)為:

dis(x,y)={d1(x,y),d2(x,y),…,dk(x,y)}

(3)

式中:di(x,y)由式(2)得到。

在算法3中可以得到節點x到所有節點v的最短的t1個距離,以及節點v到所有節點y的最短的t2個距離,其中(v,δxv)∈Lr(x),(v,δvy)∈L(y), 1≤t1,t2≤k。這些距離所經過的路徑點v滿足vmax(x,y),即w是比x和y訪問更晚的節點,可以得到新的路徑|P′|<|P|且P′=(x,…,v,…,w,…,v,…,y)。由于每個節點的自循環標記表中的路徑可以訪問所有訪問順序靠后的節點,因此在計算距離公式時需要加上C(v)來得到top-k最短距離。

2.4 復雜度分析

在索引算法中,算法2中每個節點都被訪問一遍,在碰到已經訪問過的節點時, 直接進行剪枝,剩余的每個點最多被訪問k次。因此,合計下來,每個節點最多被訪問k次,因此其時間復雜度為O(nk)。算法3中,假定距離標記表的平均長度為len,那么我們總共需要訪問n·len個節點,平均需要訪問O(n/m)條邊,判斷每條邊是否加入到隊列中所需的查詢需要O(len·klogk)的時間,因此總的時間復雜度為O(klogk(m·len+n·len2))。而在實際世界的圖中,經過剪枝后得到的索引集的長度len=n,并且需要的k也是一個相對比較小的值。至于空間復雜度,自循環標記表需要O(nk)的空間,而距離標記表需要O(nk·len)的空間,因此總的空間復雜度為O(nk·len)。

在查詢算法中,每個點的距離標記表的平均長度為len,每個點的距離存儲的數量為k,因此查詢算法的時間復雜度為O(len·klogk)。

3 動態更新算法分析

在實際世界里,數據及其關系是不斷增加的,在圖中更新基本的操作是點或邊的增加。在很多實際動態改變網絡中,移除點或邊的操作是比較少的。在文獻[16]中幾種動態更新的圖中,可以看到增加的操作要遠比刪除操作頻繁,因此本文更新算法只是針對點或者邊增加的更新。

在索引集中增加一個點的記錄很簡單:如果在圖中新增加了一個新的節點a,只需要在索引集中加入L(a)=(a,0)、Lr(a)=(a,0)、C(a)=(0,∞,…,∞)。因此更新算法著重于對邊的更新。

引理2假設G和G′為兩個圖,其中E(G)?E(G′),那么對于任意兩個節點x和y,x,y∈V(G)∩V(G′),有d(x,y)≥d′(x,y)。

由引理2中可知,如果一個圖包含另一個圖的所有邊,在這個圖中的兩點間的距離要小于等于在另一個圖中的距離。對于圖中的任意兩個節點,若它們之間有一條通過新邊更短的路徑,則這兩點間距離減小;反之,則這兩點距離不變。

3.1 索引集更新

假設增加的一條邊為(a,b),權重為eab,并且a和b兩個節點已經在圖中。我們只需要對新加入的邊影響的節點進行部分索引集的更新,即增加部分數據或更改部分數據。索引更新算法的偽代碼如下:

算法4更新索引算法

Input: index of G,(a,b)

Output: index of G′

1 procedure insert_G,(a,b)

2 UpdateCircle(L,Lr,C,a,b)

3 for all vi∈L(a) from lower vido

4 D←(d1(vi,a)+eab,d2(vi,a)+eab, …,dk(vi,a)+eab)

5 UpdateIndex(G,L,Lr,vi,b,D)

6 end for

7 return (G,L,C)

8 procedure insert_Gr,(b,a)

9 UpdateCircle(L,Lr,C,a,b)

10 for all vi∈Lr(b) from lower vido

11 D←(d1(vi,b)+eab,d2(vi,b)+eab,…,dk(vi,b)+eab)

12 UpdateIndex(Gr,Lr,L,vi,a,D)

13 end for

14 return (Gr,Lr,C)

由上述的代碼可知,索引集更新中包括自循環標記表集和距離標記表集更新,具體過程在算法5和算法6中給出。

圖1和圖2是一個示例。圖1標出了從最左邊節點出發的到所有節點的top-2最短距離;圖2標出了在加入了一條邊之后,從最左邊節點出發得到的所有節點的 top-2 最短距離。可以看到只有b和在b之后的節點會有改變,其余節點不變。

圖1 一個示例:圖中距離是與最左邊節點top-2距離

圖2 一個示例:圖增加一條邊后與最左邊節點top-2距離

引理3如果在圖更新后,節點vi到u的距離改變了,那么所有更新了的top-k最短路徑一定經過了邊(a,b)。

引理4假設一條從節點vi到u的top-k最短路徑P改變了,其中u≠a,b,那么路徑中經過了節點b之后的節點w到vi的top-k距離都改變了。

由引理4可知,如果兩個節點u到v的top-k最短距離改變了,那么改變的只是經過b之后的節點,在路徑到達a之前所經過的節點的距離不變。因此在更新中,不需要從(vi,0)開始修改,只需要從(b,d(vi,a)+g)開始。由于在索引集中得到的d(vi,a)有k個,因此初始的集為(b,d1(vi,a)+g),(b,d2(vi,a)+g),…,(b,dk(vi,a)+g)。

在圖G中增加的邊為(a,b),相應地在圖Gr中增加的邊為(b,a),因此圖G和Gr的距離標記表都需要更新。

3.2 自循環標記表更新

與算法1相對應,需要有自循環標記表和距離標記表更新的兩個過程。由引理3和引理4可知,不是所有節點的自循環標記表需要更新,只有在自循環經過的路徑中包含了邊(a,b)后才會改變該節點的自循環標記表。更新自循環標記表算法偽代碼如下:

算法5更新自循環標記表算法

Input: L,Lr,C,a,b

Output: C′

1 for all vi∈L(a)∩Lr(b) from lower viand vi≠a,b do

2 for all(vi,δvi,a)∈L(a) and (vi,δvi,b)∈Lr(b) do

3 if δvi,a+δvi,b+eab

4 α=query(L,C,vi,a)

5 β=query(Lr,C,b,vi)

6 δ=α+β+eab

7 for j=1 to k do

8 if δj

9 end for

10 end if

11 end for

12 end for

13 return C′

如果vi∈L(a)∩Lr(b),即vi→a和b→vi可達,加入邊(a,b)后,更新了vi的自循環標記表,可形成vi→a→b→vi。 因此需要更新的節點最多為max(|L(a)|,|Lr(b)|)個,其余的節點不需要更新。

3.3 距離標記表更新

如3.2所述,自循環標記表的更新只需要更新L(a)∩Lr(b)個節點。相應地,對于距離標記表的更新需要更新L(a)∪Lr(b)中的節點。如果一個節點v不屬于L(a)∪Lr(b),在計算其中一個PrunedBFS時,該節點已經被剪枝,所以增加的邊對被剪枝刪除的節點沒有影響。算法偽代碼如下:

算法6更新距離標記表算法

Input: index of G,L,Lr,vt,a,D

Output: index of L′

1 Q ← empty

2 for i=1 to |D| do

3 push(a,Di) in Q

4 end for

5 while Q is not empty do

6 pop up(u,δ) from Q

7 if TopKQuery(IN,vt,u,t)>δthen

8 add(v,δ) to L(u)

9 for all(u,w)∈E′and w≥vtdo

10 push(w,δ+e(u,w)) in Q

11 end for

12 end while

13 return L′

在算法中的TopKQuery(IN,x,y,t)定義如下:

TopKQuery(IN,x,y,t)={δxvu+δvivi+δviy}

(4)

其中:i≤t,(vi,δxvi)∈Lr(x),δvivi∈C(vi),(vi,δviy)∈L(y)。

它計算的是節點x和y在索引集IN下的第k個距離,如果節點x和y之間沒有k條路徑,那么定義TopKQuery(IN,x,y,t)=∞。

3.4 正確性證明

在經過更新算法后得到新的索引集IN′,要證明此為圖G′的索引集。對于任意兩個節點x和y,以及t(0≤t≤n),我們定義:

d1(x,y,t)=mink(d(x,vi)+d(vi,vi),d(vi,y))

如果t=0或者沒有一個vi滿足i≤t,那么d1(x,y,t)=∞。

引理5對于任意兩個節點x和y,以及i(0≤i≤n),IN是圖G的2-hop cover索引集,我們有TopKQuery(IN,x,y,i-1)=d1(x,y,i-1)。對于任意一個節點u,如果節點vi滿足d′(vi,u)

證明:由條件可知d′(vi,u)

(vi,…,z,…,z,…,a,b=w0,w1…,ws-1,u=ws)

由此可知對于所有的wi,都有d′(vi,wj)

由條件可知d′(vi,u)

同理d′(vi,wj)

因此在進行距離標記表更新算法中,對于所有的節點wj,我們有:

TopKQuery(IN,vi,wj,i)≥

得到TopKQuery(IN,vi,wj,i)≥d′(vi,wj),所以在算法6中,節點wj沒有被剪枝,該新的距離加入了L′(u)中。因此在新的索引集IN′中,一定有(vi,d(vi,u))∈L′(u)。

綜上所述,引理5成立。

引理6對于任意兩個節點x和y,以及i(0≤i≤n),更新算法得到的索引集IN′,也有TopKQuery(IN′,x,y,i)=d1(x,y,i)。

證明:如果i=0,兩個節點不可達,距離都為無窮大,成立。現在假設對于所有的i>0,對于任意兩個節點x和y,都有TopKQuery(IN′,x,y,i-1)=d1(x,y,i-1)。

只需要證明:

TopKQuery(IN′,x,y,i)=d1(x,y,i)

令δ=d1(x,y,t)。

1) 如果δ=d1(x,y,i-1),那么一定成立。

2) 如果δ

δ=d′(x,vi)+d′(vi,vi)+d′(vi,y)

只需要證明:

(vi,d′(x,vi))∈L′r(x)
(vi,d′(vi,y))∈L′(y)

在這里證明(vi,d′(vi,y))∈L′(y)。

如果d′(vi,y)=d(vi,y),從d1(x,y,t)=mink(d(x,vi)+d(vi,vi),d(vi,y))

可知,在經過節點vi到節點y的路徑中,不存在一個節點vl(l

d1(vi,y,i-1)>d1(vi,y,i)

如果(vi,d′(vi,y))?L(y),那么:

TopKQuery(IN,viy,i)=d1(vi,y,i-1)>d1(vi,y,i)存在矛盾,因此一定有(vi,d′(vi,y))∈L(y)?L′(y)。

如果d′(vi,y)

綜上所述,引理6成立。

由引理5和引理6可知,如果增加了邊(a,b),所有經過該邊的路徑距離都加入到了新的索引集中,得到的IN′為圖G′的索引集。

3.5 復雜度分析

點更新的時間消耗是O(1)。

對于邊更新,假設距離標記表的平均長度為len,即|L(v)|=len,|Lr(v)|=len。

那么在UpdateCircle(L,Lr,C,a,b)中要更新的節點數為O(len),由于需要對每個要更新節點查詢當前的k短路,每次查詢時間為O(klogk),因此在線自循環標記表更新的時間復雜度為O(len·klogk)。假設在UpdateIndex(G,L,Lr,vi,b,D)中需要訪問的節點數量為O(len),那么訪問的點對數為O(kt),在算法6中的TopKQuery需要O(len·klogk)的時間,算法4中需要訪問該過程O(len)次,因此該過程的時間復雜度為O(len2·k2tlogk)。在剪枝之后len,t?n,在平時應用中k也是比較小的數字。

所以更新算法的時間復雜度為O(len·klogk+len2k2tlogk)。

在更新算法中,使用的依舊是初始索引集的存儲空間。被修改的節點數量為O(len),每個節點最多修改了O(k)個值,因此空間復雜度為O(len·k)。

因此,動態top-k更新算法的時間復雜度為O(len·klogk+len2·k2tlogk),空間復雜度為O(len·k)。

3.6 實例分析

本節將通過實驗評估動態更新算法的時間和空間開銷。實驗用計算機的配置為Intel Core i5(2.7 GHz)處理器,8 GB內存。實驗程序運行于單個CPU核心上。數據集選用了Wikipedia vote network[17],其中包括7 115個節點,103 689條邊,在這里默認所有邊的權值為1。

實驗中,更新量定義為網絡更新時增加的邊數,每增加一條邊,運行一次更新算法。定義實驗數據的單位平均值為實測值除以更新量,即每增加一條邊的平均數據指標。

實際結果顯示,更新算法的總開銷和單位平均開銷與k和更新量有關。

表1列出了更新量為1 000時不同的k對應的更新算法的時間開銷。可以看到,隨著k的增加,單位平均更新時間也在逐漸增加;在k不大于32時,更新算法的時間開銷較小,性能較好。

表1 k對更新算法時間開銷的影響

表2列出了更新量為1 000時不同的k對應的更新算法的空間開銷和訪問節點數量。可以看到,訪問節點單位平均數量約等于索引集中每節點索引的平均長度;k不大于32時,訪問節點單位平均數量和索引集總長度單位平均增長率都較小,此時更新算法的空間開銷和訪問節點數量較小。

表2 k對更新算法空間開銷和平均訪問節點數量的影響

圖3顯示了k=16時更新量對單位平均更新時間與訪問節點單位平均數量的影響。可以看到三項數值都隨著更新量的增大而增大,但增長較為緩慢,且自循環標記表單位平均更新時間僅為微秒級。由此說明算法在大更新量下仍能保持較好的性能。

圖3 更新量對更新算法時間開銷和訪問節點數量的影響

4 結 語

本文對于動態更新的有向帶權圖中的top-k最短路徑問題設計并分析了基于2-hop cover的算法。更新算法利用預處理得到的索引集來更新,從而顯著地減少了圖改變后計算top-k的時間,只需要消耗更新算法的時間,而且有效地利用了原圖的數據集和性質。

本文在理論上證明了更新算法的正確性,分析了時間復雜度,證明其更新時間遠小于將更新后的圖看作靜態圖重新預處理的計算時間,并用實驗評估了更新算法處理大規模網絡的實際性能,證實了其實用性。動態top-k更新算法在大規模圖上可以進行快速更新和查詢。在未來研究工作中,我們將研究對于批量的更新如何提高效率、減少時間,以此來提高更新算法在實際應用中的性能。

主站蜘蛛池模板: 三级欧美在线| 激情無極限的亚洲一区免费 | 91精品啪在线观看国产91九色| 中文字幕亚洲另类天堂| 国产福利2021最新在线观看| 无码aaa视频| 嫩草国产在线| 无码电影在线观看| 国产AV毛片| 青青久久91| 日本人妻丰满熟妇区| 国产簧片免费在线播放| 熟妇丰满人妻av无码区| 91无码视频在线观看| 国产成人a在线观看视频| 国产又色又刺激高潮免费看| 亚洲IV视频免费在线光看| 欧美啪啪视频免码| 青草午夜精品视频在线观看| 久久精品欧美一区二区| 亚洲欧美日韩天堂| 手机精品福利在线观看| 免费无码AV片在线观看国产| www.99在线观看| 亚洲欧美另类久久久精品播放的| 久久综合伊人77777| 国产理论精品| 亚洲天堂2014| 国内精自线i品一区202| 国产精品不卡片视频免费观看| 高潮爽到爆的喷水女主播视频| 亚洲一区二区三区在线视频| 天堂在线www网亚洲| 中文字幕天无码久久精品视频免费 | 亚洲网综合| 青草视频在线观看国产| 免费一级无码在线网站| 日韩精品高清自在线| 九九热这里只有国产精品| 亚洲第一成网站| 91无码人妻精品一区二区蜜桃| 亚洲乱码视频| P尤物久久99国产综合精品| 国产又爽又黄无遮挡免费观看| 在线国产毛片手机小视频| 99久久性生片| 国产成人无码久久久久毛片| 国产小视频a在线观看| 国产女人爽到高潮的免费视频| 在线中文字幕网| 色综合久久久久8天国| 国产理论精品| 91精品国产麻豆国产自产在线| 国产成人综合亚洲欧美在| 亚洲精品在线影院| 99在线小视频| 国产成人亚洲无码淙合青草| 99在线观看国产| 欧美日韩成人在线观看| 中文字幕 91| 国产av剧情无码精品色午夜| 国产精品欧美亚洲韩国日本不卡| 成人一区在线| 亚洲精品福利网站| 99成人在线观看| 55夜色66夜色国产精品视频| 美女潮喷出白浆在线观看视频| 亚洲第一页在线观看| 四虎精品免费久久| 午夜一区二区三区| 欧美黄网站免费观看| 免费看一级毛片波多结衣| 免费A级毛片无码免费视频| 亚洲免费毛片| 黄色国产在线| 亚洲欧美人成电影在线观看| 国产精品真实对白精彩久久| 色婷婷亚洲综合五月| 波多野结衣无码AV在线| 久久久国产精品免费视频| 亚洲免费毛片| 亚洲一区二区三区麻豆|