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

鋼筋混凝土結構改進型分離式數值模型*

2022-07-11 23:48:50馮曉偉肖桂仲
爆炸與沖擊 2022年6期
關鍵詞:界面混凝土實驗

曾 繁,馮曉偉,黃 超,徐 權,肖桂仲,4,田 榮

(1. 中國工程物理研究院高性能數值模擬軟件中心,北京 100088;2. 北京應用物理與計算數學研究所,北京 100088;3. 中國工程物理研究院總體工程研究所,四川 綿陽 621999;4. 南京理工大學理學院,江蘇 南京 210094)

鋼筋混凝土作為結構材料廣泛應用于土木工程中,其力學性能非常復雜,與混凝土基體、鋼筋及它們之間相互作用等力學行為緊密相關。當鋼筋混凝土結構承受爆炸荷載時,鋼筋與混凝土間的粘結應力上升,當超過一定的限制值后,鋼筋與混凝土間發生滑移。盡管界面的粘結應力-滑移關系建立在微觀尺度層面,卻顯著影響結構在爆炸荷載下的力學性能。

在鋼筋混凝土結構的有限元分析方法中,圍繞鋼筋的數值建模,常用有限元模型有3 種:(1) 彌散式模型;(2) 嵌入式模型;(3) 分離式模型。

彌散式鋼筋模型常用于由均勻分布的鋼筋組成的墻體、樓板等表面型結構。該模型將鋼筋混凝土簡化為均質材料,易于使用。但是對于梁柱結構中箍筋和縱筋等復雜的鋼筋幾何問題,彌散式模型由于離散精度問題而不被采納。

在嵌入式模型中,鋼筋混凝土構件橫截面分為混凝土部分和鋼筋部分,依據平截面假定計算截面應變,再根據鋼筋和混凝土應力-應變關系和力平衡條件得到單元剛度矩陣,其中鋼筋強化效果等效為嵌入在混凝土單元中的附加軸向分量。嵌入式模型使得鋼筋的運用完全獨立于其幾何形狀,且不增加數值模型的自由度。該模型假定嵌入的鋼筋與混凝土完美粘結,即鋼筋節點位移與混凝土節點位移一致。

在分離式模型中,分別對鋼筋和混凝土精細建模,混凝土采用實體單元離散,鋼筋采用一維桿或者梁單元離散。在數值精度方面,該模型優于嵌入式模型,且通過大量的鋼筋混凝土構件爆炸實驗驗證,故而在構件一級的抗爆分析中常為工程師所青睞。由于用于離散鋼筋的單元需布置于混凝土單元的邊界,這導致混凝土的網格嚴重受限于鋼筋幾何形狀。

同時,在分離式模型中通過顯式構造界面單元,可以考慮界面的粘結滑移效應。例如,師燕超等在混凝土網格節點和鋼筋節點間構造一維的虛擬彈簧單元,通過一維滑動接觸模型考慮二者間的粘結滑移效應。Kwak 等通過混凝土網格節點與鋼筋網格節點共點方式,構造了無物理尺寸的一維連接單元。De Groot 等在混凝土與鋼筋之間定義二維粘結帶,考慮平移、旋轉等不同類型的相互作用。基于鋼筋、混凝土網格節點建立的界面模型能較好地捕捉粘結滑移效應。但是界面模型給計算帶來額外的開銷。更重要的一點是,基于節點的界面幾何建模的復雜度隨著工程結構幾何復雜度的增加而急劇增長,將耗費巨大的時間人力成本,這會極大地制約該類模型在鋼筋混凝土結構級的工程應用。因此針對大型鋼筋混凝土建筑結構,如何高效高精度地建立鋼筋混凝土數值模型依舊是一個開放性問題。

根據混合物理論,鋼筋混凝土的力學行為可由混凝土基體和強化的鋼筋兩項組分的本構關系來描述,并且無需對鋼筋進行顯式離散化。該理論具有不同組分間應變相容性的局限性假定,可以通過修正鋼筋本構關系,對混合物理論進行修正。本文構建一種基于修正混合物理論的改進型分離式數值模型,以期給出在鋼筋本構方程中如何考慮界面粘結滑移機制的簡單有效的修正方法,并避免顯式構造界面即可考慮鋼筋尺度的局部效應;對改進型分離式數值模型進行鋼筋混凝土構件級、結構級的實驗驗證,以探索該模型在宏觀大型鋼筋混凝土結構有限元分析中的實用性和可靠性。

1 改進型分離式數值模型

參照《混凝土結構設計規范(GB50010—2010)》,本文中以6 層鋼筋混凝土框架結構為例詳細闡述改進型分離式數值模型的離散步驟:(1) 混凝土結構幾何建模;(2) 鋼筋骨架線幾何建模;(3) 采用實體單元對框架結構進行網格離散;(4) 根據鋼筋骨架線幾何位置及鋼筋特征尺寸信息,對包含鋼筋骨架線的網格定義成鋼筋混凝土網格,未包含鋼筋骨架線的網格定義成素混凝土網格(這種離散方式巧妙地避免了傳統分離式建模中混凝土網格節點與鋼筋網格節點的復雜拓撲約束關系,極大縮短建模時間);(5) 定義兩類網格的材料本構模型,素混凝土網格采用經典Holmquist-Johnson-Cook (HJC)模型,鋼筋混凝土網格采用本文中提出的基于修正混合物理論的鋼筋混凝土復合材料本構模型。

圖1 為網格尺寸為25 mm 六面體單元的有限元模型,網格數為1.08×10。圖1(a)為整體框架結構有限元模型,1~3 層為混凝土部分,4~6 層為鋼筋骨架線,其中圖1(b)為局部放大圖。圖1(c)展示了構件級混凝土網格和鋼筋混凝土網格。對于大型混凝土結構而言,如采用實體單元—結構單元的傳統分離式數值模型建模,以25 mm 的網格為例,需處理千萬級鋼筋混凝土節點與混凝土節點的約束關系,前處理在網格剖分過程中將耗費大量的時間人力成本,不利于工程實用化?;诟倪M型分離式數值模型,可快速實現工程級有限元模擬,且網格分辨率達到鋼筋尺寸。

圖1 鋼筋混凝土框架結構分離式有限元模型Fig. 1 Improved discrete finite element model of reinforced concrete frame structures

1.1 修正型混合物理論模型

根據經典混合物理論,每個組分變形滿足應變協調假定,即:

式中:ε和 ε( i=1,2···)分別為應變和第個組分的應變張量。應變協調假定在復合材料理論中具有一定的局限性。對于由混凝土基體和鋼筋組成的復合材料而言,當鋼筋和混凝土間的界面發生變形或者滑移時,界面滑移將導致復合材料應變場不連續,并影響混凝土與鋼筋間的應力傳遞,從而導致鋼筋中應力下降。通過考慮界面變形或者滑移,應變協調方程(式(1))可寫成:

式中:下標f 和m 分別代表鋼筋和混凝土,ε為界面變形或者滑移的應變張量。不用顯式構造界面,而將界面的粘結滑移效應引入到鋼筋應力應變關系中。等效鋼筋彈塑性本構關系可表達為:

基于修正混合物理論,兼顧混凝土基體和具有粘結滑移效應的鋼筋力學行為,復合材料本構方程可表達為:

1.2鋼筋等效模型

Belarbi 等根據實驗結果獲取了嵌入在混凝土中鋼筋的實際屈服強度滿足如下關系式:

式中:ρ為有效配筋率,為混凝土在開裂應變為8×10時的拉伸斷裂應力。

在粘結滑移模型中,鋼筋混凝土粘結破壞有3 種形式:

(1)混凝土劈裂破壞,當由混凝土保護層或者箍筋提供的環向約束較弱時,鋼筋周邊環向拉伸應力將造成混凝土劈裂失效。

(2)拉拔式破壞,當環向約束較強時,避免了混凝土劈裂破壞模式;當界面粘結強度較小時候,界面剪應力容易達到粘結強度,從而誘發拉拔式破壞模式。

(3)鋼筋斷裂失效,當界面粘結長度較大、粘結強度較高時,鋼筋應力優先達到最大拉伸應力,從而誘發鋼筋斷裂失效。

圖2鋼筋尺度模型Fig.2Schematicofthemodelatsteelbarscale

在爆炸荷載鋼筋混凝土構件毀傷破壞過程中,通常表現為鋼筋屈服、混凝土剝落等失效形式。即鋼筋屈服時,界面剪應力滿足上升段本構關系,根據圖2 中力平衡關系給出滑移位移計算公式:

鋼筋等效彈性模量可進一步由式(5)、(6)和(7)計算。

鋼筋等效模型中,除了修正屈服強度、彈性模量等參數,同時需對塑性硬化模量進行修正。等效塑性模量為:

其中鋼筋應力可由圖2 中力平衡關系和界面最大剪應力給出,δ=δ。

因此,通過式(5)、(6)、(8),可獲得考慮界面粘結滑移的鋼筋等效應力-應變關系,其具體表達式為:

圖3 彈性模量縮放因子 /Ef 與粘結區長度l 的關系Fig. 3 Relationship between the scaling factor of the elastic modulus ( /Ef ) and the length of the bonding zone (l) under different physical parameters

圖4 塑性硬化模量縮放因子/Hf 與粘結區長度l 的關系Fig. 4 Relationship between the scaling factor of the plastic hardening modulus ( /Hf ) and the length of the bonding zone (l) under different physical parameters

從圖3 和圖4 中可以得出:鋼筋等效彈性、塑性硬化模量與粘結區長度呈正相關性,且均收斂到穩定值。塑性模量的收斂速度大于彈性模量的收斂速度。鋼筋直徑越粗、屈服強度越高、混凝土強度越高,則鋼筋等效模量越大,規律趨勢與Dehestani 等的結果一致。當粘結區長度較小時,鋼筋混凝土粘結容易發生拉拔式破壞,界面變形受鋼筋直徑、屈服強度、混凝土強度等因素影響較大。此時考慮界面粘結滑移的等效模量對鋼筋直徑、屈服強度、混凝土強度等因素較為敏感。當粘結區長度較大時,容易誘發鋼筋斷裂的破壞模式,其等效模量主要與鋼筋模型相關。

1.3 混凝土HJC 本構模型

混凝土材料在爆炸荷載作用下將處于大應變、高圍壓及高應變率的狀態。在眾多混凝土動態本構模型中,HJC 本構模型因其簡明合理的描述和計算程序的適應性,在混凝土爆炸問題中獲得了廣泛應用。因此,復合材料中混凝土基體部分采用HJC 本構模型,其應力-應變關系為:

式中: σ為混凝土應力,為混凝土抗壓強度,為內聚力參數,為損傷變量,為壓力硬化系數,為壓力,為壓力硬化指數,為應變率參數。

1.4 鋼筋混凝土復合材料本構模型

結合鋼筋等效模型及混凝土HJC 模型,復合材料本構方程(式(4))可表達為:

2 層次化驗證與確認

2.1 構件級精度驗證

改進型數值模型的功能模塊已在沖擊波結構毀傷程序中實現,該程序是一款大規模流固耦合并行程序,能夠高效高精度地模擬爆炸沖擊類問題,其中流場分析模塊是基于JASMIN 框架研發,結構動力學響應分析模塊是基于JAUMIN 框架研發,數值模擬結果采用大規??梢暦治銎脚_TeraVAP 進行后處理分析。為驗證數值模型在鋼筋混凝土結構爆炸毀傷方面的計算精度以及粘結滑移模型的有效性,選取了軸壓作用下鋼筋混凝土柱及板爆炸實驗。

2.1.1 軸壓作用柱爆炸實驗

采用18.2 kg 巖石乳化炸藥(等效11.1 kg TNT),爆高為1.5 m。鋼筋混凝土柱長2500 mm,截面為邊長200 mm 的正方形。主筋選用4 根直徑為20 mm 的HRB335 級鋼筋,對稱布置,配筋率為3.14%,鋼筋抗拉屈服強度為466.7 MPa,極限強度為601.3 MPa,彈性模量為200 GPa,塑性硬化模量為1.12 GPa。箍筋選用直徑為8 mm HRB235 級鋼筋,間距為150 mm,配筋率為0.34%,鋼筋抗拉屈服強度為483.5 MPa,極限強度為582.4 MPa,彈性模量為200 GPa,塑性硬化模量為0.94 GPa?;炷翞镃40 級商品混凝土,保護層厚度為20 mm。爆炸實驗前在柱端部施加軸向壓力12 MPa,爆炸過程中柱構件下表面通過LVDT 位移計測量位移時程,位移計布置于柱跨中位置。

基于沖擊波結構毀傷程序分別建立考慮鋼筋與混凝土間粘結滑移的改進型分離式數值模型和不考慮粘結滑移效應的分離式數值模型,網格尺寸均為10 mm。素混凝土網格采用HJC 混凝土本構模型,模型參數為:初始密度為2.4 g/cm,剪切模量為14.4 GPa,泊松比為0.20,=50.16 MPa,=0.28,=1.84,=0.84,=0.007。對于鋼筋混凝土網格,在考慮界面粘結滑移的數值模型中,參照Dehestani 等的研究,界面粘結區長度=300 mm,通過縱筋、箍筋等效材料參數(表1),及混凝土參數,再結合式(12),給出鋼筋混凝土網格中復合材料模型參數。在不考慮粘結滑移數值模型中,縱筋、箍筋材料參數為原始參數。

表1 鋼筋等效材料參數Table 1 Equivalent material parameters of rebar

對于爆炸荷載的模擬,根據藥量大小和目標距離,采用ConWep 爆炸荷載模型,通過結合反射超壓和入射超壓,并把入射角θ 的影響考慮進來,計算得到目標上壓力,直接作用在結構上,計算模型為:

在近場爆炸實驗中,通過預備性爆炸實驗,標定ConWep 模型有效性。圖5 為爆心距離1.5 m 位置處反射超壓時程曲線的ConWep 模擬結果與實驗結果對比情況。從圖5 可以看出,通過修正ConWep 模型中藥量參數,獲取與實驗測量結果相當的超壓時程曲線模擬結果,因此該模型可作為有效的構件荷載輸入。

圖5 反射超壓時程曲線對比Fig. 5 Comparison of the reflected overpressure-time history curves

圖6 為柱跨中豎向位移及破壞形貌的模擬結果與實驗結果的對比。從對比結果可以看出:考慮鋼筋與混凝土間粘結滑移效應時,數值模擬的峰值位移和殘余位移與實驗值誤差分別為6.28%和8.54%,試驗件損傷破壞形貌與模擬的損傷破壞形貌吻合較好。相比于不考慮粘結滑移效應有限元模型,改進型分離式模型的模擬結果與實驗結果更接近,這說明在爆炸荷載作用下鋼筋混凝土結構動力學響應有限元模擬中,應該考慮鋼筋與混凝土間粘結滑移。

圖6 柱豎向位移時程曲線及破壞形貌對比Fig. 6 Comparison of the vertical displacement-time history curves and failure morphology of columns

2.1.2 板爆炸實驗

實驗中藥量為10 kg TNT,爆高2.0 m。鋼筋混凝土板幾何尺寸為1200 mm×500 mm×50 mm,鋼筋型號為直徑8 mm 的HRB235 級鋼筋,采用雙層雙向配筋,主方向間距100 mm,次方向按構造配筋間距200 mm,保護層厚度10 mm,配筋率為0.96%,鋼筋參數詳見2.1.1 節。混凝土為C60 級商品混凝土,保護層厚度為10 mm。

類似地,基于沖擊波結構毀傷程序分別建立考慮粘結滑移的改進型分離式數值模型和不考慮粘結滑移效應的分離式數值模型,網格尺寸均為10 mm。素混凝土單元采用HJC 混凝土本構模型,混凝土強度為39.2 MPa,模型參數參照2.1.1 節內容。鋼筋等效材料參數為:等效強度為420.3 MPa,等效彈性模量為170.1 GPa,等效塑性硬化模量657.0 MPa,鋼筋混凝土單元中鋼筋體積分數為50%。

圖7 為C30 板跨中豎向位移的模擬結果與實驗結果的對比。從對比結果可以看出:考慮粘結滑移效應時,數值模擬的殘余位移與實驗結果更為接近。進一步說明在爆炸荷載作用下鋼筋混凝土結構動力學響應有限元模擬中,鋼筋與混凝土間粘結滑移的必要性。

圖7 板豎向位移時程曲線對比Fig. 7 Comparison of the vertical displacement-time history curves of the slab

研究有限元網格尺寸對爆炸荷載下板變形的影響,選取10、5、2.5 mm 等3 種網格尺寸,其中對應網格數分別為30000、240000、1920000。圖8 為3 種網格尺寸的數值模擬結果:數值模擬結果與網格尺寸的大小相關,如10 mm 與5 mm 的模擬結果略有差異;隨著網格尺寸逐漸減小,跨中位移時程曲線差異減少并逐步收斂(5 mm 和2.5 mm 模擬結果);進一步的細化網格尺寸可有限地提高計算精度,但是計算成本成倍增加,因此10 mm 的網格尺寸能夠兼顧計算精度及計算成本。

圖8 網格尺寸對板跨中豎向位移時程曲線的影響Fig. 8 Effect of mesh size on the vertical displacement-time history curves of the slab

2.2 結構級精度驗證

通過結構級爆炸實驗,進一步驗證改進型分離式數值模型計算精度。Woodson 等和Baylot 等開展了1/4 縮尺的兩層框架結構外場爆炸實驗。為了避免爆炸沖擊波與框架結構強耦合相互作用的干擾,即爆炸沖擊波在結構內外表面復雜入射、反射、繞射等現象,選取了文獻[26-27]中含磚墻的框架結構爆炸實驗用于驗證模型。該框架結構寬3.2 m,深1.53 m,高2 m??蚣芙Y構由6 根柱支撐,結構中柱截面尺寸為89 mm×89 mm,柱高為2 m,背面中柱截面尺寸為76 mm×76 mm,其余的柱截面尺寸為152 mm×152 mm。底層和頂層樓板厚度均為41 mm,層高0.91 m。立柱高于頂層樓板0.178 m。邊梁位于中柱和邊柱間,截面尺寸為81 mm×54 mm。四塊填充墻均由13×15 塊磚墻堆砌而成,磚墻尺寸為寬99 mm,厚48 mm,高48 mm。實驗中采用7.1 kg C4 半球形炸藥(1 kg C4 炸藥等效為1.19 kg TNT),爆點位置離地面高度為305 mm,距離底層中柱1.07 m。爆炸前,第2 層中柱頂部施加1.7 t 重物,模擬2.1 MPa 預應力荷載。圖9 為改進型分離式數值模型,網格尺寸為5 mm,總計814 萬單元,其中暴露在外部分為鋼筋混凝土網格。圖9 中插圖為爆炸前實驗設置。素混凝土網格采用HJC 混凝土本構模型,其基本參數參照2.1.1 節。鋼筋混凝土網格采用復合材料模型,參數由混凝土強度和鋼筋基本參數給出。在炸藥起爆前,在目標柱頂部施加2.1 MPa 預應力分析,當結構所有節點速度低于0.1 m/s 時,便認為結構達到新的靜力平衡。隨后開展爆炸模擬,總模擬時間為0.3 s。

圖9 改進型分離式數值有限元模型及實驗(單位:m)Fig. 9 Improved discrete finite element model and experiment (Unit: m)

圖10 為框架結構破壞形貌的模擬結果與實驗結果的定性對比,其中圖10(a)為爆炸荷載作用剛結束時框架結構損傷云圖,圖10(b)為結構最終狀態損傷破壞云圖,圖10(c)為鋼筋骨架的塑性區云圖,圖10(d)為實驗結果。模擬結果表明:(1) 迎爆面磚墻均失效;(2) 低層中柱損傷區域最早出現在柱上下兩端,隨后發展至柱跨中背爆面處,由此可見該柱失效模式為剪切為主,彎曲為輔的耦合失效(圖10(a)和10(b));(3) 低層中柱鋼筋骨架的塑性區域主要集中在柱兩端部及跨中背面位置;(4) 低層中柱失效后,失去支撐的邊梁通過其抗彎承載力提供框架的抗倒塌能力,在“梁機制”作用下框架結構未發生垮塌。

圖10 框架結構破壞形貌對比Fig. 10 Comparison of the failure morphology of the frame structure

圖11 為中柱位移時程曲線的定量對比結果。低層中柱彎曲變形,水平方向永久性位移為0.103 m,與實驗測量永久性位移0.114 m 基本一致。與不考慮粘結滑移相比,改進型數值模型計算結果與試驗結果更接近。綜合上述定性與定量驗證,可進一步驗證改進型分離式數值模型的有效性。

圖11 低層中柱位移時程曲線對比Fig. 11 Comparison of the displacement-time history curves of the lower mid column

3 結 論

(1)提出一種考慮鋼筋與混凝土間粘結滑移機制的改進型分離式數值模型。將界面粘結滑移本構模型引入鋼筋彈塑性本構關系中,對鋼筋應力-應變曲線進行了修正。修正后的彈性、塑性硬化模量與界面粘結區長度、鋼筋直徑、混凝土強度等物性參數呈正相關性。

(2)改進型數值模型通過了鋼筋混凝土構件級、結構級爆炸實驗的精度驗證,與不考慮界面粘結滑移數值模型相比,改進型模型的模擬結果與實驗結果更接近。進一步說明界面效應在爆炸荷載作用下鋼筋混凝土結構動力學響應分析中的重要性。

(3)改進型分離式數值模型由于沒有對鋼筋及其界面的顯式離散要求,使得鋼筋的運用完全獨立于其幾何形狀,同時對混凝土網格沒有約束,并且不增加計算成本,因此該模型可適用于鋼筋混凝土宏觀結構層面分析。

猜你喜歡
界面混凝土實驗
記一次有趣的實驗
混凝土試驗之家
現代裝飾(2022年5期)2022-10-13 08:48:04
關于不同聚合物對混凝土修復的研究
國企黨委前置研究的“四個界面”
當代陜西(2020年13期)2020-08-24 08:22:02
混凝土預制塊模板在堆石混凝土壩中的應用
做個怪怪長實驗
混凝土,了不起
基于FANUC PICTURE的虛擬軸坐標顯示界面開發方法研究
人機交互界面發展趨勢研究
NO與NO2相互轉化實驗的改進
主站蜘蛛池模板: 亚洲码在线中文在线观看| 亚洲成人在线网| 亚洲精品日产AⅤ| 欧美亚洲另类在线观看| 尤物在线观看乱码| 欧美成人午夜视频免看| 国内99精品激情视频精品| 国产视频一区二区在线观看| 香蕉色综合| 免费看一级毛片波多结衣| 激情视频综合网| 在线观看国产网址你懂的| 欧美色香蕉| 视频在线观看一区二区| 色欲色欲久久综合网| 99视频全部免费| 亚洲国产精品人久久电影| 黑色丝袜高跟国产在线91| 亚洲黄色网站视频| 国产1区2区在线观看| 精品国产aⅴ一区二区三区| 91久久偷偷做嫩草影院| 国产极品粉嫩小泬免费看| 精品自窥自偷在线看| 免费看a级毛片| 国产亚洲欧美在线视频| 亚洲国产成人无码AV在线影院L| 亚洲综合经典在线一区二区| 国产激情无码一区二区APP | 成人av专区精品无码国产| 青青热久麻豆精品视频在线观看| 又爽又黄又无遮挡网站| 黄色网站在线观看无码| 欧美黄网在线| 精品一区二区久久久久网站| 精品91视频| 伊人久久久久久久久久| 亚洲日产2021三区在线| 成人小视频网| 最新日本中文字幕| 国产久操视频| 日韩AV手机在线观看蜜芽| 日韩在线播放欧美字幕| 亚洲欧美国产高清va在线播放| 亚洲午夜片| 国产日韩丝袜一二三区| 精品一区二区三区波多野结衣 | 欧美成人综合在线| 日韩在线成年视频人网站观看| 亚洲国产欧美中日韩成人综合视频| 在线观看av永久| 久久精品中文无码资源站| 2021国产精品自产拍在线| 国产成人午夜福利免费无码r| 精品99在线观看| a毛片免费看| 国产在线视频欧美亚综合| 亚洲精品视频免费看| 免费一级大毛片a一观看不卡| 中文字幕亚洲专区第19页| 国产毛片久久国产| 亚洲无码精彩视频在线观看| 88av在线播放| 欧美中文字幕一区| 国产欧美日韩18| 试看120秒男女啪啪免费| 免费亚洲成人| 亚洲欧美一区二区三区麻豆| 无遮挡国产高潮视频免费观看| 爱做久久久久久| 一本大道香蕉中文日本不卡高清二区| 免费不卡视频| 国产成人无码播放| 亚洲 欧美 中文 AⅤ在线视频| 久久亚洲高清国产| 精品国产免费观看| 九九免费观看全部免费视频| 91久久性奴调教国产免费| www.国产福利| 伊人色在线视频| 精品免费在线视频| 国产永久在线观看|