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

基于集成學習的沖擊載荷識別與定位方法

2022-12-02 09:12:34陽志光郭晨宇楊帆蔣亮亮王斌
強度與環境 2022年5期
關鍵詞:結構模型

陽志光 郭晨宇 楊帆 蔣亮亮 王斌

(北京宇航系統工程研究所,北京 100076)

0 引言

結構在設計階段因為不了解真實沖擊載荷,往往會采用很高的冗余度。獲取沖擊載荷的大小,有助于優化結構設計,同時對結構可靠性分析、健康監測提供重要的信息[1]。但是作用在結構上的沖擊載荷往往難以直接測量,比如火箭發動機點火瞬間對箭體結構的沖擊[2]、星箭分離時的沖擊載荷[3]。尤其是航空航天結構的設計以重量和體積為重要指標[4,5],安裝力傳感器會帶來很大的負擔,改變結構的整體特性,甚至在一些特殊的結構部位和載荷工況下,無法安裝力傳感器來直接測量沖擊載荷。由于加速度或應變響應易于測量,且響應傳感器靈巧輕便,更易于布置,因此通過動態響應來逆向識別載荷是一種更實用的方法。

目前已有許多載荷識別方法被提出,主要可以分為頻域法[1]和時域法[4]。基于頻域法的載荷識別技術發展較早,相對更加成熟,但是要求信號樣本具有一定的長度,通常只適用于穩態載荷的識別[6]。在時域內識別載荷因結果更加直觀,且適用于瞬態載荷,近些年也得到了發展[7]。動態載荷逆向識別面臨的一個主要問題是傳遞函數矩陣秩不足或者條件數過大,導致解的不唯一或不穩定的病態性問題[8]。目前可以采用奇異值分解方法[9]、正則化方法[10]以及基函數展開式[11]等方法,通過引入額外的超參數,增加約束,來解決病態問題。其中超參數的選取,對于克服病態問題,準確逆向識別載荷至關重要[12]。實際工程結構往往因載荷量級大、極致輕質化、連接部件多等原因而在力學行為上呈現明顯的非線性,并且環境噪聲也會給結構系統帶來不確定性影響[13],導致難以根據結構建立一個準確的模型。近些年隨著人工智能的發展,以及計算能力的迅速提升,基于數據驅動的機器學習方法在載荷識別中得到應用[14]。Wang[15]通過支持向量機識別作用在圓筒結構上的多源隨機動載荷;Zhou[16]利用深度循環網絡識別非線性復合材料板上的沖擊載荷。通過結構試驗或者建模仿真生成充足的訓練數據,并結合領域知識從數據中提取合適的特征,進行監督式學習,訓練得到機器學習模型作為代理模型,可以取代沖擊載荷和響應之間難以求解的卷積映射關系;并且隨著訓練數據量增加,涵蓋的情況更加豐富,模型對噪聲、非線性情況的魯棒性也越強。

實際上在大多數情況下,沖擊載荷作用在結構上的位置是隨機的,識別沖擊載荷作用位置,也是載荷識別問題的重要組成部分。通過應力波傳播的路徑和時間,可以定位沖擊載荷的位置,但只適用于梁、板等簡單的結構[17]。另一種定位沖擊載荷的方法是通過殘差優化處理[18]來最小化模擬響應和測量響應之間的誤差函數,相較于第一種方法,在實際應用中穩定可靠,但在連續結構中計算的時間會隨著自由度數的增加而增加[19]。基于機器學習的沖擊載荷定位方法,可以快速甚至實時識別載荷大小并且定位沖擊載荷的位置。Qiu[20]基于結構應力響應的余弦相似度,提出一種結合相似性度量的模式識別方法快速定位隨機作用的沖擊載荷,并通過在平面鋼板上試驗得到了計算效率高、識別準確的結果,但方法是否適用于更復雜的結構形式,還有待進一步研究。

目前沖擊載荷識別難以在工程結構中得到大規模應用,主要是受到結構非線性力學行為和環境噪聲帶來的難以建立模型問題和逆向求解時的病態問題,并且僅通過波的傳播過程或者殘差優化方法來定位沖擊載荷,難以建立準確、適用性廣泛、計算量少的公式,因此以響應輸入和載荷輸出為數據驅動的機器學習方法更具有應用前景。

考慮到單個機器學習模型所面臨的欠擬合和過擬合的風險,而深度學習模型參數眾多,訓練時間長且難以收斂。本文采用集成學習方法識別沖擊載荷,基于數據驅動免去了建立模型和逆向求解的問題,直接建立沖擊載荷和結構響應之間的映射關系。通過試驗獲取“載荷-響應”數據集,訓練以決策樹為基礎預測器進行集成學習的自適應提升AdaBoost模型,根據結構響應對沖擊載荷做出預測,并在真實的薄壁圓筒結構上對方法進行驗證。

1 載荷識別原理

1.1 逆向識別問題

假設一個力學結構為線性定常系統,并且初始條件為 (0) 0y= 和 (0) 0y=˙ ,則載荷和響應之間的關系可以通過如下卷積關系表示:

其中h(t)代表在t時刻結構的沖擊響應函數,f(t)代表沖擊載荷在t時刻的幅值,符號*代表卷積運算,延遲時間算子τ滿足t≥τ。

結構響應y(t)是通過對應時刻的沖擊f(t)和系統傳遞函數h(t)卷積運算得到。卷積運算通過傅里葉變換,可以轉換為頻域內頻響函數H(ω)和載荷F(ω)的乘積運算。假設結構系統上響應自由度數為i,載荷自由度數為j,則頻域內離散化推導如下:

其中{Y(ω)}i×1和{F(ω)}j×1分別為包含i和j個分量的一維響應向量和載荷向量,ω為角頻率。

載荷作為激勵,引起結構系統的響應,是一個正向過程。通過響應結果來推導載荷激勵,即載荷識別,是一個逆向問題。最直接得到載荷的方法是對頻響函數矩陣求逆矩陣并與響應向量相乘:

但實際上傳遞函數矩陣求逆時往往面臨著秩虧損、條件數過大等問題,導致逆向求解載荷的結果不唯一或者對噪聲非常敏感。雖然可以通過改善響應測點布局、正則化約束方法等提升逆向計算的準確度,但這又面臨著測點布局選取和正則化基函數、參數的優化問題,并且隨著結構的復雜程度增加,優化難度也大大增加。更重要的是,許多實際的大型復雜結構由于不滿足線性或定常假設,其材料參數、幾何特征、邊界條件無法準確的知道,難以從力學基理上構建結構系統的傳遞函數。但是,一定存在函數,ijF可以描述結構上載荷輸入和響應輸出之間的關系:

函數,ijF描述了在t時刻作用于結構上j點處的沖擊力fj(t)與結構上i點處動態響應yi(t)之間的關系,n為結構的自由度數。

由前面的分析可知,從力學基理上構建系統模型,并逆向求解得到沖擊載荷和響應之間的傳遞函數是十分復雜困難的。采用機器學習方法進行監督式學習,可以直接構建沖擊載荷和動態響應之間的映射函數來替代力學模型,實現由動態響應識別沖擊載荷。

1.2 集成學習原理

目前機器學習方法已經發展出多種成熟的學習模型。線性模型、支持向量機[21]原理清晰,結構形式簡潔,但是學習性能有限,難以建立復雜的映射關系;決策樹[22]、神經網絡[23]隨著模型內部復雜度的增加,可以極大地提升學習能力,但是單個模型容易受極端樣本的影響,陷入局部最優。集成學習針對單個機器學習模型所面臨的欠擬合和過擬合的風險,同時訓練多個基礎模型,綜合給出判斷,同時保證預測結果較低的偏差與方差。其中決策樹模型[22]因其決策原理清晰、訓練迅速,適合作為集成學習模型的基礎預測器。

1.2.1 決策樹

決策樹算法整體形狀類似于樹狀結構,數據實例從根節點出發,通過對某個特征屬性判斷決策,分類到相應地子節點后,再次選取對某個特征屬性進行決策,不斷決策分類直到延伸至葉節點。葉節點根據其包含樣本數量中最多的類別,對應一個類;在回歸預測任務中葉節點對樣本數量值取平均值給出預測。

決策樹模型在分類任務的訓練過程中,通過CART算法在某個節點處,選取用于決策的特征k和進行分類的閾值tk,算法嘗試最小化分類成本函數:

其中G1、G2代表左、右子集的不純度,通常是計算基尼不純度,m、m1、m2分別是當前節點實例數量和左、右子集的實例數量。

根據動態響應來預測具體的載荷數值在機器學習中數據回歸任務,CART算法分裂訓練集的方式由最小化純度轉變為最小化均方誤差MSE:

其中節點處分裂子集的均方誤差計算公式如下:

決策樹極少對數據做出假設,因此適用于多種數據模型。但是如果不加以限制,決策樹結構將一直分裂生長直至嚴密擬合訓練數據,反而會導致在面對全新的數據進行泛化時效果不佳,出現過度擬合。為避免過擬合,通過設置正則化超參數來限制樹的生長,降低決策樹的自由度。比如超參數“最大深度”控制決策樹的最大生長深度,“葉節點最小樣本數”控制決策樹末端節點所必須有的最小樣本數量。

1.2.2 自適應提升法

提升法是指將多個較弱的預測器結合成一個強學習器的任意集成方法,目前常用的提升法有自適應提升法(AdaBoost)和梯度提升法(Gradient Boost)[22]。以決策樹為基礎預測器,AdaBoost依次訓練多棵決策樹,每棵決策樹都是對之前決策樹的預測進行糾正,針對前序決策樹擬合不足的訓練實例給予更多關注,從而使新的決策樹越來越專注于難纏的問題,預測效果逐步提升。

基于AdaBoost算法訓練模型時,每個樣本的權值w(i)最初設置為1/m,其中m為樣本實例數量。第一個預測器完成訓練后,計算其加權誤差率r1,公式如下:

第j個預測器的權重αj的計算如公式(9),其中η為學習率超參數。預測器的準確率越高,則該預測器的權重就越高。當前預測器完成訓練后,對實例的權重進行更新,更新規則如公式(10)所示,然后將左右樣本實例的權重歸一化。

在更新權重后的樣本實例上重新訓練預測器,如此重復,直至達到所需數量的預測器或者得到滿意的預測效果,算法停止。在預測的時候,AdaBoost計算所有預測器的預測結果,并使用每個預測器的權重αj對它們進行加權。

2 沖擊載荷識別與定位方法

以決策樹為基礎預測器,采用自適應提升法訓練多棵決策樹得到集成學習模型,實現沖擊載荷的逆向識別與定位。模型的訓練依賴于大量的動態響應與載荷標簽的數據集。結構受沖擊后采集的動態響應數據量龐大,并且往往受噪聲影響,如果直接輸入給機器學習模型會嚴重拖慢學習效率,需要對響應數據濾波去噪,并提取數量少但信息熵高的特征。

2.1 數據預處理與特征提取

為充分獲取結構受沖擊后的動態響應信息,通常將多個高頻應變或加速度傳感器布置在結構上。在彈性假設范圍內,結構受沖擊后的動態響應數據呈現出在基線附近上下振蕩并且幅值逐漸衰減的時間序列模式,振蕩的周期和衰減的快慢與結構的模態和阻尼相關。然而由于外界環境干擾、采集電信號不穩定等因素,實際測量的響應數據往往會出現基線漂移、噪聲混雜的情況。

采用中值濾波方法處理基線漂移問題,其基本原理是將時間響應序列中某一時刻的點用該點鄰域中各點值的中值代替。選用較大尺寸的窗口,中值濾波得到初始漂移的基線輪廓,用序列數據減去漂移的基線,就得到了以零值為基線振蕩的響應數據。

實際結構的動態響應經過中值濾波后,采用小波變換方法進行降噪。選用‘sym8’小波函數對響應信號進行小波變換后,信息被儲存在變換后的小波系數中。其中響應信號的小波系數較大,噪聲的小波系數較小。通過選取合適的閾值,認為小于此閾值的小波系數是由噪聲產生,將其置為零從而到達去噪的目的,然后還原響應信號。

工程上關注沖擊力的大小和沖擊作用的位置是否會對結構產生破壞,因此提取每個沖擊工況下沖擊力的峰值和位置坐標作為監督式學習的標簽。接著對初始響應數據進行特征提取,得到與沖擊載荷密切相關、信息熵高并且數量盡可能少的輸入特征。

根據定性的動力學理論分析,結構受沖擊載荷作用會發生振動,沖擊力越大,振動得越劇烈,因此將響應的振蕩峰值作為特征之一;沖擊通過波的形式在結構上傳播,先后被結構上的傳感器感受到,距離近的傳感器更早感受到沖擊作用,數據波動也更劇烈,因此將到達振蕩峰值的時間作為特征之一;響應初始振蕩的脈寬也與沖擊相關,將其作為提取的特征之一。至于振蕩的周期、衰減率,與結構自身的模態相關,而與外部沖擊無關,因此不作為提取的特征。在發生沖擊后的前一小段時間內,結構響應與沖擊載荷密切相關,其后主要與結構自身的模態信息相關,提取的特征包括響應傳感器記錄到的最大峰值、到達峰值時間、峰值脈寬,以及自沖擊后前5個振蕩的峰值、到達峰值時間和峰值脈寬。假設結構上共布置有nx個高頻應變傳感器,則提取的特征數nf為:

2.2 集成學習模型訓練

采集多個沖擊工況的數據,以沖擊力峰值和沖擊位置為標簽,對響應數據濾波降噪并提取特征,采用監督式學習訓練機器學習模型。數據集按照4:1的比例被分割為訓練集和測試集,訓練集用來訓練集成學習模型,測試集用于對完成訓練的模型的載荷識別效果進行評估。

以決策樹為基礎預測器的AbaBoost模型,其學習效果受模型超參數設置的影響,采用網格篩選結合交叉驗證的方法,將訓練集進一步分割為新的訓練集與驗證集,來對模型超參數:葉節點最小樣本數min_samples_leaf、集成決策樹數量n_estimators和學習率η進行篩選。計算不同參數組合下模型在驗證集上的均方根誤差RMSE,篩選得到最優參數組合,進一步提升模型的學習性能。

采用3折交叉驗證,在初步篩選得到min_samples_leaf=2,對超參數n_estimators和η進一步篩選的結果如圖 1??梢娫趎_estimators=200,η=0.2時,AdaBoost模型可以達到最佳預測效果。完成訓練的模型被部署到結構系統上,用于識別沖擊載荷。

圖1 采用3折交叉驗證網格篩選的均方根誤差Fig.1 Rmse of Grid-search with 3-fold cross-validation

2.3 載荷識別步驟

訓練用于識別沖擊載荷的Adaboost模型的主要步驟如下:

Step 1:數據采集。在機械結構上作用隨機沖擊載荷,記錄沖擊力峰值和沖擊位置,采用高頻傳感器采集響應數據;

Step 2:數據預處理。對響應數據中值濾波和小波變換降噪,提取響應數據中多個振蕩峰值、到達峰值時間、峰值脈寬,作為特征;

Step 3:模型訓練。提取特征作為輸入,沖擊力和沖擊位置作為輸出,采用網格篩選結合交叉驗證的方法優化模型超參數,訓練得到AdaBoost模型;

Step 4:模型部署。將數據降噪、特征提取、AdaBoost模型集成為系統,實現輸入原始響應數據,識別輸出沖擊載荷的峰值與位置。

將完成訓練的AdaBoost模型部署到結構上進行載荷識別的流程如下:

圖2 沖擊載荷識別與定位流程Fig.2 Process of impact load identification and localization

3 試驗系統搭建

鋁制薄壁圓筒具有密度低、承壓能力強的優秀性質,在航天運載器中經常作為艙段結構,在該結構的縮比試驗件上驗證沖擊載荷識別方法。圓筒結構高度為600mm,直徑為1000mm,壁厚1.5mm,端框厚度為4mm。圓筒具有空間曲面構型和上下突出的端框,并且采用4個吊點的懸掛邊界條件,具有明顯的幾何非線性,在承受沖擊時結構響應復雜。

在圓筒內部各個象限均勻粘貼應變片,來記錄結構應變響應。另外通過一個外表面粘貼的加速度傳感器作為觸發傳感器。當加速度超過量程10%時,記錄共計1.5s時長的沖擊力和應變響應,其中前100ms為觸發前的信號,目的是完整地記錄沖擊力過程。上述沖擊力、應變和加速度均是通過DH5923高速動態數據采集系統進行采集,設定采樣頻率為5kHz。

考慮到圓筒為軸對稱結構,通過軸向高度和周向角度兩個坐標定位敲擊的位置。沿著軸向均勻劃分16行,周向劃分64列,如圖3所示。力錘敲擊這些網格點,記錄敲擊點的高度和角度作為位置標簽,用以監督式學習。最終總共采集了1684條敲擊數據,每條數據包含沖擊力的歷程、沖擊作用在圓筒表面的高度和角度,以及應變傳感器的響應數據。其中記錄到某個樣本的沖擊力和某個應變片的應變響應如圖4所示,可見沖擊力的形狀比較復雜,應變響應變化迅速,數據量龐大,并且受噪聲影響。

圖3 圓筒外表面網格劃分及內部應變片分布Fig.3 Distribution of grids and strain gaudes on the surface of the cylinder

圖4 沖擊載荷和動態響應Fig.4 Impact load and dynamic strain response

4 試驗結果與討論

采用平均絕對誤差MAE和相對誤差RE作為載荷識別效果的評價指標,計算公式如下:

4.1 沖擊載荷識別與定位

將采集到動態響應數據進行中值濾波,并選用選用‘sym8’小波函數進行去噪。下圖是在某個沖擊工況下的應變響應數據以及濾波后的響應數據。可見經過濾波降噪,應變響應數據中基線漂移和噪聲問題被很好地解決。

圖5 動態應變響應濾波和去噪Fig.5 Filtering and denoising on original strain response

從濾波去噪后的應變響應數據中提取峰值、到達峰值時間和脈寬作為特征,訓練AbaBoost模型分別用于識別沖擊力的峰值和沖擊位置的高度和角度。

訓練用于識別角度的模型時,0°和360°在數值上差異最大,但實際上是在同一位置,誤差梯度為零,導致在0°識別角度偏大,在360°識別角度偏小。解決辦法是對角度標簽進行三角函數轉換,從而得到在0°和360°處連續一致的標簽。但由于單個三角函數在2π周期內的映射不唯一,將角度標簽分別進行sin和cos三角函數轉換,訓練兩個模型,根據兩個識別結果判斷角度范圍,從而逆向識別出唯一的角度。角度所在范圍與sin、cos函數轉換后正負號如下表所示:

表1 三角函數轉換區間Table1 positive and negative signs of the range of angles

下圖為在測試樣本上,AdaBoost模型對沖擊峰值、位置高度和角度的識別結果對實際值的擬合效果。由力錘敲擊產生的沖擊力分布在0~550N的范圍內,從圖6 (a)可見識別沖擊力峰值對實際峰值的擬合性較好,但是對于峰值較大和較小的樣本,擬合稀疏發散,識別誤差較大。這主要是因為在手動敲擊圓筒時,沖擊力樣本分布類似于正態分布,導致機器學習在訓練時更關注數量上占優勢的中段沖擊力樣本。在圖6 (b)中可見,模型對沖擊力作用點的擬合效果良好。

圖6 AbaBoost模型的識別結果Fig.6 Identification and loaclisation of impact load using AdaBoost model

AdaBoost模型識別沖擊載荷的平均絕對誤差MAE和相對誤差RE及95%置信區間半徑如表 2所示,識別沖擊力峰值的MAE=20.56±2.49,MRE=13.44%±1.30%;識別高度的絕對誤差MAE=47.05±1.72,以圓筒高度 600mm 為基數的相對誤差MRE=7.84%±0.79%;識別角度的絕對誤差MAE=1.45±0.49,以圓筒環向一周360°為基數的相對誤差MRE=0.40%±0.14%??梢姴捎肁daboost模型可以較好地識別沖擊力大小,并準確定位沖擊的位置。

表2 AdaBoost模型識別沖擊載荷的誤差Table2 Impact load identification error of AdaBoost Model

4.2 邊界條件影響

圓筒為空間曲面結構,并且上下兩端邊界為突出的端框,具有幾何非線性。當圓筒承受沖擊載荷時應力波從沖擊區域四散傳播,到達端框邊界時發生反射交匯,導致此區域處應力應變響應復雜。從上一節的載荷識別結果可見,在靠近上下端框的邊界區域識別誤差較大。探討去除靠近邊界的樣本,是否能改善載荷識別的準確度。將下端框0~100mm和上端框500~600mm區域定義為邊界區域,表3是從測 試集中剔除邊界區域內樣本后,在中段區域內各個模型識別沖擊載荷峰值以及位置與全筒段區域上識別相對誤差的對比。

表3 非邊界區域沖擊載荷識別誤差Table3 Impact load identification error in the middle area

在遠離端框邊界的中段區域,沖擊峰值識別的相對誤差從13.44%降低到了11.88%,沖擊位置高度和角度識別的相對誤差分別從7.84%、0.40%降低到了7.03%、0.26%。從圖7可見,在遠離邊界的區域,模型對沖擊載荷的識別效果與是實際載荷的擬合效果也更好。邊界區域應力環境復雜,非線性程度高,會導致較大的識別誤差,因此使載荷作用在遠離邊界的區域可以改善載荷識別與定位的準確度。

圖7 在遠離邊界區域沖擊載荷的識別擬合效果Fig.7 Fitting effect of predicted impact load to actual impact load

5 結論

本文提出基于集成學習識別復雜結構上沖擊載荷并定位載荷作用位置的方法。以決策樹 為基礎預測器,采用自適應提升方法進行集成,通過沖擊試驗獲取“載荷-響應”數據集,進行監督式學習,訓練得到AdaBoost模型,實現由動態響應逆向識別沖擊載荷,并在薄壁圓筒結構上進行驗證,結果表明:

1)集成學習方法可以根據基礎預測器的學習效果給予相應的權重,利用群體的智慧提升載荷識別的效果,并有效避免陷入局部最優的過擬合風險;

2)模型的訓練依賴于數據的支持,從原始動態響應中去除噪聲干擾,并提取數量少、信息熵高的特征,可以提升訓練效率和載荷識別的準確度;

3)在靠近結構邊界的區域,應力環境復雜,容易產生較大的識別誤差,因此設計結構使載荷作用在遠離邊界的區域可以改善逆向識別的效果。

本文提出的集成學習方法可以識別并定位復雜結構上的沖擊載荷,其接近于10%的識別誤差滿足工程應用需求。但從動態響應中提取合適的特征依賴于工程經驗和大量的試驗,后續將針對動態響應特征提取和提升載荷識別準確度開展進一步研究。

猜你喜歡
結構模型
一半模型
《形而上學》△卷的結構和位置
哲學評論(2021年2期)2021-08-22 01:53:34
重要模型『一線三等角』
重尾非線性自回歸模型自加權M-估計的漸近分布
論結構
中華詩詞(2019年7期)2019-11-25 01:43:04
新型平衡塊結構的應用
模具制造(2019年3期)2019-06-06 02:10:54
論《日出》的結構
3D打印中的模型分割與打包
FLUKA幾何模型到CAD幾何模型轉換方法初步研究
創新治理結構促進中小企業持續成長
現代企業(2015年9期)2015-02-28 18:56:50
主站蜘蛛池模板: 91精品人妻互换| 成人福利免费在线观看| 国产区人妖精品人妖精品视频| 日韩第一页在线| 伊人91在线| 国产麻豆aⅴ精品无码| 成人精品午夜福利在线播放| 波多野结衣一区二区三区AV| 在线无码九区| 91极品美女高潮叫床在线观看| 54pao国产成人免费视频| 区国产精品搜索视频| 亚洲精选高清无码| 国产人前露出系列视频| 免费av一区二区三区在线| 尤物视频一区| 东京热一区二区三区无码视频| 欧美日本在线| 91久久国产热精品免费| 久久久久久久久18禁秘| 米奇精品一区二区三区| 热99精品视频| 99re在线免费视频| 九九线精品视频在线观看| 国产精品成人观看视频国产 | 国产制服丝袜无码视频| 十八禁美女裸体网站| 亚洲精品视频免费| 日本久久久久久免费网络| 亚洲精品另类| 精品国产中文一级毛片在线看| 青青草国产免费国产| 国产一在线| 亚洲综合激情另类专区| 女人天堂av免费| 亚洲人成日本在线观看| 国产精品综合久久久 | AⅤ色综合久久天堂AV色综合| 亚洲中文字幕无码爆乳| 香蕉久人久人青草青草| 国产精品自在线天天看片| 亚洲第一在线播放| 澳门av无码| 亚洲精品中文字幕无乱码| 91在线国内在线播放老师| 欧美激情综合| 亚洲最新在线| 亚洲午夜福利在线| 人人看人人鲁狠狠高清| 天天色天天综合| 欧美色视频在线| 久久国产精品娇妻素人| 88国产经典欧美一区二区三区| 国产视频久久久久| 日韩黄色大片免费看| 国产又粗又猛又爽视频| 人人妻人人澡人人爽欧美一区| 国产流白浆视频| 亚洲一区二区约美女探花| 国产精品久久久久久久伊一| 亚亚洲乱码一二三四区| 小说区 亚洲 自拍 另类| 亚洲AV免费一区二区三区| 国产成人啪视频一区二区三区| 国产97公开成人免费视频| 毛片最新网址| 国产无码在线调教| 国产一级小视频| 日本妇乱子伦视频| 久久中文无码精品| 亚洲一区二区成人| 国产精品综合久久久| 亚洲精品成人片在线观看| 国产波多野结衣中文在线播放| 欧美激情二区三区| 亚洲AV无码不卡无码| 97av视频在线观看| 在线观看国产黄色| 狠狠色噜噜狠狠狠狠色综合久| 亚洲AV无码乱码在线观看裸奔| 久久五月天国产自| 免费女人18毛片a级毛片视频|