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

基于D-S證據理論的航天設備壽命預測方法

2016-11-09 09:43:30丁瑞陸寧云程月華姜斌邢琰
中國空間科學技術 2016年4期
關鍵詞:分配理論融合

丁瑞,陸寧云,*,程月華,姜斌,邢琰

1.南京航空航天大學 自動化學院,南京 211016 2.南京航空航天大學 航天學院,南京 211016 3.北京控制工程研究所,北京 100191 4.空間智能控制技術國家重點實驗室,北京 100191

?

基于D-S證據理論的航天設備壽命預測方法

丁瑞1,陸寧云1,*,程月華2,姜斌1,邢琰3,4

1.南京航空航天大學 自動化學院,南京 211016 2.南京航空航天大學 航天學院,南京 211016 3.北京控制工程研究所,北京 100191 4.空間智能控制技術國家重點實驗室,北京 100191

以動量輪為研究對象,提出一種基于D-S證據理論和Bayes理論的信息融合及壽命預測方法。首先,挖掘多個壽命信息源的內、外部信息作為Bayes多源驗前信息,并構建證據集合;其次,利用合成規則對證據集合進行合成,得到合理的驗前融合權重的分配結果;然后,利用Bayes方法求解壽命參數的融合驗后分布,并計算壽命參數的估計值;最后,根據參數估計值預測設備的剩余壽命。仿真結果表明,基于多源壽命信息融合的預測結果與單一來源的壽命預測結果相比,更加接近于設備的真實壽命信息。

壽命預測;多源信息融合;D-S證據理論;Bayes方法;動量輪

航天設備運行機理復雜,運行環境惡劣,設備性能退化過程中容易導致惡性事故發生。為保障航天器安全可靠運行,確保航天任務順利實施,準確的設備剩余壽命預測方法是航天工程中的迫切需求。航天設備具有典型的小子樣特點,不適合使用傳統的基于大樣本壽命數據的預測方法[1-2]。然而,航天設備在設計、研制、試驗和使用等階段均存在不同來源、不同環境及不同層次的壽命信息,對這些多源壽命信息進行合理融合,將有助于得到更準確的設備剩余壽命預測結果。

信息融合的最初定義是指,對從單個或多個傳感器信息源獲取的數據和信息進行關聯、相關和綜合,從而對態勢和威脅及其重要程度進行全面及時評估的信息處理過程[3]。廣義上,融合就是將來自多信息源的信息和數據進行綜合處理,從而得出更為準確可靠的結論[4]。D-S證據理論、Bayes方法、模糊測度積分和專家系統等均是基于信息融合的可靠性評估技術中的常見方法。其中,D-S證據理論因其可以靈活處理不同可靠性數據和信息的形式差異,并具有直接表達不確定未知信息的能力,受到越來越多的重視和應用。

D-S證據理論是一種基于“證據”和“組合”來處理不確定性推理問題的數學方法[5],已被嘗試應用于故障診斷和健康管理的相關領域[6-7]。譬如,文獻[8]研究了基于D-S證據理論的旋轉機械綜合故障診斷方法;文獻[9]將證據理論與神經網絡相結合,旨在解決航空發動機試車臺試驗中發動機磨損的故障診斷問題;文獻[10]將D-S數據融合技術應用于多Agent的衛星故障診斷中,提高了診斷的準確性和多Agent系統的智能性;文獻[11]提出一種基于粗糙集理論和D-S證據理論的設備技術狀態評估模型。

針對航天設備的多源壽命信息,本文將D-S證據理論和Bayes理論相互結合,提出一種基于多源信息融合的壽命預測方法。首先,針對所研究的對象收集多種來源的壽命信息,包括專家經驗、性能退化數據、相關(或相似)設備的壽命信息等;其次,將這些可靠性信息作為壽命預測的驗前信息,挖掘其內外部信息,構建證據集合;最后,運用D-S證據合成規則對多證據的基本概率分配(Basic Probability Assignment,BPA)進行合成,將合成后的BPA作為各驗前分布的融合權重分配結果,最終計算研究對象剩余壽命的驗后分布曲線。仿真試驗結果表明了本文方法的可行性和合理性。

1 D-S證據理論及Bayes可靠性評估

本文的主要目的是,將設備的多源壽命信息視作Bayes驗前信息,且驗前分布形式為威布爾分布,對分布中的尺度參數br進行基于D-S證據理論的融合估計求解。首先運用合成規則確定各驗前信息的融合權重,得到尺度參數br的融合驗前分布;再根據Bayes理論計算參數的驗后分布估計值,并用估計值確定設備壽命的威布爾分布函數,最終得到設備壽命的預測值。本節將簡要介紹D-S證據合成規則和Bayes驗前信息融合的基本思想。

1.1證據合成的基本概念

D-S合成規則反映的是多組證據聯合作用的結果。給定幾組同一識別框架下的基于不同證據的基本概率分配函數,如果這幾個證據不是完全沖突的,那么就可以利用D-S合成規則,計算出多證據聯合作用下的基本概率分配函數[13]。

(1)

在基于信息融合的航天設備壽命預測問題中,焦元A1,A2,…,Ai和B1,B2,…,Bj代表多個壽命信息來源,m(C)是指經多個證據合成后分配到各信息源上的基本可信數。

1.2Bayes驗前信息融合

本文選用威布爾分布作為航天設備剩余壽命的分布函數,因為它是航空航天設備常用的一種失效分布,尤其適用于機電類產品磨損累計失效的分布描述[14]。

假設設備的剩余壽命T服從威布爾分布

(2)

(3)

式中:Γ(?)為伽馬函數;Q′(br)為Q(br)的一階導數;αi和βi為參數,i表示第i個信息源。

假定現有k個關于尺度參數br的驗前信息源,驗前分布分別為πi(br)(i=1,2,…,k),若得到各信息源的權重wi,則融合后br的驗前分布為π(br)=wiπi(br);另外假定通過現場試驗已獲得該設備的n個壽命數據,即真實壽命信息的集合,記為T={T1,T2,…,Tn};由指數分布的共軛性質可以計算尺度參數br的驗后分布

(4)

進而,根據最小方差估計原理可以求出尺度參數br的最小方差估計值

(5)

后文中將應用D-S證據理論中的基本概率分配和證據合成,確定各信息源的權重wi;與其他權重確定算法相比,本文方法可以統一處理形式各異的各類壽命信息,使得權重確定過程更加合理、準確。

2 基于D-S證據理論的壽命信息  融合方法

由于航天設備很難獲得大樣本數據,本文將充分挖掘專家經驗、遙測性能參數、相似系統壽命數據等驗前信息,形成有關設備壽命的證據集合,運用D-S證據合成規則進行多個證據的BPA合成,作為各驗前分布的融合權重分配,并利用權重得出驗前分布融合結果,結合現場真實分布信息,最終計算得出研究對象壽命的Bayes驗后分布曲線,具體實現流程如圖1所示。

圖1 D-S方法流程Fig.1 Flowchart of D-S method

2.1驗前分布的信息挖掘及證據構成

本文將設備壽命信息劃分為內部信息和外部信息[14]。內部信息是指設備壽命數據本身提供的信息,如數據容量、分布特征度;外部信息是壽命數據采集過程中的一些影響因素,如數據采集方法和水平、數據處理精度及信息來源的可信度。

(1)內部信息證據

不同來源壽命信息所包含的數據量大小存在明顯差異,而數據量大小是影響數據驅動預測方法性能的主要因素,因此第一個D-S證據由數據容量特性構成。本文采用計算各信息源數據量的百分比作為數據容量基本概率分配(BPA)函數m1的構造方法。譬如,某設備有4種壽命信息源S1~S4,數據容量分別為130、400、255和140,則通過計算各容量占總容量的百分比可得數據容量證據E1的BPA結果為

E1={m1(S1)=0.141,m1(S2)=0.432,

m1(S3)=0.276,m1(S4)=0.151}

除了數據容量,不同壽命信息源的驗前分布參數與現場真實分布參數之間的差異,也可度量該信息源的融合權重大小,因此第二個D-S證據為參數分布證據E2某設備4種信息源S1~S4的壽命分布參數θ的驗前分布分別為π1(θ)=guass(4.9,1),π2(θ)=guass(5.2,1),π3(θ)=guass(4.2,1),π4(θ)=guass(4.5,1),均服從高斯分布;假定由現場試驗數據得到的θ真實分布為π(θ)=guass(5,1),下面來推導E2的BPA函數m2。

(6)

(2)外部信息證據

不同信息來源的壽命數據一般有著不同的采集方法和技術水平、數據處理的精度也不盡相同。此外,不同的信息來源也具有不同的可信度,這些外部信息也應當作為證據納入D-S證據理論的證據合成中。

例如,文獻[14]在分析某型柴油機的累積失效率和故障里程數之間的關系時,4種來源的故障里程數信息在數據質量、預處理方法和傳感器精度等方面存在差異,經過數位專家的評估分析后,對信息源S1~S4賦予了不同的可信度,分別為0.995、0.987、0.926和0.912,該論文中得到的基于外部信息的可信度證據E3的BPA函數m3結果為

E3={m3(S1)=0.260 5,m3(S2)=0.258 4,

m3(S3)=0.242 4,m3(S4)=0.238 7}

同理,若獲得了其他可量化的外部信息,也可以根據實際情況轉化為融合權重分配這同一框架下的相應證據集合,加入影響信息源驗前分布融合權重的考慮因素中。

2.2D-S證據理論融合過程模型

利用證據理論的合成規則(即式(1)),對第2.1節所提出的多組證據集合進行融合,可得到各信息源驗前分布的權重分配集合。基于內部信息的各組證據的BPA融合形式為Wj=E1⊕E2,j,其中E2是分k個區間依次求解得到的,故E2有k個值。

由第1.1節可知,證據合成的次序對最終結果沒有影響。若仍需繼續考慮其他內外部信息的證據集合,權重計算公式可以擴展為通用形式

(7)

并利用得到的權重分配矩陣Wj(k×n),其中,k表示求取參數分布證據E2時所劃分的區間數,n表示驗前信息來源的個數。結合式(3),可得到k組融合驗前分布πj(br)的集合

(8)

總結上述步驟,基于D-S證據理論的壽命信息融合過程如圖2所示。其中,n表示證據數目,t表示信息源數目,mj(Si)表示在同一證據Ej的基本概率分配下的第i個信息源的基本可信數,i=1,2,…,t,j=1,2,…,n。通過D-S合成規則,各證據下的BPA函數mj進行合成,得到證據組合后的基本可信數m(Si),即為組合后的各信息源置信度大小。進一步將置信度作為融合權重,至此,最終得到各信息源的驗前融合權重Wi。

圖2 基于D-S證據理論的融合過程模型Fig.2 Fusion process model based on D-S evidence theory

3 算法驗證及分析

動量輪是衛星三軸穩定控制系統中的主要執行部件,動量輪的壽命直接影響衛星的使用壽命,實現動量輪剩余壽命的準確預測,對提高衛星可靠性進而延長衛星使用壽命具有重要意義[15]。

本文以XX型動量輪為研究對象,采用文獻[16]中挖掘到的3種信息源S1~S3的驗前分布信息和一組現場試驗壽命數據分布信息,進行基于D-S證據理論的壽命信息融合方法的仿真驗證。這里選擇威布爾分布作為動量輪剩余壽命的分布形式,對尺度參數br的估計值進行融合求解。其中,3種來源的驗前信息分別為相似型號動量輪的壽命數據、基于在軌溫度性能的遙測數據推導得出的偽壽命數據、基于專家經驗的壽命分布參數的區間估計。

3.1驗前壽命信息挖掘

表1 現場真實分布及各驗前分布參數表

圖3 尺度參數各驗前分布的概率密度曲線Fig.3 Prior distributions curves of scale parameter

3.2證據構成及融合

已知3種來源S1~S3的驗前信息的原始數據量分別為13、5和10,通過計算各信息源數據量的百分比,得到數據容量證據E1的BPA,結果如表2示??傻?/p>

E1={m1(S1)=0.4643,

8m1(S2)=0.1786,m1(S3)=0.3571}

表2 信息源數據容量特征表

為了分析各信息源和現場試驗尺度參數的分布情況,對πi(br)取積分,得到尺度參數br的累積分布函數,為了獲得精確的參數分布證據E2的BPA,本文將時間軸0~800(月)劃分成20個區間,確定各信息源在每個區間內的特征樣本點及分布情況,如圖4所示。

圖4 尺度參數樣本分布Fig.4 Samples distribution of scale parameter

這樣,從尺度參數br累積概率的分布情況中統計得到特征樣本點的數值大小,并進行接下來的分析。各信息源S1~S3和現場試驗S0的特征樣本點的數值和區間劃分如表3所示。

分別計算第i個信息源在第j個區間內的樣本點Pi,j與現場真實分布的樣本點P0,j之間的距離L=|Pi,j- P0,j|,接著根據式(6)構建出各區間內的參數分布證據BPA,構建的20組BPA集合,如表4所示。

根據第2.2節所示方法及D-S證據理論合成規則,可以對數據容量證據E1和參數分布證據E2,j進行融合(j=1,2,…,20)。

表3 尺度參數樣本點累積概率表

由于目前沒有挖掘到更多的外部信息量,假設所有的信息源數據都采用統一的方法進行收集,并且在收集過程中有相同的技術水平,數據預處理精度一致,收集的數據有相同的可信度,這表明除了驗前分布本身的內部信息之外并沒有其他的外部信息應該在權重分配過程中考慮。所以,E1和E2,j經過證據合成后的BPA即為各信息源驗前分布的權重集合,即Wj=E1⊕E2,j。權重分配結果如表5所示。

進一步地,可以根據表6結合各驗前分布πi(br)計算得到融合后的尺度參數br的累積概率分布曲線,如圖5所示。由圖觀察可知,融合后的參數累積概率分布與其他各單一信息來源的驗前分布相比,更接近于現場真實分布。在認定現場數據最能反映產品真實信息的前提下,認為融合后的尺度參數br的累積概率分布曲線是合理的。

表4 參數分布證據BPA表

為了驗證經過D-S證據組合后的權重分配結果更加合理,分別使用數據容量證據E1和參數分布證據E2這兩類單一證據BPA作為驗前分布的融合權重的分配結果,即分別用m1和m2直接表示Wj。計算由這兩種單一證據所決定的權重分配,在各區間內融合后的特征樣本點與現場數據樣本點之間的距離,進行比較分析,結果如圖6所示。

表5 驗前分布權重分配向量

圖5 融合后的尺度參數累積概率分布曲線Fig.5 Cumulative probability distribution curve of scale parameter after fusion

圖6 三種融合樣本點分布距離比較Fig.6 Distances comparison of three types fusiondistribution samples

證據ddmax證據10.0103550.0309證據20.008230.0636證據組合0.007470.0412

3.3基于D-S證據合成的動量輪剩余壽命估計

根據第3.2節得到的各信息源驗前分布的融合權重向量集合Wj,利用式(4)和式(7),可以得到動量輪剩余壽命的尺度參數br的驗后概率密度

經判定,該評估結果和工程實際相符合。由此可見,本文提出的基于D-S證據理論的驗前分布融合權重分配的方法可行,可以解決工程實際問題。該方法改進了文獻[16]考量驗前分布相似性的基于JS距離的權重確定方法,充分挖掘了各驗前信息源的內外部信息,并建立相應的證據和證據的基本概率分配函數,加入到證據合成的過程中,使得各信息源的權重分配更加合理準確。

圖7 動量輪剩余壽命分布曲線Fig.7 Residual life distribution curve of momentum wheel

4 結束語

本文采用D-S證據合成方法來研究動量輪剩余壽命的驗前分布融合問題,充分挖掘了各驗前信息源的內、外部信息,構建了基于各類信息的D-S證據基本概率分配函數,并利用合成規則對多證據BPA進行合成,最終得到合理的各信息源的權重分配結果。

與現有研究方法相比,本文主要創新性在于,融合權重的確定不僅僅依據各信息源數據本身的特性,而是通過挖掘多種信息,構成D-S證據集合,再利用證據合成規則,來確定的權重分配結果。在對XX型動量輪的可靠性評估和剩余壽命預測過程中,針對多信息源所挖掘到的信息,構建了數據容量證據E1和參數分布證據E2的基本概率分配函數,通過兩類證據的BPA的合成,得到各驗前分布的權重分配并計算得出尺度參數br的最小方差估計值及估計區間。結果表明,通過本文方法獲得的動量輪剩余壽命的威布爾分布符合實際工程經驗,可以解決工程實際問題。

由于沒有挖掘到XX型動量輪3種信息源的更多的外部信息量,本文假設所有的信息源數據都采用統一的收集方法進行收集,并具有相同的技術水平,數據預處理方式一致,信息來源具有相同的可信度。實際上,若存在可量化的外部信息,也應作為證據納入對權重分配的考慮中。因此,本文后續工作還應探究外部信息的挖掘和證據的構建方法,使得權重的分配更加有據合理,剩余壽命的評估結果更加真實可信。

References)

[1]YE X, MA Y, MENG H, et al. Degradation failure model of electromagnetic relay[C]∥Proceeding of ICEC-ICREPEC 2012.Beijing, 2012: 116-123.

[2]邵瑞芝, 范本堯. 長壽命通信衛星的可靠性研究[J]. 中國空間科學技術, 1996, 16(4):24-33.

SHAO R Z, FAN BY. Reliability study of long-life communication satellites[J]. Chinese Space Science and Technology, 1996, 16(4):24-33(in Chinese).

[3]方艮海. 產品可靠性評估中的多源信息融合技術研究[D].合肥:合肥工業大學,2006.

FANG G H. Research on the multi-source information fusion techniques in the process of reliability assessment[D]. Hefei:Hefei University of Technology, 2006(in Chinese).

[4]莊釗文, 郁文賢, 王浩,等. 信息融合技術在可靠性評估中的應用[J]. 系統工程與電子技術, 1999, 22(4):6-9.

ZHUANG Z W, YU W X, WANG H, et al. Information fusion and application in reliability assessment[J]. Systems Engineering and Electronics, 1999, 22(4):6-9(in Chinese).

[5]SHAFER G. Perspectives on the theory and practice of belief functions[J]. International Journal of Approximate Reasoning, 1990, 4(5-6): 323-362.

[6]BAE H R, GRANDHI R V, CANFIELD R A. An approximation approach for uncertainty quantification using evidence theory[J]. Reliability Engineering and System Safety, 2004, 85(1-3): 215-225.

[7]AGARWAL H, RENAUD J E, PRESTON E L, et al. Uncertainty quantification using evidence theory in multidisciplinary design optimization[J]. Reliability Engineering and System Safety, 2004, 85(1-3): 281-294.

[8]高洪濤, 王敏. 證據理論在旋轉機械綜合故障診斷中的應用[J]. 大連理工大學學報. 2001, 41(4): 459-462.

GAO H T, WANG M . Application of D-S evidential reasoning to the comprehensive fault diagnosis of rotating machine[J]. Journal of Dalian University of Technology, 2001, 41(4): 459-462(in Chinese).

[9]陳果. 基于神經網絡和D-S證據理論的發動機磨損故障融合診斷[J]. 航空動力學報, 2005, 20(2): 303-308.

CHEN G. Fusion diagnosis of engine wearing fault based on neural network and D-S evidence theory[J]. Journal of Aerospace Power, 2005, 20(2): 303-308(in Chinese).

[10]范顯峰, 姜興渭. 基于多Agent的衛星故障診斷融合技術研究[J]. 中國空間科學技術, 2003, 23(2):39-44.

FAN X F, JIANG X W. Research of multi-agent based satellite fault diagnosis and fusion technology[J]. Chinese Space Science and Technology, 2003, 23(2): 303-308(in Chinese).

[11]耿俊豹, 邱瑋, 孔祥純,等. 基于粗糙集和D-S證據理論的設備技術狀態評估[J]. 系統工程與電子技術, 2008, 30(1): 112-115.

GENG J B, QIU W, KONG X C, et al. Techical condition evaluation for devices based on rough set theory and D-S evidence theory[J]. Systems Engineering and Electronics, 2008, 30(1): 112-115(in Chinese).

[12]ZHU P, XIONG W, QIN N, et al. D-S theory based on an improved PSO for data fusion[J]. Journal of Networks, 2012, 7(2):369-376.

[13]董彥民, 賀小帆, 劉文珽. 基于不同壽命分布的DFR值換算關系[J]. 北京航空航天大學學報, 2011, 37(12): 1524-1528.

DONG Y M, HE X F, LIU W T. Conversion relation of detail fatigue rating based on different fatigue life distribution[J]. Journal of Beijing University of Aeronautics and Astronautics, 2011, 37(12): 1524-1528(in Chinese).

[14]孫銳. 基于D-S證據理論的信息融合及在可靠性數據處理中的應用研究[D]. 成都:電子科技大學, 2012.

SUN R. Research on D-S evidence theory based information fusion and its application in reliability data processing[D]. Chengdu:University of Electronic Science and Technology of China, 2012(in Chinese).

[15]JIN G, MATTHEWS D, FAN Y, et al. Physics of failure-based degradation modeling and lifetime prediction of the momentum wheel in a dynamic covariate environment[J]. Engineering Failure Analysis, 2013, 28(3):222-240.

[16]劉強. 基于失效物理的性能可靠性技術及應用研究[D].長沙:國防科學技術大學, 2011.

LIU Q. Research on the performance reliability technology and the application based on physics of failure[D].Changsha:Nation University of Defense Technology, 2011(in Chinese).

(編輯:高珍)

Lifetime prediction of aerospace equipment based on D-S evidence theory

DING Rui1, LU Ningyun1,*, CHENG Yuehua2, JIANG Bin1, XING Yan3,4

1.College of Automation Engineering, Nanjing University of Aeronautics & Astronautics, Nanjing 211016, China 2.College of Aerospace Engineering, Nanjing University of Aeronautics & Astronautics, Nanjing 211016, China 3.Beijing Institute of Control Engineering, Beijing 100191, China 4.State Key Laboratory of Space Intelligent Control, Beijing 100191, China

An information fusion method was proposed based on D-S evidence theory and Bayes theory for lifetime prediction of important aerospace equipment, momentum wheel. Firstly, multi-source life information were collected and mined to obtain the prior distribution of the momentum wheel′s lifetime, in order to build D-S evidence collections.Secondly, D-S combination rule was used to obtain reasonable weight allocation for prior distributions. After that, fusion posterior distribution was figured out and the values of life parameters were also estimated. Finally, according to parameters′ estimation, the lifetime prediction for momentum wheel was derived. Simulation result shows that the prediction using the proposed method is closer to the real lifetime measurements.

lifetime prediction; multi-source information fusion;D-S evidence theory; Bayes method; Momentum wheel

10.16708/j.cnki.1000-758X.2016.0044

2015-11-26;

2015-12-30;錄用日期:2016-05-11;

時間:2016-07-1213:26:34

http:∥www.cnki.net/kcms/detail/11.1859.V.2016.0712.1326.001.html

國家自然科學基金(61374141,61203091)

丁瑞(1991-),女,碩士研究生,dingrui504@163.com

陸寧云(1977-),女,教授,博士生導師,luningyun@nuaa.edu.cn,主要研究方向為復雜系統數據驅動建模、故障診斷與預測的理論和應用

TB114

A

http:∥zgkj.cast.cn

引用格式:丁瑞,陸寧云,程月華,等.基于D-S證據理論的航天設備壽命預測方法[J].中國空間科學技術, 2016,36(4):58-66.DINGR,LUNY,CHENGYH,etal.LifetimepredictionofaerospaceequipmentbasedonD-Sevidencetheory, 2016,36(4):58-66(inChinese).

猜你喜歡
分配理論融合
堅持理論創新
當代陜西(2022年5期)2022-04-19 12:10:18
村企黨建聯建融合共贏
今日農業(2021年19期)2022-01-12 06:16:36
神秘的混沌理論
融合菜
理論創新 引領百年
從創新出發,與高考數列相遇、融合
相關于撓理論的Baer模
《融合》
現代出版(2020年3期)2020-06-20 07:10:34
應答器THR和TFFR分配及SIL等級探討
遺產的分配
主站蜘蛛池模板: 国产成人综合日韩精品无码首页| 免费人欧美成又黄又爽的视频| 成年片色大黄全免费网站久久| 亚洲第一香蕉视频| 日韩av高清无码一区二区三区| 四虎影视8848永久精品| 亚洲狼网站狼狼鲁亚洲下载| 91久久青青草原精品国产| 亚洲开心婷婷中文字幕| 无码有码中文字幕| 91精品专区| 国产无码在线调教| 国产午夜精品鲁丝片| 97人人模人人爽人人喊小说| 99热国产这里只有精品无卡顿"| 国产亚洲精品无码专| 先锋资源久久| 欧美精品亚洲日韩a| 小蝌蚪亚洲精品国产| 日本亚洲成高清一区二区三区| 亚洲精品波多野结衣| 影音先锋亚洲无码| 沈阳少妇高潮在线| 拍国产真实乱人偷精品| 亚洲欧洲日韩久久狠狠爱 | 久草中文网| 美女被躁出白浆视频播放| 欧美亚洲国产一区| 一级看片免费视频| 欧美69视频在线| 欧美啪啪视频免码| 99精品国产自在现线观看| 激情视频综合网| 露脸一二三区国语对白| 中文字幕乱码二三区免费| 亚洲无码免费黄色网址| 中文字幕免费视频| 日韩中文字幕亚洲无线码| 国产午夜看片| 国产一级精品毛片基地| 国产女人18水真多毛片18精品| 国产高清在线观看91精品| 亚洲无码日韩一区| 三上悠亚在线精品二区| 亚洲精品视频在线观看视频| 亚洲第一在线播放| 久久国产成人精品国产成人亚洲 | 国产成人久久综合一区| 欧美区一区| 99这里只有精品免费视频| 国产一级无码不卡视频| 国产精品妖精视频| 高清码无在线看| 欧美国产精品不卡在线观看| 老汉色老汉首页a亚洲| 亚洲人成亚洲精品| 波多野结衣久久精品| 中文字幕啪啪| 3D动漫精品啪啪一区二区下载| 日本伊人色综合网| 特级毛片免费视频| 91小视频在线观看免费版高清| 亚洲人妖在线| 欧美天天干| 久久国产精品麻豆系列| 久久免费视频6| 国产精品无码翘臀在线看纯欲| 亚洲欧美成aⅴ人在线观看| 18禁影院亚洲专区| 亚洲精品自产拍在线观看APP| 99国产精品一区二区| 国产日韩欧美视频| 国产精品永久免费嫩草研究院| 九色视频一区| 亚洲一区二区在线无码| 中国特黄美女一级视频| 成人精品视频一区二区在线| 波多野结衣视频网站| 亚洲免费三区| 亚洲av片在线免费观看| 国产一区二区人大臿蕉香蕉| 国产免费黄|