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

基于AquaCrop模型的夏玉米生長模擬及灌溉制度優化

2019-07-18 07:09:36陳超飛柳雙環郭大辛徐芳平馬孝義
干旱地區農業研究 2019年3期
關鍵詞:產量模型

陳超飛,柳雙環,郭大辛,徐芳平,馬孝義

(1.西北農林科技大學水利與建筑工程學院,陜西 楊凌 712100;2.西北農林科技大學旱區農業水土工程教育部重點實驗室,陜西 楊凌 712100)

玉米是中國第一大糧食作物。2015年中國玉米種植面積達3 811.9萬hm2,總產量達2.246億t[1]。關中地區是陜西糧食主產區,主要實行冬小麥-夏玉米輪作制[2]。該區夏季容易出現極端天氣,降雨分配不均導致伏旱時有發生[3],這對夏玉米的生產極為不利,氣候變化、水資源短缺和田間管理不當更加劇了生產的不穩定性[4]。因此,研究夏玉米節水灌溉對保障關中地區糧食穩產高產、提高水資源利用率有重要意義。

田間試驗往往耗時耗力,且易受氣候變化和病蟲害等環境因素影響,結合作物模型進行農田或區域尺度的作物生產模擬、預測和評估成為國際上的研究熱點[5]。作物模型克服了傳統試驗研究對象單一、試驗周期長、時間和空間有限制的缺點,有助于田間管理決策[6]。相較其他作物模型而言,AquaCrop模型以良好的適用性、較少的參數和較好的模擬效果吸引了大批科研工作者的關注和應用[7]。世界各地學者對該模型適用性展開的大量研究表明,利用該模型能較好模擬不同水分處理多種農作物的產量形成機制[8-10]。國內學者應用AquaCrop模型對冬小麥[11-12],玉米[13-14],棉花[15],大豆[16],大蔥[17]和西北胡麻[18]等多種作物進行模型校準和驗證,評價了生物量、產量和水分利用效率等方面的精度和適用性。目前,AquaCrop模型主要用于模擬不同作物產量和水分利用效率方面的研究。國內較缺乏對夏玉米不同水氮組合情境下的模擬研究,且AquaCrop作物模型在關中地區針對夏玉米模擬的研究更是少見報道。

于2016年和2017年6-10月在陜西武功縣高新農業示范園針對夏玉米設置不同梯度灌水和施氮處理的大田試驗。利用試驗數據對AquaCrop模型進行校驗,評價其適用性,最后模擬不同降水年型下不同灌溉策略對夏玉米產量和水分利用效率的影響,以期為模型在關中地區的應用和夏玉米不同典型年的灌溉制度提供參考。

1 材料與方法

1.1 研究區概況

田間試驗分別于2016年和2017年6-10月在陜西武功縣高新農業示范園(34°21′N,108°03′E)進行,冬小麥-夏玉米輪作是該區主要的糧食種植制度。試驗區位于關中平原中西部,地處渭河一、二階地。該區多年平均降水量596.8 mm,夏玉米生長季6-9月的多年平均降水量為360.2 mm,有效降水量316.5 mm;平均氣溫為13.3℃,平均日照時數2 095 h,年輻射總量480.90 kJ,無霜期221 d;土壤為典型土,0~20 cm土層容重為1.32 g·cm-3,pH值為8.14,有機質含量10.20 g·kg-1,全氮含量0.96 g·kg-1,堿解氮質量分數為60.3 mg·kg-1,速效磷質量分數為10.28 mg·kg-1,速效鉀質量分數為128.5 mg·kg-1。

1.2 試驗方案設計

2016-2017年田間試驗在同一地塊進行,設置水、氮2個因素。水分設4個水平:I0,雨養(不灌水);I1,低水30 mm;I2,適水60 mm;I3高水90 mm。氮素4個水平:N0,不施氮;N1,為推薦施氮量的67%,120 kg·hm-2;N2,推薦施氮量,180 kg·hm-2;N3,推薦施氮量的133%,240 kg·hm-2,共16個處理,3次重復。試驗小區為東西走向,寬5.7 m,長28 m,面積159.6 m2。裂區設計,水分為主區,肥料為副區,小區間由壟隔開,壟寬0.5 m,高0.2 m。灌水方式為畦灌,灌水時間依據旱情決定,2016年灌水日期為8月13日,2017年灌水日期為7月11日。供試氮肥為尿素(含N≥46%),磷肥為過磷酸鈣(含P2O5≥16%),鉀肥為硫酸鉀(含K2O≥50%)。各處理施磷肥(60 kg·hm-2)和鉀肥(60 kg·hm-2),均作為基肥在播前一次性施入。玉米供試品種為鄭單958,行距65 cm,株距30 cm,密度為51 300株·hm-2。2016年于6月14日播種,10月1日成熟;2017年于6月15日播種,10月7日成熟。其他管理措施如雜草控制、病蟲害管理等按照一般大田管理規范進行。

1.3 觀測指標及處理

1.3.1 葉面積指數 每隔10~15 d測量標記選擇供試植株的所有綠葉葉片的長度(從葉領到葉尖)及寬度(葉片最寬處),葉面積指數(leaf area index,LAI)采用下式計算[19]:

(1)

式中,ρ為種植密度(株·hm-2);n為第j株玉米的總葉片數;m為測定株數;Lij為葉片葉領到葉尖的長度(m);Bij為葉片最寬處寬度(m);0.75為與葉形有關的葉面積回歸系數。

1.3.2 冠層覆蓋度 冠層覆蓋度(canopy cover,CC)是指土壤表面綠色冠層覆蓋的面積百分比,由LAI計算獲得,其計算公式如下[20]:

CC=1.005[1-exp(-0.6×LAI)]1.2

(2)

1.3.3 土壤含水量 本試驗使用德國IMKO公司制造的TRIME-TDR測定土壤含水量(soil water content,SWC),測量深度為1m,間隔為20 cm。每次播種前和收獲后,用土鉆取土測量土壤含水量,并校準TRIME-TDR儀器。

1.3.4 地上部生物量及產量 玉米進入拔節期后,每隔10~15 d測定地上部生物量(aboveground biomass,B),每小區取3株能夠代表作物平均長勢的植株,剪去地下部分后于105℃下殺青0.5 h,75℃烘至恒重,測其干重,乘以種植密度即為地上部生物量。玉米成熟后,選取能夠代表作物整體長勢的20株成熟玉米,果穗風干后經人工脫粒,75℃烘至恒重后測算產量。

1.3.5 水分利用效率 水分利用效率(water use efficiency,WUE)是指單位水分消耗獲得的生物量或產量,產量水分利用效率WUE可按下式計算:

(3)

式中,Y為籽粒產量(kg·hm-2);ET為作物耗水量,(mm),由AquaCrop模型模擬計算獲得。

1.4 構建AquaCrop模型數據庫

AquaCrop模型需要輸入的數據包括氣象、作物、土壤、田間管理等。

1.4.1 氣象數據 武功氣象站1961-2017年的逐日降水、氣溫等氣象數據來源于國家氣象科學數據共享服務平臺(http://data.cma.cn/site/index.html)中國地面氣候資料日值數據集(V3.0)。參考作物蒸發蒸騰量ET0采用Penman-Montieth方程[21]進行計算。2016年和2017年夏玉米生育期內逐日最高和最低氣溫、太陽輻射、參考蒸發蒸騰量、最小相對濕度及降水量見圖1。

1.4.2 土壤數據 模型輸入的土壤數據包括土層深度及層數、土壤質地、容重、各土層的凋萎含水量、田間持水量和飽和含水量等參數。試驗區土壤為典型土,用土鉆分別在0~20、20~40、40~60、60~80 cm和80~100 cm 5個深度取土,風干過2 mm篩后采用英產Mastersizer 2000型激光粒度儀測定土壤質地。各土層容重采用環刀法測定,采用離心法測定土壤水分特征曲線,用RETC軟件擬合得到凋萎含水量、田間持水量和飽和含水量,輸入AquaCrop模型建立土壤數據庫文件,詳細土壤參數如表1所示。

圖1 2016-2017年武功縣夏玉米生育期氣象數據Fig.1 Meteorological data of summer maize growth period in Wugong County from 2016 to 2017

1.4.3 作物參數數據 作物參數中的植株密度和物候期等參數根據試驗記錄直接輸入模型。歸一化水分生產力、參考收獲指數,冠層增長衰減系數和水分脅迫響應系數等重要參數根據模型手冊和文獻提供的取值范圍采用“試錯法”進行修正,部分默認參數如基底溫度、上限溫度等采用模型手冊上的推薦值。參照Vanuytrecht等[22]介紹的校準順序:第一步,校準冠層覆蓋度;第二步,校準地上部生物量;第三步,校準產量。利用2017年田間試驗數據進行模型調試,使模型模擬結果與實測結果相符合,經過校準,得到‘鄭單958’的品種參數,部分參數詳見表2。

AquaCrop模型沒有明確考慮養分的循環與平衡,而是采用半定量的評估方法來描述養分脅迫對生物量的影響[23]。模型引入相對于最大地上部生物量的比例系數Brel,公式如下:

(4)

式中,Bstress為無肥料脅迫和水分脅迫條件下得到的地上部生物量(t·hm-2);Bref為有肥料脅迫但無水分脅迫條件下得到的地上部生物量(t·hm-2);Brel的范圍為0%~100%,0%表示有養分脅迫時作物生長失敗,100%表示沒有養分脅迫,Bstress和Bref參照當地生產現狀、前人研究成果[24-26]以及本次試驗數據確定。模型內嵌自動校準程序,使用4個脅迫參數來表現肥料脅迫的響應:冠層擴張(Ksexp,f),最大冠層覆蓋度(KsCCx),生物量水分生產力(KsWP)和冠層衰減(fCDecline)[27]。通過校準,結果如表3。

表1 試驗區土壤數據

表2 AquaCrop模型玉米校準參數

1.5 模型評價指標

在本研究中,采用以下統計指標來評價模型校準和驗證的精度:決定系數R2,均方根誤差RMSE,標準均方根誤差NRMSE,一致性指數d和納什效率系數ENS。

(5)

(6)

(7)

(8)

(9)

表3 土壤肥料脅迫參數校準表

1.6 不同降水年型灌溉情景模擬

利用Pearson-Ⅲ型分布模型對1961-2017年夏玉米生育期降水量進行分析,得到濕潤年1990年、平水年1989年和干旱年1999年降水量分別為377.5、283.1 mm和207.3 mm。選擇苗期、拔節期、抽雄期和灌漿期4個關鍵生育期制定模擬灌溉方案,土壤初始含水量設置為田間持水量;肥料脅迫設置為100%,即無肥料脅迫;灌溉模擬情景設置為雨養(P0)、灌1水(P1~P4)和灌2水(P5~P10),灌水定額分別為60 mm和120 mm,按全排列共11種處理,各灌溉方案見表4。

2 結果與分析

2.1 模型校準結果

利用2017年的試驗數據進行校準,試驗共有4個水分處理,限于篇幅,主要分析充分灌水處理I3與雨養處理I0。不同施氮水平下夏玉米冠層覆蓋度校準結果如圖2,地上部生物量校準結果如圖3,土壤含水量校準結果如圖4。

2017年各處理的冠層覆蓋度和地上部生物量的校準誤差統計指標見表5。所有處理下冠層覆蓋度的R2、RMSE和ENS分別為0.868、11.83%、0.743,表明AquaCrop模型能夠較為準確地模擬出不同水氮組合情景下冠層覆蓋度的動態變化過程。所有處理地上部生物量的R2、RMSE和ENS分別為0.975、0.970 t·hm-2、0.968,表明AquaCrop模型能夠準確地模擬出不同水氮組合情景地上生物量的動態變化過程。所有處理下土壤含水量的R2、RMSE和ENS分別為0.743、20.1 mm、0.490,表明AquaCrop模型能夠較為準確地模擬出不同水氮組合情景下土壤含水量的動態變化過程。

表4 模擬灌溉方案/mm

注:N0,N1,N2,N3:0%,67%,100%和133%的推薦施氮量;I0,I1,I2,I3:雨養,虧水,適水,充分灌水。下同。Note: N0, N1, N2, N3: Application of 0%, 67%, 100% and 133% of the recommended nitrogen rate; I0, I1, I2, I3:Rainfed, deficit, suitable and adequate irrigation. The same below.圖2 充分灌水與雨養條件下夏玉米冠層覆蓋度校準結果Fig.2 Calibration results of summer maize canopy cover under full irrigation and rainfed conditions.

圖3 充分灌水與雨養條件下夏玉米地上部生物量動態模擬Fig.3 Calibration results of summer maize aboveground biomass under full irrigation and rainfed conditions

圖4 充分灌水與雨養條件下夏玉米土壤含水量校準結果Fig.4 Calibration results of summer maize soil water content under full irrigation and rainfed conditions

如圖5(a)所示,所有試驗處理下夏玉米產量的R2、RMSE、NRMSE、d和ENS分別為0.879、0.297 t·hm-2、4.390、0.959和0.812。如圖5(b)所示,最終生物量的R2、RMSE、NRMSE、d和ENS分別為0.873、1.108 t·hm-2、6.977、0.904和0.502。各處理的相對誤差見表6,產量的相對誤差最小為0.60%,對應處理為I3N2,最大為7.63%,對應處理為I1N3;最終生物量的相對誤差最小為1.27%,對應處理為I1N0,最大為11.04%,對應處理為I0N2。

2.2 模型驗證結果

采用2016年的試驗數據進行模型驗證,土壤含水量的驗證結果如圖6所示,誤差統計見表7,所有試驗處理下土壤含水量的R2、RMSE、ENS分別為0.691、15.7 mm和0.671。如圖7(a)所示,所有試驗處理下夏玉米產量的R2、RMSE、NRMSE、d和ENS分別為0.919、0.249 t·hm-2、4.112、0.977和0.915。如圖7(b),最終生物量的R2、RMSE、NRMSE、d和ENS分別為0.860、0.977 t·hm-2、6.407、0.933和0.694。各處理的相對誤差見表8,產量的相對誤差最小為1.51%,對應處理為I3N0,最大為6.65%,對應處理為I3N1;最終生物量的相對誤差最小為1.24%,對應處理為I3N2,最大為12.05%,對應處理為I0N3。

表5 2017年不同施氮及灌溉水平下冠層覆蓋度、地上部生物量和土壤含水量校準誤差統計

圖5 2017年所有灌溉施氮處理下產量(a)和最終生物量(b)校準結果Fig.5 Calibration results of yield (a) and final biomass (b) under all differentnitrogen and irrigation levels in 2017

表6 2017年所有灌溉施氮處理下產量、最終生物量校準結果

圖6 充分灌水與雨養條件下夏玉米土壤含水量驗證結果Fig.6 Validation results of summer maize soil water content under adequate irrigation and rainfed conditions

表7 2016年不同施氮及灌溉水平下土壤含水量驗證誤差統計

2.3 不同降水年型灌溉模擬結果

利用校準后的模型對3種不同降水年型進行不同灌溉方案的模擬,夏玉米產量與水分利用效率的模擬結果如圖8、表9所示。不同灌溉方案下,濕潤年產量在7.899~7.998 t·hm-2,平水年產量在7.724~8.066 t·hm-2,干旱年產量在1.984~7.349 t·hm-2。濕潤年由于降水量豐沛,雨養和灌溉條件下的水分利用效率基本不變,其模擬值為2.19 kg·m-3或2.21 kg·m-3;平水年水分利用效率模擬值為2.50~2.72 kg·m-3;干旱年水分利用效率模擬值為0.86~2.33 kg·m-3。

3 討 論

3.1 AquaCrop模型的適用性

本研究利用AquaCrop模型對關中地區夏玉米栽培進行校準和驗證,夏玉米產量和地上部生物量的模擬值與實測值一致性高,模擬誤差在可接受范圍內,且能準確模擬土壤含水量的變化趨勢,與國內外的研究成果基本吻合。如倪玲等[6]評價了AquaCrop模型在黃土塬區的適用性,結果顯示夏玉米模擬產量與實測產量間的決定系數R2為0.927,相對誤差在-2.479%~11.182%之間;模擬地上部生物量與實測地上部生物量間的R2為0.784,模型對產量的模擬效果優于對生物量的模擬,也與本研究結果相似。Abedinpour等[27]模擬了不同水氮組合情景下的玉米產量、生物量及水分生產力的變化特征,其中產量R2為0.96,RMSE為0.1 t·hm-2,ENS為0.98;生物量的R2為0.9,RMSE為0.75 t·hm-2,ENS為0.95。Ran等[14]模擬不同虧缺條件下制種玉米冠層覆蓋度的模擬值與實測值的R2為0.818,RMSE為12.9%,ENS為0.811;生物量R2為0.929,RMSE為1.972 t·hm-2,ENS為0.903;土壤含水量R2為0.736,RMSE為33.1%,ENS為0.364。這些都表明AquaCrop模型可以準確模擬出不同田間管理條件下玉米產量和地上部生物量,說明AquaCrop模型可以為關中地區夏玉米田間灌溉和決策提供指導。

圖7 2016年不同灌溉施氮處理下產量(a)和最終生物量(b)驗證結果Fig.7 Validation results of yield (a) and final biomass (b) under all different nitrogen and irrigation levels in 2016

表8 2016年不同灌溉施氮產量、最終生物量驗證結果

圖8 不同降水年型不同灌溉方案的產量(a)和水分利用效率(b)模擬結果Fig.8 Simulated yield (a) and WUE (b) in different irrigation schemes under wet, normal and dry year

表9 不同降水年型不同灌溉方案的產量與水分利用效率模擬結果

與此同時,本文也反映出模型存在一定不足,模型輕微高估了雨養I0條件下的冠層覆蓋度和地上部生物量(如圖2(e)~2(h)與圖3(e)~3(h)),即水分脅迫越嚴重,模型的模擬誤差越大,這與李玥等[18]研究胡麻的結果類似;模型整體高估了土壤含水量,尤其是2017年的所有處理(如圖4(a)~4(h)),這與Ran等[14]的研究結果相似;同時模型低估了衰落期的冠層覆蓋度,這是由于2017年玉米生育后期降水量過多,日照時數較少,延緩了玉米的衰老,這可能是模型對于衰老期的水分脅迫響應過于敏感[19],模擬的衰老速度過快造成。因此模型對于不同生育期的水分脅迫響應仍需改進[14]。

3.2 針對不同降水年型灌溉方案的制定

不同生育期受旱會對夏玉米的生長發育及產量造成不同程度的影響[28]。生長初期地上部生長緩慢,比較耐旱,一般無需灌水,但若長時間遭遇伏旱,則容易造成弱苗,最終導致大幅減產。玉米進入拔節期后,植株生長旺盛,此時缺水會引起營養體生長不良,明顯制約夏玉米株高和葉面積的增長及干物質的積累[29]。抽雄期受旱會導致穗粒數減少,灌漿期則明顯減少百粒重,均會造成減產[19]。因此,本研究按照夏玉米的生育期實際劃分情況,分別統計3種降水年型在苗期、拔節期、抽雄期和灌漿期的有效降水量。濕潤年4個生育期有效降水量分別為206.8、57.6、86.9 mm和26.2 mm,各生育期無明顯水分脅迫;平水年4個生育期有效降水量分別為105.7、14.4、131 mm和32 mm,拔節期出現水分脅迫;干旱年4個生育期有效降水量分別為129.7、6.5、22.9 mm和48.2 mm,拔節期和抽雄期出現水分脅迫。

不同降水年型灌溉模擬結果顯示:雨養條件下,濕潤年產量最高,但水分利用效率低于平水年,這是由于濕潤年的降水在各生育期分布相對較均勻,各生育期未發生明顯的水分脅迫,所以產量最高;濕潤年的耗水量為363 mm,大于平水年297.2 mm,而產量相差不大,因此導致水分利用效率較低。

灌1水情況下,濕潤年和平水年均在拔節期灌水60 mm時產量達到最大,分別為7.998 t·hm-2和8.066 t·hm-2,對應的水分利用效率也最大,分別為2.21 kg·m-3和2.72 kg·m-3,這是因為平水年在拔節期出現水分脅迫,此時灌水可以提高產量與水分利用效率。對于干旱年,抽雄期灌水60 mm產量最大,為6.618 t·hm-2,對應水分利用效率為2.23 kg·m-3;拔節期灌水60 mm的產量略低,為5.749 t·hm-2,對應水分利用效率為2.05 kg·m-3。這是因為干旱年在拔節期和抽雄期缺水嚴重,在這2個生育期灌水均可提高產量,模型模擬結果顯示抽雄水獲得的產量要高于拔節水,表明模型對于營養生長階段的水分脅迫響應較為敏感。

灌2水情況下,苗期和拔節期各灌水60 mm,即方案P5,此時濕潤年和平水年的產量最大,分別為7.998、8.170 t·hm-2,且灌水日期越延后產量越低,但如果存在拔節期灌水的組合,如方案P8和P9,產量均較其他生育期的灌水組合大。干旱年在拔節期和抽雄期各灌水60 mm時(方案P8)產量最大,為7.349 t·hm-2,且水分利用效率也達到最大,為2.33 kg·m-3,這些模擬結果都表明在玉米拔節期或抽雄期灌水可提高夏玉米產量,與前人研究成果極為吻合[28],基本說明了不同灌溉時間和灌溉次數對夏玉米產量的影響。

濕潤年灌水后的產量和水分利用效率增加幅度較小,如方案P2的產量較雨養條件下的產量僅提高1.2%,水分利用效率提高不足1%,因此若能保證出苗整齊,濕潤年不灌水也能夠獲得較高產量和水分利用效率,同時也達到節約水資源的目的。平水年雖然方案P5可獲得最高產量8.170 t·hm-2,但相較方案P2的增幅較小,為1.3%,且水分利用效率低于P2,增加灌水量并沒有顯著提升產量,反而會降低水分利用效率,造成水資源浪費,因此平水年可在拔節期灌水60 mm,推薦灌溉方案P2,可節約50%的灌溉水資源,產量同樣能維持在較高水平。干旱年采用灌溉方案P8能夠獲得最高的產量和水分利用效率,相較于方案P3,產量和水分利用效率增幅分別為11%和4.5%,考慮到干旱年整體降水較少,推薦灌溉方案P8,即在拔節期和抽雄期各灌水60 mm,可實現水資源充分利用和作物高產。

4 結 論

AquaCrop模型可以較好地模擬關中地區不同水氮條件下夏玉米冠層覆蓋度、生物量、產量和土壤含水率。通過不同降水年型不同灌溉方案的模擬,依據模擬結果及夏玉米生理需水規律,同時為了實現高產和高水分利用效率,得出不同年型的優化灌溉制度為:若能保證出苗整齊,濕潤年可不灌水,同樣能夠獲得較高產量和水分利用效率,同時也達到節約水資源的目的;平水年推薦在拔節期灌水60 mm,能節約50%水資源,亦能保證穩產;干旱年推薦在拔節期和抽雄期各灌水60 mm,可獲得高產和高水分利用效率。

猜你喜歡
產量模型
一半模型
2022年11月份我國鋅產量同比增長2.9% 鉛產量同比增長5.6%
提高玉米產量 膜下滴灌有效
今日農業(2021年14期)2021-11-25 23:57:29
世界致密油產量發展趨勢
重要模型『一線三等角』
海水稻產量測評平均產量逐年遞增
今日農業(2020年20期)2020-11-26 06:09:10
重尾非線性自回歸模型自加權M-估計的漸近分布
2018年我國主要水果產量按省(區、市)分布
2018上半年我國PVC產量數據
聚氯乙烯(2018年9期)2018-02-18 01:11:34
3D打印中的模型分割與打包
主站蜘蛛池模板: 青青草原国产av福利网站| 国产国拍精品视频免费看 | 在线国产你懂的| 永久免费无码成人网站| 国产精品福利一区二区久久| 久久久久九九精品影院| 亚洲欧美另类视频| 国产在线精品人成导航| 欧美日韩中文字幕在线| 久久综合九九亚洲一区| 亚洲第一成年人网站| 日韩美女福利视频| 在线毛片网站| 欧美成人午夜在线全部免费| 日本午夜精品一本在线观看| 精品国产免费人成在线观看| 国产乱人视频免费观看| 波多野结衣一二三| 成人福利在线视频| 伊人久久福利中文字幕| 亚洲区视频在线观看| 日韩午夜片| 欧美自拍另类欧美综合图区| 国产熟女一级毛片| 91美女视频在线| 好吊妞欧美视频免费| 2021国产精品自产拍在线观看| 天天视频在线91频| 中文字幕乱码二三区免费| 国产成人精品亚洲77美色| 92午夜福利影院一区二区三区| 欧美精品啪啪| 国产91丝袜在线播放动漫 | 国产欧美日韩免费| AV不卡在线永久免费观看| 国产欧美日韩免费| 日韩在线视频网| 囯产av无码片毛片一级| 又大又硬又爽免费视频| 综合成人国产| 人人爽人人爽人人片| 在线国产欧美| 国产成人av一区二区三区| 亚洲日韩国产精品综合在线观看| 欧美第一页在线| 91年精品国产福利线观看久久| 日本尹人综合香蕉在线观看| …亚洲 欧洲 另类 春色| 成人免费视频一区二区三区| 亚洲激情99| 国产成人区在线观看视频| 成年免费在线观看| 亚洲首页在线观看| 欧美午夜在线观看| 色国产视频| 久夜色精品国产噜噜| 久久久久国产精品嫩草影院| 亚洲免费毛片| 国产精品任我爽爆在线播放6080 | 香蕉久久国产精品免| 男女性午夜福利网站| 成人免费网站久久久| 青草娱乐极品免费视频| 青青草一区| 91国内视频在线观看| 国产成人夜色91| 国产91丝袜在线播放动漫 | 欧美全免费aaaaaa特黄在线| 久久国产热| 亚洲成人精品| 国产欧美日韩精品第二区| 久久国产热| 四虎精品国产永久在线观看| 色成人亚洲| 97无码免费人妻超级碰碰碰| 久久6免费视频| 国产一区免费在线观看| 亚洲Aⅴ无码专区在线观看q| 免费无码AV片在线观看中文| 激情亚洲天堂| 久久亚洲国产视频| 欧美另类精品一区二区三区|