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

裝備平行仿真中演化建模框架研究*

2018-06-13 08:20:16葛承壟朱元昌邸彥強胡志偉
火力與指揮控制 2018年5期
關鍵詞:模型

葛承壟,朱元昌,邸彥強,胡志偉

(陸軍工程大學石家莊校區,石家莊 050003)

0 引言

裝備平行仿真技術是近年來軍用仿真領域的新興技術,不少國內仿真專家和學者提及裝備平行仿真技術的基本思想[1]并提出具體應用設想,已經成為軍用仿真領域的研究熱點[2],然而少有文獻明確給出裝備平行仿真的技術框架,可參考文獻等已有成果寥寥無幾。從裝備作戰指揮決策領域和裝備維修保障領域對軍用仿真技術的新需求出發介紹裝備平行仿真技術的概念。

在裝備作戰指揮決策領域中,作戰輔助決策系統采用固定的計算模型,不能根據戰場態勢進行模型自適應調整,影響指揮決策的準確性,利用仿真系統預測戰場態勢并為指揮員提供實時決策支持成為指揮控制系統的發展趨勢[3]。在裝備維修保障領域中,故障預測與健康管理(Prognostic and Health Management,PHM)技術是實現裝備視情維修的主要技術途徑,其中裝備剩余壽命(Remaining Useful Life,RUL)預測是 PHM[4]的主要研究內容之一。現有RUL預測方法大都是離線預測方法,且預測模型無法根據裝備信息進行模型更新,不具備自適應能力,使得RUL預測的時效性差、準確度不高,利用實時退化數據演化預測模型,實現具有自更新能力的RUL預測成為解決這些問題的有效途徑[5]。

以上兩個裝備運用領域中遇到的突出問題可以歸納為領域模型固定,不具備自適應演化能力,要求研究具有模型演化能力的軍用仿真技術即裝備平行仿真技術。文獻[6-7]提出了裝備平行仿真技術的概念和技術內涵,裝備平行仿真技術是指實際武器裝備和平行仿真系統通過特定的接口設備連接在一起,平行仿真系統通過傳感器等設備以在線的方式實時獲取裝備信息,用于演化裝備參考模型,提高仿真結果的可信度,仿真結果通過執行器反饋給實際武器裝備,從而提高武器裝備的運用效能和保障效能。裝備平行仿真示意圖如圖1所示。在空間上,實際武器裝備和與之對應的平行仿真系統構成裝備平行系統,二者之間是雙向交互、互利共生的關系。裝備平行仿真可用于裝備指揮決策和維修保障等領域中,解決領域模型固定、不具備自適應演化能力這一問題。裝備平行仿真的理論緣起與平行系統理論[8]、共生仿真[9]、在線仿真[10]和動態數據驅動應用系統(Dynamic Data Driven Application System,DDDAS)[11]等理論范式密切相關,文獻[7]已作詳細評述。作為新興仿真技術,裝備平行仿真具有虛實共生、雙向交互、模型演化、高效運行等特點,其中模型演化是其在建模方法上區別于以往仿真技術的主要所在,但是至今仍缺乏完整的技術框架。

圖1 裝備平行仿真示意圖

1 裝備平行仿真中的演化建模問題

1.1 裝備平行仿真中的仿真模型與模型演化

利用不同的建模方式可以得到刻畫裝備不同特性的模型,設集合M為所有可用于描述裝備特性的模型集,Θ為當前時刻的仿真模型。隨著裝備狀態、行為的變化,為提高平行仿真系統的仿真可信度,需要根據當前裝備狀態和行為信息演化當前模型Θ,以適應裝備狀態和行為的變化。對于特定的仿真應用和裝備狀態來說,存在適宜模型集以高逼真度仿真實際武器裝備的特定狀態和行為,稱此模型集為參考模型集MR,有,是否成立取決于時間,隨著裝備運行狀態和行為的變化,平行仿真系統需要對Θ進行演化修正甚至模型更替,當前模型的演化修正、新模型加入和舊模型剔除使得MR是時變模型集。模型演化并沒有嚴格定義,參考已成熟的軟件動態演化的概念[12],本文認為裝備平行仿真中的模型演化是指平行仿真系統為適應動態變化的裝備狀態、滿足不斷變化的仿真需求而對仿真模型進行的自適應調整,并從仿真模型的一般組成對模型演化進行細化。一般來講,完整的仿真模型Θ由模型輸入input、模型參數模型輸出output構成,特定建模方式下三者的組合決定了模型結構,如圖2所示。故模型演化可以分為模型輸入演化、模型參數演化、模型輸出演化和模型結構演化,模型輸入演化是指減少模型的冗余輸入[13],以優化模型結構;模型參數演化是指模型參數的更新/優化[14];模型輸出演化是指模型輸出的自適應校正[15];模型結構演化是指模型的適宜性選擇/動態更替[16]。由于在特定仿真應用中,模型輸入的數據源往往是固定的,因此,這里僅考慮模型參數演化、模型輸出校正和模型適宜性選擇。本文主要關注面向裝備RUL預測的平行仿真中演化建模問題。

圖2 仿真模型的一般表示示意圖

1.2 面向裝備RUL預測的平行仿真中演化建模框架

裝備平行仿真中的仿真模型是能夠反映裝備某一特性的模型,不同種類的裝備平行仿真中模型的內涵不同,在面向裝備RUL預測的平行仿真中首先應該解決仿真模型定位即基礎預測模型的問題。

1.2.1 基礎預測模型

仿真模型是仿真執行的基礎,面向裝備RUL預測的平行仿真中仿真模型屬于預測模型。受裝備復雜度和裝備運行環境等的影響,實際武器裝備狀態具有明顯的不確定性、時變性和動態性,因此,實際武器裝備可以看作動態系統。實際武器裝備由正常工作到發生功能性故障的過程是性能逐漸退化的過程,且描述動態系統常用的方法是建立系統的狀態空間模型(State Space Model,SSM),所以基于SSM的裝備性能退化建模是面向裝備RUL預測的平行仿真中的建模方向。SSM的一般形式包括狀態方程和觀測方程,即

其中,xk為狀態變量,zk為外界觀測,和分別是均值為0、(協)方差為Qk-1和Rk且相互獨立的過程噪聲和觀測噪聲,為有界映射。

1.2.2 演化建模框架

面向裝備RUL預測的平行仿真中演化建模框架如圖3所示,模型演化主要包括3個階段:首先,平行仿真系統根據裝備實時退化狀態和模型在線選擇算法從模型庫中在線選擇適宜裝備退化SSM,作為RUL基礎預測模型;其次,平行仿真系統利用數據同化算法對裝備退化狀態進行跟蹤,實現對裝備退化SSM的輸出校正,得到更為準確的裝備退化后驗狀態;平行仿真系統根據裝備當前和歷史退化狀態,利用參數在線估計算法對預測模型中的未知參數(包括超參數)進行準確實時估計;最后,平行仿真系統利用演化后的裝備退化SSM,結合裝備失效閾值,外推裝備退化狀態到達失效閾值的時間,并經過多次Monte-Carlo仿真得到裝備RUL及其概率密度函數。

圖3 面向裝備RUL預測的平行仿真中演化建模框架

面向裝備RUL預測的平行仿真中演化建模是一種裝備退化數據驅動的演化建模,實現裝備退化狀態的在線、實時、準確感知是完成演化建模的前提。裝備退化狀態感知涉及傳感器技術、數據預處理技術、虛實接口技術等,在此不作展開討論,重點研究討論模型在線適宜性選擇技術、數據同化算法和參數在線估計算法。

2 數據驅動的模型在線適宜性選擇

與以往仿真技術不同,在裝備平行仿真中將模型的驗證過程與模型構建過程一同處理,模型的構建確保了仿真模型的可信度。裝備平行仿真中的模型構建屬于數據驅動的模型構建過程,平行仿真系統能根據裝備當前和歷史狀態數據以在線形式從模型庫中調用適宜模型,其本質屬于模型適宜性選擇過程。裝備性能退化過程受復雜工況、人為因素和其他隨機因素等的影響,具有時變性強、波動大、多階段退化等特點,在面向裝備RUL預測的平行仿真中,根據裝備性能退化軌跡的特征,可以將裝備性能退化SSM劃分為線性退化SSM和非線性退化SSM,模型適宜性選擇指的平行仿真系統根據裝備當前和歷史性能退化數據從兩類退化SSM中動態選擇適宜模型。根據平行仿真系統模型輸入維數(特征退化量的數量)的不同,提出兩種模型適宜性選擇方法。

2.1 單維退化數據情況下模型在線適宜性選擇

部分裝備的性能退化可用單一退化特征量表征,如用漂移特征量表征導彈裝備中陀螺儀的性能退化過程[17]、用振動信號的均方根值描述軸承的性能退化過程[18]等。與DDDAS類似,裝備平行仿真能將模型構建與模型驗證一并處理[19],因此,可將模型適宜性選擇過程看作多仿真模型的在線驗證,并提出時頻域聯合分析方法對裝備性能退化SSM進行適宜性選擇,即在時域、頻域中分別采用灰色關聯度分析、最大熵譜估計對仿真模型進行選擇。

灰色關聯度分析[20]是模型驗證領域中的經典時域分析方法,屬于幾何分析范疇,它通過分析比較仿真模型輸出數據和參考數據變化趨勢的相似程度即關聯度來定量分析二者的一致性,關聯度越大,模型越為適宜,反之模型越不適宜。設參考數據序列即裝備性能退化真實數據為,并設有m個用于仿真裝備性能退化的SSM,第i個仿真模型輸出序列為,其中。仿真模型輸出序列xi與裝備退化數據x0在時刻k處的差異可以表示為

稱為關聯系數,ζ為分辨系數且,ζ可減小序列極值對二者關聯度的影響。xi和x0的關聯度定義為

對式(3)取最大值有,關聯度最大者對應的模型即為當前退化狀態下的適宜SSM。

最大熵譜估計是常用的現代譜估計方法,克服了古典加窗譜分析精度低、能量泄露等不足,能在頻域中對仿真模型進行有效驗證,具體算法參見文獻[21]。利用灰色關聯度分析和最大熵譜估計對仿真模型進行適宜性選擇,是分別在時域和頻域兩個不同角度得到的選擇結果,存在選擇結果沖突的情況,可利用D-S證據理論中沖突理論合成方法對二者的選擇結果進行融合[22],得到最終的仿真模型選擇結果。

2.2 多維退化數據情況下模型在線適宜性選擇

武器裝備性能退化過程的復雜性和不確定性決定了存在表征性能退化過程的多特征量[23]。多特征量問題集中表現在特征量多、特征量間相關,使得特征量之間存在信息重疊、數據維數高。為降低RUL預測的計算復雜度,需對多退化特征量進行數據降維,得到反映原有多退化特征量信息的綜合退化特征量。核主成分分析(Kernel Principle Component Analysis,KPCA)是常用的數據降維方法,相較于主成分分析(Principle Component Analysis,PCA)方法,KPCA引入核方法,能夠保留更多的原始信息,更適宜處理多維退化數據的特征提取[24]。設高維退化特征量為,定義非線性映射函數用于實現從輸入空間向特征空間(再生核Hilbert空間)的變換,即,輸入空間樣本xb、xc在特征空間中的距離用點積表示,定義核函數K,使得

其中,高斯核函數是常用核函數,即

σ為核帶寬參數,非線性變換完畢后,在特征空間執行PCA算法得到主成分。模型在線適宜性選擇的本質是模式識別問題,即多維退化數據對退化SSM的選擇。誤差反向傳播神經網絡(Back Propagation Neural Network,BPNN)是常用的模式識別方法。針對多維退化數據情況下模型在線適宜性選擇問題,提出一種基于KPCA和BPNN的模型適宜性選擇方法,如圖4所示。

首先對M維退化特征量采用KPCA算法得到維相互獨立的主成分;L維主成分經歸一化處理后,得到BP神經網絡的輸入,經模式識別后,得到適宜的裝備性能退化SSM。模式識別前,需要利用裝備性能退化歷史數據對BP神經網絡進行訓練,得到輸入層qj和隱藏層yi之間的權值矩陣Wji、隱藏層yi和輸出層ol之間的權值矩陣Til,并將最新的識別過程輸入輸出數據動態加入到神經網絡訓練數據庫中,提高模式識別的準確度。

圖4 多維退化數據情況下基于KPCA和BP神經網絡的模型適宜性選擇

3 基于數據同化的模型輸出校正

數據同化(Data Assimilation,DA)緣起于數值天氣預報(Numerical Weather Prediction,NWP)領域,指模型輸出數據與觀測數據變得相近的過程[25],是聯系模型預測和觀測數據的核心部分。在面向裝備RUL預測的平行仿真中,借鑒數據同化的思想,將最新的裝備性能退化數據引入裝備退化SSM中,利用順序數據同化算法校正預測模型輸出,提高RUL預測精度。順序DA算法又稱濾波算法,常用的算法包括卡爾曼濾波(Kalman Filter,KF)及其改進算法[26]、粒子濾波(Particle Filter,PF)及其改進算法[27]、層次Bayesian方法等,這里結合RUL預測需求重點研究KF和PF算法。

3.1 卡爾曼濾波算法

KF是典型的最小方差估計方法,最優KF問題就是給定裝備性能退化數據觀測序列,找到性能退化特征量即系統狀態 xk+1的最優線性估計,使得估計誤差的方差最小。將式(1)具體化,裝備性能退化離散狀態空間模型可表示為

其中,Fk和Hk是系統矩陣和觀測矩陣;狀態初值一般取,狀態協方差P0可任意假定。

以式(6)為基礎預測模型,應用KF算法對裝備性能退化過程進行跟蹤預測。KF算法分為預測和更新兩個步驟:預測階段利用當前時刻k的狀態及其協方差估計Pk|k得到下一時刻k+1二者的先驗估計;更新階段利用k+1時刻最新觀測數據zk+1和先驗估計得到狀態xk+1的后驗估計。KF算法可參見文獻[26]。

由于KF算法要求SSM是線性模型且w和v為相互獨立的高斯白噪聲,因此,KF算法適宜于退化過程屬于線性高斯的武器裝備。為提高KF算法的適用性,國內外學者分別提出了擴展卡爾曼濾波(Extended Kalman Filter,EKF)、集合卡爾曼濾波(Ensemble Kalman Filter,EnKF)、無跡卡爾曼濾波(Unscented Kalman Filter,UKF)等改進方法,不再展開研究。

3.2 粒子濾波算法

實際武器裝備退化狀態的變化過程大都是非線性的,不適宜應用KF算法,對應的退化SSM是非線性非高斯模型。考慮非線性SSM可表示為式(1),此時和分別是均值為0、協方差為Qk+1和Rk且相互獨立的過程噪聲和觀測噪聲,為有界非線性映射。對于非線性非高斯的裝備性能退化狀態估計問題,以Monte-Carlo方法和序貫重要性采樣(Sequential Importance Sampling,SIS)為基礎的粒子濾波算法具有較大優勢。PF算法通過對來自于裝備退化狀態概率密度函數的采樣集進行預測和更新得到裝備退化狀態的后驗Bayesian估計。PF算法可參見文獻[27]。

4 基于參數在線估計的模型參數演化

在式(1)的裝備退化SSM中,存在未知參數Q、R,當基于隨機過程如Wiener過程、Gamma過程等進行SSM建模時還存在與這些分布相關的超參數,這些參數構成了SSM的未知參數集。平行仿真系統利用模型參數在線估計算法對未知參數集Θ進行在線估計,使得參數值隨著裝備退化狀態的變化進行自適應調整,從而提高裝備退化狀態預測的準確性,可選的算法包括期望最大化(Expectation Maximization,EM)算法[28-29]、馬爾可夫鏈蒙特卡羅(Markov Chain Monte Carlo,MCMC)方法[30]等。

4.1 期望最大化算法

EM算法是一種迭代算法,主要用于求取后驗分布的極大似然估計,能有效估計具有隱含狀態SSM中的未知參數,它的每一次迭代由兩步組成,即E步(求期望)和M步(極大化)。一般地,以表示Θ的基于裝備退化觀測數據z的后驗分布密度函數,稱為觀測后驗分布,表示添加隱含退化狀態x后得到的關于Θ的后驗分布密度函數,稱為添加后驗分布,表示在給定Θ和退化觀測數據z下隱含退化狀態x的條件分布密度函數。為計算觀測后驗分布的極大似然估計值,記為第i次迭代開始時的極大似然估計值,則第i+1次迭代的兩步為:

E步:將記為關于x的條件分布并求數學期望,從而把隱含退化狀態x利用積分運算消掉,即

M步:將極大化,即找到新參數,使

如此形成了一次迭代。將上述E步和M步進行迭代直至或時停止迭代。EM算法能與數據同化算法聯合使用,對裝備退化狀態、參數進行聯合估計,即在E步可利用不斷更新的裝備退化觀測數據z和DA算法更新裝備退化狀態,在M步得到SSM中未知參數的估計值,進而將經模型演化獲得的退化狀態和SSM未知參數用于RUL預測。

4.2MCMC方法

EM得到是參數后驗分布的極大似然估計值,對于Θ的后驗均值、后驗方差等的獲取需要利用MCMC方法,尤其是在多退化特征量情況下,參數集Θ中參數的聯合分布是復雜高維且非標準形式的分布,此時EM算法難以實施。

MCMC方法屬于統計試驗方法,它從一個目標概率分布中抽樣,得到收斂、平穩的一定容量樣本,由于各個參數的樣本值分布近似于其邊緣分布,故把各個參數樣本值帶入定積分近似計算公式中就得到該參數的數字特征。通過MCMC方法可得到參數后驗聯合分布,在Bayesian框架下,目標概率分布就是參數的后驗聯合分布。MCMC算法可參見文獻[30]。

5 結論

作為新興軍用仿真技術,裝備平行仿真能有效滿足裝備指揮決策領域和維修保障領域中對仿真技術的需求,其中數據驅動的模型演化是其主要技術內涵,通過討論裝備平行仿真中的演化建模問題,明確了演化建模內涵,并建立了面向裝備RUL預測的平行仿真中演化建模框架,能實現對裝備RUL的實時、在線、自適應預測。

[1]胡曉峰.大數據時代對建模仿真的挑戰與思考[J].軍事運籌與系統工程,2013,27(4):5-12.

[2]李曉婉,竇林濤,程健慶.平行仿真技術在指控系統中的應用構想[C]//2015全國仿真技術學術會議論文集,2015:36-39.

[3]王飛躍.指控5.0:平行時代的智能指揮與控制體系[J].指揮與控制學報,2015,1(1):107-120.

[4]SCANFF E,FELDMAN K L,GHELAM S,et al.Life cycle cost impact of using prognostic health management(PHM)for helicopter avionics [J].Microelectronics Reliability,2007(47):1857-1864.

[5]張仕新,昝翔,李浩,等.狀態維修理論及剩余壽命預測的研究現狀與展望[J].兵工自動化,2014,33(9):15-20.

[6]葛承壟,朱元昌,邸彥強,等.裝備精確維修平行仿真系統及關鍵技術研究[J].現代防御技術,2016,44(6):162-168.

[7]GE C L,ZHU Y C,DI Y Q,et al.Equipment residual useful life prediction oriented parallel simulation framework[C]//Proceedings of AsiaSim/SCS AutumnSim 2016,Part I,CCIS 643:377-386.

[8]王飛躍.平行控制:數據驅動的計算控制方法[J].自動化學報,2013,39(4):293-302.

[9]AYDT H,CAI W,TURNER S J,et al.Symbiotic simulation for optimization oftooloperations in semiconductor manufacturing[C]//Proceedings of the Winter Simulation Conference,2011:2093-2104.

[10]吳金平,陸銘華.指揮控制系統中仿真決策系統的嵌入研究[J].火力與指揮控制,2015,40(6):84-87.

[11]DAREMA F.Dynamic data driven applications systems:new capabilities for application simulations and measurements[C]//Proceedings of the International Conference on Computational Science2005, Berlin,Springer-Verlag,2005,LNCS 3515:610-615.

[12]李青山,王璐,褚華,等.一種基于智能體技術的軟件自適應動態演化機制[J].軟件學報,2015,26(4):760-777.

[13]郭俊,周建中,王浩,等.系統理論水文模型結構與參數多目標優化[J].水力發電學報,2014,33(2):1-7.

[14]侯景偉,孔云峰,孫九林.蟻群算法在需水預測模型參數優化中的應用[J].計算機應用,2012,32(10):2952-2955.

[15]馬建文,秦思嫻.數據同化算法研究現狀綜述[J].地球科學進展,2012,27(7):747-757.

[16]楊廣斌,劉鵬舉,唐小明.動態數據驅動的林火蔓延模型適宜性選擇[J].林業科學,2011,47(1):107-112.

[17]馮磊,王宏力,周志杰,等.基于狀態空間的慣性測量組合剩余壽命在線預測[J].清華大學學報(自然科學版),2014,54(4):508-514.

[18]MORTADA M.Diagnosis of rotor bearings using logical analysis of data [J].Journal of Quality in Maintenance Engineering,2011,17(4):371-397.

[19]周云,黃柯棣,胡德文.動態數據驅動應用系統的概念研究[J].系統仿真學報,2009,21(8):2138-2141.

[20]寧小磊,吳穎霞,于天明,等.基于改進灰色關聯分析的仿真模型綜合驗證方法[J].兵工學報,2016,37(2):338-347.

[21]王建華,符文星,董敏周,等.最大熵譜估計在空空導彈仿真模型驗證中的應用[J].彈箭與制導學報,2005,25(4):848-850.

[22]張鑫,牟龍華.基于局部沖突消除的證據合成法則[J].系統工程與電子技術,2016,38(7):1594-1599.

[23]王小林,郭波,程志君.基于非線性漂移Wiener過程的產品實時可靠性評估[J].中南大學學報,2013,44(8):3203-3209.

[24]SCHOLKOPF B,SMOLA A J,MULLER K R.Nonlinear component analysis as a kernel eigenvalue problem[J].Neural Computation,1998,10(5):1299-1319.

[25]RODELL M,HOUSER P R,JAMBOR U,et al.The global land data assimilation system [J].Bulletin of the American Meteorological Society,2004,85(3):381-394.

[26]REICHLE R H,WALKER J P,KOSTER R D,et al.Extended versus ensemble Kalman filtering for land data assimilation[J].Journal of Hydrometeorology,2002,3(6):728-740.

[27]畢海蕓,馬建文.粒子濾波算法在數據同化中的應用研究進展[J].遙感技術與應用,2014,29(5):701-710.

[28]茆詩松,王靜龍,濮曉龍.高等數理統計[M].北京:高等教育出版社,1998.

[29]DEMPSTER A P,LAIRD N M,RUBIN D B.Maximum likelihood from incomplete data via the EM algorithm[J].Journal of the Royal Statistical Society:Series B,1977,39(1):1-38.

[30]王浩偉,徐廷學,劉勇.基于隨機參數Gamma過程的剩余壽命預測方法[J].浙江大學學報(工學版),2015,49(4):699-704,762.

猜你喜歡
模型
一半模型
一種去中心化的域名服務本地化模型
適用于BDS-3 PPP的隨機模型
提煉模型 突破難點
函數模型及應用
p150Glued在帕金森病模型中的表達及分布
函數模型及應用
重要模型『一線三等角』
重尾非線性自回歸模型自加權M-估計的漸近分布
3D打印中的模型分割與打包
主站蜘蛛池模板: 国产精品三级专区| 尤物国产在线| 欧美激情视频一区二区三区免费| 亚洲天堂免费| 国产一区二区三区在线观看免费| 久久久久人妻精品一区三寸蜜桃| 91www在线观看| 亚洲欧美成aⅴ人在线观看| 少妇精品久久久一区二区三区| 欧美亚洲第一页| 国产午夜看片| 成人另类稀缺在线观看| 久久青草视频| 国产一级二级在线观看| 亚洲制服中文字幕一区二区| 欧美激情第一区| 四虎影视国产精品| 久久亚洲中文字幕精品一区| 福利小视频在线播放| 亚洲欧美综合另类图片小说区| 偷拍久久网| 成人国产精品网站在线看| 欧美综合成人| 久久久久亚洲精品成人网| 欧美一级特黄aaaaaa在线看片| 黄色网址免费在线| 91毛片网| 国内精品小视频福利网址| 色视频国产| h网址在线观看| 国产欧美精品专区一区二区| 国产全黄a一级毛片| 精品国产美女福到在线直播| 婷婷色中文| 99热亚洲精品6码| 在线观看国产黄色| 欧美国产综合色视频| 成人国产精品2021| 亚洲性一区| 国产精品第三页在线看| 久久伊人操| 亚洲最大福利视频网| 久久精品人妻中文系列| 久一在线视频| 亚洲品质国产精品无码| 毛片一级在线| 首页亚洲国产丝袜长腿综合| 亚洲综合经典在线一区二区| 青青网在线国产| 91午夜福利在线观看精品| 国产福利小视频在线播放观看| a天堂视频在线| 无码视频国产精品一区二区 | 精品成人一区二区| 免费观看无遮挡www的小视频| 久久国产毛片| 91精品视频网站| 啪啪啪亚洲无码| 成人午夜精品一级毛片| 97国内精品久久久久不卡| 亚洲国产理论片在线播放| 毛片网站免费在线观看| 久久久久无码国产精品不卡 | 欧美影院久久| 五月天在线网站| 国产精品思思热在线| 亚洲日本一本dvd高清| 美女内射视频WWW网站午夜 | 国产又爽又黄无遮挡免费观看| 亚洲日本精品一区二区| 欧美成人一级| www.亚洲天堂| 六月婷婷激情综合| 亚洲女同一区二区| 992Tv视频国产精品| 91无码网站| 98超碰在线观看| 亚洲无线国产观看| 51国产偷自视频区视频手机观看| 性做久久久久久久免费看| 青草视频网站在线观看| 欧美三级视频网站|