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

循環加載下壓實黃土的邊界面塑性本構模型

2024-04-29 00:44:03代倩廖紅建康孝森孫玉軍周恒
西北大學學報(自然科學版) 2024年1期

代倩 廖紅建 康孝森 孫玉軍 周恒

收稿日期:2023-10-15

基金項目:陜西省教育廳科研計劃項目(22JK0445);陜西省秦創原“科學家+工程師”隊伍建設項目(2023KXJ-178)。

第一作者:代倩,女,博士生,從事黃土動力學試驗與數值模擬等研究,daiqian68@stu.xjtu.edu.cn。

通信作者:廖紅建,女,博士,教授,從事巖土本構關系、巖土工程減災與防災等研究,hjliao@mail.xjtu.edu.cn。

摘要? 壓實黃土的動應力應變關系預測是黃土高原平山填方工程長期運營涉及的關鍵問題之一,已有模型預測壓實黃土的動本構關系難度較大。該文在邊界面塑性框架下,引入動態映射法則,提出了一個反映應變累積的循環加載塑性模量表達式,構建了循環加載條件下壓實黃土邊界面塑性本構模型,分析了循環應力比與圍壓對應力應變關系與模量的影響規律,采用壓實粉質黃土動三軸試驗結果進行了驗證。結果表明,該循環加載塑性模量基本反映了壓實黃土的應變累積行為,所構建的本構模型考慮了循環應力比對動應力應變關系與模量的影響,可為填方壓實黃土長期變形分析提供參考。

關鍵詞? 壓實黃土;循環加載;應力應變;本構模型;塑性模量

中圖分類號: TU473? DOI:10.16152/j.cnki.xdxbzr.2024-01-004

Bounding surface plasticity model of compactedloess under cyclic loading

DAI Qian1,2,? LIAO Hongjian1,2, KANG Xiaosen3, SUN Yujun4, ZHOU Heng5

(1.School of Human Settlements and? Civil Engineering, Xian Jiaotong University, Xian 710049, China;

2.School of Civil Engineering,?? Xian Jiaotong University City College, Xian 710018, China;

3.School of Geological Engineering and Geomatics, Changan University, Xian 710064, China;

4.State Power Investment Corporation, Yellow River Upper Reaches Hydropower Development Co., Ltd., Xining 810001,? China;

5.Powerchina Northwest Engineering Co., Ltd., Xian 710065, China)

Abstract? The prediction of the dynamic stress-strain relationship of compacted loess is essential in the projects of gully reclamation for farming and cutting mountains for city building in Loess Plateau. It is still difficult to predict the behaviors. To fill the gap, this paper formulates a bounding surface plasticity model for compacted loess subjected to cyclic loading. Specifically, a plastic modulus for cyclic loading is proposed to reflect the behaviors of? accumulated strain, and an updated mapping rule is introduced. Furthermore, the influence of confining pressure and cyclic stress ratio on stress-strain relationship and modulus is shown by several case studies, which is verified by experimental tests of compacted loess. The results show that the plastic modulus reflects? accumulated strain of compacted loess.The constitutive model considers the effect of cyclic stress ratio on dynamic stress-strain relationship and modulus of compacted loess because of the proposed plastic modulus and updated mapping rule. The results can give reference for analysis of long-term deformation of loess-filled foundation.

Keywords? compacted loess; cyclic loading; stress-strain; constitutive model; plastic modulus

黃土地區平山造城、 治溝造地、 高速鐵路等工程涉及不同規模的填方壓實黃土(見圖1), 長期循環荷載下壓實黃土產生累積變形,準確預測壓實黃土的循環加載累積變形對于陜西、 山西及甘肅等地區的黃土填方工程運營具有重要的實際意義。

目前,循環加載條件下砂土與黏土本構模型研究較為廣泛和深入。砂土模型一般采用錐面屈服面,基于記憶面[1]或邊界面等[2]構建砂土循環加載本構模型[3-4],考慮相變狀態參量[5]、組構演化[6]及非共軸特性[7]等。黏土循環加載本構模型一般采用修正劍橋模型屈服面,基于邊界面塑性,采用映射法則與塑性模量反映循環加載過程中的塑性應變累積[8-9],求解復雜動力邊值問題[10]。還將非飽和土水力特性引入邊界面塑性框架,建立了熱水力耦合模型[11],以及考慮非飽和狀態參量和顆粒破碎參數[12]及損傷特性[13]。而黃土循環加載本構模擬研究較少,已有研究主要是單調加載下黃土的損傷特性[14-15]、各向異性[16]、結構性[17-18]、臨界狀態線不唯一[19]等。可借鑒以往黏土本構模型構建黃土的循環加載本構模型,崔廣芹在邊界面塑性框架下研究了黃土的動本構模型[20]。在土的動本構模型中,映射中心與塑性模量關系到加載過程中塑性應變累積,映射中心可能沿應力加載方向運動[21]。塑性模量影響循環加載條件下黃土偏應變累積穩定與非穩定狀態,合適的循環加載塑性模量對預測壓實黃土動力特性至關重要[10,22]。

雖然上述學者開展了黃土動本構模型的研究,但模型尚未考慮到黃土的臨界狀態線不唯一的基本特性,尤其是粉質黃土中粉粒組含量高,屬于間隔級配土一類,且大孔隙架空結構導致其初始孔隙比范圍較大,臨界狀態線明顯依賴于初始孔隙比[19,23]。而臨界狀態線對于循環加載下黃土的動力行為預測至關重要,循環加載下黃土的本構關系預測仍然是一個挑戰。本文在邊界面塑性理論框架下,基于Kang等建立的邊界面塑性本構模型[24],引入動態映射法則,提出一個循環加載條件下壓實黃土塑性模量表達式,推導壓實黃土動本構模型,模擬不同循環應力比時壓實黃土動本構關系。研究結果可為填方黃土變形分析提供參考。

1? 本構模型

1.1? 考慮黃土臨界壓力比的邊界面

傳統各向異性邊界面不能反映黃土的臨界狀態線不唯一性[19,23],可通過臨界壓力比這一指標反映(見圖1),故采用考慮臨界壓力比的邊界面[24],

[AKf-]=([AKq-D]-[AKp-D]′α)2+[AKM-]2[AKp-D]′([AKp-D]′-[AKp-D]′0)=0? [JY](1)

式中:([AKp-D]′,[AKq-D])為各向異性結構性邊界面上的像應力;α為各向異性不變量;p′[TX-]0為等向硬化參數;[AKM-]由下式確定。

[AKM-]=(+(1-)[AKp-D]′[AKp-D]′0-1)[KF(]M2-α2[KF)]? [JY](2)

式中:M為不同洛德角時臨界狀態線斜率,M=gθMc;由下式確定[24]。

=[KF(][SX(]M-αc[]M+αc[SX)][KF)][SX(]r[](r-1)32[SX)]-[SX(]1[]r-1[SX)]? [JY](3)

式中:αc為臨界狀態時α的值,臨界壓力比r反映臨界狀態線與等向壓縮線的距離。

1.2? 加載面與塑性勢面

對于超固結土,初始時加載面在邊界面內,破壞時加載面與邊界面一致,加載面為[24]

[AKf-]=(q-p′α)2+[AKM-]2p′(p′-p′0)=0? [JY](4)

式中:(p′,q)為當前應力,位于加載面上。循環加載過程加載面與邊界面相似,相似中心為反轉應力點。

塑性勢函數[24]

g=([AKq-D]-[AKp-D]′α)2+M2[AKp-D]′([AKp-D]′-[AKp-D]′g)? [JY](5)

式中:[AKp-D]′g反映塑性勢面大小。

1.3? 等向硬化與旋轉硬化

等向硬化采用修正劍橋模型硬化法則,

d[AKp-D]′0=〈L〉[SX(][AKp-D]′0[]Cp[SX)][SX(]g[][AKp-D]′[SX)]? [JY](6)

式中:Cp為λ-κ1+e0,λ與κ分別為壓縮線與等向加卸載線的斜率。

初始各向異性與誘發各向異性采用Dafalias與Taiebat提出的各向異性演化表達式[25],

dα=〈L〉crpa(αb-α)[SX(][AKp-D]′[][AKp-D]′0[SX)]? [JY](7)

式中: pa為標準大氣壓; cr為常數, 控制各向異性結構性邊界面q軸方向的旋轉率。 〈·〉為Maca-uley算子。αb為各向異性不變量α的上限值,表達式為[25]

αb=±[SX(]M[]ζ[SX)](1-eζ|η|M)? [JY](8)

式中:ζ為控制參數。

1.4? 循環加載映射法則

邊界面塑性框架下,預測循環加載條件下應變累積,需采用更新映射中心法則[26],循環加載條件下,映射中心從上一個應力反轉點更新至下一個應力反轉點(見圖2)。壓實黃土首次加載與后續循環加卸載所產生的累積應變規律不同,首次加載時累積應變顯著大于后續循環加載一次時應變累積。

首次加載時映射中心設為原點,單調加載本構模型亦設原點為映射中心,后續循環加載時映射中心設定為應力反轉點,已廣泛應用于循環加載本構模型中。采用Seidalinov與Taiebat所提出的更新映射法則[27]如下,

[AKp-D]′=p′mc+ρ(p′-p′mc)

[AKq-D]=qmc+ρ(q-qmc)[JY](9)

式中:ρ為相似比;(p′mc,qmc)為映射中心的應力狀態。將像應力代入邊界面表達式中即可求得相似比ρ,算得當前應力與像應力在應力空間的距離,確定循環加卸載中塑性模量值。

1.5? 循環加載塑性模量

塑性模量是邊界面塑性理論最為重要的一部分,直接關系到塑性變形累積計算的準確性。在邊界面塑性框架下,塑性模量Kbδ=Kb+Kδ。Kb為各向異性結構性邊界面上像應力相關的塑性模量分量,由各向異性結構性邊界面一致性條件確定,而經典塑性理論只有Kb這一項,只能反映應力狀態達到屈服面時塑性應變,不能反映壓實黃土初始塑性應變累積與循環加載塑性應變累積效應。Kδ為當前應力狀態相關的塑性模量分量,由當前應力與像應力在應力空間中坐標決定。

基于Seidalinov與Taiebat[26]的塑性模量公式,考慮循環加載下壓實黃土變形規律,提出一個改進循環加載塑性模量Kδ,

Kδ=[AKp-D]′04(ρ-1)2e-ad∫[SX(]1[]ε0[SX)]|Δεpq |? [JY](10)

式中:ad為塑性模量衰減參數;ε0為參考應變值0.15;Δεpq為塑性偏應變增量。

循環加載條件下,當ad越大,Kδ衰減越快,反映隨著塑性應變累積過程中壓實黃土更容易破壞。

(ρ-1)2反映當前應力狀態相關的塑性模量Kb大小與當前應力和像應力在應力空間中的距離密切相關。若當前應力遠離像應力,如當前應力接近映射中心,則(ρ-1)2達到較大值,則Kδ較大,塑性應變累積極小,接近彈性狀態。若當前應力逐漸接近像應力,例如當前應力在映射中心與像應力中間時,則(ρ-1)2逐漸接近0,Kδ逐漸接近0,塑性應變累積率逐漸增大,表現為材料由彈性狀態向彈塑性狀態過渡。當前應力與像應力相等時,則(ρ-1)2為0,Kδ為0,塑性應變累積率達到最大,表現為材料彈塑性狀態,這一狀態與經典塑性力學彈塑性狀態一致。

由邊界面式(1)的一致性條件,推導得塑性模量Kb表達式為

Kb=cpa(α-αb)[SX(][AKp-D]′[][AKp-D]′0[SX)][SX(][AKf-][]α[SX)]-[SX(][AKp-D]′0[]Cp[SX)][SX(]g[][AKp-D]′[SX)][SX(][AKf-][][AKp-D]′0[SX)]? [JY](11)

1.6? 應力應變關系

根據非關聯流動法則,

dεp=〈L〉[SX(]g[][AKσ-D][SX)]? [JY](12)

式中:塑性乘子表達為

L=[SX(]1[]Kbδ[SX)][JB((][SX(]f[][AKσ-D][SX)]:d[AKσ-D][JB))]? [JY](13)

根據廣義胡克定律,

dσ=Dedεe? [JY](14)

式中:為De彈性矩陣。

由應變可加性計算總應變增量,

dε=dεe+dεp? [JY](15)

由上述表達式及一致性條件,得增量型應力應變關系為

dσ=Depdεe? [JY](16)

式中:彈塑性剛度矩陣為

Dep=De-(Kbδ+nDem)-1DemnDe [JY](17)

式中:m和n分別為

m=(g/[AKσ-D])/‖g/[AKσ-D]‖? [JY](18)

n=([AKf-]/[AKσ-D])/‖[AKf-]/[AKσ-D]‖? [JY](19)

2? 算例分析

2.1? 模型參數

模型含10個材料參數,通過常規三軸試驗、循環加載試驗及試算確定。①強度參數:黏聚力c,斜率Mc為常規三軸壓縮路徑下臨界狀態線的斜率,Me為常規三軸拉伸路徑下臨界狀態線的斜率,通過常規三軸試驗即可確定。②壓縮參數與回彈參數:λ和κ分別為v-ln p′平面上等向壓縮線斜率與卸載加載線斜率,通過等向壓縮試驗確定。③彈性參數:泊松比μ設為0.3。④臨界壓力比:r反映的是臨界狀態時有效平均主應力p′與硬化參數p′0之間的關系,在v-ln p′平面上分析臨界狀態線與等向壓縮線確定(見圖1)。⑤旋轉硬化參數:cr控制邊界面旋轉,ζ反映了應力路徑對旋轉硬化參數上限值的影響,通過試算確定[25]。⑥塑性模量參數:ad控制塑性模量增大或減小,根據循環加卸載試驗應變累積確定,若偏應變穩定于某一值且小于破壞應變,則ad取正值;若偏應變大于破壞應變且不穩定,則ad取負值。

2.2? 算例

通過循環加載算例分析不同循環應力比條件下應力應變滯回圈、循環加載非線性、動楊氏模量的變化規律。材料參數值根據一般土性參數設定(見表1),結果見圖3。

由圖3(a)可知,循環加載過程中應力應變曲線非線性,應力應變滯回圈逐漸趨于穩定,循環應力比越大,滯回圈越顯著,反映了循環加載下黃土的彈塑性行為。進一步,隨著循環次數的增加,累積塑性應變增大,但變化量越來越小,尤其是循環應力比CSR(σd/(2σc))越大,塑性應變累積越快。

由圖3(b)可知,循環加載過程中動楊氏模量逐步衰減,當循環應力比為常數,動楊氏模量隨著塑性應變累積而減小,當累積塑性應變為某一特定值,動楊氏模量隨著循環應力比增大而減小。該規律反映了循環加載過程中土結構損傷演化與應變累積和循環應力比的相關性,循環應力比較小引起結構損傷較小,循環應力比較大引起結構損傷較大,符合循環加載下黃土變形規律。

3? 驗證與討論

應用該模型預測循環加載條件下飽和壓實黃土的動應力應變關系、塑性應變累積及動楊氏模量衰減規律[27]。材料參數取值根據文獻[24]確定(見表2)。

由圖4、5可知,模型模擬了飽和壓實黃土的循環加載非線性及應力應變滯回效應。首先,當循環應力比較大時,飽和黃土應力狀態接近邊界面,應力應變表現出強非線性,軸向應變累積顯著,且首次加載應力應變曲線與后續加載應力應變曲線顯著不同。模型通過引入更新映射中心法則,首次加載時映射中心為原點,后續加載時映射中心為應力反轉點,較好模擬了飽和壓實黃土首次加載與后續加載應力應變關系。其次,飽和壓實黃土應力應變滯回效應顯著,循環應力比越大滯回圈愈大,且隨著加載次數增大逐漸趨于穩定,模型較好地模擬了應力應變滯回圈疏密變化,經典彈塑性模型難以實現。

由圖6可知,模型基本預測了飽和壓實黃土的塑性應變累積行為。隨著循環次數的增大,累積軸向應變增大,軸向應變增加速率逐漸減小。循環應力比越大,累積軸向應變越大,土樣更容易達到破壞條件,如CSR=0.43對應的軸向應變遠大于CSR=0.35對應的軸向應變(見圖6)。一方面,循環應力比越大,應力狀態越靠近邊界面,導致塑性模量較小且應變增量較大,更易產生塑性應變累積。另一方面,循環應力比越大,應力狀態更快到達臨界狀態線,更快達到破壞條件。

由圖7可知,模型基本預測了飽和壓實黃土動楊氏模量衰減規律。動楊氏模量計算值小于動楊氏模量實驗值,但趨勢一致,且在一個數量級。動楊氏模量隨著圍壓的增大而增大,反映了圍壓對模量的影響。隨著循環加載次數的增大,動楊氏模量衰減,在102次循環次數內衰減顯著,反映了循環加載下黃土結構劣化規律。

4? 結語

為預測壓實黃土的動本構關系,本文修正了既有邊界面塑性本構模型,預測分析了飽和壓實粉質黃土動三軸試驗結果,得出以下主要結論。

在邊界面塑性理論框架下,引入了一個動態映射法則,提出了一個反映應變累積的循環加載塑性模量表達式,構建了循環加載下壓實黃土邊界面塑性本構模型,分析了循環應力比對應力應變關系與動楊氏模量的影響規律。

采用該模型預測分析了循環加載下飽和壓實黃土動本構關系,模型計算結果與試驗結果分析表明所提出的循環加載塑性模量能夠反映壓實黃土的循環加載非線性、塑性應變累積、應力應變滯回及動楊氏模量衰減規律。研究方法與結果可為填方壓實黃土長期變形分析提供參考。

參考文獻

[1]CORTI R, DIAMBRA A, WOOD D M, et al. Memory surface hardening model for granular soils under repeated loading conditions[J]. Journal of Engineering Mechanics, 2016, 142(12):04016102.

[2]WANG Z L, DAFALIAS Y F, SHEN C K. Bounding surface hypoplasticity model for sand[J].Journal of Engineering Mechanics, 1990, 116(5):983-1001.

[3]方火浪, 沈揚, 鄭浩,等. 砂土三維多重機構邊界面模型[J].巖土工程學報, 2017, 39 (7): 1189-1195.

FANG H L, SHEN Y, ZHENG H, et al. Three-dimensional multi-mechanism bounding surface model for sands[J].Chinese Journal of Geotechnical Engineering, 2017, 39 (7): 1189-1195.

[4]董建勛, 劉海笑, 李洲. 適用于砂土循環加載分析的邊界面塑性模型[J]. 巖土力學, 2019, 40(2): 684-692.

DONG J X, LIU H X, LI Z. A bounding surface plasticity model of sand for cyclic loading analysis[J]. Rock and Soil Mechanics, 2019, 40(2): 684-692.

[5]遲明杰, 趙成剛, 李小軍. 剪脹性砂土邊界面模型的研究[J]. 工程地質學報, 2008, 16(3): 128-134.

CHI M J, ZHAO C G, LI X J. A bounding surface plasticity based constitutive model for dilatant sand [J]. Journal of Engineering Geology, 2008, 16(3): 128-134.

[6]WOO S I, SALGADO R. Bounding surface modeling of sand with consideration of fabric and its evolution during monotonic shearing[J]. International Journal of Solids and Structures, 2015, 63:277-288.

[7]FANG H L, ZHENG H, ZHENG J. Micromechanics-based multimechanism bounding surface model for sands[J]. International Journal of Plasticity, 2017, 90:242-266.

[8]XIAO Y, LIU H L, LIU H, et al. Unified plastic modulus in the bounding surface plasticity model[J]. Science China Technological Sciences, 2016, 59(6):932-940.

[9]HU C, LIU H X, HUANG W. Anisotropic bounding-surface plasticity model for the cyclic shakedown and degradation of saturated clay[J]. Computers and Geotechnics, 2012, 44: 34-47.

[10]KIM T, JUNG Y H. A new perspective on bounding surface plasticity: The moving projection origin[J]. KSCE Journal of Civil Engineering, 2017, 21(3):652-658.

[11]ZHOU C, NG C W W. Simulating the cyclic behaviour of unsaturated soil at various temperatures using a bounding surface model[J]. Géotechnique, 2016, 66(4):344-350.

[12]RUSSELL A R, KHALILI N. A unified bounding surface plasticity model for unsaturated soils[J]. International Journal for Numerical and Analytical Methods in Geomechanics, 2006, 30(3):181-212.

[13]YANG C, CUI Y J, PEREIRA J M, et al. A constitutive model for unsaturated cemented soils under cyclic loading[J]. Computers and Geotechnics, 2008, 35(6):853-859.

[14]FU Y K, GAO Z W, HONG Y, et al. Destructuration of saturated natural loess: From experiments to constitutive modeling[J]. International Journal of Damage Mechanics, 2021, 30(4):575-594.

[15]褚峰, 張宏剛, 邵生俊. 隴東Q3結構性黃土壓剪損傷本構模型試驗研究[J].巖土力學, 2019, 40(10): 3855-3870.

CHU F, ZHANG H G, SHAO S J. Experimental study of constitutive model of Longdong Q3 structural loess with compressive and shearing damage [J]. Rock and Soil Mechanics, 2019, 40(10): 3855-3870.

[16]侯樂樂, 翁效林, 崔藝鋮, 等.? K0固結結構性黃土的各向異性本構模型[J].巖石力學與工程學報, 2022, 41(10): 2124-2134.

HOU L L, WENG X L, CUI Y C, et al. Critical state constitutive model of K0 consolidated structural loess[J].Chinese Journal of Rock Mechanics and Engineering, 2022, 41(10): 2124-2134.

[17]姚志華, 陳正漢, 方祥位, 等. 非飽和原狀黃土彈塑性損傷流固耦合模型及其初步應用[J].巖土力學, 2019, 40(1): 216-226.

YAO Z H, CHEN Z H, FANG X W, et al. Elastoplastic damage seepage-consolidation coupled model of unsaturated intact loess and its application [J]. Rock and Soil Mechanics, 2019, 40(1): 216-226.

[18]蔣明鏡. 現代土力學研究的新視野: 宏微觀土力學[J].巖土工程學報, 2019, 41(2): 195-254.

JIANG M J. New paradigm for modern soil mechanics:Geomechanics from micro to macro [J]. Chinese Journal of Geotechnical Engineering, 2019, 41(2): 195-254.

[19]康孝森, 廖紅建, 黃強兵, 等. 結構性黃土的臨界壓力比及其邊界面塑性本構模型預測[J].巖土力學, 2022, 43(6): 1469-1480.

KANG X S, LIAO H J, HUANG Q B, et al. Spacing ratio of structural loess and its prediction using bounding surface plasticity model [J]. Rock and Soil Mechanics, 2022, 43(6): 1469-1480.

[20]崔廣芹. 飽和黃土動力本構模型及地鐵隧道周圍土層變形分析[D]. 西安:西安建筑科技大學, 2014.

[21]JUNG Y H, LEE J H. Experimental validation of the radial mapping rule in bounding surface plasticity model[J]. Journal of the Korean Geotechnical Society, 2013, 29(1): 171-181.

[22]KAN M E, TAIEBAT H A, KHALILI N. Simplified mapping rule for bounding surface simulation of complex loading paths in granular materials[J]. International Journal of Geomechanics, 2014, 14(2):239-253.

[23]XU L, COOP M R. The mechanics of a saturated silty loess with a transitional mode[J]. Géotechnique, 2017, 67(7): 581-596.

[24]KANG X S, LIAO H J, HUANG Q B, et al. Enhanced anisotropic bounding surface plasticity model considering modified spacing ratio of anisotropically consolidated clay[J]. Acta Geotechnica, 2022, 17(6): 2213-2233.

[25]DAFALIAS Y F, TAIEBAT M. Anatomy of rotational hardening in clay plasticity[J]. Géotechnique, 2013, 63(16): 1406-1418.

[26]SEIDALINOV G, TAIEBAT M. Bounding surface SANICLAY plasticity model for cyclic clay behavior[J]. International Journal for Numerical and Analytical Methods in Geomechanics, 2014, 38(7):702-724.

[27]代倩, 廖紅建, 康孝森, 等. 動荷載下填方體壓實黃土動應變與動孔壓變化規律研究[J].巖土工程學報, 2021, 43(S1): 235-240.

DAI Q, LIAO H J, KANG X S, et al. Behaviors of dynamic strain and pore pressure of compacted loess in loess-filled foundation induced by dynamic loading[J]. Chinese Journal of Geotechnical Engineering, 2021, 43(S1): 235-240.

(編? 輯? 張? 歡)

主站蜘蛛池模板: 久久久久国产一区二区| 人妻无码中文字幕第一区| 91色老久久精品偷偷蜜臀| 国产精品女同一区三区五区| 久久久噜噜噜久久中文字幕色伊伊 | 欧美日韩激情| 综合天天色| 国产毛片久久国产| 67194在线午夜亚洲| 怡红院美国分院一区二区| 58av国产精品| 色亚洲激情综合精品无码视频| 欧美激情视频一区| 五月天婷婷网亚洲综合在线| 91精品视频播放| 亚洲精品无码av中文字幕| 日韩精品亚洲人旧成在线| 国内丰满少妇猛烈精品播 | 国产一区亚洲一区| 91视频精品| 欧美人在线一区二区三区| 98超碰在线观看| 狠狠色丁香婷婷| 日韩在线网址| 一级爱做片免费观看久久| 精品无码专区亚洲| 亚洲国产精品日韩av专区| 欧美精品一区在线看| 97久久人人超碰国产精品| 无码日韩视频| 麻豆精品在线| 日本一本正道综合久久dvd| 一区二区欧美日韩高清免费 | 中文无码影院| 日韩国产黄色网站| 国产精品自在在线午夜区app| 色天堂无毒不卡| 大香网伊人久久综合网2020| 亚洲天堂精品视频| 久久99热这里只有精品免费看| 鲁鲁鲁爽爽爽在线视频观看| 99热这里只有精品5| 尤物精品视频一区二区三区| av一区二区三区高清久久| 尤物精品视频一区二区三区| 久久男人资源站| 国产精品va免费视频| 欧美一级高清免费a| 九九视频免费在线观看| 成人亚洲国产| 激情爆乳一区二区| 国产成人永久免费视频| 欧美一级在线看| 自拍偷拍欧美| 在线欧美国产| 青青青国产精品国产精品美女| 日韩精品无码免费专网站| 欧美精品在线免费| 91在线激情在线观看| 日韩午夜伦| 久久夜色精品| 婷婷亚洲综合五月天在线| 国产jizz| 57pao国产成视频免费播放| 五月婷婷丁香综合| 亚洲一区色| 亚洲成av人无码综合在线观看| 久久精品只有这里有| 五月婷婷丁香色| 人妻少妇久久久久久97人妻| 精品无码日韩国产不卡av| 国产你懂得| 一级毛片在线播放免费| 成年免费在线观看| 国产成人在线无码免费视频| 欧美在线黄| 蜜芽一区二区国产精品| 精品少妇人妻一区二区| 无码丝袜人妻| 2021国产在线视频| 一区二区三区国产精品视频| 久久婷婷五月综合色一区二区|