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

冬小麥葉片和冠層尺度的光合-蒸騰模擬

2023-10-10 08:52:12張傳偉沈彥俊張玉翠
干旱地區農業研究 2023年5期
關鍵詞:模型

趙 陽,張傳偉,沈彥俊,張玉翠

(1.中國科學院遺傳與發育生物學研究所農業資源研究中心/中國科學院農業水資源重點實驗室/河北省節水農業重點實驗室,河北 石家莊 050022;2.中國科學院大學現代農業科學學院,北京 100049)

植物光合是生態系統有機物的重要來源,蒸騰作用是與之相伴的水分散失過程,兩者之間的耦合模擬是評價全球和區域初級生產力、模擬作物生長、研究陸面過程與氣候相互作用和預測生態環境變化的重要手段,對研究地球表層能量物質遷移、地球生物化學循環等有重要意義[1]。

植物光合和蒸騰作用耦合的常用方法是把氣孔導度和蒸騰作用耦合到一個光合作用模型中[2-4]。氣孔作為控制植物和外界進行水分與氣體交換的通道,通過影響CO2和水汽的擴散對光合作用和蒸騰作用實現控制。因此,氣孔導度是水碳耦合模型的關鍵變量。氣孔導度受光合有效輻射、相對濕度、土壤含水量、溫度、CO2濃度5種環境因子的影響,并且隨著外界條件的改變而改變,相比將其設為一個常數代入模型中應用,引入氣孔導度模型不僅可以提高模擬精度[5],而且對預測植物應對氣候條件的改變更具有價值。目前,氣孔導度模型的形式眾多,可以廣義地劃分為3種類型——乘和或者多重因素限制模型(Jarvis,1976[6]);基于氣孔導度和光合速率線性關系的半經驗模型(Ball-Berry類[7-8])以及基于脫落酸(ABA)信號、葉片水勢或水力控制的水分脅迫響應模型[9]。其中由于Ball-Berry類模型(以下簡稱BB模型)所需要的參數較少且易獲取,被廣泛用于葉片[2,10-11]、冠層[12-13]以及區域和全球系統模型(ESMs)中,用于模擬水和碳的循環過程。

在進行模型的耦合時,植物的功能類型、研究尺度及研究目的不同,模型耦合的側重點就不同,所采用的模型表達形式也不同,這造成模型復雜程度的較大差異,對實際應用影響較大。Yu等[14]基于氣孔控制光合與蒸騰的生理機制,引入了CO2內部導度的概念,開發了基于氣孔行為的光合-蒸騰耦合模型(SMPT-SB模型),模型所需參數易獲取,并且在玉米、大豆及樹木的葉片尺度應用,獲得了較好的結果。任傳友等[15]將該模型進行尺度提升,在長白山溫帶闊葉林冠層尺度研究中獲得了良好的適用性。以上研究均將光合速率和蒸騰速率獨立進行研究,尚未涵蓋水碳耦合方面的研究。張寶忠等[16]利用該模型計算了夏玉米葉片及冠層的光合、蒸騰及水分利用效率,與實測值一致性較好,在不同尺度研究水碳耦合關系。然而,利用該模型模擬冬小麥的光合和蒸騰特性還未有定論。此外,冠層尺度的水碳耦合是衡量農田蒸散和作物生產力以及農田水分利用狀況的重要內容,同時也是區域或全球蒸散模型應用和驗證的基礎。

因此,本文的研究目的主要有:(1)探究冬小麥葉片和冠層尺度的光合及蒸騰作用日變化規律,驗證SMPT-SB模型對于冬小麥的適用性;(2)借助該模型實現作物水分利用研究從葉片水平到冠層尺度的擴展,明確冬小麥葉片和冠層水碳交換的關系,為提高冬小麥多尺度水分利用效率提供科學依據。

1 材料與方法

1.1 研究區概況

試驗在中國科學院欒城農業生態系統試驗站(37°53′N,114°41′E,海拔50.1 m)進行,該站位于太行山山前平原,屬于半干旱半濕潤季風氣候,種植模式以冬小麥-夏玉米一年兩熟為主。年平均溫度為11.5℃,1月月均溫度最低(-5.6℃),7月月均溫度最高(26.7℃);年均降雨量487 mm,70%的降水出現在7—9月[17]。2011—2015年均蒸散量690 mm[18],冬小麥品種為‘科農199’,種植期間通常進行2~3次灌溉。小麥季關鍵生育期的劃分見表1。

表1 冬小麥生育時期的劃分Table 1 Growth periods of winter wheat

1.2 田間觀測試驗

在冬小麥生長發育期間選取4個晴日進行葉片尺度氣體交換數據的測量,9∶00—18∶00每小時觀測1組數據,冠層尺度數據由渦動相關系統測得,從中選出相對應的4個晴日數據。

葉片氣體交換參數的獲取:利用LI-6400便攜式光合作用系統(LICOR Inc.,Lincoln,NE,USA),測定小麥第一片完全展開葉片的凈光合速率 (AL,μmol·m-2·s-1)、氣孔導度(gs,mol·m-2·s-1)、蒸騰速率(TL,mmol·m-2·s-1)等葉片氣體交換參數以及葉溫(T1,℃)、相對濕度(hs)等環境因子。觀測時采用透明葉室進行自然光觀測,隨機選取6個小區,每個小區隨機選擇3株小麥,測量時每株計數3次,最后求均值作為每個小區的代表值。

冠層尺度數據的獲取:大田中建有通量觀測塔,渦度相關系統的安裝高度為3.5 m,由三維超聲風速儀(CSAT3,Campbell,USA)和開路式CO2/H2O分析儀(Li-7500A,LICOR,USA)組成,采樣頻率為10 Hz,數據采集器為CR3000,同時輸出30 min平均值,測量數據主要包括生態系統的水碳通量、三維風速等。光合有效輻射、溫濕度、降雨量等由塔上安裝的光合有效輻射表(LI-190SL,LICOR,USA)、溫濕度計(HMP155,Vaisala,Finland)、雨量筒(TE525MM,Texas Electronics,USA)等觀測獲得。

葉面積指數和株高的觀測:在葉片氣體交換參數測量的同時進行葉面積和株高的測量,每個試驗小區隨機選擇10株小麥,用直尺量取每株植株所有有效葉片的長和寬,并記錄莖數和作物種植密度,最后計算得到整個田塊的葉面積指數。

1.3 不同尺度光合-蒸騰耦合模型

1.3.1 葉片尺度水碳耦合模型 忽略水分子進出氣孔對CO2吸收造成的影響,葉片碳吸收過程受到葉片氣孔、邊界層及內部生物化學及光化學過程的影響,葉片尺度的光合速率[16]可以表示為:

(1)

其中,gsc、gbc、gic分別為葉片氣孔、邊界層、葉片內部對CO2的導度(mol·m-2·s-1);Ca為環境CO2濃度(μmol·mol-1);Γ*為葉片無暗呼吸時的CO2補償點(μmol·mol-1)。

葉片的蒸騰受到邊界層阻力和氣孔阻力的共同作用,因此,葉片尺度的蒸騰速率[16]可以表示為:

(2)

其中,Wi、Wa分別為葉片內外水汽濃度(mmol·mol-1);gsw、gbw分別為葉片氣孔和邊界層對水汽的導度(mol·m-2·s-1)。

葉片氣孔導度利用Luening修正的Ball-Woodrow-Berry模型進行估算:

(3)

其中,m為參數;f(Ds)為水汽響應函數;Cs為葉片周圍CO2濃度(μmol·mol-1);Γ為有暗呼吸時的CO2補償點(μmol·mol-1)。g0表示當凈光合速率趨向于0時的最小氣孔導度(mol·m-2·s-1),當考慮氣孔對CO2擴散的影響時可忽略;而Γ和Γ*在田間試驗條件下相差很小,對模擬結果造成的誤差可忽略,因此,假設Γ=Γ*,同時滿足關系gsw=1.56gbw和gbw=1.37gbc,將公式(3)代入公式(1)、(2)中,即可得到葉片尺度的光合和蒸騰速率模型:

(4)

(5)

1.3.2 冠層尺度水碳耦合模型 與葉片尺度不同的是,冠層尺度考慮水碳傳輸阻力時,需要考慮到空氣動力學阻力的影響。因此,冠層尺度的光合速率模型[16]為:

(6)

式中,AC為冠層光合速率(μmol·m-2·s-1);Ca-C為參考高度處冠層大氣CO2濃度(μmol·mol-1);Γ*-C為冠層無暗呼吸時的CO2補償點(μmol·mol-1);gbw-C為冠層對水汽的邊界層導度(mol ·m-2·s-1);gic-C為冠層內部導度(mol·m-2·s-1);ga-C為空氣動力學導度(mol·m-2·s-1)。

冠層尺度的蒸騰速率模型為:

(7)

式中,TC為冠層蒸騰速率(mmol·m-2·s-1);Wi-C、Wa-C為冠層內和冠層參考高度處的水汽濃度(mmol·mol-1);g0-C為冠層凈光合速率趨向于0時的氣孔導度(mol·m-2·s-1)。

1.4 模型參數和水汽響應函數f(Ds)的確定

1.4.1 葉片和冠層阻力參數的確定

(1)內部阻力的計算。凈同化速率不僅取決于氣體的擴散,同時也受到葉肉中生化和光化學過程的控制,如光反應、暗反應等,Jones等[19]提出葉片內部導度:

(8)

其中,Ci為胞間CO2濃度(μmol·mol-1)。

研究表明,對于C3植物而言,胞內CO2濃度與外界CO2濃度之間的比值通常保持為一定值,因此,冠層內部CO2濃度可近似采用下面的關系式求解[19]:

Ci=0.7Ca

(9)

那么冠層內部導度為:

(10)

而Γ*可用下式[20]計算:

Γ*=42.7+1.68(T-25)+0.012(T-25)2

(11)

(2)邊界層阻力的計算。葉片尺度的邊界層阻力(rbw,m2·s·mol-1)可由LI-6400直接測量得到,冠層尺度的邊界層阻力可以由下式表示[21]:

(12)

式中,α為冠層內風速衰減系數(α=3);W為葉寬幅(cm),根據試驗中測量的葉片寬度取1.5 cm;Uh為冠層高度h處的風速(m·s-1);rbw-C為冠層對水汽的邊界層阻力(m2·s·mol-1);2.24×10-2為單位由mol·m-2·s-1向m2·s·mol-1轉換的轉換系數。不同時期的葉面積指數(LAI)見表2。

表2 不同生育期冬小麥葉面積指數及株高Table 2 Leaf area index and plant height of winter wheat at different growth periods

(3)空氣動力學阻力的計算。空氣動力學阻力由冠層內和冠層上方空氣流動決定,冬小麥冠層在中性層結條件下,動力學阻力可由下式計算:

(13)

式中,UZ為參考高度處的風速(m·s-1);k為卡曼常數,其值取0.41;Z為觀測濕度和風速的高度(m);d為零平面位移(m);Z0為表面粗糙長度(m)。

對于小麥而言,d和Z0的值可由下式進行估算[22]:

d=0.63h

(14)

Z0=0.13h

(15)

式中,h為植株的高度(m)。

1.4.2 水汽響應函數 Ball-Woodrow-Berry氣孔導度模型中的參數m代表了氣孔導度對光合速率、CO2濃度、相對濕度以及溫度的綜合敏感性,其值由于水分脅迫程度、物種以及作物功能類型的不同有所差異。氣孔導度對水汽的響應函數具有多種形式,本研究中對比了常見的3種形式(圖1),當f(Ds)=RH時,斜率m=9.05,這也與一些C3植物模型中的m=9.0[22-23]相符,氣孔導度模擬值和實測值的決定系數(R2=0.75)也在3種形式中最高,因此,選擇此函數作為水碳耦合時的水汽響應函數。

1.5 模型驗證

將同時期的葉片和冠層尺度儀器觀測的光合和蒸騰數據與模擬結果進行對比,建立模擬和觀測值之間的相關性分析。通過計算模擬值與實測值的回歸斜率b、決定系數R2、均方根誤差RMSE對模型模擬效果進行評價。

2 結果與討論

2.1 葉片光合及蒸騰速率日變化實測與模擬

在外界環境條件影響下,光合與蒸騰的日變化可以反映出冬小麥通過氣孔開閉對光合與蒸騰的調節作用。從實測結果來看(圖2),冬小麥的光合速率呈先增后減的趨勢,通常在11∶00左右達到峰值。上午隨著光合有效輻射的增加以及相對濕度的減小,植物葉片氣孔導度增大,光合速率逐漸增大;而14:00之后光合有效輻射快速下降,相對濕度也逐漸增大,光合速率迅速下降。4月17日以及30日12∶00時光合和蒸騰速率快速下降,其原因是溫度較高且輻射強,部分氣孔關閉以減少水分的散失,存在明顯的“光合午休”現象。由于5月9日14∶00—15∶00不時有云層遮蔽陽光,平均光合速率下降較快。實測值與模擬值的日變化趨勢一致,模型完整地模擬出了小麥葉片光合在一天中的波動情況。實測值顯示,從整個生育期來看,4月17日的日平均光合速率為17.94 μmol·m-2·s-1,4月30日光合速率日均值達到最大20.95 μmol·m-2·s-1,這兩個時期內小麥生長速度較快,5月9日光合速率日均值相對抽穗期下降了20%左右,5月17日達到最低13.85 μmol·m-2·s-1。模擬和實測的日均值最大相對誤差為0.05%,每小時模擬結果的相對誤差在0.14~3.85 μmol·m-2·s-1之間。

冬小麥葉片在4個典型日的蒸騰速率日變化實測和模擬值如圖3所示,與光合速率類似,小麥葉片的蒸騰速率也呈現先增后減的趨勢。4個典型日中,由于4月17日當天晴朗無風,葉片周圍水汽凝聚,因此當天蒸騰速率變化波動較小。4月17日及30日蒸騰速率相差不大,但從4月30日到5月9日,蒸騰速率明顯增大,達到了8.05 mmol·m-2·s-1,主要因為小麥葉面積隨著發育進程而增大,水分消耗速率加快。模型模擬的小麥蒸騰速率日變化趨勢與實測結果一致,并且實測值與模擬值相差的絕對值較小,最大差值僅為0.41 mmol·m-2·s-1,分時段的模擬結果中也不超過2 mmol·m-2·s-1。因此,模型在誤差控制和趨勢模擬方面表現良好,可以較準確地反映出小麥在一天之內的水分消耗和CO2吸收情況。

圖3 葉片蒸騰速率日變化實測和模擬值對比Fig.3 Comparison of measured and simulated diurnal variation of leaf transpiration rate

2.2 冠層光合和蒸騰速率日變化實測與模擬

冠層尺度的光合和蒸騰速率變化更能反映農田整體情況。如圖4冠層尺度光合速率實測和模擬值日變化對比所示,與葉片尺度相同的是,冠層尺度光合速率也呈現上午逐漸增加、14∶00以后快速下降、至18∶00下降至0附近的趨勢。與葉片尺度的光合速率變化不同,冠層尺度光合速率日變化沒有明顯的“光合午休”現象,冠層尺度的光合速率最大日平均值出現在5月9日(23.81 μmol·m-2·s-1),因此可以推斷,影響冠層光合速率變化的環境因子與葉片尺度有所不同,冠層光合速率并不是葉片光合速率的簡單加和。從模型的模擬效果上來看,模擬值與實測值也有著相同的變化趨勢,模型的誤差也較小,最大誤差僅為0.84 μmol·m-2·s-1。

圖4 冠層光合速率日變化實測與模擬值對比Fig.4 Comparison of measured and simulated diurnal variation of canopy photosynthetic rate

相比于冠層光合速率的模擬結果,冠層蒸騰速率效果略低,但仍然保持了與實測結果變化趨勢的一致性(圖5)。需要說明的是,模型中將冬小麥冠層視為“大葉模型”,將蒸騰視為唯一涌源,但渦動相關系統所觀測的“蒸騰”速率實際為“蒸騰+蒸發”速率,雖然這種誤差會隨著小麥葉面積指數的增加逐漸減小,但仍然不可忽略。這也是冠層蒸騰速率模擬精度低于光合速率的原因之一。但從實際的模擬情況來看,模型仍然適用于冠層蒸騰速率的模擬,分時段模擬誤差在-1.63~2.21 mmol·m-2·s-1之間,日均誤差不超過0.15 mmol·m-2·s-1。

圖5 冠層蒸騰速率日變化實測與模擬值對比Fig.5 Comparison of measured and simulated diurnal variation of canopy transpiration rate

2.3 不同尺度光合和蒸騰速率模擬效果評價

利用LI-6400便攜式光合儀觀測的葉片尺度和渦度相關系統觀測的冠層尺度的實測數據,與模型模擬的數據進行相關性分析表明(圖1,圖6),不論是葉片還是冠層尺度,各散點均勻、緊密地分布在1∶1線的兩側,并且沒有明顯偏離的點。葉片和冠層尺度光合、蒸騰速率模擬的回歸斜率在0.8~1.0之間,葉片和冠層光合速率模擬的RMSE分別為2.14 μmol·m-2·s-1和2.48 μmol·m-2·s-1(表3),而蒸騰速率模擬的RMSE分別為0.86 mmol·m-2·s-1

圖6 葉片和冠層尺度光合及蒸騰速率模擬值及實測值對比Fig.6 Comparison of simulated and measured photosynthetic and transpiration rates at leaf and canopy scales

表3 葉片、冠層尺度光合和蒸騰速率模擬統計參數Table 3 Simulation statistical parameters of photosynthetic and transpiration rates at leaf and canopy scales

和0.30 mmol·m-2·s-1,R2除冠層蒸騰速率為0.88外,其他均在0.9以上,相關性顯著。但從葉片和冠層兩個尺度比較來看,葉片尺度的模擬效果要優于冠層尺度,這是因為將冠層假設為“一片葉子”,忽略了土壤蒸發,即使小麥的冠層高度小,但是仍然存在輻射和風速在冠層中衰減等因素,整個冠層并不是完全的均一狀態;且葉片尺度向冠層尺度的提升,需要考慮輻射、濕度、溫度、土壤水分、風速等各方面的綜合影響,這是一個復雜的過程。即便如此,通過對模擬結果與實測值的對比分析,該模型的模擬效果仍然處于較高水平,至于上面所提到的缺點,需要盡可能提高儀器測量精度,同時也需要將更多的影響因素逐漸歸納至模型中。

2.4 不同m值確定方式對模擬結果的影響

氣孔導度對水汽響應函數的斜率值m作為模型的重要參數,通常會對模擬蒸騰[24]和估算碳通量[25]產生顯著影響,其值通常在物種間或者同一物種內都有較大差異,目前尚難以量化這種差異。因此,在明確“大葉模型”適用于冬小麥的前提下,模擬冠層光合和蒸騰速率時忽略了土壤蒸發以及冠層中不同部位葉片光合能力不同的影響,采取了不同時期擬合不同m值的方式,將其模擬結果與各個時期采用同一個m值的模擬結果進行了對比。如圖7所示,在4月17日和5月17日,采用一個恒定的m=9.6模擬時,模擬效果仍然較好,其模擬的趨勢與實測結果一致,且模擬值與實測值之間絕對誤差小。這兩天的相對濕度日均值(9∶00—18∶00)較大并且相差很小,4月17日為54%,5月17日為51%;而4月30日和5月9日兩天的相對濕度較小分別為22%和31%,導致該因子對光合速率模擬的影響作用突出,其m值恒定為9.6時冠層光合速率AC-mm模型模擬值相對于實測值明顯低估。雖然仍能反映出基本相同的變化趨勢,但與實測值之間的絕對誤差較大,不能反映真實光合情況。

注:AC為實測值,AC-m和AC-mm分別為分時期擬合m值和恒定m值擬合結果。Note: AC is the measured value, AC-m and AC-mm are the fitting results of period-specific fitting m value and constant m value, respectively.圖7 不同m值率定方式下光合速率日變化實測與模擬趨勢對比Fig.7 Comparison of measured and simulated diurnal variation trends of photosynthetic rate under different m value calibration methods

m值的確定對于冠層尺度蒸騰速率的模擬效果也有影響(圖8),當采用恒定的m值時,模型亦能大致反映出冠層尺度蒸騰速率的日變化情況,但其模擬值相對于實測值出現明顯波動(5月17日),模擬效果并不理想。冠層尺度的光合和蒸騰速率如果用恒定的m值進行模擬,那么“m×RH”這個因子在相對濕度變化較大時,不能真實準確反映出其實際情況,從而降低模型的模擬精度。因此,冠層尺度的模擬分時期模擬參數m值會大大提高模型模擬精度。

注:TC為實測值,TC-m和TC-mm分別為分時期擬合m值和恒定m值擬合結果。Note: TC is the measured value, TC-m and TC-mm are the fitting results of period-specific fitting m value and constant m value, respectively.圖8 不同m值率定方式下蒸騰速率日變化實測與模擬趨勢對比Fig.8 Comparison of measured and simulated diurnal variation trends of transpiration rate under different m value calibration methods

3 結 論

本研究基于SMPT-SB模型,將其與改進的Ball氣孔導度模型耦合,建立了葉片和冠層尺度的水碳耦合模型。系統分析和模擬了冬小麥葉片和冠層尺度光合、蒸騰速率的日變化趨勢,討論了模型在兩個尺度的適用性。冬小麥葉片和冠層光合、蒸騰速率存在“單峰”和“雙峰”兩種形式,這是當輻射大、溫度高時作物自身氣孔調節造成的。耦合在氣孔導度模型中的參數m,在葉片尺度模擬時其值通常在8~10之間,但將其提升到冠層尺度后,在相對濕度差異較大的時期,其值存在較大誤差,m值無法精確反映氣孔導度的變化,從而降低模擬效果,而分時段率定m值則會大幅提高模型模擬精度。葉片尺度模擬結果優于冠層尺度。模型對于小麥在日尺度上的水分消耗和碳吸收變化模擬效果良好,模擬值與實測值具有相同的變化趨勢,并且絕對誤差控制合理,因此該模型適用于冬小麥葉片及冠層尺度的土壤-植物-大氣之間水碳交換的模擬。

猜你喜歡
模型
一半模型
一種去中心化的域名服務本地化模型
適用于BDS-3 PPP的隨機模型
提煉模型 突破難點
函數模型及應用
p150Glued在帕金森病模型中的表達及分布
函數模型及應用
重要模型『一線三等角』
重尾非線性自回歸模型自加權M-估計的漸近分布
3D打印中的模型分割與打包
主站蜘蛛池模板: 国产SUV精品一区二区6| 日本三级精品| 尤物国产在线| 国产人碰人摸人爱免费视频| 欧美色伊人| 亚洲女同欧美在线| 国产浮力第一页永久地址| AV不卡在线永久免费观看| 99久久国产综合精品女同| 天天综合网在线| 99re66精品视频在线观看| 久久精品一品道久久精品| 麻豆精品在线| 首页亚洲国产丝袜长腿综合| 国产成人午夜福利免费无码r| 日本精品中文字幕在线不卡| 91区国产福利在线观看午夜| 中文字幕久久精品波多野结| 色精品视频| 日韩人妻无码制服丝袜视频| 欧美亚洲激情| 九九热视频精品在线| 怡红院美国分院一区二区| 日本一本在线视频| 亚洲精品另类| 91小视频在线观看免费版高清| 国产成人1024精品下载| 91精品网站| 国产成人一区| 午夜国产精品视频| 亚洲欧洲日韩国产综合在线二区| 老色鬼久久亚洲AV综合| 亚洲女同一区二区| 99中文字幕亚洲一区二区| 亚洲男人的天堂在线观看| 国产啪在线91| 在线亚洲精品福利网址导航| 亚洲av日韩综合一区尤物| 久久精品国产国语对白| 性做久久久久久久免费看| 亚洲无码四虎黄色网站| 亚洲va视频| 日韩精品免费一线在线观看| 色婷婷啪啪| 看国产一级毛片| 97精品国产高清久久久久蜜芽| 在线观看无码av五月花| 亚洲va视频| 精品视频91| 亚洲国产综合第一精品小说| 亚洲自偷自拍另类小说| 亚洲专区一区二区在线观看| 色AV色 综合网站| 日韩在线2020专区| 在线精品欧美日韩| 日本道中文字幕久久一区| 亚洲九九视频| 国产成+人+综合+亚洲欧美| 久久99蜜桃精品久久久久小说| 国产美女在线免费观看| 亚洲午夜天堂| 久久男人视频| 免费福利视频网站| 日韩在线网址| 亚洲欧美日韩视频一区| 亚洲美女一级毛片| 99精品在线看| 精品自窥自偷在线看| 亚洲综合色区在线播放2019| 亚洲福利片无码最新在线播放| 欧美成人第一页| 九色视频最新网址| 成人亚洲视频| 欧美激情综合一区二区| 国产好痛疼轻点好爽的视频| 四虎国产在线观看| 玖玖精品在线| 日本一区二区不卡视频| 久996视频精品免费观看| 国产福利小视频在线播放观看| 欧美.成人.综合在线| 国产人人乐人人爱|