亚洲免费av电影一区二区三区,日韩爱爱视频,51精品视频一区二区三区,91视频爱爱,日韩欧美在线播放视频,中文字幕少妇AV,亚洲电影中文字幕,久久久久亚洲av成人网址,久久综合视频网站,国产在线不卡免费播放

        ?

        基于材料參數更新的RC框架結構時變易損性分析*

        2021-08-26 01:37:50楊思昭王憲杰
        建筑結構 2021年15期
        關鍵詞:經時易損性服役

        楊思昭,王憲杰,2

        (1 云南大學建筑與規(guī)劃學院,昆明 650500;2 浙江大學建筑工程學院,杭州 310058)

        0 概述

        地震易損性分析能夠預測在不同強度地震作用下,建筑結構達到或超越不同性能水準的概率[1-2]。但傳統的易損性分析主要是靜態(tài)的,未考慮材料性能衰退對結構抗震性能的影響。事實上,結構在全壽命周期內,隨著服役齡期的增加,由于受力構件有效截面面積的減小和材料力學性能的退化,結構在同一地面運動加速度作用下發(fā)生破壞的概率增加[3-4]。因此,在進行地震易損性評估的研究中,考慮材料時變特性的影響更能真實反映不同服役期建筑結構的地震易損性。

        但由于耐久性研究存在時間跨度長、影響因素復雜、不確定性大等問題,目前尚未有統一、可靠的理論模型來刻畫材料、構件和結構層面的長周期、隨機衰變特性。與此同時,在實際工程的檢測、監(jiān)測和檢修等作業(yè)中,已取得了服役結構的諸多性能參數。但實際工程中的實測數據(特別是有損檢測)由于樣本數量少而缺乏足夠的代表性,僅依靠小樣本檢測數據評估建筑結構的實際工作性能缺乏足夠的完備性[5]。結合現有的實測數據,引入信息更新理論對初始設計模型不斷修正和完善,可實現對具體結構性能退化的合理預測。針對實時信息更新與系統吸收,經典貝葉斯方法可以有效實現對先驗分布的相關參數進行重新估計和檢驗,得到系統后驗分布數據[6-7]。因此,采用貝葉斯方法進行結構耐久性分析既能利用多年來耐久性基礎研究對結構內在規(guī)律的理性認識,又能考慮結構的實際情況[8],從而準確刻化結構實時動態(tài)性能。

        退化結構的時變結構抗力是非平穩(wěn)隨機過程,要求必須以隨機過程理論為基礎,研究隨機荷載效應超過隨機抗力的概率[9]。即通過理論分析研究結構的時變易損性時,同時需要考慮結構自身參數的不確定性和外界激勵的隨機性。此外,在傳統的地震易損性研究中,為了便于問題的求解,通常需要基于以下幾點假定進行分析[2-4,10]:結構地震需求參數和給定地震動水平作用之間服從冪指數回歸關系;在給定地震動水平作用下,結構地震需求參數服從對數正態(tài)分布;結構的概率抗震能力服從對數正態(tài)分布。顯然,基于此求解得到的地震易損性結果只是經驗估計值,不能反映結構響應的真實概率分布。李杰[11]提出的概率密度演化理論不僅為復合隨機問題的求解提供了理論依據,而且可以分析得到滿足真實分布的結構響應概率密度分布函數,為精細化評估既有結構的地震易損性提供了新思路。

        本文通過總結建筑材料耐久性的相關研究成果,建立材料力學性能退化模型;引入不同服役齡期材料力學性能的實測數據,并采用貝葉斯方法進行數據更新。基于更新的材料性能退化模型與概率密度演化理論對某7層帶內廊式對稱三跨RC框架結構的地震易損性進行分析,得到RC框架結構抗震性能隨服役齡期的一般變化規(guī)律。

        1 材料性能時變特性分析與貝葉斯方法

        1.1 隨機參數更新貝葉斯方法

        已有研究[12-13]表明,一般大氣環(huán)境中既有建筑物混凝土的抗壓強度服從正態(tài)分布,其均值和標準差是服役時間t的函數。采用共軛先驗分布確定均值μ的后驗分布。當方差σ2已知時,均值μ的先驗分布可視為服從正態(tài)分布N(θ,τ2),即:

        (1)

        式中:θ,τ為超參數,可根據先驗預測公式進行確定。

        先驗分布反映對研究對象內在規(guī)律的理性認識,不可避免的存在主觀判斷,尚需考慮其實際情況。設t時刻材料物理性能參數的實測樣本值{x1,x2,x3,…,xn}取自總體樣本N(μ,σ2),其中方差σ2已知。則此樣本x的似然函數可表示為:

        (2)

        進而,由貝葉斯方法可將結構物理參數μ的后驗分布表示為:

        (3)

        將式(1),(2)帶入式(3),經過一定的數學推導可得μ的后驗分布為:

        (4)

        (5)

        圖1給出了采用貝葉斯方法對材料性能參數進行信息更新的一般流程。

        圖1 材料性能參數貝葉斯動態(tài)信息更新過程

        1.2 混凝土經時抗壓強度預測模型

        一般而言,混凝土經時抗壓強度的均值在初期隨時間增大而增大,但增長速率逐漸減緩,在后期則隨時間增大而降低;標準差隨時間增大而增大,且不受均值影響。牛荻濤等[12]通過分析長期暴露試驗和實測的服役混凝土強度,建立了一般大氣環(huán)境中混凝土強度均值μf(t)和標準差σf(t)的經時變換模型:

        μf(t)=1.452 9exp[-0.024 6(lnt-

        1.715 4)2]·μ′0

        (6)

        σf(t)=(0.030 5t+1.236 8)·σ′0

        (7)

        式中μ′0和σ′0分別為混凝土28d強度均值和標準差。

        本文依據《混凝土結構設計規(guī)范》(GB 50010—2010)給出的混凝土抗壓強度均值和標準值之間的關系,將混凝土強度均值μf(t)的經時變化模型轉換為混凝土強度標準值fcu(t)的經時變化模型:

        fcu(t)=μf(t)-1.645σf(t)

        (8)

        1.3 鋼筋材料性能時變模型

        相較于混凝土經時抗壓強度的研究,鋼筋經時銹蝕深度已有較為成熟的研究成果[14],且已在規(guī)范中進行體現[15]。處于一般大氣環(huán)境中的RC框架結構,鋼筋銹蝕開始時間ti可按下式進行計算:

        (9)

        式中:C為混凝土保護層厚度;k為碳化系數;x0為碳化殘量。具體參數取值詳見規(guī)范[15]。

        鋼筋銹蝕是一個電化學腐蝕過程,其銹蝕速率與鋼筋表面的含氧量有關,而含氧量取決于混凝土質量、保護層厚度和環(huán)境條件。在溫度為20℃,相對濕度為75%的典型環(huán)境中,鋼筋初始銹蝕電流密度icoor(1)為[14]:

        (10)

        式中w/c為水灰比,當已知混凝土抗壓強度標準值fcu(t)時,水灰比可按下式計算:

        (11)

        進而,鋼筋開始銹蝕之后的某時刻tc,鋼筋銹蝕電流密度icoor(tc)為:

        (12)

        通常情況下,當銹蝕電流密度icoor為1μA/cm2時,相應的鋼筋銹蝕深度速率為11.6μm/年。

        2 經時材料本構關系定義

        2.1 經時混凝土的本構關系

        相較于大量已有的時不變或素混凝土本構關系模型,考慮銹蝕箍筋約束影響建立的經時混凝土本構關系更符合實際服役情況。在考慮了體積配箍率、箍筋屈服強度等影響因素建立的Kent-Scott-Park模型基礎上,通過引入峰值應力和應變軟化修正系數[16]建立了銹蝕矩形箍筋約束混凝土本構關系模型,如式(13)所示:

        (13)

        式中:σc和εc分別為銹蝕箍筋約束混凝土的應力和應變;εcc為銹蝕箍筋約束混凝土的峰值點應變,εcc=0.002c2K;εcu為銹蝕箍筋約束混凝土的極限應變;c1,c2和c3分別為非銹蝕箍筋約束混凝土的峰值應力、峰值點應變和應變軟化調整系數,按式(14)確定;k1,k2分別為考慮箍筋銹蝕影響的峰值應力和應變軟化修正系數,按式(15)確定;K為考慮箍筋約束影響的混凝土強度增強系數,按式(16)確定;Zm為應變軟化段斜率,按式(17)確定;fc(t)為混凝土的軸心抗壓強度。

        c1=0.157 0λt+0.963 4

        (14-1)

        c2=-1.556 8λt+1.232 0

        (14-2)

        c3=17.744 0λt+0.974 2

        (14-3)

        k1=(0.291 4λt-0.179 1)ω+1.0

        (15-1)

        k2=exp[(-1.585 2λt+0.960 8)ω]

        (15-2)

        (16)

        (17)

        式中:λt為箍筋特征值;ω為銹脹裂縫寬度;ρs為體積配箍率;fyh為箍筋屈服強度;h″為從箍筋外邊緣算起的核心混凝土寬度;sh為箍筋間距。

        2.2 銹蝕鋼筋的本構關系

        實際工程中的銹蝕鋼筋,當其銹蝕率ηs<80%時,本構關系可按式(18)進行確定[17]:

        (18)

        式中:σs和εs分別為銹蝕鋼筋的應力和應變;Es0為未銹蝕鋼筋的彈性模量;fyc,fuc分別為銹蝕鋼筋的屈服強度和極限強度,按式(19)確定;εsyc,εshc和εsuc分別為銹蝕鋼筋的屈服應變、強化應變和極限應變,按式(20)確定。

        (19-1)

        (19-2)

        (20-1)

        (20-2)

        εsuc=e-2.501ηsεsu0

        (20-3)

        式中:fy0,fu0分別為未銹蝕鋼筋的屈服強度和極限強度;εsh0,εsu0分別為未銹蝕鋼筋的強化應變和極限應變;ηs,cr為鋼筋屈服平臺消失時的臨界銹蝕率,帶肋鋼筋取20%,光圓鋼筋取10%。

        3 基于PDEM的結構時變易損性分析

        時變易損性表征既有建筑結構在不同服役時期遭遇可能強度地震作用時,結構達到或超越各級性能水準的概率??杀硎緸椋?/p>

        FR(a,t)=Pf[PL|SA=a,T=t]

        (21)

        式中:FR(a,t)為結構的時變易損性;Pf為失效概率;PL為結構的性能水平;SA為地震動強度指標,SA=a;T為結構服役期,T=t。

        顯然,結構響應θ(SA,T)達到或超過目標性能水準θc的概率可表示為:

        Pf[SA,T]=Pr[θc≤θ(SA,T)]

        (22)

        通過以上分析可知,劃分合理可信的目標性能水準和較為精細化求解結構在不同地震動強度和不同服役年限條件下的結構響應是時變易損性分析的基礎,以下分別就這兩方面進行討論。

        3.1 概率抗震能力分析

        《建筑抗震設計規(guī)范》(GB 50011—2010)蘊含豐富的基于性能抗震設計思想,其目標性能的選定依托于對性能水平的合理劃分。針對RC框架結構,可采用結構層間位移角來定義結構的性能水平,將RC框架結構不同破壞程度對應的最大層間位移角限值列于表1。在此指出,為了便于問題的考慮,本文將各目標破壞狀態(tài)限值視為確定量進行考慮。

        破壞等級與最大層間位移角限值的關系 表1

        3.2 基于PDEM的概率地震需求分析

        傳統的易損性分析方法通常假定概率地震需求參數服從對數正態(tài)分布,是基于經驗的一種近似分析方法。區(qū)別于傳統方法,基于概率守恒原理提出的概率密度演化理論[11,18]可以準確求取結構動力響應的概率密度函數及其演化過程,為結構的概率地震需求分析提供了新思路。

        假設結構反應需求參數Z(t)為所考察的狀態(tài)量,對于某一物理解答存在、唯一且連續(xù)依賴于初始條件的概率保守系統(Z(t),Ψ),其中Ψ=(Ψ1,Ψ2,…,Ψs)為隨機參數向量,s為隨機變量的個數,其聯合概率密度函數為PZΨ(z,ψ,t),根據概率守恒原理的隨機事件描述,則有:

        (23)

        式中:D/Dt(·)為全導數;Ωt×ΩΨ為初始狀態(tài)空間在t時刻對應的解區(qū)域,對其經過一系列的數學推導,即求得廣義概率密度演化方程。

        (24)

        設初始條件為:

        pZΨ(z,ψ,t0)=δ(z-z0)pΨ(ψ)

        (25)

        式中:pΨ(ψ)為Ψ的聯合概率密度函數;δ(·)為Dirac函數;z0為確定性初始值。

        進而,采用TVD格式的有限差分法可求得Z(t)的時變概率密度函數:

        (26)

        式中ΩΨ為Ψ的分布空間。

        顯然,所求的時變概率密度函數pZ(z,t)即反映了給定地震動強度與結構響應之間的概率分布關系。

        3.3 概率密度演化方程求解的TVD格式

        為方便討論,取式(24)中的任一偏微分方程,記為:

        (27)

        針對形如式(27)的一維偏微分方程,采用有限差分法求解是行之有效的思路之一[19],且目前已發(fā)展了單邊差分格式、Lax-Wendroff格式和TVD格式三種求解格式。其中,單邊差分格式精度較低,耗散效應明顯;Lax-Wendroff格式不能保證結果的非負性,色散效應突出;基于上述兩種基本格式,并施加適當形式的通量限制器可求得具有總變差不增的TVD格式為:

        (28)

        (29)

        式中:φ0(r)=max[0,min(2r,1),min(r,2)];

        3.4 基本求解流程

        目前,國內外已發(fā)展了多種結構易損性分析方法[20]。為了較為精細地考慮既有建筑結構的抗震性能,本文給出了采用貝葉斯方法進行材料參數更新,并基于概率密度演化理論研究信息更新結構的時變易損性分析方法,具體步驟包括如下。

        (1)材料性能參數的貝葉斯更新:1)選取合理的材料物理性能經時變化模型作為先驗預測信息;2)實測不同服役齡期建筑材料的實時物理性能;3)采用貝葉斯方法對先驗模型進行信息更新,獲得更符合實際情況的材料性能衰退模型。

        (2)隨機非線性力學模型的建立:1)建立不同服役齡期條件下結構的非線性力學模型;2)選取一系列滿足場地條件的不同震級、不同斷層距的地震動時程記錄;3)分析結構參數、外部激勵存在的隨機性,確定隨機變量的個數,進而采用數論選點法進行離散代表點的選取和賦得概率的計算。

        (3)時變易損性分析:1)采用TVD格式的有限差分法求解廣義概率密度演化方程,獲得不同服役齡期、不同地震動強度作用下結構需求參數的概率密度分布函數;2)定義結構的破壞等級及相應的量化指標;3)計算不同服役齡期、不同強度地震作用下結構響應超過不同破壞等級限值的條件概率,進而繪制以地震動參數為變量的時變易損性曲線(面)。

        4 算例分析

        4.1 RC框架結構設計

        本文以7層帶內廊式對稱三跨鋼筋混凝土框架結構辦公樓為例。結構設計使用年限為50年,建筑抗震設防類別為丙類;抗震設防烈度為8度(0.2g),Ⅱ類場地,設計地震分組第二組;地面粗糙度類別C類,基本風壓0.3kN/m2。

        框架結構平面布置如圖2所示。結構底層層高3.9m,其余層層高3.3m。柱截面尺寸為0.55m×0.55m,梁截面尺寸為0.50m×0.20m,板厚為100mm。梁板柱混凝土強度等級均為C30。梁柱主筋采用HRB400級鋼筋,梁柱箍筋、板采用HRB335級鋼筋。梁柱混凝土保護層厚度20mm。填充墻為190mm厚的混凝土空心砌塊,活荷載標準值按《建筑結構荷載規(guī)范》(GB 50009—2012)規(guī)定取值。

        圖2 框架結構平面布置圖

        采用YJK軟件進行荷載組合、內力分析和截面配筋等初步設計工作。其中,結構的柱配筋見圖3。

        圖3 結構的柱配筋示意圖

        4.2 混凝土經時強度更新的貝葉斯方法

        文獻[13]采用回彈法和鉆芯法實測一般大氣環(huán)境中不同服役年限民用建筑物的混凝土抗壓強度。由于鉆芯法更加準確,本文選取基于此檢測得到的服役混凝土相對強度樣本,如表2所示。當已知混凝土的立方體抗壓強度標準值時,乘以相應的相對強度即得到混凝土時變強度標準值。

        實測服役混凝土相對強度樣本值 表2

        將式(6),(7)視為先驗預測模型,并通過式(8)建立混凝土經時抗壓強度平均值和標準值的換算關系。圖4為先驗混凝土經時相對強度標準值的預測模型和實測樣本值的對比。顯然,采用牛荻濤模型[12]預測的混凝土經時強度標準值高于實測樣本值,即該預測模型不能很好地預測混凝土經時強度的真實情況。若直接采用該模型進行結構非線性響應分析,預測結果將比真實情況更大,不利于對實際服役結構進行檢修決策和優(yōu)化加固設計等作業(yè)。

        圖4 先驗混凝土經時強度預測值與試驗值對比

        根據貝葉斯方法,將表2中的實測數據樣本值視為似然函數,并用其將先驗預測模型進行多次修正,得到信息更新后混凝土經時強度的預測值,如圖5所示。分析可知:貝葉斯方法可有效實現先驗分布和實測數據的綜合。一方面,更新后的模型曲線與先驗預測曲線的發(fā)展趨勢存在相似的變化規(guī)律,即貝葉斯方法繼承了對先驗分布理性規(guī)律的認識;另一方面,貝葉斯模型曲線與先驗分布曲線存在一定的偏移,經過若干次貝葉斯更新即可使預測結果改善,使更新值更加接近實測樣本值,即貝葉斯方法能充分考慮結構的實際服役情況。

        將混凝土的經時抗壓強度fcu(t)視為滿足正態(tài)分布的隨機變量,為了充分考慮材料的實際服役情況,采用貝葉斯方法對其均值進行參數更新(圖5);而將鋼筋的時變銹蝕深度視為確定性變量。故本文考察的7層RC框架結構中,將每一層混凝土的時變強度視為一個隨機變量,共包含7個隨機變量。表3列出了服役齡期條件下材料參數的基本時變信息。鋼筋和混凝土采用第2節(jié)定義的經時材料本構關系模型。

        圖5 混凝土經時相對強度均值的貝葉斯更新

        不同服役齡期下的材料性能 表3

        4.3 地震波的選取

        地震波的合理選取直接影響結構易損性的分析效果。本文基于PEER地震動數據庫,取場地剪切波速Vs,30為260~510m/s,對應中國的Ⅱ類場地[21]。選取的20條地震動記錄均勻分布于以下5個地震動條帶中[22]:1)SMSR(小震級、短距離):5.5≤M≤6.5,15km≤R≤30km;2)SMLR(小震級、長距離):5.5≤M≤6.5,30km≤R≤50km;3)LMSR(大震級、短距離):6.5≤M≤7.5,15km≤R≤30km;4)LMLR(大震級、長距離):6.5≤M≤7.5,30km≤R≤50km;5)NEAR(近場地震):6.5≤M≤7.5,0km≤R≤15km。其中M為震級,R為斷層距,SA1為結構第一自振周期對應的加速度反應譜值。所選的地震波信息見表4。

        地震動記錄 表4

        4.4 時變易損性分析

        基于最大層間位移角的概率密度分布函數求得結構在不同服役年限、不同地震動水平、不同性能等級作用下的失效概率,采用分段三次Hermite函數擬合得到服役結構在不同服役齡期條件下地震易損性曲線如圖6所示。為了更加直觀地反映服役期對結構失效概率的影響,進而擬合得到三維時變易損性曲面如圖7所示。

        圖7 時變易損性曲面

        從圖6可知,隨著結構破壞程度的增加,結構易損性曲線在不同服役齡期之間的差異越大。如:當選用結構最大層間位移角為1/550時(完好),易損性曲線在不同服役齡期之間差異很?。欢x用結構最大層間位移角為1/60時(不嚴重破壞),易損性曲線在不同服役齡期之間差異性較大。即當處于不同服役齡期的既有建筑結構遭遇“小震”作用時,其易損性曲線不會存在過大的偏差;而當其遭遇“大震”作用時,材料的實時力學性能很大程度上決定了結構的安全性。另一方面,結合表3可知,結構服役齡期分別為14年和20年時,鋼筋均未發(fā)生銹蝕,而混凝土強度均值處于相當水平;而結構服役齡期分別為1年和50年時,50年混凝土的相對強度略大于1年時的相對強度,且1年時鋼筋未發(fā)生銹蝕,而50年時鋼筋則發(fā)生了程度不大的銹蝕。故而,在本文考察的6個服役齡期內,結構服役1年和50年、14年和20年生成的地震易損性曲線是大致一致的,即結構抗震性能相似。

        圖6 不同服役齡期條件的地震易損性曲線

        綜合分析圖6和圖7的結構時變易損性變化規(guī)律,可以發(fā)現:結構在服役初期,由于混凝土徐變、硬化等因素的影響,其強度處于增強的過程,在同一地面運動加速度作用下,其抗震性能隨著服役齡期的增加而提高。結構服役齡期達10年左右,混凝土強度發(fā)展到較為穩(wěn)定的階段,并在今后10~20年之間處于平穩(wěn)過程,該時期結構的抗震性能最好。結構服役齡期達25~30年時,混凝土強度開始降低,且鋼筋開始銹蝕,其抗震性能開始呈現衰減的趨勢。結構服役齡期為50~60年時(超過結構設計使用年限),混凝土強度呈現較為明顯的降低趨勢,且鋼筋銹蝕仍持續(xù)發(fā)展,結構的抗震性能亦有一定程度的降低。

        此外,可以注意到,傳統的易損性分析方法通常假定結構地震需求參數服從對數正態(tài)分布,所得的分析結果只是經驗近似值。而采用概率密度演化理論可直接分析求得地震需求參數的真實分布模式,其計算結果相較于傳統方法應更加可靠。限于篇幅,圖8僅給出了結構在服役齡期1年時分別采用本文選用方法與傳統分析方法所得的地震易損性曲線。顯然,當地震作用較小時,兩種分析方法的計算結果差異不大;而當地震作用較大時,選用本文方法所得的失效概率略大于傳統方法。即當地震作用較小時,可選用傳統方法進行結構抗震性能分析,而當地震作用較大時,傳統分析方法所求得的計算結果偏于不安全。

        圖8 服役齡期1年時地震易損性計算結果對比

        5 結論

        本文考慮材料性能的時變規(guī)律對結構抗震性能的影響,采用實測數據進行材料信息更新,對一考慮隨機參數的7層RC框架結構進行時變易損性分析,結論如下:

        (1)貝葉斯方法可有效實現先驗分布和實測數據的綜合,經過更新的貝葉斯數據既可保持已有理性認識的變化規(guī)律,又能更加接近實際服役情況。

        (2)RC框架結構的抗震性能與材料的經時性能有關,前期抗震性能隨著服役齡期的增加而提高,之后一段時期內保持在較為穩(wěn)定的階段,當服役齡期達到30年前后,由于材料性能的退化,結構的抗震性能隨服役齡期的增加呈現降低的趨勢。

        (3)結構的時變易損性曲線在不同服役期是不一致的,且與劃分的性能水平有關。隨著結構破壞等級的增加,結構的易損性曲線在不同服役年限之間存在的差異越大。

        (4)當地震作用較小時,采用本文方法與傳統方法求得的結構失效概率差異不大;而當地震作用較大時,本文方法所得的結構失效概率略大于傳統方法。

        利用本文所提方法,可以求得滿足實際概率分布的服役結構時變易損性,研究成果對既有建筑結構的檢修決策和優(yōu)化加固提供參考,具有良好的工程實用性。

        猜你喜歡
        經時易損性服役
        基于IDA的預應力混凝土連續(xù)梁橋易損性分析
        工程與建設(2019年5期)2020-01-19 06:22:48
        材料服役行為研究:助力國家名片——材料服役行為分論壇側記
        張語涵(一首)
        文史雜志(2019年1期)2019-01-21 02:02:04
        盧進海
        基于PSDM和IDA法的深水隔震橋梁地震易損性分析比較
        減水劑對混凝土經時損失的影響
        四川水泥(2016年2期)2016-08-16 10:01:04
        基于溫度的預拌水泥混凝土坍落度經時損失的試驗研究
        2015年中考熱身單項選擇填空精練100題
        基于性能的FRP加固RC框架結構地震易損性分析
        潮州市湘橋區(qū)洪澇災害承災體易損性及其變化
        亚洲天堂av福利在线| 国产麻豆一精品一AV一免费软件 | 91综合久久婷婷久久| 色婷婷亚洲一区二区三区在线| 亚洲中文字幕无码av永久 | 美腿丝袜中文字幕在线观看| 国产亚洲精品国产精品| 国产精品成人久久电影| 午夜亚洲国产理论片亚洲2020| 91国内偷拍一区二区三区| 国产小视频在线看不卡| 亚洲aⅴ在线无码播放毛片一线天| 国产成人精品三级在线影院| 亚洲一区二区三区自拍麻豆| 后入内射国产一区二区| 中文字幕日韩一区二区三区不卡| 久久久亚洲精品午夜福利| 午夜视频一区二区三区四区| 天天摸夜夜摸夜夜狠狠摸| A午夜精品福利在线| 一区二区三区视频在线免费观看 | 国产精品免费无遮挡无码永久视频| 国产zzjjzzjj视频全免费| 欧美一级视频在线| 日本一二三四区在线观看| 久久不见久久见中文字幕免费 | 日产精品高潮一区二区三区5月| 又粗又粗又黄又硬又深色的| 久久精品免视看国产明星 | 熟女一区二区中文字幕| a级毛片免费完整视频| 动漫在线无码一区| 国产精品美女主播在线| 人成午夜免费视频无码| 国产精品多人P群无码| 亚洲免费av第一区第二区| 免费无遮挡无码永久在线观看视频 | 国产欧美一区二区精品仙草咪| 久久久高清免费视频| 国产av一区二区亚洲精品| 国产一区二区三精品久久久无广告 |