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

基于時序Petri網的機器人柔性作業車間無死鎖調度優化算法

2024-12-31 00:00:00陳海波張超隆董建明
浙江理工大學學報 2024年11期

摘 要: 研究了一類源于芯片生產等智能制造領域的機器人柔性作業車間調度問題。針對該類問題在算法設計過程中,需要同時考慮工件加工、機器人運輸和死鎖求解等帶來的問題,提出了一種基于時序Petri網的無死鎖調度優化算法。首先,對問題進行時序Petri網建模,提出了一種基于Petri網變遷串的解表示形式,以方便死鎖求解和算法尋優;其次,通過對問題死鎖的結構分析和分類,提出了一種死鎖判斷和求解算法,并證明了算法可在多項式時間內求解任意類型死鎖,同時也能在一定程度上保證解的優良結構;最后,提出了問題的一種離散蜂群算法求解方案,在算法尋優過程中進行實時死鎖求解以保證解的可行性和較快的收斂速度。不同規模實例的數值實驗和與問題最優解及最優解下界的比較分析表明,提出的算法對不同類型實例都體現了較好的性能和較低的時間復雜度。該研究為復雜作業車間實時調度問題的算法研究提供了新思路和方法。

關鍵詞: 作業車間調度;離散蜂群算法;單機器人;Petri網;阻塞和死鎖

中圖分類號: TP301.6

文獻標志碼: A

文章編號: 1673-3851 (2024)11-0839-12

引文格式:陳海波,張超隆,董建明. 基于時序Petri網的機器人柔性作業車間無死鎖調度優化算法[J]. 浙江理工大學學報(自然科學),2024,51(6):839-850.

Reference Format:" CHEN Haibo,ZHANG" Chaolong,DONG Jianming. A deadlock-free scheduling algorithm of robotic flexible job shops based on temporal Petri nets[J]. Journal of Zhejiang Sci-Tech University,2024,51(6):839-850.

A deadlock-free scheduling algorithm of robotic flexible job shops based on temporal Petri nets

CHEN Haibo1, ZHANG Chaolong1, DONG Jianming2

(1.School of Computer Science and Technology (School of Artificial Intelligence), Zhejiang Sci-Tech University, Hangzhou 310018, China;

2.School of Management and E-Business, Zhejiang Gongshang University, Hangzhou 310018, China)

Abstract:" A robotic flexible job shop scheduling problem originating from intelligent manufacturing fields such as chip manufacturing was studied. To solve this problem, it is necessary to consider the difficulties caused by job processing, robot transportation and deadlock solving. As a result, a deadlock-free scheduling optimization algorithm based on temporal Petri net was proposed. Firstly, the problem was modeled by temporal Petri net, and a solution representation based on Petri net transition string was proposed to facilitate deadlock solving and algorithm optimization. Secondly, through the analysis and classification of the deadlock structure, a deadlock detection and resolution algorithm was proposed, and it was proved that the algorithm could solve any type of deadlock in polynomial time, and at the same time could ensure the good structure of the solution to a certain extent. Finally, a discrete swarm algorithm was proposed to solve the problem. In the optimization process of the algorithm, deadlock was solved to ensure the feasibility of the solution and fast convergence speed. The numerical experiments of different scale instances and the comparison with the optimal solution and the lower bound of the optimal solution showed that the proposed algorithm had good performance and lower time complexity for different types of instances. The study provides a new idea and method for the algorithm research of complex job shop real-time scheduling problems.

Key words: job shop scheduling; discrete bee colony algorithm; single robot; Petri nets; blocking and deadlock

0 引 言

隨著工業互聯網及工業4.0的逐步深入,制造業的升級朝著無人化和智能化推進,滿足“多樣化、小規模、周期可控”的柔性化生產制造將成為趨勢。大量無人輔助制造工具,例如自動導引機器人(Automated guided vehicle, AGV)、碼垛機器人等設備的使用,大大提高了生產制造的效率,同時也給柔性生產調度提出了更高的要求。調度理論主要研究如何將有限的資源進行優化并合理安排任務的執行順序使得設定的目標(例如時間最短、成本最低或利潤最高等)達到最優,已被大量用于生產調度與排程等諸多領域[1]。柔性作業車間調度問題(Flexible job-shop scheduling psroblem, FJSP)研究復雜柔性車間生產任務的優化調度,是運籌學、調度理論和工業工程等研究的重要問題之一[2]。在傳統作業車間調度問題中,工件具有若干道工序,需要在不同的機器上進行加工,且不同工件的工序數量及其在機器上的加工順序可以不同,同一工件的不同工序不能同時加工。若同一臺機器可以加工不同類型的工序,從而產生了一類更為一般的柔性作業車間調度問題。柔性車間作業調度要求機器具有多功能屬性,可以用于加工不同類型的工序,所以更加符合小批量定制等新型生產模式的要求。

經典的柔性作業車間調度問題通常不考慮工件(工序)在機器間的移動,并且假設工件(工序)在機器上加工完后即被移走,或具有一個無限大的緩沖區可以存放加工完的工件。在現代智能制造及無人工廠場景下,除了需要優化生產環節,產品的運輸及路徑的優化調度同樣至關重要[3]。特別是在一些典型柔性制造場景下,通常工件(工序)加工完后需要由無人運輸工具(例如AGV)進行運輸至下一道工序所需的機器上;由于加工區域狹小,機器通常并沒有緩沖空間可以存放加工完的工件,因此工件可能要滯留在被加工機器上,等待運輸機器人空閑時將完工工件移出機器并將其運至后續工序機器等待可用時進行加工,這種緩沖約束被稱為“阻塞”。當運輸工具數量有限且調度不當時,機器的“阻塞”可能會導致“死鎖”現象的發生,此時整個加工過程將會停滯。所以,在實時作業車間調度環境下,如何保證機器柔性的情況下,將工件加工和運輸工具調度進行協同優化,并實現無死鎖調度與算法尋優協同至關重要。

為了適應個性化定制等新型制造模式,在傳統的作業車間調度問題(JSP)模型基礎上,柔性作業車間調度問題得到了廣泛的研究[2]。例如,將傳統的每臺機器只能加工一類工序,擴展成每臺機器可以加工多種類型工序,或者每道工序可以選擇多臺平行機器中的一臺進行加工,甚至考慮加工時間等參數不確定情況下的動態柔性作業車間調度問題[4-5]。典型(柔性)作業車間調度問題假設機器緩沖區是無限的,即工件在機器上加工完成后機器即可釋放。事實上,在大量實際制造場景下,由于制造空間狹小等原因,通常機器緩沖區有限甚至無緩沖區,所以已有研究人員開始研究緩沖區的不同設置對算法設計的影響。例如,Liu等[6]研究了具有4種緩沖區約束(即無等待、無緩沖、有限緩沖和無限緩沖約束)情況下的廣義JSP問題,給出了問題的混合整數規劃模型和多個啟發式算法。Bektur等[7]研究了工件具有到達時間且機器具有有限緩沖區情況下的不同類機調度問題,給出了問題的混合整數規劃模型、禁忌搜索算法(TS)和模擬退火(SA)算法求解算法。考慮緩沖區有限情況下的模型,需要同時考慮工件的加工和緩沖區的合理使用,所以算法設計變得更加困難。

當緩沖區有限情況下,工件加工完后通常需要借助運輸工具(機器人)將完工工件運出緩沖區,并同時運至下一階段的機器,此時加工與運輸協同調度效率顯得異常重要。同時,當緩沖區容量為0時,即機器無緩沖區時,考慮阻塞和死鎖情況下的調度,不僅符合現實而且是解決緩沖區有限且考慮工件運輸情況下調度問題的關鍵[8]。無死鎖調度決策在實時調度系統中較為常見。例如,Li等[9]研究了一類機器無緩沖區的柔性裝配系統(AFSSP)調度問題,工件在裝配過程中會存在死鎖約束。Xing等[10]研究了柔性裝配系統(FAS)的死鎖控制問題,使用Petri網(Petri Net)分析了FAS問題的活性結構,并設計了Petri網控制器以控制FAS的無死鎖運行。李浚等[11]研究了機器人柔性流水車間調度問題,提出了一種時序Petri網監督學習的啟發式算法,并且設計了時序Petri網的全連接神經網絡模型。Luo等[12]研究柔性制造系統(FMS)中的死鎖避免策略,提出了一種類似銀行家算法的死鎖避免策略,并分析了所設計策略的時間復雜度。這些研究主要集中在裝配線調度問題,且研究的核心集中在問題的混合整數規劃模型和基于Petri網的死鎖控制策略兩個方面。然而,采用基于Petri網的死鎖仿真往往會產生狀態爆炸等現實問題,因而在實際應用中受到一定的限制。

對于考慮機器人運輸的調度問題研究,現實調度系統的機器環境主要涉及車間作業環境。例如,Hurink等[13]將傳統的作業車間調度問題分為兩個階段研究,先確定工件調度,再確定工件的運輸路線規劃并將機器設置成無限緩沖區,同時給出了問題的一種禁忌搜索算法。Nouri等[14]研究了一類具有運輸時間和多機器人的作業車間調度問題(JSPT-MR),該問題將機器設置成無限緩沖區,并給出了一種基于混合多智能體的啟發式算法。晏曉凡[15]研究了針對物料機器人指派和作業車間的聯合調度問題,該調度環境為無限緩沖區,并且設計了一種改進灰狼優化算法進行求解。上述問題都假設機器的緩沖區無限,此時機器人的運輸與工件的加工不會形成死鎖。Sun等[16]則研究了一類機器緩沖區容量為0的作業車間調度問題,且問題中同時考慮機器人對工件的運輸和死鎖情形,并給出了該類問題的一個混合整數規劃模型,對小規模問題實例進行了規劃模型的求解,但未給出問題的有效算法。

本文將進一步解決文獻[16]中提出的問題,即給出該類問題的無死鎖調度算法。該類問題在芯片生產等領域具有廣泛應用背景。例如,芯片生產分為清洗、光刻、蝕刻等若干步驟單元,其中清洗單元是其中一個至關重要的步驟。在清洗單元,晶圓需要經過浸泡去除有機物、氧化層去除、顆粒去除、金屬離子去除等多道工序,且工序間具有嚴格的先后順序要求,并在多個清洗站由專門的機器負責完成,如槽式清洗機、DI水沖洗機、超聲波清洗機、干燥設備等。清洗單元對環境潔凈度要求極高,例如在英飛凌晶圓處理環境中,28 L空氣中顆粒度需小于1。因此,清洗單元內部部署機器人,能在保持空氣高潔凈度的同時,精確地將晶圓移動至不同的清洗機器進行處理。如何優化機器人的調度和晶圓清洗流程的排序是清洗單元高效運行的關鍵。本文針對該類機器緩沖區容量為0且考慮機器人運輸的柔性作業車間調度問題(Robotic flexible job-shop scheduling problem, RFJSP),給出了該問題的首個無死鎖調度的高效算法。首先對問題進行時序Petri網建模,將調度問題的解轉化成相應的Petri網變遷串;其次,對問題解的結構進行分析,將死鎖類型進行分類,并基于Petri網變遷串的解結構設計了相應的高效死鎖檢測和求解算法,并證明算法在多項式時間內可以求解任意死鎖類型;最后,設計了該問題的人工蜂群算法解決方案。在蜂群算法設計中,解結構采用Petri網變遷串,且當算法運行時遇到死鎖導致解不可行時,則采用設計死鎖判定和求解算法進行死鎖的求解,使算法的尋優過程既能保證解結構的優良特性又始終在可行解中進行。

1 RFJSP問題的形式化描述及時序Petri網建模

首先給出RFJSP問題的形式化描述:n個工件的集合記為J={j1,j2,…,ji,jn},需要在m臺機器上進行加工,機器集合記為M={M1,M2,M3,…,Mm}。對于每個工件ji(i=1,2,3,…,n),由θi個工序組成,所以每個工件也可以表示為ji={Oi,1,Oi,2,…,Oθi},且工序的加工順序給定。對于每道工序Oi,j,需要在對應機器M(Oi,j)∈{M1,M2,M3,…,Mm}上加工,加工時間記為ρi,j。工件在機器上的加工需要滿足作業車間環境(Job shop)要求,即每個工件工序的加工順序可以不同,且同一臺機器可以加工同一工件的多道工序(機器柔性),但是同一時刻同一工件的不同工序不能同時加工。

有一臺移動機器人記為R,負責工件在不同機器上的運輸,且機器人每次只能運輸一個工件,例如將工件j從當前工序Oi,j所在的機器M(Oi,j)運輸到下一道工序Oi,j+1所在的機器M(Oi,j+1)。每臺機器沒有緩沖區,即加工完相應工序的工件需要一直占用機器(機器阻塞),直到機器人將工件移出機器。所以,機器人的移動方式有兩種,即空載移動和負載移動。空載移動是指機器人在沒有攜帶任何工件的情況下移動,而負載移動是指機器人在攜帶工件的情況下移動至指定機器。假定,每個工件的第1道工序都需要通過移動機器人R從初始站S中運出,且最后一道工序完工后由移動機器人運輸至存儲站D。機器人在相鄰機器間的運輸時間記為單位時間,且以跨機器數來計算運輸總時間。

問題的優化目標是合理安排工件(工序)在機器上的加工順序與機器人運輸順序協同,使得最后完工工件的完工時刻(makespan)最早。本文中對RFJSP問題的解進行Petri網建模,首先介紹了Petri網的基本概念,然后設計了RFJSP問題的時序Petri網模型。

1.1 Petri網

Petri網(下稱“PN”)是一個具有兩種節點類型的有向二分圖,分別稱為“Place”(庫所)和“Transition”(變遷)。其中節點通過有向弧連接,且有向弧從一個庫所指向一個變遷,或從一個變遷指向一個庫所。在圖表示中,庫所、變遷和弧分別用圓圈、方框和箭頭表示[17]。所以,一個PN可以表示為一個三元組N=(P,T,F),其中:P是一組有限的庫所集合;T是一組有限變遷的集合,滿足P≠,T≠,P∩T=;F(P×T)∪(T×P)表示有向弧的集合[1]。每條弧上可以添加權重,用來表示變遷觸發的條件。弧權重為1的PN稱為普通Petri網。對于x∈P∪T,x的前集用·x={y∈P∪T|(y,x)∈F}表示,同樣x的后集可以表示為x·={y∈P∪T|(x,y)∈F}。給定集合Z={0,1,2,3,…}以及,標Zk={1,2,3,…,k}記ψ是一個映射ψ:(PZ),給定一個庫所p∈P和一個標記ψ,ψ(p)表示在標記ψ時庫所p中的令牌(token)的數量。帶有初始標記ψ0的PN稱為有標記的PN,表示為(N,ψ0)。

如果存在一個變遷t(t∈T),p∈·t,ψ(p)gt;0,那么變遷t稱為使能狀態(可以被觸發),用ψ[tgt;來表示。對于普通PN,在標記ψ狀態下,使能狀態下的變遷t被觸發從而產生新的標記ψ′,用ψ[tgt;ψ′表示。由于是普通PN,即弧權重為1,每次變遷需要1個令牌,所以ψ′(p)會出現如下3種情況。第1種情況:ψ′(p)=ψ(p)-1,p∈·t/t表示變遷t的前集的所有庫所中的令牌(Token)的數量減少1個單位;第2種情況:ψ′(p)=ψ(p)+1,p∈·t/t,表示變遷t的后集的所有庫所中的令牌(Token)的數量增加1個單位;第3種情況:ψ′(p)=ψ(p),pt·∪t·,表示未發生變遷的t前后集的庫所中的令牌數量不變。如果變遷串λ=t1,t2,ti,…,tk,變遷ti被觸發的狀態時刻的標記為ψi[ti,當滿足條件ψi[tigt;ψ′i+1,i∈Zk,ψ1=ψ,那么就稱變遷串λ在標識ψ下一定是可達的。

時序Petri網在普通Petri網的基礎上引入了時間概念,以更精確地模擬和分析系統中各個事件的時間行為。允許令牌在庫所中滯留一段時間以便滿足變遷的條件,因此,時序Petri網能夠更準確地反映系統中的動態行為。

1.2 RFJSP問題的時序PN模型

本文使用時序PN對RFJSP問題進行建模。首先需要將RFJSP問題的元素映射為PN的基本元素,即將庫所、變遷和弧與RFJSP問題中的元素一一對應,并且給定一個RFJSP問題的實例,則可以產生一個時序Petri網模型與之對應:

a)對于給定問題實例,每道工序對應2個狀態庫所,分別表示工序正在被機器人運輸至指定機器狀態和工序正在指定機器上進行加工的狀態,且這兩個庫所都帶時間窗口,時間窗口包括運輸時間和加工時間。

b)每道工序對應2個變遷,分別表示工件(工序)開始運輸和工件(工序)開始加工。

c)對于每個工件的第1道工序和最后一道工序,分別設置1個表示開始狀態的初始庫所和1個表示終止狀態的存儲庫所。

d)對于每個工件的最后一道工序,設置2個虛擬變遷,分別表示最后一道工序完工后開始運往存儲站和最后一道工序運至存儲站并釋放機器人。

e)整個模型中有1個機器人資源庫所和m個機器資源庫所。

f)每個工件的初始庫所、機器人資源庫所和m個機器資源庫所初始各有1個令牌。所以,對于實例中的每個工件ji,(ji∈J),可以將其每道工序在PN模型中的庫所和變遷按照規定的加工順序組成如下“庫所-變遷”序列ξi:

ξi=pi,S,ti,1,pi,1,ti,2,pi,2,tt,3,…,ti,l,pi,l,…,ti,L,pi,L,ti,(L+1),pi,D(1)

其中:pi,S表示工件ji的初始庫所,ti,1表示工件ji的第1道工序開始運輸變遷,pi,1表示第1道工序正在運輸狀態庫所,ti,2表示第1道工序開始加工變遷,pi,2表示第1道工序正在指定機器上加工狀態庫所。依此類推。ti,L表示最后一道工序開始運輸變遷,pi,L表示最后一道工序正在運輸狀態庫所,ti,(L+1)表示最后一道工序運至存儲站變遷,pi,D表示工件ji完工且已到達存儲站狀態庫所。

表1給出了RFJSP問題的一個具體建模實例,其中包含3臺機器M1,M2,M3,1臺移動機器人R,4個工件j1,j2,j3,j4,每個工件都有5道工序。

以工件j1為例,由式(1)得到其“庫所-變遷”序列如下:

ξ1=p1,S,t1,1,p1,1,t1,2,p1,2,t1,3,…,t1,l,p1,l,…,t1,11,p1,11,t1,12,p1,D。

其具體的Petri網示意圖如圖1所示。

同j1和j4,工件j2和j3同樣可以得到相應的“庫所-變遷”序列和相應的PN模型表示。同樣以工件j1的第1道工序為例,其各庫所狀態變化和變遷觸發規則可描述如下:工件j1初始化在S站中,經過變遷t1,1請求到機器人R(如果R庫所存在令牌),當變遷t1,1觸發之后令牌轉移至p1,2狀態庫所,且工件j1在庫所p1,2狀態中停留一定時間(機器人移動到機器所需時間);當機器M3可用時(機器M3庫所存在令牌),則變遷t1,2觸發并轉移至庫所p1,3進入在機器M3上加工狀態,同時歸還令牌給機器人R庫所。令牌在庫所p1,3需停留一定時間(至少是第1道工序的加工時間)來完成工序。工件j1的后續工序可按照上述規則執行。工件j2、j3、j4同j1。

1.3 解的表示形式

本文根據RFJSP問題的時序PN模型表示定義RFJSP問題解的形式。確定RFJSP問題解主要有兩個方面:一是確定工件(工序)之間的加工順序,二是確定機器人運輸工件的順序(運輸路徑)。但是,在單機器人和機器零緩沖區情況下,任意可行解中工件的加工順序與機器的運輸順序是一致的。所以,考慮使用PN模型中的變遷串來表示RFJSP問題的解(可行或不可行)。例如,針對表1實例,給定工序的加工順序的局部為(O1,1,O4,1,O1,2,O1,3,O4,2,…),可以構造如下變遷串λ來表示問題的解:

λ=t1,1,t1,2,t4,1,t4,2,t1,3,t1,4,t1,5,t1,6,t4,3,t4,4,…。

上述變遷編號如圖1所示,并依此規則可以將完整解進行表示。為了保證解的基本可行性,即滿足每個工件工序的加工順序必須滿足指定的順序,在構造變遷串時,同一工件工序的變遷串標號保持遞增順序以保證解滿足工序先后順序的硬約束。同時,可以根據變遷串對于初始標記是否可達來判斷當前變遷串所對應的解是否可行。例如,若初始標記ψ=p1,S+p2,S+p3,S+p4,S+pM1+pM2+pM3+pR,則上述變遷串λ所表示的局部解是可行的。

同時,對于任一工序都包含兩個變遷,即工件開始運輸變遷和工序開始加工變遷(對于每個工件的最后一道工序另外分別增加到虛擬存儲站的兩個變遷)。不難發現,當單機器人情況下,工序的運輸變遷和加工變遷在變遷串中必須相鄰且都先后被觸發;否則,該變遷串所表示的解不可行且一定會發生死鎖(關于死鎖的描述見2.1),因為運輸變遷觸發后若不能及時觸發加工變遷,則機器人無法釋放且產生死鎖。所以,本文將解的變遷串表示成運輸變遷和加工變遷成對出現的形式,以避免死鎖發生。例如上述變遷串λ就是成對出現,所以可以進一步表示成如下形式:

λ=(t1,1,t1,2),(t4,1,t4,2),(t1,3,t1,4),(t1,5,t1,6),(t4,3,t4,4),…。

2 死鎖的判定與求解

2.1 阻塞與死鎖

在RFJSP問題中,由于機器上沒有緩沖區,所以當工件加工完成后無法及時將工件運離機器,則機器將會出現阻塞。若機器阻塞導致后續工序無法繼續被加工則產生死鎖。下面根據不同的死鎖結構和求解死鎖的不同方案將RFJSP問題可能存在的死鎖分成兩類如下。

第1類死鎖:以表1為例,假設工序加工順序的局部為(O1,1,O2,1,O1,2,…),初始標記ψ0=p1,S+p2,S+p3,S+pM1+pM2+pM3+pR,變遷串λ1=t1,1,t1,2,t2,1,t2,2,…可以表示該加工順序所表示的解。對于λ1的子串(t1,1,t1,2,t2,1)對于標記ψ0是可達的,即ψ0[λ1gt;ψ1,ψ1=p1,3+p2,2+p3,S+p4,S+pM1+pM2而對于λ1的子串(t1,1,t1,2,t2,1,t2,2),對于標記ψ0不可達,所以該解不可行且在變遷t2,2處發生死鎖。

因為工序O1,1已經在機器M3上加工完成,此時等待機器人將其卸載,并運送至下一臺機器M1上繼續加工第2道工序O1,2。但是此時,工件j2的工序O2,1需要緊跟O1,1之后加工并已經請求到機器人(負載狀態)。這導致了工件j1阻塞在機器M3上(等待機器人將其運走),工件j2占用機器人等待機器M3空閑,從而導致死鎖發生,如圖2(a)所示。該死鎖發生的原因是O1,1和O2,1需要使用同一臺機器且工序O2,1的加工順序不應排在O1,1之后,同時可以簡單通過將工序O2,1后移至后面相應位置進行加工即可求解死鎖。本文記符合該特征類別的死鎖為“死鎖I”。

第2類死鎖:以表1實例為例,假設工序加工順序的局部為(O1,1,O1,2,O2,1,O4,1,…),初始標記ψ0=p1,S+p2,S+p3,S+p4,S+pM1+pM2+pM3+pR。

變遷串λ1=t1,1,t1,2,t1,3,t1,4,t2,1,t2,2,t4,1,t4,2,t2,3,t2,4,…可以表示該加工順序所表示的解。對于λ1的子串(t1,1,t1,2,t1,3,t1,4,t2,1,t2,2,t4,1,t4,2,t2,3)對于標記ψ0是可達的,即ψ0[λ1gt;ψ1,且ψ1=p1,5+p2,4+p4,3+pM3。而對于λ1的子串(t1,1,t1,2,t1,3,t1,4,t2,1,t2,2,t4,1,t4,2,t2,3,t2,4)則對于標記ψ0不可達,所以該解不可行且在變遷t2,4處發生死鎖,如圖2(b)所示。

因為,當變遷t1,1,t1,2,t1,3,t1,4,t2,1,t2,2,t4,1,t4,2,t2,3被觸發后,此時機器M1和M3有工件在加工或被阻塞,M2空閑且工序O2,2對應工件j2已申請到機器人的運輸。此時,由于工序O2,2需要M1,而M1被阻塞,所以此時表示O2,2開始加工的變遷t2,4無法被觸發,從而產生死鎖。事實上,此時機器人去裝載任意一個工件都將會造成死鎖,且該類死鎖無法通過簡單移動t2,4及其后面的變遷位置來求解,需要調整O2,2之前已經加工的工序的順序。記該類特征的死鎖類型為“死鎖Ⅱ”。

2.2 死鎖類型的判定及求解算法

給定RFJSP問題解的一個變遷串,需要首先判斷是否有死鎖且分析死鎖類型,然后求解死鎖。假設給定問題實例解的一個變遷串表示形式為λ=(t1,1,t1,2),(t4,1,t4,2),(t1,3,t1,4),(t1,5,t1,6),(t4,3,t4,4),…。由于每道工序的運輸變遷和加工變遷成對出現,為了算法描述方便,將解的形式簡化表示成λ=α1,α2,α3,…,αθ+n-1,αθ+n,其中αi表示第i對變遷串,即表示第i順序被加工的工序,n和θ分別為工件總數和工序總數。本文設計死鎖判定和求解算法如下,記為算法DF(λ)。

算法 DF(λ)

輸入:RFJSP問題實例解的變遷串表示λ,給定初始標記ψ0。

Step1. 對變遷串λ從左到右進行遍歷,依據初始標記ψ0,判斷當前變遷串是否可被觸發。

Step2. 若當前變遷對αi無法被觸發,此時根據如下規則判斷死鎖類型及求解當前死鎖。

Step2.1. 若從變遷對αi+1開始,依次往后遍歷,當遍歷到的變遷對記為αl時,將αl插入到αi之前,即第i個變遷串位置時,可使αl被觸發,則將其插入并保持αi及之后的變遷對的相對順序不變,并對變遷對下標重新編號。繼續從αi+1位置開始遍歷,進行Step2操作尋找下一個死鎖變遷對。αi處發生了死鎖Ⅰ。

Step2.2. 若從變遷對αi+1開始,依次往后遍歷,不存在某個變遷對αl,使得將其插入至第i個變遷串位置時可被觸發,則說明當前αi處發生了死鎖Ⅱ,繼續按照如下規則解鎖:當該情況發生時,說明當前從αi開始之后的變遷串所對應工序所需的機器都處于阻塞狀態,此時記這些被阻塞的機器中,最早被阻塞的機器所對應的被觸發變遷對為αk,其對應的工序為Os,t。

Step2.2.1. 若Os,t是工件js的最后一道工序,則將完工工件返回存儲站D的變遷對插入到αk之后,使得機器M(Os,t)被釋放,αk之后的其他變遷對相對位置不變,并對變遷對重新編號。從變遷對αk+1開始回到Step2繼續往后遍歷。

Step2.2.2. 若Os,t不是工件Os,t+1的最后一道工序,則將工序Os,t+1所對應的變遷對插入到αk之后,使得機器M(Os,t)被釋放,αk之后的其他變遷對相對位置不變,并對變遷對重新編號。從變遷對αk+1開始回到Step2繼續往后遍歷。

Step3. 若遍歷到最后一個變遷對αθ+n,則遍歷結束,此時得到的變遷串即為無死鎖變遷串,并將其輸出,算法結束。

引理1 對于RFJSP問題任意解的變遷串λ,在初始標記ψ0下,若λ中存在死鎖,算法DF(λ)一定可以在多項式時間內求解死鎖,并將λ轉化成無死鎖的可行解。

證明 算法DF(λ)的基本思想是,從左到右遍歷變遷串λ,若當前遍歷位置的變遷串發生死鎖Ⅰ,則將當前位置之后的某個變遷對插入到當前位置之前并被觸發,然后繼續往后遍歷。顯然,若遍歷過程中,每次遇到的死鎖都是死鎖Ⅰ,即每次都能在當前死鎖變遷串位置之后找到變遷對插入到當前位置之前并被觸發,當遍歷結束即可獲得一個無死鎖的可行解且遍歷過程可在多項式時間內完成。

若當遍歷過程中,當前位置的變遷串無法被觸發,且將當前位置之后的任意一個變遷對插入到當前位置之前都無法被觸發,即此時發生死鎖Ⅱ;按照算法DF(λ),此時將需要在當前位置之前的某個位置之后插入一個變遷對,即Step2.2中的操作。為了證明Step2.2可求解死鎖Ⅱ,且可在多項式時間內完成,只需證明Step2.2可行,且每次求解死鎖Ⅱ都可以使得無死鎖變遷對增加即可。

首先,證明Step2.2的可行性。當死鎖Ⅱ發生時,當前變遷串αi及之后所有變遷串所對應的工序所需的機器都被阻塞,記這些機器的集合為Mz,即當前變遷串之前變遷串所對應的工序在這些機器上加工,并且之后并沒有同一工件的工序被加工而使這些機器被釋放。顯然按照Step2.2,在Mz中找到一臺最早被阻塞的機器(記為Mk)和相應被觸發的變遷串(記為αk),若在Mk上被觸發變遷串Mk所對應工序是對應工件的最后一道工序,則可將對應工件的返回存儲站D變遷對插入到該變遷串之后即可,此時插入變遷串一定可被觸發且同時釋放機器Mk,即Step2.2.1可行。若變遷串αk所對應工序不是對應工件的最后一道工序,則必然可以找到其下一道工序所對應的變遷串(記為αx),并將其插入到該變遷串之后即可,此時αx一定可被觸發且同時釋放機器Mk。因為αx在原位置一定未被觸發,即i≤x,且Mk是Mz中最早被阻塞的機器且αx所需機器(記為Mc)也在Mz中且在αk前未被阻塞,所以此時αx必可被觸發,即Step2.2.2可行。所以Step2.2可行。

其次,當每執行一次Step2.2,插入的變遷串αx及之前變遷串都可被觸發。若繼續遍歷過程中遇到死鎖I,則通過Step2.1使得可被觸發的變遷串增加。若繼續遇到死鎖Ⅱ,則假設原變遷串中機器Mc在變遷對αc處被阻塞,如圖3所示。顯然,除了αc以外,原串變遷中不可被觸發的變遷依然在變遷對αi之后,且被阻塞的第1臺機器所在的位置至少在αx的新位置,即第k+1變遷對位置。所以繼續執行Step2.2,新的插入位置至少會在變遷串αk+1之后,即αx的新位置之后,所以每執行一次Step2.2,無死鎖變遷對的長度至少增加1。顯然,Step2.2執行的次數是多項式時間的,所以綜上可得算法可在多項式時間內完成。

綜上引理1得證。

2.3 目標函數計算

本文考慮RFJSP問題的目標函數是極小化最后完工工件的完工時間,也就是需要首先計算最后一個變遷串所對應的工序的開工時間。按照1.3中對變遷串的簡化定義,假設當前解的變遷串為λ=α1,α2,α3,…,αθ+n-1,αθ+n,其中:αθ+n-1為最后加工的工序所對應的運輸變遷和加工變遷對,αθ+n則表示最后加工的工序所對應的工件運輸至存儲站所對應的變遷對。

假設變遷對αk所對應的工序為Oi,j,并令f(Oi,j)表示αk中加工變遷的觸發時間,即工序Oi,j在機器上的開工時間。令Oi,j為αk-1所對應的工序,αk′為工序Oi,j-1所對應的變遷對,則由變遷觸發的順序規則和同一工件工序加工的順序約束,可得f(Oi′,j′)≤f(Oi,j),f(Oi,j-1)+ρi,j-1≤f(Oi,j)。

所以Oi,j的開工時間f(Oi,j)可計算如下:

f(Oi,j)=

max{f(Oi,j)+rt3+rt4,rt1+rt2}, "j=1;

max{f(Oi,j)+rt3+rt4,f(Oi,j-1)+ρi,j-1+rt4}, "jgt;1

(2)

其中:rt1表示機器人從機器M(Oi′,j′)到初始站S的運輸時間,rt2表示機器人從初始站S到機器M(Oi,j)的運輸時間,rt3表示機器人從機器M(Oi′,j′)到M(Oi,j-1)的運輸時間,rt4表示機器人從機器M(Oi,j-1)到M(Oi,j)的運輸時間。當Oi,j是解中第1個被加工的工序時,rt1=0,即機器人初始位置在初始站S。

根據式(2),對解的變遷串從左到右進行遍歷并計算每個工序的開工時間,則計算得到最后一個開工工序的開工時間時,即可得到解的目標函數值。

2.4 問題最優解的下界估計

記RFJSP問題的最優目標函數值為C*max,下面給出最優解的一個下界估計如下。

引理2 問題RFJSP的最優目標函數值滿足

LBRFJSP≤C*max,且

LBRFJSP=max{max1≤k≤m(∑ni=1μki)+min_rtS+min_rtD,max1≤i≤n(∑mk=1μki+∑θij=1rtij)}(3)

其中:n表示工件總數,m表示機器總數,μki表示工件ji在機器Mk上加工時間的總和,min_rtS表示從初始站到最近機器的運輸時間,min_rtD表示從存儲站到最近機器的運輸時間,∑θij=1rtij則表示工件ji的所有相鄰工序之間的運輸時間之和。

證明 因為表達式max1≤k≤m(∑ni=1μki)表示每臺機器上需要加工的工序總時間的最大值,顯然在最優解中也至少需要在機器上加工這些時間。每臺機器上的第1道工序和最后一道工序需要分別從來自初始站和運輸至存儲,需要至少min_rtS+min_rtD時間。所以最優目標函數值max1≤k≤m(∑ni=1μki)+min_rtS+min_rtD≤C*max。

由于問題RFJSP的機器環境是作業車間環境,即同一工件工序的加工不能重疊,且相鄰工序的加工之間需要通過機器人運輸,所以對于每個工件ji,其所有工件的加工時間之和為∑mk=1μki,其所有涉及的運輸時間為∑θij=1rtij,在最優解中,至少需要∑mk=1μki+∑θij=1rtij時間才能將工件加工完并運送至存儲站。所以最優目標函數值max1≤i≤n(∑mk=1μki+∑θij=1rtij)≤C*max。

綜上引理2得證。

3 離散人工蜂群算法設計

本文在對問題RFJSP進行時序PN建模后,將運用離散人工蜂群算法(Discrete artificial bee colony, DABC)[18]對問題進行求解。在求解過程中,將設計算法對產生的每個解進行死鎖判定,并對產生死鎖的不可行解進行求解,以使得解變成無死鎖的可行解。離散人工蜂群算法是一種用于在離散參數空間中進行全局尋優的群體智能優化算法,其靈感來自工蜂的覓食行為。DABC算法中用食物源的位置表示解,且DABC算法中存在3種類型的個體:雇傭蜂、觀察蜂和偵察蜂。雇傭蜂在其鄰域內尋找更好的解。基于雇傭蜂的狀態,觀察蜂選擇若干解并進行進一步搜索,以提高解群體的質量。偵察蜂為陷入局部最優解的個體生成新的解。DABC算法中有3個控制參數,即蜜源的數量SN(初始解的數量)、迭代次數、偵察峰出動的標準。本文設計的基于RFJSP問題的DABC算法如下。

3.1 種群初始化

種群初始化需要設定一些參數,主要包括算法迭代次數,解空間(蜜源(SN))的數量,雇傭蜂、觀察蜂、偵察蜂的數量,以及解被搜索閾值(Limit)。在DABC算法中,每一個解都是一組由PN生成的變遷序。在RFJSP問題中,一個解就是一個編碼。例如假設有兩個工件分別為ji和jj,ji和jj均只有兩道工序,那么根據1.3中形式化描述,ji在Petri網的變遷串為

λ1=(ti,1,ti,2),(ti,3,ti,4),(ti,5,ti,6)。

同樣,工件jj在Petri網的變遷串為

λj=(tj,1,tj,2),(tj,3,tj,4),(tj,5,tj,6)。

那么該編碼方式為ω=λi∪+λj,其中使用符號“∪+”表示集合的隨機交錯合并,λi中的變遷對與λj中的變遷對隨機合并,但是屬于λi中的變遷對必須保持在λi中的相對順序。例如

ω1=(ti,1,ti,2),(tj,1,tj,2),(ti,3,ti,4),(tj,3,tj,4),(ti,5,ti,6),(tj,5,tj,6),

其中:ω1表示一個解(蜜源)。采用隨機集合交錯合并方式可以一定程度保證解空間的多樣性。

3.2 雇傭蜂階段

雇傭蜂通過局部搜索策略在當前蜜源附近探索新蜜源,以搜索更好的解。根據鄰近蜜源的優先級雇傭蜂將舊蜜源更新為新蜜源,具體操作如下:

a)在ω中任意交叉兩個變遷對i,j,其他元素相對位置不發生變化得到ωnew。

b)ωnew在經過交叉操作之后可能存在死鎖情況,調用DF算法將ωnew修正為ω′new。

c)求解ω′new的目標函數,將新解ω′new與舊解ω比較,如果新解更優,則ω′new取代ω的位置,否則保持不變。

3.3 觀察蜂階段

雇傭蜂將蜜源信息與觀察蜂共享之后,觀察蜂將對蜜源區域進行搜索,觀察蜂的出動有助于提高搜索精度和全局收斂性,使用可變鄰近搜索的搜索方法,在解ω′new基礎上進行交換和隨機插入。

3.4 偵察蜂階段

隨著算法迭代次數增加,蜜源發生變化的概率逐漸減小。為了避免陷入局部最優,必須對蜜蜂所攜帶的蜜源限定次數,同時提高種群的多樣性。如果蜜源在限定次數內未發生任何變化,那么該蜜源的攜帶者就會轉化為偵察蜂,蜜源會被重新初始化。

4 數值實驗及結果分析

本文針對不同問題規模和輸入參數的問題實例進行大規模數值實驗,以驗證算法在不同情況下的效果和時間復雜度。數值實驗在配備了具有3.2 GiHz頻率中央處理器(CPU)的個人PC機上進行,且該機配置了8 GiB RAM。算法實現采用Python語言編程,并運用了Python的threading模塊多線程技術進行并行處理。由于該問題現有文獻沒有給出基準實例且尚無算法結果,因此本文采用相關問題[16]類似方法,隨機生成不同類型的實例以測試DABC算法的性能。生成8種不同類型的實例如表2所示,其中:第1列表示類型索引;第2列用一個三元組(n,θ,m)表示該類型實例的規模,分別表示工件的數量、工件工序數量的上界和機器數量。同時,工序加工時間也分兩種類型,分別取[1,10]與[10,100]之間的隨機整數。機器人在任意相鄰兩臺機器之間移動的時間是1,且移動時間依據相隔的機器數遞增。針對每種類型,本文隨機構建10個不同的實例,其中每個實例的工序數量取(1,θ)之間的隨機數。上述8種不同類型分別針對機器數量、工件數量、工序數量和加工時間大小這4個體現問題規模的參數出發,進行例舉和組合,從而分析算法在不同規模特征實例下的效果。

同時,本文DABC算法的相關參數設置如下:迭代次數Niter=100,種群數目SN=50,limit數目=10。為了防止執行算法出現棧溢出,設置深度優先搜索樹層數最大為20,CPU限制時間為3600 s。為了驗證算法的性能,比較對象是RFJSP問題最優解的下界。通過算法解CDABC與下界LBRFJSP之間的相對誤差Gap=CDABC-LBRFJSPLBRFJSP×100%來表征算法DABC的算法解相對問題最優解的性能。分別用Max_Gap、Min_Gap、Gap、Av_time分別表示在同一類實例中算法解目標函數的最大相對誤差、最小相對誤差、平均相對誤差、平均運行時間。算法求解不同類型實例的結果如表3—表5和圖4所示。

首先,對于小規模實例類型EX1,本文對隨機產生的10個實例,首先通過運用CPLEX求解器對問題RFJSP的整數規劃模型[16]進行求解,得到這些實例的最優解。同時,計算這些實例的最優解下界。分別計算本文算法解與最優解之間的相對誤差Gap1和與下界之間的相對誤差Gap2,具體結果如表3所示。從表3可知,對于這些實例,每個實例運行10次后取平均后的算法解呈現非常好的性能。相對于最優解而言,平均相對誤差都小于10%。事實上,對于這些實例,在每個實例的10次算法執行中都曾得到最優解,所以算法具有較好的全局尋優性。同時,從表3中可以看到與下界比較得到的相對誤差Gap2與Gap1比較,大了近30%,原因是問題的復雜性導致下界的估計相對最優解而言偏離較大。所以,下面對于規模較大實例的分析中,由于無法在短時間內得到最優解,所以本文只分析算法解與下界之間的相對誤差一定程度上也可以較好刻畫解的質量。

由表4中可知,算法對于每一類實例在工件加工時間較小情況下(且與運輸時間差距較小)的實例體現了較好的效果,與最優解下界比較得到的相對誤差都在50%左右。即使對于規模最大的兩個類型EX7和EX8,算法的性能較規模較小實例并未下降太多。事實上,由于算法解的比較對象是問題的下界,而問題下界比最優解小很多,所以算法解相對實際最優解的相對誤差將會遠小于50%(可以參考算法在相對規模較小實例中的表現,見表3)。對于算法執行的時間效率,本文算法尤為出色,對于小規模與中規模實例,算法都在10 min內得到滿意結果,對于大規模實例,算法也可以在1 h內得到滿意解。

從表5中可知,對于每一類實例在工件加工時間較大情況下,此時運輸時間相對于加工時間較小,算法同樣也體現了較好的求解效果。無論實例規模大小,工件加工時間的變大對算法的影響較小,與表4中體現的性能差別不大。表4和表5中的不同實例不僅體現了工件加工時間和問題規模對算法的影響,同時對加工時間與運輸時間差距不同情況下的算法分析,一定程度上也說明本文算法的魯棒性。當工件加工時間與運輸時間差距較小時,機器人運輸策略除了影響死鎖的產生以外,對最終的目標函數也會產生較大影響。當運輸時間遠小于工件加工時間時,機器人運輸策略對于目標函數值的影響較小,主要起到對工件加工順序的支持和死鎖的控制作用。從表4和表5中可以看到,本文算法受運輸時間大小的影響較小。

進一步對本文算法的收斂性能進行分析。通過對每一類實例算法的迭代次數與解的質量進行分析,如圖4所示。圖4(a)圖像針對的是工件加工時間ρ∈[1,10]情況下的實例類型EX1—EX4進行求解的結果,圖4(b)圖像針對的是工件加工時間ρ∈[10,100]情況下實例類型EX5—EX8的計算結果。從圖4(a)的迭代過程曲線可知,對于實例類型EX1—EX4,當算法迭代次數在45次的時候,算法即收斂到較好解。從圖4(b)的迭代過程曲線可知,對于較大規模實例類型EX5—EX8,算法在迭代80次左右時也收斂到較好滿意解。所以本文算法具有滿意的收斂速度,并在算法的尋優效果及收斂速度上達到了較好的平衡。

5 結 論

本文研究了一類依賴機器人運輸的柔性作業車間調度問題,針對該類在芯片制造等領域具有廣泛應用背景的實時調度問題,現有文獻只給出了問題的整數規劃模型,本文研究的目的是給出該問題的一個高效無死鎖調度算法。該研究的難點在于最優調度方案的搜索過程中,如何保證解的可行性即如何檢測并求解死鎖,并保證死鎖求解與算法尋優具有較好的收斂性。本文算法的基本思路是首先對問題進行時序Petri網建模,目的不是簡單使用傳統的Petri網進行死鎖仿真,而是給出解的一種適合死鎖判斷并求解的形式,同時方便進行解的解碼與迭代。其次,通過對問題死鎖產生的原因和類型進行分析,給出通用的死鎖判定和求解算法,并證明算法對于任意可能發生的死鎖都能在多項式時間內進行求解。最后,給出問題的蜂群算法求解方案。在求解過程中,對于產生的死鎖解,通過快速檢測和求解,在保證解的可行性情況下盡可能保留原有解的結構。通過理論分析和大規模數值實驗分析,得到本文的死鎖求解算法不僅能夠有效求解任意死鎖,同時極大地保留了工序之間的相對順序,使得蜂群算法在尋優過程中能夠更多保留優質解的特性,從而加快算法收斂。通過不同規模實例的數值實驗和與問題最優解及其下界的比較,也驗證了本文算法在不同輸入參數變化情況下都體現了較好的性能及收斂速度,特別是對于較大規模實例的求解,在收斂速度和解的質量上體系了較好的平衡,具有較好的實際應用前景。

本文結果是該類問題首個有效的算法結果。未來可以優化算法,進一步降低死鎖求解的時間復雜度,提高啟發式算法尋優的效率。此外,研究多機器人情形下問題的模型,并設計高效算法,將使本文研究更具現實意義和實用價值。

參考文獻:

[1]胡覺亮,李志林,董建明.工件加工需要額外資源的平行機調度問題[J].浙江理工大學學報(自然科學版),2022,47(5):791-797.

[2]Ali Chaudhry I, Ali Khan A. A research survey: review of flexible job shop scheduling techniques[J]. International Transactions in Operational Research, 2016, 23(3): 551-591.

[3]曹珍,韓曙光.考慮充電調度的電動無人車配送路徑規劃問題研究[J].浙江理工大學學報(自然科學),2023,49(6):784-794.

[4]Yan P Y, Liu S Q, Sun T F, et al. A dynamic scheduling approach for optimizing the material handling operations in a robotic cell[J]. Computers amp; Operations Research, 2018, 99: 166-177.

[5]Wei L X, He J X, Guo Z Y, et al. A multi-objective migrating birds optimization algorithm based on game theory for dynamic flexible job shop scheduling problem[J]. Expert Systems with Applications, 2023, 227: 120268.

[6]Liu S Q, Kozan E, Masoud M, et al. Job shop scheduling with a combination of four buffering constraints[J]. International Journal of Production Research, 2018, 56(9): 3274-3293.

[7]Bektur G, Sara T. A mathematical model and heuristic algorithms for an unrelated parallel machine scheduling problem with sequence-dependent setup times, machine eligibility restrictions and a common server[J]. Computers amp; Operations Research, 2019, 103: 46-63.

[8]Drobouchevitch I G, Neil Geismar H, Sriskandarajah C. Throughput optimization in robotic cells with input and output machine buffers: A comparative study of two key models[J]. European Journal of Operational Research, 2010, 206(3): 623-633.

[9]Li X L, Xing K Y, Lu Q C. Hybrid particle swarm optimization algorithm for scheduling flexible assembly systems with blocking and deadlock constraints[J]. Engineering Applications of Artificial Intelligence, 2021, 105: 104411.

[10]Xing K Y, Wang F, Zhou M C, et al. Deadlock characterization and control of flexible assembly systems with Petri nets[J]. Automatica, 2018, 87: 358-364.

[11]李浚,羅繼亮,李旭航,等.基于Petri網和監督學習的機器人柔性流水車間調度方法[J/OL].控制理論應用:1-9(2024-04-18)[2024-05-01]. http:∥kns.cnki.net.elib.zstu.edu.cn:80/kcms/detail/44.1240.tp.20240415.2248.004.html.

[12]Luo J C, Liu Z Q, Zhou M C. A Petri net based deadlock avoidance policy for flexible manufacturing systems with assembly operations and multiple resource acquisition[J]. IEEE Transactions on Industrial Informatics, 2019, 15(6): 3379-3387.

[13]Hurink J, Knust S. Tabu search algorithms for job-shop problems with a single transport robot[J]. European Journal of Operational Research, 2005, 162(1): 99-111.

[14]Nouri H E, Driss O B, Ghdira K. Hybrid metaheuristics for scheduling of machines and transport robots in job shop environment[J]. Applied Intelligence, 2016, 45(3): 808-828.

[15]晏曉凡.改進灰狼優化算法求解機器人作業車間調度問題[J].計算技術與自動化,2023,42(2):158-163.

[16]Sun Y G, Chung S H, Wen X, et al. Novel robotic job-shop scheduling models with deadlock and robot movement considerations[J]. Transportation Research Part E: Logistics and Transportation Review, 2021, 149: 102273.

[17]李大成, 羅繼亮, 孫莎莎, 等. 基于平行 Petri 網的制造系統調度與控制一體化方法[J]. 自動化學報, 2023, 49(4): 845-856.

[18]Bai D Y, Bai X Y, Li H R, et al. Blocking flowshop scheduling problems with release dates[J]. Swarm and Evolutionary Computation, 2022, 74: 101140.

(責任編輯:康 鋒)

主站蜘蛛池模板: 久久人体视频| 国产一区成人| 亚洲av日韩综合一区尤物| 中文字幕色在线| 亚洲欧洲日韩久久狠狠爱| 欧美日本激情| 精品国产aⅴ一区二区三区| 中文一区二区视频| 亚洲午夜福利精品无码不卡 | 欧美伊人色综合久久天天| 国产一级在线播放| 永久在线播放| 亚洲天堂日韩在线| 国内精品一区二区在线观看| www.youjizz.com久久| …亚洲 欧洲 另类 春色| 日韩国产综合精选| 人妻一区二区三区无码精品一区| 国产高清不卡| 国产一区在线视频观看| 香蕉国产精品视频| 日本久久久久久免费网络| 天堂成人av| 国产女人在线观看| yjizz视频最新网站在线| 欧美日本中文| 国产精品视频导航| 国产成人精品男人的天堂| 青青操国产视频| 国内精品久久九九国产精品| 亚洲国产第一区二区香蕉| 色综合日本| 国产乱人伦精品一区二区| 成人国产精品一级毛片天堂 | 区国产精品搜索视频| 国产视频只有无码精品| 亚洲第一香蕉视频| 国产chinese男男gay视频网| 一本色道久久88| 久久综合色播五月男人的天堂| 玖玖精品在线| 精品无码一区二区三区在线视频| 免费国产在线精品一区| 中文字幕永久在线看| 91青青草视频| 五月婷婷伊人网| 亚洲精品制服丝袜二区| 日韩性网站| 国产男人的天堂| 久久99国产视频| 成人福利一区二区视频在线| 国产福利拍拍拍| 国产九九精品视频| 无码'专区第一页| 欧美第二区| 精品亚洲欧美中文字幕在线看| 色哟哟色院91精品网站| 一级毛片不卡片免费观看| 中国特黄美女一级视频| 亚洲人成网址| 国产精品亚洲日韩AⅤ在线观看| 国产白丝av| 国产在线98福利播放视频免费| 国产精品视频观看裸模| 国产福利小视频高清在线观看| 欧美一区二区人人喊爽| 精品国产欧美精品v| 国产精品美乳| 亚洲一区无码在线| 午夜小视频在线| 日韩视频福利| 成年人福利视频| 伊人中文网| 色综合激情网| 国产精品永久久久久| 久久网欧美| 天天干天天色综合网| 无码AV日韩一二三区| 成人av手机在线观看| 精品黑人一区二区三区| 亚洲精品自在线拍| 国产欧美日韩在线在线不卡视频|