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

內蒙古河套灌區不同土層有機質空間變異的分形

2017-12-18 06:17:52張紅玲屈忠義張棟良
干旱地區農業研究 2017年5期
關鍵詞:研究

張 娜,張紅玲,屈忠義,張棟良

(1.內蒙古農業大學水利與土木建筑工程學院,內蒙古呼和浩特010018;2.寧夏回族自治區水利科學研究院,寧夏銀川750021)

內蒙古河套灌區不同土層有機質空間變異的分形

張 娜1,2,張紅玲2,屈忠義1,張棟良1

(1.內蒙古農業大學水利與土木建筑工程學院,內蒙古呼和浩特010018;2.寧夏回族自治區水利科學研究院,寧夏銀川750021)

以內蒙古河套灌區解放閘灌域為研究對象,以區域內間隔4 km等間距112個取樣點共448個樣本為基礎,通過多重分形及聯合多重分形方法,對有機質含量在不同土層(0~20 cm、20~40 cm、40~70 cm與70~100 cm)的空間變異特性及其相互關系進行研究。結果顯示,有機質含量0~20 cm至70~100 cm土層的多重分形譜寬度依次為 0.933、0.795、0.824、0.805,即在0~20 cm土層空間變異性最大,而20~40 cm土層空間變異性最小。不同土層多重分形譜曲線均表現為左鉤狀(Δf>0),表明在灌域內土壤有機質的空間分布中,數值較大的數據占主導地位。同時從土壤表層到埋深 1m處 D(1)/D(0)依次為 0.983、0.980、0.973、0.965,其值逐漸減小,即空間分布逐漸由密集區域向稀疏區域轉化。聯合多重分形分析表明0~20 cm土層與20~40 cm、40~70 cm土層空間變異相關程度較高,而與70~100 cm土層空間變異相關程度較低。

有機質;不同土層;空間異質性;多重分形;聯合多重分形;內蒙古河套灌區

由于受到自然因素和人為因素的共同作用,土壤特性在空間分布上表現出不同程度的變異性[1-2],且研究區域和尺度不同時,土壤特性空間分布的復雜性和異質性程度并不相同[3-4]。此外,土壤是三維自然空間實體,在水平方向和垂直方向上都存在空間變異性[5]。垂直方向上不同土層土壤特性的空間變異特征及其主要影響因素并不完全相同[6],且垂直方向的不同土層土壤特性空間變異性之間存在一定的相互關系。這種相互關系的研究有助于呈現土壤特性在三維空間上的變異特征,使土壤過程的模擬更貼近真實情況,也有助于實現通過某一土層土壤特性的空間變異特征反映其它土層土壤特性的空間變異特征。因此,定量分析不同土層土壤特性的空間變異特征及其相互關系已成為目前水土科學的研究熱點之一。

近年來,各國學者對不同尺度土壤特性參數的空間變異進行了探討,從不同角度提出了一些值得借鑒的研究方法,如地質統計學、方差分析、自相關分析、譜分析、小波分析和分形理論等[7]。分形理論作為一種能夠解釋研究對象復雜結構并定量表征均一程度的工具,已被廣泛應用于土壤特性空間分布研究中[8-11]。多重分形方法將分形集聯合起來形成多重分形譜用以描述土壤特性的空間變異特性。該方法不需要數據符合某一種分布,且其結果不受研究對象取樣尺度的影響,是以研究對象尺度為基礎劃分多尺度特征的結果。同時聯合多重分形理論可以定量分析和確定同一幾何支撐下不同研究對象在多尺度上的相互關系[12-13]。

有機質含量是反映土壤地力的重要指標,其空間變異性的研究具有重要意義。然而,目前應用分形理論進行土壤有機質的空間分布特征研究,主要集中在田間尺度及某一土層上,而大區域尺度上不同土層空間變異性的相互關系研究較少。為此,本文以內蒙古河套灌區解放閘灌域內間隔4 km等間距112個取樣點共448個樣本為基礎,通過多重分形及聯合多重分形方法,研究土壤有機質不同土層空間變異及其相互關系,對于全面理解土壤有機質空間分布具有重要意義,并為農業技術研究中野外采樣系統設計提供一定的參考。

1 研究區概況及研究方法

1.1 研究區概況

內蒙古河套灌區位于黃河流域上游,是我國第三大灌區且為亞洲最大的灌區,總土地面積1.12×104km2,其中灌溉地面積達到5.70×103km2。本文以河套灌區解放閘灌域作為研究對象,解放閘灌域是河套灌區第二大灌域,地處東經 106°43′~107°27′,北緯 40°34′~41°14′,東南緊鄰黃河,西北地處狼山腳下。西北以總排干溝為界,東南則以總干渠為界,大部分位于內蒙古自治區杭錦后旗所轄境內。解放閘灌域總控制面積2.1568×105hm2,灌溉面積1.421×105hm2,其中耕地 1.2193×105hm2,林果地面積 5.507×103hm2,牧草地面積 1.461×104hm2。灌域南北長約87 km,東西寬約78 km,呈三角形狀,解放閘灌域屬中溫帶高原型氣候,風大雨少,氣候干燥。地形較為平緩,地勢由東南向西北微度傾斜,坡度約為0.02%。海拔在1 030~1 046m之間。

1.2 樣點布設及土樣采集

于2011年8月15日至9月10日期間,根據解放閘灌域實際控制面積,通過ArcGIS在數字化圖(1∶12.5萬)上按照4 km×4 km進行網格設計,按0~20、20~40、40~70 cm與 70~100 cm共4層進行取樣。取樣區域北部以總排干溝為界,東部以四排干溝為界。采樣點位置均用GPS定位,遍布整個灌域,基本屬于平均分布。經統計,共取有效點112個。取樣點位置如圖1所示。

1.3 測試項目及方法

按照土壤學實驗規范,土樣經過風干,研磨,過篩(2mm)后,參考《土壤、水、植物理化分析教程》[14]中的濃硫酸-重鉻酸鉀容量法測定。

1.4 研究方法

多重分形方法主要用于分析單一變量的空間變異性,描述多重分形的參量有兩種,一種為 D(q)和q,另一種為α(q)和 f(q),計算過程可參考文獻[12-13]。

式中,D(q)為廣義維數;q為質量概率;pi(δ)為統計矩;δ為尺度;pi(δ)α(q)為研究對象μ的奇異指數;f(q)為研究對象奇異指數為奇異指數的維數分布函數;奇異指數α(q)和奇異指數維數分布函數f(q)的關系曲線稱為多重分形譜。

其中研究對象的多重分形譜的寬度越大,表示研究對象的空間變異性越強;若研究對象的多重分形譜左偏趨勢比較明顯,表示研究對象的高值信息對研究對象的空間分布特征影響顯著,若研究對象的多重分形譜右偏趨勢比較明顯,表示研究對象的低值信息對研究對象的空間分布特征影響顯著。

將多重分形拓展為聯合多重分形后,可利用聯合多重分形方法確定不同變量在多尺度上的相關性,計算公式可參考文獻[15-17]。

圖1 取樣點分布示意圖Fig.1 Distributionmap of sampling sites

聯合 多 重 分 形 參 數 α1(q1,q2)、α2(q1,q2)、f(α1,α2)的三維圖稱為研究對象1和研究對象2的聯合多重分形譜。當確定了聯合多重分形譜參數α1(q1,q2)、α2(q1,q2)、f(α1,α2),可以通過以下兩種方法分析研究對象1和研究對象2之間的相互關系:一種方法是將聯合多重分形譜三維圖轉化為二維圖,用灰度值表示 f(α1,α2),然后將 f(α1,α2)投影到α1(q1,q2)、α2(q1,q2)的平面上,得到研究對象1和研究對象2聯合多重分形譜的灰度圖,如果研究對象1和研究對象2聯合多重分形譜的灰度圖相對集中且沿對角線方向延伸,則可定性的說明研究對象1和研究對象2具有某些相同的分布或者相關程度較高,如果研究對象1和研究對象2聯合多重分形譜的灰度圖分布比較離散,則反之;另一種方法是計算α1(q1,q2)、α2(q1,q2)的相關系數,如果α1(q1,q2)、α2(q1,q2)的相關性顯著,則可定量的說明在多尺度上研究對象1與研究對象2的相關程度高,反之,可定量的說明研究對象1與研究對象2的的相關程度低。

2 結果與分析

2.1 有機質含量的統計特征值

不同土層土壤有機質的統計特征結果(表1)顯示,土壤有機質含量的均值隨著土層深度的增加而減少,即從土壤表層到埋深1米處土壤有機質含量遞減。這一規律產生的原因為:對于原始土壤來說,微生物是土壤有機質的最早來源,而隨著生物的進化和成土過程的發展,動植物殘體就成為土壤有機質的基本來源。在通常的自然植被條件下,土壤有機質絕大部分直接來源于土壤上生長的植物殘體和根系分泌物,故土壤表層的有機質含量較高,且土壤表層晝夜溫差較大,微生物分解活動較弱,土壤有機質更易于積累。

表1 土壤有機質含量描述性統計特征值Table 1 Descriptive statistics of soil organicmatter concentration

變異系數的大小反映總體內部各樣本之間的變異程度,一般可將樣本的變異程度分為3級:CV<10%為弱變異性,10%≤CV≤100%為中等變異性,CV>100%為強變異性[13]。由表1可知,不同土層有機質含量的變異系數介于37.52%~60.40%,均屬中等變異性。變異系數隨著土層深度的增加而增加。

2.2 基于多重分形配分函數的分形特征判定

根據多重分形配分函數的性質[17],當 log(ε)~log(χq(ε))關系為線性時,表明研究區域內變量分布具有多重分形特征[18]。擬合不同土層下土壤有機質 log(ε)與 log(χq(ε))的關系,圖 2給出了雙對數坐標下不同 q值對應的配分函數值。由圖2可知,在 -10≤q≤10時,log(ε)~log(χq(ε))之間成一定的線性關系,表明研究區域內土壤有機質含量具有多重分形特征。表2給出了不同土層的有機質含量 log(ε)-log(χq(ε))擬合的斜率 K和決定性系數R2,可知 log(ε)-log(χq(ε))之間均有良好的線性關系,決定性系數 R2位于0.81~1.00之間,在整個空間尺度上表現出良好的標度不變性,同時,當 q取不同值時,log(ε)與 log(χq(ε))之間擬合直線的斜率均不同,這表明了有機質含量在不同土層的空間分布上具有多重分形特征。

2.3 基于廣義維數譜的多重分形分析

根據多重分形廣義維數譜算法對不同土層有機質含量進行多重分形分析,在-10≤q≤10的變化范圍內得到土壤有機質分布的廣義維數譜 D(q),有機質不同土層廣義維數譜曲線 q-D(q)見圖3。

根據廣義分形維數隨質量概率統計矩階數變化的趨勢對有機質含量在研究區域上的多重分形特征進行進一步判定。當質量概率統計矩階數大于0時,如果隨質量概率統計矩階數的增加,廣義分形維數依次減小,則可判定研究對象的多重分形特征比較明顯[19-20]。不同土層的廣義分形維數均隨質量概率統計矩階數的增大而減小,這進一步說明有機質在研究區域內的空間分布格局具有多重分形特征,可以利用多重分形方法研究其空間異質性。αmax-αmin的差異來獲取研究對象的空間分布及其不規則特征是多重分形分析的一個重要方面。Δα越大表示分布越不均勻[21],多重分形譜α-f(α)曲線呈左鉤狀時,即Δf>0,說明數據中較大值的數量占主導地位,反之,較小數值的數量占主導地位[21]。

表 2 log(ε)—log(χq(ε))擬合的斜率 K和決定性系數 R2Table 2 The slope K and determination coefficients R2 of the fitting line between log(ε)and log(χq(ε))

圖2 不同土層土壤有機質的log(ε)與log(χ(qε))關系Fig.2 The relationship of log(ε)and log(χ(qε))of soil organicmatter in different layers

圖3 不同土層有機質含量的廣義分形維數Fig.3 Generalizedmultifractal dimension of soil organicmatter concentration at different soil layer

圖4給出土壤有機質在10≥q≥-10情況下不同土層的多重分形譜α-f(α)曲線,表3給出土壤有機質含量的多重分形譜參數。

定量描述土壤有機質分布的離散程度可以通過信息熵維數和盒維數的比值 D(1)/D(0)反映。D(1)/D(0)越接近于 1時表明土壤特征分布主要集中在密集區,接近于0時說明分布集中于稀疏區域[1]。可知從土壤表層到埋深 1 m處 D(1)/D(0)逐漸較小,認為有機質的空間分布逐漸由密集區域向稀疏區域轉化。

2.4 土壤有機質空間變異性的多重分形分析

根據多重分形譜形狀和多重分形譜譜寬Δα=

圖4 不同土層有機質含量的多重分形譜曲線Fig.4 Multifractal spectrums of soil organic matter concentration at different soil layer

表3 有機質含量空間分布多重分形譜參數Table 3 Multifractal parameters of soil organicmatter concentration

由圖4可知土壤有機質的多重分形譜均為“鐘形”或類似開口向下的二次拋物線,但不同土層的多重分形譜的“頂點”、開口大小存在差異,說明在不同土層其空間分布及其不規則性是有差別的。由圖4及表3可知不同土層間多重分形譜的譜寬Δα具有一定差異:0~20、20~40、40~70 cm及 70~100 cm土層的多重分形譜寬度依次為 0.933、0.795、0.824、0.805,可見有機質含量不同土層的空間變異性由大到小依次為:0~20、40~70、70~100 cm及 20~40 cm土層。

由表3可知,不同土層的Δf均大于0,其在多重分形譜曲線上均表現為左鉤狀,表明在內蒙古河套灌區解放閘灌域土壤有機質的空間分布中,數值較大的數據占主導地位,其概率分布較大。

2.5 不同土層的聯合多重分形

為了探討有機質含量在不同土層空間變異性的相互關系,同時便于用淺層易于取樣測定的數據,以減少深層取樣測試工作量。繪制了0~20 cm土層與20~40、40~70 cm及70~100 cm土層的聯合多重分 形 譜,分 別 用 αOM0~20cm、αOM20~40cm及 f(αOM0~20cm、αOM20~40cm)表示α1(q1,q2)、α2(q1,q2)、f(α1,α2),其它土層依次類推。為了便于分析將其聯合多重分形譜投影到平面上得到聯合多重分形譜的灰度圖(圖6)。從圖6可知,0~20 cm土層與20~40、40~70 cm的聯合分形譜灰度圖相對集中且沿對角線方向延伸,表明有機質含量在0~20 cm土層的空間分布與20~40 cm及40~70 cm土層空間變異的相似性較高。因此,0~20 cm土層有機質含量的空間變異性可較好反映其在20~40 cm及40~70 cm土層的空間變異性狀況;而有機質含量在0~20 cm土層與70~100 cm土層的聯合多重分形灰度圖分布較為分散,說明二者的相關性很小,因此其空間變異不能相互表達。

為進一步定量化分析土壤有機質含量0~20 cm土層與20~40、40~70 cm及70~100 cm土層空間變異性的相互關系,求解了其對應奇異指數之間的相關性。0~20 cm土層與20~40、40~70 cm及70~100 cm土層對應奇異指數的相關系數分別為0.654、0.362、0.032,且 0~20 cm土層與 20~40、40~70 cm土層在0.01水平上顯著,而與70~100 cm土層不顯著。進一步說明有機質含量空間變異性在0~20 cm土層與20~40 cm及40~70 cm土層相關程度較高。

通過一元函數建立0~20 cm土層分別與20~40、40~70 cm及70~100 cm土層的轉換模型,見表4。

圖5 0~20 cm土層有機質含量與其它土層的聯合多重分形譜Fig.5 Jointmultifractal spectrums between organicmatter concentration in 0~20 cm soil layer and other layers inmid-scale

可知:0~20 cm僅與20~40 cm土層的預測效果較好,而與40~70 cm及70~100 cm土層的回歸效果較差。

圖6 0~20 cm土層有機質含量與其它土層的聯合多重分形譜的灰度圖Fig.6 Gray scale image of jointmultifractal spectrums between organicmatter concentration in 0~20 cm soil layer and other layers inmid-scale

表4 0~20 cm土層有機質與其它土層轉化模型Table 4 The transformationmodel of organicmatter in 0~20 cm and other soil layers

這一結果與聯合多重分形分析的結果相差較大,首先聯合多重分形研究兩個具有多重分形特征的變量在幾何支撐上的關系,其涵蓋了研究對象的非均勻性和各向異性等信息,體現的是變量空間分布的信息,而一元函數僅僅表達了一個變量隨另一個變量的線性變化趨勢,并不能很好地體現變量空間分布特征。可見,一元函數雖可簡便地建立不同土層變量之間的轉換函數,但并不能很好地體現變量空間變化的信息。基于聯合多重分形結果建立兩個變量考慮空間變異情況的轉換函數有待于進一步深入研究。

3 結 論

通過多重分形及聯合多重分形方法分析了內蒙古河套灌區解放閘灌域有機質含量在不同土層的空間變異特性及其相關性,得到如下結論:

1)有機質含量在不同土層均為中等變異性,且變異系數隨著土層深度的增加而增加。

2)內蒙古河套灌區解放閘灌域有機質含量在空間分布上具有多重分形特征,其在不同土層的空間變異性由大到小依次為:0~20、40~70、70~100 cm及20~40 cm土層。同時在有機質的空間分布中,數值較大的數據占主導地位。

3)聯合多重分形表明0~20 cm土層與20~40、40~70 cm土層空間變異相關程度較高,而與70~100 cm土層空間變異的相關程度較低。

4)通過一元函數建立0~20 cm土層有機質含量分別與20~40、40~70 cm及70~100 cm土層的轉換模型,其回歸效果0~20 cm土層僅與20~40 cm土層較好,而與40~70 cm及70~100 cm土層較差。

5)基于聯合多重分形結果建立兩個變量考慮空間變異情況的轉換函數有待于進一步深入研究。

[1] 管孝艷,楊培嶺,任樹梅,等.基于多重分形理論的壤土粒徑分布非均勻性分析[J].應用基礎與工程科學學報,2009,17(2):196-205.

[2] Jacinthe PA,LaiR.Spatial variability of soil propertiesand trace gas fluxes in reclaimed mine land of southeastern Ohio[J].Geoderma,2006,136(3-4):598-608.

[3] Eghball B,Schepers JS,Negahban M,et al.Spatial and temporal variability of soil nitrate and corn yield:multifractalanalysis[J].A-gronomy Journal,2003,95(2):339-346.

[4] Caniego F J,Espejo R,Martin M A,etal.Multifractal scaling of soil spatial variability[J].Ecological Modelling,2005,182:291-303.

[5] 邵明安,王全九,黃明斌.土壤物理學[M].北京:高等教育出版社,2006.

[6] 劉繼龍,馬孝義,付 強,等.不同土層土壤特性空間變異性關系的聯合多重分形研究[J].農業機械學報,2012,43(5):37-42.

[7] 王 軍,邱 揚.土地質量的空間變異與尺度效應研究進展[J].地理科學進展,2005,24(4):28-35.

[8] Kravchenko A N,Bullock DG,BoastCW.Jointmultifractal analysis of soil yield and terrain slope[J].Agronomy Journal,2000,92:1279-1290.

[9] Lovejoy S,Schertzer D.Multifractal,cloud radiances and rain[J].Journal of Hydrology,2006,322:59-88.

[10] Cheng Q.Multifractal interpolation[C]//Proceedings of the Fifth Annual Conference of the International Association for Mathematical Geology.Trondheim,Norway:[s.n.],1999:245-250.

[11] Li Y,LiM,Robert H.Single and jointmultifractal analysis of soil particle size distributions[J].Pedosphere,2011,21(1):75-83.

[12] Zeleke T B,Si B C.Characterizing scale-dependent spatial relationships between soil properties using multifractal techniques[J].Geoderma,2006,134(3-4):440-452.

[13] Zeleke TB,SiBC.Scaling relationshipsbetween saturated hydraulic conductivity and soil physical properties[J].Soil Sci Soc Am J,2005,69(6):1691-1702.

[14] 張 韞.土壤、水、植物理化分析教程[M].北京:中國林業出版社,2011.

[15] Perfect E,Centry RW,Sukop M C,et al.Multifractal Sierpinski carpets:theory and application to upscaling effective saturated hydraulic conductivity[J].Geoderma,2006,134:240-252.

[16] Scheuring I,RiediRH.Application ofmultifractal to the analysisof vegetation pattern[J].Journal of Vegetable Science,1994,5(4):489-496.

[17] Wang D,Fu B J,ZhaoW W,et al.Multifractal characteristics of soil particle size distribution under different land-use types on the Loess Plateau,China[J].Catena,2008,72(1):29-36.

[18] Lee C K.Multifractal characteristics in air pollutant concentration time series[J].The Air Soil Pollut,2002,135(1-4):389-409.

[19] Eghball B,Schepers JS,Negahban M,et al.Spatial andtemporal variability of soil nitrate and corn yield:multifractal analysis[J].A-gron J,2003,95(2):339-346.

[20] Peitgen H O,Jurgens H,Saupe D.Chaos and Fractals:New Frontier of Science[M].NewYork:Springer-Verlag,1992:737.

[21] Kravchenko A N,CharlesW B,Donald GB.Multifractal analysisof soilspatial variability[J].Agronomy Journal,1999,91(6):1033-1041.

[22] Zeleke T B,Si B C.Scaling properties of to-pographic indices and crop yield:multifractal and jointmultifractal approaches[J].Agronomy Journal,2004,96(4):1082-1090.

Fractal study on organicmatter spatial heterogeneity of different soil layers in Inner Mongolia Hetao Irrigation District

ZHANG Na1,2,ZHANG Hong-ling2,QU Zhong-yi1,ZHANG Dong-liang1
(1.College ofWɑter Resourcesɑnd Civil Engineering,Inner MongoliɑAgriculturɑl University,Hohhot,Inner Mongoliɑ010018,Chinɑ;2.Scientific Reseɑrch Institute ofWɑter Conservɑncy of Ningxiɑ,Yinchuɑn,Ningxiɑ750021,Chinɑ)

The Jiefangzha irrigation district in the Inner Mongolia Hetao irrigation area was taken as research subject.Based on 448 samples of112 sampling points at4 km of interval distance in the region,and usingmulti-fractal and jointmulti-fractalmethod,we studied the spatial variability of organicmatter concentration in different soil layers(0~20 cm,20~40 cm,40~70 cm,70~100 cm)and their characteristics in themultiscale correlation.Results showed that the organicmatter concentration in 20~40 cm soil layer had the smallest spatial variability,and in the surface was the strongest(Δαof 0~20 cm to 70~100 cm soil layerwas 0.933、0.795、0.824 and 0.805).Themulti-fractal spectrum curve of different layers showed a left hook(Δf>0),indicating that the larger data dominated.Also,the spatial distribution gradually transformed from the dense region to the sparse region from the surface to the 1m depth(D(1)/D(0)from the surface to the 1m depth was0.983、0.980、0.973 and 0.965).The jointmulti-fractal contourmap and the corresponding singular index of different soil layers showed that the correlation between 0~20 cm and 20~40 cm、40~70 cm were relatively high,butwith 70~100 cm soil layer was relatively low.

organicmatter;different soil layers;spatial variation;multifractal;jointmultifractal;Inner Mongolia hetao irrigation district

S154.4

A

1000-7601(2017)05-0157-07

10.7606/j.issn.1000-7601.2017.05.23

2016-08-04

2016-10-27

國家自然科學基金項目(51469020)

張 娜(1988—),內蒙古呼和浩特人,博士,研究方向為非充分灌溉理論。E-mail:710840895@qq.com。

屈忠義(1969—),內蒙古巴彥淖爾人,博士,教授,主要從事節水灌溉理論與區域水鹽監測方面的科學研究。E-mail:

猜你喜歡
研究
FMS與YBT相關性的實證研究
2020年國內翻譯研究述評
遼代千人邑研究述論
視錯覺在平面設計中的應用與研究
科技傳播(2019年22期)2020-01-14 03:06:54
關于遼朝“一國兩制”研究的回顧與思考
EMA伺服控制系統研究
基于聲、光、磁、觸摸多功能控制的研究
電子制作(2018年11期)2018-08-04 03:26:04
新版C-NCAP側面碰撞假人損傷研究
關于反傾銷會計研究的思考
焊接膜層脫落的攻關研究
電子制作(2017年23期)2017-02-02 07:17:19
主站蜘蛛池模板: 在线观看国产小视频| 国产福利拍拍拍| 欧美亚洲国产精品第一页| 啪啪啪亚洲无码| 高清不卡毛片| 97无码免费人妻超级碰碰碰| 老司机午夜精品视频你懂的| 国产Av无码精品色午夜| 国模沟沟一区二区三区| 77777亚洲午夜久久多人| 一本大道在线一本久道| 精品中文字幕一区在线| 无码专区第一页| 国产综合精品一区二区| 色婷婷综合激情视频免费看| 996免费视频国产在线播放| 亚洲区视频在线观看| аⅴ资源中文在线天堂| 亚洲精品你懂的| 亚洲国产成人麻豆精品| 国产99热| 十八禁美女裸体网站| 日韩不卡高清视频| 鲁鲁鲁爽爽爽在线视频观看 | 日韩中文欧美| 国产精品私拍在线爆乳| 久久久久久尹人网香蕉| 波多野结衣一区二区三区AV| 天天躁日日躁狠狠躁中文字幕| 中文字幕免费播放| 狠狠做深爱婷婷综合一区| 亚洲aⅴ天堂| 99性视频| 精品国产美女福到在线不卡f| 亚洲精品无码AV电影在线播放| 国产福利小视频高清在线观看| 怡春院欧美一区二区三区免费| 国产一级毛片在线| 超清无码一区二区三区| 22sihu国产精品视频影视资讯| 国产97视频在线| 一级毛片在线播放免费| 人妻无码AⅤ中文字| 亚洲精品第五页| 又黄又湿又爽的视频| 国产乱子伦精品视频| 色综合成人| 国产精品第一区| 最新日本中文字幕| 91福利一区二区三区| 91精品国产自产91精品资源| 91成人在线观看| aa级毛片毛片免费观看久| 国产在线观看人成激情视频| 欧美无专区| 国产免费自拍视频| 日本影院一区| 五月天丁香婷婷综合久久| 免费在线视频a| 亚洲最大情网站在线观看| 综合色在线| 国产爽妇精品| 热久久这里是精品6免费观看| 国产精品青青| 亚洲精品欧美日本中文字幕| 国产精品专区第一页在线观看| 国产第一页免费浮力影院| 自偷自拍三级全三级视频| 欧美啪啪视频免码| 欧美伦理一区| 日本久久网站| 人妖无码第一页| 日本免费高清一区| 99热亚洲精品6码| 精品少妇人妻无码久久| 免费毛片在线| 国产又粗又爽视频| 亚洲最大看欧美片网站地址| 在线精品视频成人网| 国产女人18水真多毛片18精品| 久久99国产精品成人欧美| 欧美亚洲一二三区|