• <tr id="yyy80"></tr>
  • <sup id="yyy80"></sup>
  • <tfoot id="yyy80"><noscript id="yyy80"></noscript></tfoot>
  • 99热精品在线国产_美女午夜性视频免费_国产精品国产高清国产av_av欧美777_自拍偷自拍亚洲精品老妇_亚洲熟女精品中文字幕_www日本黄色视频网_国产精品野战在线观看 ?

    采用廣義 Pareto 分布進行車輛荷載效應極值估計的研究

    2016-01-15 03:33:19馮海月,伊廷華,陳斌
    振動與沖擊 2015年15期
    關(guān)鍵詞:閾值

    采用廣義Pareto分布進行車輛荷載效應極值估計的研究

    馮海月1,伊廷華1,2,陳斌3

    (1. 大連理工大學土木工程學院,遼寧大連116023;2. 同濟大學土木工程防災國家重點實驗室,上海200092;3.浙江工業(yè)大學建筑工程學院,杭州310003)

    摘要:為了解決車輛荷載效應數(shù)據(jù)間的相關(guān)性問題并能夠充分利用樣本數(shù)據(jù),提出了一種改進的獨立風暴法。首先提出了一種三次類閾值法進行數(shù)據(jù)的初步分析來獲得閾值;然后考慮到簡支梁橋的車輛荷載效應主要受單輛重車影響,提出了改進的獨立風暴法進行獨立同分布樣本數(shù)據(jù)的提??;最后采用廣義Pareto分布進行車輛荷載效應極值估計;利用該方法,對動態(tài)稱重系統(tǒng)記錄的國內(nèi)某大橋?qū)崪y數(shù)據(jù)進行了車輛荷載效應極值估計,并與超越閾值法、獨立風暴法進行了對比分析,結(jié)果表明:在較短評估期(T<20年)內(nèi),三種方法均可以較好地預測荷載效應極值;而在中長評估期內(nèi),改進的獨立風暴法預測值較高,預測結(jié)果偏于安全。

    關(guān)鍵詞:車輛荷載效應;廣義Pareto分布;閾值;獨立風暴法;動態(tài)稱重系統(tǒng)

    中圖分類號:TU312.1文獻標志碼:A

    基金項目:國家自然科學基金(10632040, 11302058)

    收稿日期:2014-06-13修改稿收到日期:2014-07-16

    基金項目:汽車仿真與控制國家重點實驗室開放基金(20121117);中央高校基本科研業(yè)務費專項資金(長安大學基礎研究支持計劃專項基金CHD2012JC021)

    收稿日期:2014-07-09修改稿收到日期:2014-09-03

    Extreme estimation for vehicle load effect based on generalized Pareto distribution

    FENGHai-yue1,YITing-hua1,2,CHENBin3(1. School of Civil Engineering, Dalian University of Technology, Dalian 116023, China;2. State Key Laboratory for Disaster Reduction in Civil Engineering, Tongji University, Shanghai 200092, China;3. College of Civil Engineering and Architecture, Zhejiang University of Technology, Hangzhou 310003, China)

    Abstract:In order to solve data correlation problems for vehicle load effect and make full use of data samples, the modified method of independent storms (MMIS) was proposed. Firstly, a triple-class threshold method was proposed to obtain the threshold with the primary data analysis. Then, considering that the vehicle load effect of simply supported beam bridges was mainly affected by a single heavy vehicle, the MMIS was adopted to extract the independent and identically distributed (IID) sample data. Finally, the extreme vehicle load effect was estimated with the generalized Pareto distribution. In the end, the estimation of the extreme vehicle load effect was performed for the measured data of a bridge recorded with the weighinmotion (WIM). The results were compared with those of the peak-over-threshold method and those of the method of independent storms. The results showed that within the shorter estimation period (T<20 years), all the three methods can be used to better predict the extreme load effect; while within the middle or longer estimation period, the estimation of the modified method of independent storms is higher and safer.

    Key words:vehicle load effect; generalized Pareto distribution; threshold; method of independent storms; weighinmotion (WIM)

    車輛荷載是影響橋梁結(jié)構(gòu)安全性和耐久性的主要因素之一。近年來,隨著我國交通貨運量呈現(xiàn)爆發(fā)式的增長,重車和超重車已逐漸成為最重要的運輸工具之一。因此,針對橋梁所處位置的交通狀況,建立精確的車輛荷載效應模型,具有重要的現(xiàn)實意義。車輛荷載效應是一個隨機過程,它受車重、軸重、車間距、軸距、多車同時作用以及車速等多因素耦合作用的影響。精確的車輛荷載效應模型不但能夠?qū)煌ㄜ囕v的限載提供基礎依據(jù),而且能夠為橋梁安全評估、壽命預測及養(yǎng)護維修提供有力的支持。

    動態(tài)稱重系統(tǒng)(Weigh In Motion,WIM)是一組安裝的傳感器和含有軟件的電子儀器,用以測量動態(tài)輪胎力和車輛通過時間并提供計算輪重、軸重、總重(如車速、軸距等)的數(shù)據(jù)[1-2]。由于它可以準確獲取實際車流數(shù)據(jù),因此通過這些數(shù)據(jù)可以進行車輛荷載效應的計算,進而利用極值理論來預測未來荷載效應的情況[3]。傳統(tǒng)基于極值的方法均是基于加載事件服從獨立同分布這一基本假設[4],而實際的車輛荷載效應屬于平穩(wěn)時間序列,數(shù)據(jù)之間具有一定的相關(guān)性,因此簡單地將荷載效應作為獨立同分布序列并不十分合理。Cook[5]通過分析一段時間內(nèi)的風暴最大值,提出了獨立風暴理論(Method of Independent Storms,MIS),并采用Gumbel分布進行擬合,得到了不同評估期的荷載極值。Bhattacharya[6]認為平穩(wěn)荷載為點過程序列,他通過分析相關(guān)獨立性,并結(jié)合閾值超過數(shù)的分簇研究,提出了極值指標的計算方法,得到了平穩(wěn)時間序列分布函數(shù)。Naess等[7]基于極值理論是漸進分布,在實際應用時無法判斷在多大程度上適用于實際數(shù)據(jù)的問題,提出了針對Gumbel 分布的最大值分布模型。Mazas等[8]針對超越閾值法(Peak Over Threshold,POT)中閾值的合理選擇問題,提出了一種雙閾值方法,并進行了不同評估期荷載取值研究。李春前[9]提出了類BM離散法,得到近似獨立荷載數(shù)據(jù),采用廣義Pareto分布(Generalized Pareto Distribution,GPD)來擬合車輛荷載效應的尾部分布,認為一天荷載效應值只代表一天最大值截口分布,從而得到了更加貼近實際的車輛荷載效應概率模型。Wang等[10]針對數(shù)據(jù)的非獨立性,提出將原始數(shù)據(jù)序列在足夠大的時間間隔中提取出來,得到漸進獨立數(shù)據(jù),并結(jié)合POT分簇算法得到最優(yōu)閾值,然后采用GPD分布進行擬合,得到了短期荷載極值分布。

    本文采用GPD模型對車輛荷載效應尾部值進行了估計研究,提出了一種三次類閾值法進行數(shù)據(jù)的初步選擇,并針對數(shù)據(jù)間的相關(guān)性問題,提出了一種改進的獨立風暴法(Modified Method of Independent Storms,MMIS);最后利用WIM記錄的某橋一個月實際車流數(shù)據(jù),對不同設計基準期的車輛荷載效應極值進行了預測,并與POT模型的方法進行對比,得到了一些有益的結(jié)論。

    1極值理論和GPD模型

    1.1極值理論

    根據(jù)極值理論[11],設獨立同分布的隨機樣本{X}={X1,X2…Xn}服從某個未知的分布函數(shù)F(x)。假設Mn=max{Xi},則存在數(shù)列an>0和bn∈R,使得當n→∞時,概率Pr{(Mn-bn)/an

    這三種分布代表了三種不同的極值行為,可以用統(tǒng)一的形式進行表示:

    (2)

    式中:μ,k∈R,σ>0;μ是位置參數(shù);σ為尺度參數(shù);k為形狀參數(shù)。H(x)稱為廣義極值分布,當k=0時,H(x)為極值Ⅰ型,分布的尾部呈指數(shù)狀;k>0時,H(x)為極值Ⅱ型,分布的尾部較長;k<0時,H(x)為極值Ⅲ型,分布具有有限的上端點。

    對于式(2),當有足夠大的閾值u,且X>u時,有超出量Y=X-u樣本服從GPD分布:

    (3)

    1.2車輛荷載效應的極值理論

    概率密度函數(shù)可以分為低尾、中間、和高尾三部分[12]。由于傳統(tǒng)的單峰或者多峰分布擬合方法可以滿足車輛荷載效應分布中低尾和中間部分的精度要求,因此本文主要進行高尾部分的研究。

    根據(jù)GPD理論,當u值選定時,樣本變量X的累積概率分布為:

    (4)

    2獨立樣本數(shù)據(jù)的獲取

    傳統(tǒng)計算車輛荷載效應極值的方法是由WIM系統(tǒng)獲得實際車輛的相關(guān)信息,然后根據(jù)科學計算軟件來模擬產(chǎn)生荷載效應,并通過一定的時間間隔將效應值劃分為若干段,通過提取每段的最大值作為分析的樣本。若選擇的時間間隔足夠長,則認為可以采用極值理論。通常選擇的時間間隔為一天,然后假設一天的荷載效應最大值為全年最大值,然后應用極值理論進行不同基準期的荷載效應最大值評估。這種方法計算簡單,應用較為普遍。但是一天的荷載效應只使用一個最大值,其他數(shù)據(jù)則不考慮會造成數(shù)據(jù)的極大浪費;或者選擇某個閾值,提取大于這個閾值的所有樣本數(shù)據(jù),并用GPD方法進行擬合,從而得到車輛荷載效應分布函數(shù)。閾值的選擇會直接影響到結(jié)果的準確性。過大的閾值會使超出量太少,造成估計量的方差較大;而過小的閾值,則會使得估計量成為有偏估計[11]。

    2.1三次類閾值法

    本文提出的三次類閾值法的基本思想是剔除荷載效應樣本值中的低尾和中間部分,從而得到MMIS法分析所需的樣本,具體步驟為:

    (1)根據(jù)樣本數(shù)據(jù)經(jīng)驗分布函數(shù),得到類閾值u1。類閾值u1應取在樣本經(jīng)驗分布函數(shù)的中間部分,這樣不僅可以將原始數(shù)據(jù)中低尾的部分剔除掉,而且還可以減少類閾值u2的選擇循環(huán)次數(shù),從而節(jié)省類閾值的選擇時間。提取大于u1的樣本{Y}={Y1,Y2,…,YN},這里N為原始樣本中大于u1的個數(shù)。

    (2)由平均超出函數(shù)法,初步得到一個類閾值u2。類閾值u2的值宜取樣本經(jīng)驗分布函數(shù)的中間部分,這樣可以保證得到的數(shù)據(jù)樣本不會丟失風暴值。平均超出函數(shù)[13]是類閾值u2的線性函數(shù),其表達式為:

    (5)

    式中:μ,k∈R,σ>0;μ,σ,k分別對應極值分布中的位置參數(shù)、尺度參數(shù)和形狀參數(shù)。

    對于樣本{Y}={Y1,Y2…YN},平均超出函數(shù)可由經(jīng)驗平均超出函數(shù)(Empirical mean excess function,EMEF)en(u2)來估計,表達式為:

    (6)

    (3)根據(jù)原始樣本{X}得到累計分布函數(shù)值為0.99所對應的樣本值,此值即為類閾值u3。在這里將類閾值u3作為三次類閾值法的閾值選擇上限。

    (4)分別將類閾值u2和u3作為MIS法的閾值進行分析。如果MIS法得到的兩個類閾值所對應的獨立風暴個數(shù)相差一個量級,則重新進行類閾值u2的選擇,直到前后選擇的類閾值u2所對應的MIS樣本數(shù)相差不超過一個量級。這樣選擇類閾值u2的原因是,MIS法對閾值的穩(wěn)定性較好,閾值適當?shù)奶岣呋驕p小對MIS分析樣本不會產(chǎn)生太大的影響[14]。

    2.2MMIS法

    MIS法首先需要確定出閾值u。Cook[5]將閾值點稱為平靜點,兩個平靜點之間為獨立風暴。即若是連續(xù)數(shù)據(jù)超過閾值u,就將其看作為是一個獨立風暴;將每個獨立風暴的最大值作為獨立樣本,忽略掉該數(shù)據(jù)段中較小的其它值。

    由于MIS法在獲取數(shù)據(jù)樣本時,僅關(guān)注兩個連續(xù)閾值之間的最大值,并不考慮連續(xù)閾值區(qū)間次大值是否大于鄰近區(qū)間的最大值,從而會造成數(shù)據(jù)的浪費;且同時考慮到簡支梁橋?qū)嶋H車輛荷載效應最大值主要是由重車控制這一特點,本文提出了一種MMIS法。具體步驟為:

    (1)由三次類閾值法得到閾值u。

    (2)將連續(xù)超過閾值u的車輛荷載效應作為一個獨立的風暴段,記錄風暴段中的最大值和次大值的位置。假設其中的兩個風暴段為{Zi}={Zi-j,Zi-j+1,…,Zi}和{Zk}={Zk-j,Zk-j+1,…,Zk};令Mi=max{Zi},Mk=max{Zk},Mk-1為序列{Zk}中次最大值。

    (3)比較Mk-1和Mi。如果Mk-1>Mi,而且Mk-1和Mi之間的距離大于某個固定值s時,則將Mk-1加入最后選取的獨立風暴樣本。s取為根據(jù)實測車輛荷載數(shù)據(jù)統(tǒng)計特性得到的單輛重車長度。獨立風暴值的分布函數(shù)為:

    (7)

    (4)假設R年中有m個獨立風暴樣本,年平均風暴發(fā)生概率為r,且r=m/R。假設一年內(nèi)車輛運行情況變化不會太大,則設計基準期T年發(fā)生的風暴總數(shù)概率分布為:

    (8)

    3工程算例

    由于實際的車輛類型多樣,各生產(chǎn)廠家生產(chǎn)的車輛軸距也不盡相同,直接進行實測車輛的分析計算比較困難。根據(jù)WIM記錄的某橋一個月實際車流數(shù)據(jù),并參考計算橋梁疲勞時的車輛荷載類型分類,將車輛類型分為七類,各類型車輛所占總數(shù)量的比例及車輛簡圖見表1。各類型車輛軸重所占車重的比例以及軸距直接在車輛簡圖中表示。由實測的車輛荷載統(tǒng)計數(shù)據(jù)通過蒙特卡洛法模擬生成隨機車輛,并根據(jù)影響線加載[9]分別得到車輛在經(jīng)過30 ~ 60 m橋時的跨中彎矩值。

    表1 車輛分類及各車型比例統(tǒng)計

    以50 m跨徑簡支梁橋分析為例。由隨機車輛產(chǎn)生的簡支梁跨中彎矩值經(jīng)驗分布函數(shù)見圖1。由不同閾值得到的平均超出函數(shù)見圖2。由經(jīng)驗分布函數(shù)選定類閾值u1為6 000 kN·m,但根據(jù)平均超出函數(shù)圖得到的初步類閾值u2為10 000 kN·m,故初步選定類閾值為u2。累計分布函數(shù)值為0.99所對應的類閾值u3為11 046 kN·m。初步類閾值u2對應的MIS法樣本個數(shù)為140個,類閾值u3對應的MIS樣本數(shù)為48個。因為樣本個數(shù)相差很多,故重新選擇類閾值u2為10 500 kN·m,此時類閾值u2所對應的樣本數(shù)為85個。再次選擇類閾值u2為10 600 kN·m,對應的MIS樣本個數(shù)為80個。因為前后選擇的類閾值u2所對應的MIS樣本個數(shù)相差較小,故選擇閾值為10 500 kN·m所對應的樣本為分析樣本。

    圖1 簡支梁橋50 m跨中彎矩頻率直方圖 Fig.1 Bending moment frequency histograms of 50 m midspan of simply supported beam bridge

    圖2 經(jīng)驗平均超出函數(shù)圖 Fig.2 Empirical mean excess function

    簡支梁橋最大車輛荷載效應一般由重車控制。根據(jù)WIM記錄的車輛荷載統(tǒng)計特點,四軸車車重較大,而且相比五、六軸車出現(xiàn)的概率較大,因此本文選取四軸車所對應的橋梁影響長度為近似s值,s取為10 m。由MIS方法得到的獨立樣本數(shù)據(jù)共85個。GPD分布擬合參數(shù)為k=-0.129 7。由于k<0是個截尾分布,參考李植淮等[16]的建議,將樣本數(shù)據(jù)順序排列,取k=0,利用式(3)可得到σ值為847.531。由MMIS法得樣本數(shù)據(jù)116個。GPD分布擬合參數(shù)為k=-0.225 64。取k=0,得到σ為1 000.815。則MIS法最大值截口分布為:

    x>1.05×104

    (9)

    MMIS法得到的最大值截口分布為:

    x>1.05×104

    (10)

    由POT方法得到的最大值截口分布為:

    x>1.05×104

    (11)

    其它跨徑橋梁類閾值、MIS法和MMIS法樣本個數(shù)見表2。

    表2 簡支梁橋各跨徑樣本個數(shù)

    采用規(guī)范[15]給出的一天荷載效應最大值分布代表一年最大值分布,取分位值滿足FT(x)=0.95的x作為彎矩設計值,計算上述三種不同方法彎矩設計值與設計基準值之間的關(guān)系如圖3所示。

    圖3 三種分析方法的彎矩設計值與基準期之間的關(guān)系 Fig.3 Bending moment v.s. service life of three methods

    由圖3可知:①三種方法都可以反映出荷載效應極值隨基準期增加而增大的基本變化規(guī)律。②POT方法由于數(shù)據(jù)具有一定的相關(guān)性,其結(jié)果明顯偏??;而MIS方法雖然選擇的數(shù)據(jù)近似獨立,但是由于選擇的數(shù)據(jù)樣本較少,其結(jié)果也偏小,但較POT估計結(jié)果準確。③隨著橋梁跨徑的增加,POT法和MIS法評估得到的車輛荷載效應極值相差越來越小。造成這種現(xiàn)象的原因是,隨著橋梁跨徑的增加,車輛荷載效應最大值不再主要是由單輛重車控制引起,而是要考慮多車同時作用的影響。④這三種方法都可以較好地反映樣本的尾部變化規(guī)律,在較短的周期(T<20年),三種方法估計結(jié)果相差不多。因此,在進行舊橋評估時,三種方法都可使用,因為舊橋的剩余壽命一般少于20年。⑤在中長期(T>20年)荷載效應極值估計時,POT和MIS法估計結(jié)果明顯偏低。因此,在對新建橋梁或擬建橋梁荷載效應計算時,本文提出的MMIS法估計結(jié)果更為保守,也偏于安全。

    4結(jié)論

    由于車輛荷載效應樣本是平穩(wěn)時間序列,具有一定的相關(guān)性。為了解決數(shù)據(jù)間的相關(guān)性問題,本文對MIS法進行了改進,提出了MMIS法。同時提出了一種三次類閾值法,通過該閾值法來得到MMIS法的樣本。并以WIM記錄的某橋一個月實際樣本數(shù)據(jù)進行了分析,得到了以下結(jié)論:

    (1)本文提出的三次類閾值法,有效地克服了POT方法中閾值需要精確選擇的問題。根據(jù)樣本經(jīng)驗分布函數(shù)進行類閾值u1的選擇;通過經(jīng)驗平均超出函數(shù)法進行類閾值u2的初步選擇,并結(jié)合MIS法具有對閾值的穩(wěn)定性較好的特點進行類閾值u2的最終確定;并將類閾值u3作為閾值選擇上限,取樣本累積分布函數(shù)值為0.99所對應的值。該方法程序?qū)崿F(xiàn)簡單,收斂效率較高。

    (2)結(jié)合簡支梁橋車輛荷載效應主要受單輛重車影響的特點,對MIS法不能充分利用原始樣本數(shù)據(jù)的缺陷進行了改進,提出了MMIS法。由于MMIS法的分析樣本不僅相互獨立,而且相比MIS法數(shù)據(jù)利用更加完全,因此得到的分析結(jié)果也更為可靠。

    (3)通過對某橋WIM記錄的一個月實際車流量分析得到的車輛荷載效應值進行分析,得到以下結(jié)果:①若評估期較短,在對舊橋承載力評估或者確定舊橋加固方案時,理論上三種方法的模擬結(jié)果相近,都可以使用。但是相比POT方法閾值的選擇、MIS方法中的平靜點確定,本文提出的用三次類閾值法進行樣本數(shù)據(jù)的選擇更為簡單,而且也較為準確。②若評估期較

    長,在對擬建橋梁或新建橋梁服役期間進行最大荷載效應估計時,POT方法和MIS方法可能會低估車輛荷載效應,隱藏著車輛設計荷載偏小的危險;而利用MMIS方法估計未來設計期限內(nèi)的荷載效應值則較為保守,偏于安全。

    參考文獻

    [1]Sivakumar B, Ghosn M, Moses F. Protocols for collecting and using traffic data in bridge design[M]. Washington: Transportation Research Board, 2011.

    [2]Kozikowski M. WIM based live load model for bridge reliability[D]. Nebraska: Dissertation to the University of Nebraska, 2009.

    [3]Nowak A S. WIM based lived load model for bridges[R]. Florida: Published in Transportation Board, 2011.

    [4]Galambos J. The asymptotic theory of extreme order statistics[M]. Malabar: Krieger, 1987.

    [5]Cook N J. Towards better estimation of extreme winds[J]. Journal of Wind Engineering and Industrial Aerodynamics, 1982, 9(3): 295-323.

    [6]Bhattacharya B. The extremal index and the maximum of a dependent stationary pulse load process observed above a high threshold[J]. Structural Safety, 2008, 30(1): 34-48.

    [8]Mazas F, Hamm L. A multi-distribution approach to POT methods for determining extreme wave heights[J]. Coastal Engineering, 2011, 58(5): 385-394.

    [9]李春前. 基于GPD模型的車輛荷載效應極值估計[D]. 北京: 清華大學,2012.

    [10]Wang Y, Xia Y, Liu X. Establishing robust short-term distributions of load extremes of offshore wind turbines[J]. Renewable Energy, 2013, 57: 606-619.

    [11]史道濟. 實用極值統(tǒng)計方法[M]. 天津: 天津科學技術(shù)出版社, 2006.

    [12]Caers J, Maes M A. Identifying tails, bounds and end-points of random variables[J]. Structural Safety, 1998, 20(1): 1-23.

    [13]趙旭. 廣義Pareto分布的統(tǒng)計推斷[D]. 北京: 北京工業(yè)大學, 2012.

    [14]An Y, Pandey M D. A comparison of methods of extreme wind speed estimation[J]. Journal of Wind Engineering and Industrial Aerodynamics, 2005, 93(7): 535-545.

    [15]中華人民共和國住房和城鄉(xiāng)建設部. 工程結(jié)構(gòu)可靠性設計統(tǒng)一標準(GB 50153-2008) [S]. 北京: 中國建筑工業(yè)出版社.

    [16]李植淮,李春前,孫健康,等.基于GPD模型的車輛荷載效應極值估計[J].工程力學,2012,29(增刊1):166-171.

    LI Zhi-huai, LI Chun-qian, SUN Jian-kang, et al. Estimation of extreme vehicle load effect based on GPD model[J]. Engineering Mechanics, 2012,29(Sup1): 166-171.

    第一作者符毅強男,碩士生,1991年生

    通信作者陳予恕男,教授,中國工程院院士,1931年生

    郵箱:yschen@hit.edu.cn

    第一作者張平男,博士,講師,1977年生

    通信作者馬建男,博士,教授,博士生導師,1957年生

    猜你喜歡
    閾值
    采用紅細胞沉降率和C-反應蛋白作為假體周圍感染的閾值
    小波閾值去噪在深小孔鉆削聲發(fā)射信號處理中的應用
    基于CS-TWR的動態(tài)閾值貪婪算法成像研究
    基于自適應閾值和連通域的隧道裂縫提取
    比值遙感蝕變信息提取及閾值確定(插圖)
    河北遙感(2017年2期)2017-08-07 14:49:00
    室內(nèi)表面平均氡析出率閾值探討
    基于遲滯比較器的雙閾值穩(wěn)壓供電控制電路
    基于改進閾值去噪的諧波檢測
    基于新閾值函數(shù)的小波閾值降噪方法
    基于小波閾值去噪的改進方案
    av在线老鸭窝| 日韩在线高清观看一区二区三区| 叶爱在线成人免费视频播放| 天天躁夜夜躁狠狠久久av| 久久久久精品久久久久真实原创| 欧美国产精品一级二级三级| av女优亚洲男人天堂| 在线亚洲精品国产二区图片欧美| 日本黄色日本黄色录像| 国产一区二区三区av在线| www.熟女人妻精品国产| 曰老女人黄片| 久热这里只有精品99| 精品国产露脸久久av麻豆| 五月伊人婷婷丁香| 日日摸夜夜添夜夜爱| 狂野欧美激情性bbbbbb| 日韩精品免费视频一区二区三区| 交换朋友夫妻互换小说| 日韩中文字幕视频在线看片| 女性被躁到高潮视频| 女人高潮潮喷娇喘18禁视频| 国产精品一区二区在线观看99| 69精品国产乱码久久久| 日日摸夜夜添夜夜爱| 美女视频免费永久观看网站| 免费少妇av软件| 欧美激情 高清一区二区三区| 人人妻人人澡人人爽人人夜夜| 日本免费在线观看一区| 久久久国产一区二区| 日本欧美国产在线视频| 亚洲精品久久成人aⅴ小说| 寂寞人妻少妇视频99o| 中国三级夫妇交换| www日本在线高清视频| 伊人久久国产一区二区| 国产av码专区亚洲av| 精品久久久久久电影网| 9191精品国产免费久久| av在线app专区| 久久女婷五月综合色啪小说| av免费在线看不卡| www日本在线高清视频| 免费少妇av软件| av.在线天堂| 久久久久久伊人网av| 嫩草影院入口| 亚洲第一av免费看| 91午夜精品亚洲一区二区三区| 国产激情久久老熟女| 五月天丁香电影| 99久久人妻综合| 国产成人午夜福利电影在线观看| 不卡视频在线观看欧美| 久久久久久久久久久久大奶| 欧美bdsm另类| 欧美精品国产亚洲| 精品少妇久久久久久888优播| 精品久久久精品久久久| 国产精品久久久久成人av| 中国三级夫妇交换| 久久99热这里只频精品6学生| av女优亚洲男人天堂| 韩国高清视频一区二区三区| 秋霞伦理黄片| 国产高清不卡午夜福利| 久久99精品国语久久久| 日韩,欧美,国产一区二区三区| 人妻一区二区av| 欧美xxⅹ黑人| 精品一区二区三区四区五区乱码 | 亚洲国产欧美网| 欧美成人午夜精品| 久久这里有精品视频免费| 人成视频在线观看免费观看| 成年女人毛片免费观看观看9 | 晚上一个人看的免费电影| 亚洲精品久久成人aⅴ小说| 国产成人午夜福利电影在线观看| 黄频高清免费视频| 亚洲欧洲国产日韩| 精品一区在线观看国产| 最近手机中文字幕大全| 天天躁夜夜躁狠狠躁躁| 91精品三级在线观看| 中文欧美无线码| 国产激情久久老熟女| av不卡在线播放| 国产亚洲av片在线观看秒播厂| 欧美日本中文国产一区发布| 亚洲av中文av极速乱| 欧美精品国产亚洲| 热99国产精品久久久久久7| 18禁观看日本| 99re6热这里在线精品视频| 亚洲av欧美aⅴ国产| 亚洲av成人精品一二三区| 亚洲av男天堂| 国产精品一区二区在线不卡| 国产精品久久久久成人av| 激情视频va一区二区三区| 极品少妇高潮喷水抽搐| 久久毛片免费看一区二区三区| av有码第一页| freevideosex欧美| 久久久久国产一级毛片高清牌| 亚洲 欧美一区二区三区| 不卡av一区二区三区| 精品人妻一区二区三区麻豆| 国产伦理片在线播放av一区| 人成视频在线观看免费观看| 一级毛片 在线播放| 观看美女的网站| a级毛片黄视频| 亚洲少妇的诱惑av| 一本一本久久a久久精品综合妖精 国产伦在线观看视频一区 | 美女主播在线视频| 欧美另类一区| www日本在线高清视频| 亚洲少妇的诱惑av| 国产黄色免费在线视频| 人人妻人人添人人爽欧美一区卜| 韩国精品一区二区三区| 天天操日日干夜夜撸| 色哟哟·www| 日韩中文字幕欧美一区二区 | 午夜福利乱码中文字幕| 丝袜人妻中文字幕| 9191精品国产免费久久| 国产xxxxx性猛交| 一级,二级,三级黄色视频| 青青草视频在线视频观看| 欧美日韩亚洲高清精品| 久久女婷五月综合色啪小说| 国产片内射在线| 国产无遮挡羞羞视频在线观看| 亚洲中文av在线| 亚洲av男天堂| 日韩一区二区三区影片| 黄色配什么色好看| 一级片'在线观看视频| 欧美成人午夜精品| 香蕉丝袜av| 国产男女内射视频| 麻豆精品久久久久久蜜桃| 国产一区二区激情短视频 | 欧美97在线视频| 亚洲精品自拍成人| 午夜免费鲁丝| 欧美在线黄色| 久久97久久精品| 国产爽快片一区二区三区| 国产有黄有色有爽视频| 高清av免费在线| 妹子高潮喷水视频| 大话2 男鬼变身卡| 人人妻人人添人人爽欧美一区卜| 国产xxxxx性猛交| 久久精品亚洲av国产电影网| 国产精品一区二区在线不卡| 国产精品久久久久久av不卡| 久久久精品国产亚洲av高清涩受| 少妇人妻久久综合中文| 人人妻人人添人人爽欧美一区卜| 久久久久精品性色| 亚洲国产毛片av蜜桃av| 1024香蕉在线观看| 精品人妻偷拍中文字幕| 国产97色在线日韩免费| 卡戴珊不雅视频在线播放| 婷婷成人精品国产| 韩国精品一区二区三区| 十八禁网站网址无遮挡| 黄片无遮挡物在线观看| 99热网站在线观看| 99九九在线精品视频| 高清在线视频一区二区三区| 永久免费av网站大全| 天堂8中文在线网| 亚洲伊人色综图| 亚洲av欧美aⅴ国产| 日本爱情动作片www.在线观看| 精品亚洲成a人片在线观看| 欧美精品人与动牲交sv欧美| 中文字幕av电影在线播放| 一级爰片在线观看| 超色免费av| 在线观看免费高清a一片| 国产xxxxx性猛交| 香蕉国产在线看| 国产在视频线精品| 午夜免费男女啪啪视频观看| 国产精品秋霞免费鲁丝片| 中文乱码字字幕精品一区二区三区| 国产成人一区二区在线| 精品一区二区免费观看| 久久久久网色| 国产男女超爽视频在线观看| 校园人妻丝袜中文字幕| 丁香六月天网| 午夜av观看不卡| 精品国产国语对白av| 国产又色又爽无遮挡免| 国产极品粉嫩免费观看在线| 三级国产精品片| 成人亚洲欧美一区二区av| 夫妻午夜视频| 国产野战对白在线观看| 最新的欧美精品一区二区| 9191精品国产免费久久| 热99久久久久精品小说推荐| 精品久久蜜臀av无| 亚洲三区欧美一区| 久久国产精品大桥未久av| 男女国产视频网站| 久久精品国产a三级三级三级| 麻豆av在线久日| 亚洲成人手机| 尾随美女入室| www.精华液| 熟妇人妻不卡中文字幕| 永久免费av网站大全| 免费观看无遮挡的男女| 国产精品久久久av美女十八| 又大又黄又爽视频免费| 国产成人aa在线观看| 午夜福利在线免费观看网站| 伦理电影免费视频| av线在线观看网站| 这个男人来自地球电影免费观看 | 男人爽女人下面视频在线观看| 日本黄色日本黄色录像| av网站在线播放免费| 在线亚洲精品国产二区图片欧美| 在线观看免费日韩欧美大片| 午夜激情av网站| 国产精品一区二区在线观看99| 久久ye,这里只有精品| 自线自在国产av| 天天躁夜夜躁狠狠久久av| 午夜日韩欧美国产| 熟女av电影| 一本色道久久久久久精品综合| 亚洲av综合色区一区| 最近最新中文字幕大全免费视频 | 欧美精品av麻豆av| 新久久久久国产一级毛片| 欧美精品一区二区免费开放| 少妇 在线观看| 国产精品av久久久久免费| 男女高潮啪啪啪动态图| 国产日韩欧美亚洲二区| 日本av免费视频播放| 美女国产高潮福利片在线看| 免费观看在线日韩| 日本爱情动作片www.在线观看| 亚洲精品久久成人aⅴ小说| 亚洲国产av影院在线观看| 韩国av在线不卡| av.在线天堂| 久久97久久精品| 国产又色又爽无遮挡免| 日韩三级伦理在线观看| 久久久久久人人人人人| 在线看a的网站| 啦啦啦视频在线资源免费观看| 另类亚洲欧美激情| 视频区图区小说| 一级黄片播放器| 免费观看在线日韩| 91精品国产国语对白视频| 国产午夜精品一二区理论片| 久久久久视频综合| 天天躁夜夜躁狠狠躁躁| 大香蕉久久成人网| 亚洲精品日本国产第一区| 欧美最新免费一区二区三区| 美女福利国产在线| 亚洲第一青青草原| 久久久国产一区二区| 亚洲国产精品国产精品| av在线观看视频网站免费| 中文乱码字字幕精品一区二区三区| 久久女婷五月综合色啪小说| 一二三四中文在线观看免费高清| 黄片小视频在线播放| 啦啦啦在线观看免费高清www| 久久精品国产鲁丝片午夜精品| 欧美精品一区二区免费开放| 国产精品一二三区在线看| 人成视频在线观看免费观看| 国产老妇伦熟女老妇高清| 人人妻人人澡人人爽人人夜夜| 在线观看美女被高潮喷水网站| 一区二区日韩欧美中文字幕| 色哟哟·www| 欧美精品高潮呻吟av久久| 26uuu在线亚洲综合色| 日韩一区二区三区影片| 成人毛片60女人毛片免费| 天天躁日日躁夜夜躁夜夜| 九草在线视频观看| 黑人欧美特级aaaaaa片| 日本欧美视频一区| 久久久久久伊人网av| 999久久久国产精品视频| 伊人亚洲综合成人网| 一本大道久久a久久精品| 国产白丝娇喘喷水9色精品| 综合色丁香网| 99国产精品免费福利视频| 91久久精品国产一区二区三区| 日韩视频在线欧美| 桃花免费在线播放| 久久久久久人妻| 日韩人妻精品一区2区三区| 亚洲,一卡二卡三卡| av国产久精品久网站免费入址| 亚洲色图综合在线观看| 欧美日韩亚洲高清精品| 国产熟女欧美一区二区| 亚洲国产欧美在线一区| 日韩一区二区三区影片| 高清视频免费观看一区二区| 亚洲美女视频黄频| 亚洲欧美中文字幕日韩二区| 91精品三级在线观看| 日韩三级伦理在线观看| 三级国产精品片| 国产精品.久久久| 欧美精品国产亚洲| 亚洲成人一二三区av| 欧美激情高清一区二区三区 | 亚洲精品乱久久久久久| 建设人人有责人人尽责人人享有的| 亚洲经典国产精华液单| 国产亚洲午夜精品一区二区久久| 飞空精品影院首页| 最近最新中文字幕大全免费视频 | 成年美女黄网站色视频大全免费| 日韩免费高清中文字幕av| 亚洲人成77777在线视频| 国产成人免费无遮挡视频| 亚洲国产日韩一区二区| 国产精品女同一区二区软件| 国产不卡av网站在线观看| 国产精品欧美亚洲77777| 天天影视国产精品| 亚洲成国产人片在线观看| 十八禁高潮呻吟视频| 久久久久精品久久久久真实原创| 亚洲图色成人| 十分钟在线观看高清视频www| 亚洲图色成人| 99久久中文字幕三级久久日本| 久久精品国产亚洲av天美| 日韩伦理黄色片| 女人被躁到高潮嗷嗷叫费观| 少妇人妻 视频| 亚洲国产成人一精品久久久| 亚洲美女视频黄频| 老熟女久久久| 国产麻豆69| 久久人人97超碰香蕉20202| 久久久久精品久久久久真实原创| www.自偷自拍.com| 亚洲av成人精品一二三区| 欧美成人午夜免费资源| 日本wwww免费看| 久久人人97超碰香蕉20202| 又大又黄又爽视频免费| 国产爽快片一区二区三区| 免费在线观看黄色视频的| 在线亚洲精品国产二区图片欧美| 国产麻豆69| 精品亚洲乱码少妇综合久久| 中文字幕色久视频| 97在线人人人人妻| 国产成人免费观看mmmm| 精品亚洲成国产av| 桃花免费在线播放| tube8黄色片| 国产成人免费观看mmmm| 丁香六月天网| 亚洲内射少妇av| 日韩一区二区视频免费看| 有码 亚洲区| 日韩一区二区视频免费看| 久久久精品94久久精品| 日韩在线高清观看一区二区三区| 黄色怎么调成土黄色| 99re6热这里在线精品视频| 看免费av毛片| 成人亚洲精品一区在线观看| 丰满少妇做爰视频| 亚洲第一区二区三区不卡| 最近手机中文字幕大全| 男女免费视频国产| 亚洲精品日本国产第一区| av网站免费在线观看视频| 日韩 亚洲 欧美在线| 午夜日本视频在线| 国产日韩欧美视频二区| 日产精品乱码卡一卡2卡三| 啦啦啦啦在线视频资源| 日韩中字成人| videos熟女内射| 爱豆传媒免费全集在线观看| 美女中出高潮动态图| 色哟哟·www| 1024香蕉在线观看| 午夜免费鲁丝| 久久国产精品大桥未久av| 亚洲精品国产av成人精品| 男女免费视频国产| 狂野欧美激情性bbbbbb| 国产黄频视频在线观看| 飞空精品影院首页| 精品一区二区三区四区五区乱码 | 又黄又粗又硬又大视频| 国产福利在线免费观看视频| 少妇 在线观看| 91aial.com中文字幕在线观看| 如何舔出高潮| av免费在线看不卡| 日韩电影二区| 精品第一国产精品| 日本黄色日本黄色录像| 亚洲婷婷狠狠爱综合网| 中文字幕另类日韩欧美亚洲嫩草| 水蜜桃什么品种好| 国产精品成人在线| 晚上一个人看的免费电影| 久久这里只有精品19| 波多野结衣av一区二区av| 久久久国产一区二区| 国产精品国产av在线观看| 观看美女的网站| 看十八女毛片水多多多| 亚洲精品国产av成人精品| 中文天堂在线官网| 久久久精品区二区三区| 日韩精品有码人妻一区| 欧美av亚洲av综合av国产av | 97精品久久久久久久久久精品| 久久99蜜桃精品久久| 午夜精品国产一区二区电影| 国产一区二区在线观看av| 欧美av亚洲av综合av国产av | 亚洲精品第二区| 免费观看a级毛片全部| 国产深夜福利视频在线观看| 一级片'在线观看视频| 18+在线观看网站| 观看美女的网站| 国产成人精品在线电影| 老司机亚洲免费影院| 亚洲精品日韩在线中文字幕| 狠狠精品人妻久久久久久综合| 亚洲伊人久久精品综合| 国产在线视频一区二区| 男男h啪啪无遮挡| 国产熟女欧美一区二区| 一级a爱视频在线免费观看| 成人手机av| 18禁观看日本| 亚洲综合色网址| 另类亚洲欧美激情| 少妇猛男粗大的猛烈进出视频| 免费观看av网站的网址| 一边摸一边做爽爽视频免费| 欧美日韩一级在线毛片| 久久久久久久亚洲中文字幕| 成人18禁高潮啪啪吃奶动态图| 女人久久www免费人成看片| 久久国内精品自在自线图片| 麻豆av在线久日| 久久久久久久大尺度免费视频| 成人国产av品久久久| 亚洲,欧美,日韩| 国产成人a∨麻豆精品| 欧美xxⅹ黑人| 激情视频va一区二区三区| 国产精品秋霞免费鲁丝片| 亚洲三级黄色毛片| 午夜免费鲁丝| av.在线天堂| 亚洲精品成人av观看孕妇| 免费在线观看黄色视频的| 国产成人精品一,二区| 国产人伦9x9x在线观看 | 国产老妇伦熟女老妇高清| 免费人妻精品一区二区三区视频| 亚洲图色成人| 一级毛片电影观看| 老司机亚洲免费影院| 亚洲精品,欧美精品| 国产在线一区二区三区精| 日本猛色少妇xxxxx猛交久久| 亚洲精品美女久久久久99蜜臀 | 久久久久久久精品精品| av.在线天堂| 青春草国产在线视频| 日韩制服骚丝袜av| 亚洲精品久久成人aⅴ小说| 看免费av毛片| 在线看a的网站| 免费看av在线观看网站| 日韩大片免费观看网站| 亚洲精品一二三| 久久午夜福利片| 老鸭窝网址在线观看| 免费黄色在线免费观看| 日韩精品有码人妻一区| h视频一区二区三区| 日本爱情动作片www.在线观看| 免费观看av网站的网址| 另类精品久久| 一本—道久久a久久精品蜜桃钙片| 日韩在线高清观看一区二区三区| 久久精品国产综合久久久| 久久久久久久久久久久大奶| 成年av动漫网址| 欧美日韩亚洲国产一区二区在线观看 | 91国产中文字幕| 国产日韩欧美亚洲二区| 精品第一国产精品| 你懂的网址亚洲精品在线观看| 建设人人有责人人尽责人人享有的| 亚洲综合色惰| 亚洲,欧美精品.| 精品久久久精品久久久| 91国产中文字幕| 丝袜喷水一区| 国产综合精华液| 久久人妻熟女aⅴ| 超碰97精品在线观看| 久久久久视频综合| 久久久久精品人妻al黑| av女优亚洲男人天堂| 中文字幕人妻熟女乱码| 在线观看免费视频网站a站| 久久国产亚洲av麻豆专区| 国产成人精品无人区| 久久久久久人人人人人| 男女边吃奶边做爰视频| 日韩熟女老妇一区二区性免费视频| 九色亚洲精品在线播放| 亚洲欧洲国产日韩| 国产成人午夜福利电影在线观看| 精品久久久精品久久久| 两性夫妻黄色片| 波多野结衣av一区二区av| 观看美女的网站| 亚洲av日韩在线播放| 999精品在线视频| 午夜激情久久久久久久| 欧美激情极品国产一区二区三区| 777久久人妻少妇嫩草av网站| 视频区图区小说| 国产午夜精品一二区理论片| 午夜老司机福利剧场| 黄片无遮挡物在线观看| 国产xxxxx性猛交| 高清不卡的av网站| 男女边摸边吃奶| 国产在线一区二区三区精| 久久午夜综合久久蜜桃| 亚洲第一区二区三区不卡| 如何舔出高潮| 国产精品欧美亚洲77777| 国产男女超爽视频在线观看| 久久久精品区二区三区| 中文精品一卡2卡3卡4更新| 亚洲国产日韩一区二区| 国产精品99久久99久久久不卡 | 午夜老司机福利剧场| 欧美精品一区二区免费开放| 美女中出高潮动态图| 亚洲成人手机| 一区二区三区精品91| 亚洲欧美色中文字幕在线| 亚洲色图综合在线观看| 亚洲激情五月婷婷啪啪| 日本91视频免费播放| 久久青草综合色| 日韩av不卡免费在线播放| videos熟女内射| 亚洲精品久久午夜乱码| 欧美精品人与动牲交sv欧美| 国产乱人偷精品视频| 国产综合精华液| 亚洲精品aⅴ在线观看| 亚洲av在线观看美女高潮| 久久久久久久久免费视频了| 亚洲精品aⅴ在线观看| 久久精品国产鲁丝片午夜精品| 久久人人爽av亚洲精品天堂| 人成视频在线观看免费观看| 日日啪夜夜爽| 日本av手机在线免费观看| 欧美精品人与动牲交sv欧美| 极品少妇高潮喷水抽搐| 色视频在线一区二区三区| 天天操日日干夜夜撸| 日日啪夜夜爽| 最近中文字幕高清免费大全6| 男女边摸边吃奶| 国产精品久久久av美女十八| 在线观看免费高清a一片| 大香蕉久久成人网| 精品视频人人做人人爽| 亚洲国产最新在线播放| 汤姆久久久久久久影院中文字幕|