本發(fā)明屬于土壤水分運(yùn)移模擬，具體涉及一種煤礦開(kāi)采裂縫地區(qū)土壤水分運(yùn)移模擬方法。

背景技術(shù)：

1、黃土高原是我國(guó)煤炭的主要產(chǎn)區(qū)，隨著煤炭科技的進(jìn)步，黃土高原礦區(qū)已基本形成高強(qiáng)度、高效率的地下開(kāi)采模式。然而，大規(guī)模地下采煤已引起嚴(yán)重的地表裂縫，對(duì)黃土高原礦區(qū)造成廣泛而劇烈的地表?yè)p害。地下采煤引起的地表裂縫對(duì)于土壤水分及其運(yùn)移過(guò)程產(chǎn)生了直接影響，導(dǎo)致該區(qū)域生態(tài)環(huán)境損害問(wèn)題日趨突出且嚴(yán)重。然而，裂縫發(fā)育過(guò)程中土壤水分運(yùn)移難以獲取，現(xiàn)有的土壤水分運(yùn)移模擬方法主要針對(duì)自然條件下的均質(zhì)土壤，而開(kāi)采裂縫改變了土壤結(jié)構(gòu)、土壤水分蒸散發(fā)速率及運(yùn)移路徑，增加了水分運(yùn)移的復(fù)雜性，導(dǎo)致裂縫區(qū)的水力特性與自然狀態(tài)下的土壤差異較大，由于本身的運(yùn)移路徑的改變以及不均勻的水力特性存在，現(xiàn)有的土壤水分運(yùn)移模擬方法無(wú)法直接適用于開(kāi)采裂縫地區(qū)。

2、因此，現(xiàn)如今缺少一種設(shè)計(jì)合理的煤礦開(kāi)采裂縫地區(qū)土壤水分運(yùn)移模擬方法，建立的煤礦開(kāi)采裂縫地區(qū)模型分割為裂縫域和基質(zhì)域，并對(duì)裂縫域和基質(zhì)域分別建立水分運(yùn)移方程進(jìn)行模擬，獲取煤礦開(kāi)采裂縫地區(qū)的土壤含水量，有效地適應(yīng)裂縫地區(qū)土壤水分運(yùn)移情況，提高了水分運(yùn)移模擬準(zhǔn)確性。

技術(shù)實(shí)現(xiàn)思路

1、本發(fā)明所要解決的技術(shù)問(wèn)題在于針對(duì)上述現(xiàn)有技術(shù)中的不足，提供一種煤礦開(kāi)采裂縫地區(qū)土壤水分運(yùn)移模擬方法，其方法步驟簡(jiǎn)單、設(shè)計(jì)合理，建立的煤礦開(kāi)采裂縫地區(qū)模型分割為裂縫域和基質(zhì)域，并對(duì)裂縫域和基質(zhì)域分別建立水分運(yùn)移方程進(jìn)行模擬，獲取煤礦開(kāi)采裂縫地區(qū)的土壤含水量，有效地適應(yīng)裂縫地區(qū)土壤水分運(yùn)移情況，提高了水分運(yùn)移模擬準(zhǔn)確性。

2、為解決上述技術(shù)問(wèn)題，本發(fā)明采用的技術(shù)方案是：一種煤礦開(kāi)采裂縫地區(qū)土壤水分運(yùn)移模擬方法，其特征在于，該方法包括以下步驟：

3、步驟一、測(cè)量數(shù)據(jù)的采集：

4、步驟101、在設(shè)定的總監(jiān)測(cè)時(shí)間中，采用裂縫監(jiān)測(cè)模塊獲取待研究的煤礦開(kāi)采裂縫地區(qū)中裂縫深度、裂縫下寬度和裂縫上寬度；采用土壤監(jiān)測(cè)模塊獲取待研究的煤礦開(kāi)采裂縫地區(qū)的土壤含水量，土壤水吸力，土壤孔隙度；采用氣象模塊獲取待研究的煤礦開(kāi)采裂縫地區(qū)的降雨強(qiáng)度、氣溫、太陽(yáng)總輻射；

5、步驟102、對(duì)待研究的煤礦開(kāi)采裂縫地區(qū)中的土壤質(zhì)地進(jìn)行測(cè)定，獲取待研究的煤礦開(kāi)采裂縫地區(qū)的土壤質(zhì)地?cái)?shù)據(jù)；

6、步驟103、將裂縫深度記作第1類測(cè)量數(shù)據(jù)，裂縫下寬度記作第2類測(cè)量數(shù)據(jù)，裂縫上寬度記作第3類測(cè)量數(shù)據(jù)，土壤含水量記作第4類測(cè)量數(shù)據(jù)，土壤水吸力記作第5類測(cè)量數(shù)據(jù)，土壤孔隙度記作第6類測(cè)量數(shù)據(jù)，降雨強(qiáng)度記作第7類測(cè)量數(shù)據(jù)，氣溫記作第8類測(cè)量數(shù)據(jù)，太陽(yáng)總輻射記作第9類測(cè)量數(shù)據(jù)；

7、步驟二、采集測(cè)量數(shù)據(jù)的預(yù)處理：

8、步驟201、采用計(jì)算機(jī)對(duì)第1類測(cè)量數(shù)據(jù)至第9類測(cè)量數(shù)據(jù)進(jìn)行預(yù)處理，得到預(yù)處理后的第1類測(cè)量數(shù)據(jù)至第9類測(cè)量數(shù)據(jù)；

9、步驟202、對(duì)預(yù)處理后的第1類測(cè)量數(shù)據(jù)至第9類測(cè)量數(shù)據(jù)分別進(jìn)行均值處理，得到測(cè)量到的裂縫深度均值、裂縫下寬度均值、裂縫上寬度均值、土壤含水量均值、土壤水吸力均值、土壤孔隙度均值、降雨強(qiáng)度均值、氣溫均值和太陽(yáng)總輻射均值；其中，將裂縫深度均值、裂縫下寬度均值、裂縫上寬度均值記作裂縫數(shù)據(jù)，將土壤含水量均值、土壤水吸力均值、土壤孔隙度均值、降雨強(qiáng)度均值、氣溫均值和太陽(yáng)總輻射均值均記作模擬初始參數(shù)；

10、步驟三、獲取煤礦開(kāi)采裂縫地區(qū)的土壤水力特性參數(shù)：

11、步驟301、采用計(jì)算機(jī)以土壤水吸力為橫坐標(biāo)，以土壤含水量為縱坐標(biāo)，繪制出預(yù)處理后的第5類測(cè)量數(shù)據(jù)和預(yù)處理后的第4類測(cè)量數(shù)據(jù)的各個(gè)點(diǎn)；

12、步驟302、采用計(jì)算機(jī)基于繪制出預(yù)處理后的第5類測(cè)量數(shù)據(jù)和預(yù)處理后的第4類測(cè)量數(shù)據(jù)的各個(gè)點(diǎn)進(jìn)行曲線擬合，得到殘余含水量θr、飽和含水量θs、土壤孔隙第一相關(guān)參數(shù)α、土壤孔隙第二相關(guān)參數(shù)n和土壤孔隙第三相關(guān)參數(shù)m；

13、步驟303、采用計(jì)算機(jī)根據(jù)公式得到非飽和導(dǎo)水率函數(shù)k(θ)；其中，ks為飽和導(dǎo)水率；θ表示土壤含水量自變量；

14、步驟四、建立煤礦開(kāi)采裂縫地區(qū)模型，并根據(jù)建立的基質(zhì)域水分運(yùn)移方程和裂縫域水分運(yùn)移方程進(jìn)行模擬：

15、步驟401、采用計(jì)算機(jī)利用hydrus軟件根據(jù)步驟202中的裂縫數(shù)據(jù)，建立待研究的煤礦開(kāi)采裂縫地區(qū)斷面圖，記作煤礦開(kāi)采裂縫地區(qū)模型0bga，并設(shè)置模擬初始參數(shù)及邊界初始條件；其中，設(shè)定煤礦開(kāi)采裂縫地區(qū)模型0bga中距離裂縫區(qū)域邊界10cm的土壤區(qū)域?yàn)榱芽p域，剩余的土壤區(qū)域?yàn)榛|(zhì)域；

16、步驟402、采用計(jì)算機(jī)利用hydrus軟件采用等間距網(wǎng)格對(duì)煤礦開(kāi)采裂縫地區(qū)模型0bga進(jìn)行網(wǎng)格化，得到煤礦開(kāi)采裂縫地區(qū)網(wǎng)格模型；

17、步驟403、采用計(jì)算機(jī)建立基質(zhì)域水分運(yùn)移方程，如下式：

18、其中，ca為基質(zhì)域的比水容重，θa為基質(zhì)域的土壤含水量函數(shù)，t為時(shí)間，表示關(guān)于時(shí)間t的偏導(dǎo)，表示關(guān)于x方向的偏導(dǎo)，ha為基質(zhì)域的水勢(shì)，表示關(guān)于z方向的偏導(dǎo)，k(θa)表示土壤含水量自變量取值θa時(shí)的非飽和導(dǎo)水率；

19、采用計(jì)算機(jī)建立裂縫域水分運(yùn)移方程，如下式：

20、其中，cc為裂縫域的比水容重，θc為裂縫域的土壤含水量函數(shù)，hc為裂縫域的水勢(shì)，k(θc)表示土壤含水量自變量取值θc時(shí)的非飽和導(dǎo)水率，γw為雙域水分交換項(xiàng)，γl為層流交換項(xiàng)，γi為水平吸滲交換項(xiàng)，ξa為基質(zhì)域體積比，ξc為裂縫域的體積比，αw為水轉(zhuǎn)換率，kc-a為雙域交界處的飽和導(dǎo)水率，且w為裂縫的寬度均值，g為重力加速度，v為土壤壓縮系數(shù)；lf為裂縫左邊壁濕潤(rùn)長(zhǎng)度，lr為裂縫右邊壁濕潤(rùn)長(zhǎng)度，表示關(guān)于lf的偏導(dǎo)，表示關(guān)于lr的偏導(dǎo)；

21、步驟404、采用計(jì)算機(jī)將建立的基質(zhì)域水分運(yùn)移方程和裂縫域水分運(yùn)移方程編譯生成dll文件并加載至hydrus軟件中；

22、步驟405、采用計(jì)算機(jī)利用hydrus軟件設(shè)定時(shí)間步長(zhǎng)，總模擬時(shí)長(zhǎng)，根據(jù)建立的基質(zhì)域水分運(yùn)移方程和裂縫域水分運(yùn)移方程對(duì)煤礦開(kāi)采裂縫地區(qū)網(wǎng)格模型進(jìn)行模擬運(yùn)行，得到煤礦開(kāi)采裂縫地區(qū)網(wǎng)格模型中各個(gè)節(jié)點(diǎn)的土壤含水量。

23、上述的一種煤礦開(kāi)采裂縫地區(qū)土壤水分運(yùn)移模擬方法，其特征在于：步驟401，具體過(guò)程如下：

24、步驟4011、采用計(jì)算機(jī)利用hydrus軟件中幾何模塊geometry，建立煤礦開(kāi)采裂縫地區(qū)模型0bga；其中，煤礦開(kāi)采裂縫地區(qū)模型0bga的x方向?yàn)榱芽p寬度方向，煤礦開(kāi)采裂縫地區(qū)模型0bga的z方向?yàn)榱芽p深度向上；煤礦開(kāi)采裂縫地區(qū)模型0bga包括裂縫區(qū)域cdbf和土壤區(qū)域，距離裂縫區(qū)域邊界10cm的土壤區(qū)域?yàn)榱芽p域，剩余的土壤區(qū)域?yàn)榛|(zhì)域；

25、步驟4012、采用計(jì)算機(jī)獲取煤礦開(kāi)采裂縫地區(qū)模型0bga的邊界0b、邊界bc、邊界cd、邊界de、邊界ef、邊界fg、邊界ga和邊界ao；

26、步驟4013、采用計(jì)算機(jī)利用hydrus軟件輸入土壤質(zhì)地?cái)?shù)據(jù)、土壤孔隙度均值、殘余含水量、飽和含水量、飽和導(dǎo)水率ks、土壤孔隙第一相關(guān)參數(shù)α、土壤孔隙第二相關(guān)參數(shù)n、土壤孔隙第三相關(guān)參數(shù)m；

27、步驟4014、采用計(jì)算機(jī)利用hydrus軟件設(shè)定初始邊界條件；其中，邊界bc、邊界de、邊界fg為定通量邊界，并給邊界bc和邊界fg賦值步驟202中測(cè)量到的土壤含水量均值和土壤水吸力均值；設(shè)定邊界ao為自由排水邊界，邊界ob、邊界ga為零通量邊界，裂縫區(qū)域兩側(cè)邊界cd和邊界ef不參與水分運(yùn)移；

28、步驟4015、采用計(jì)算機(jī)利用hydrus軟件輸入步驟202中降雨強(qiáng)度均值、氣溫均值和太陽(yáng)總輻射均值。

29、上述的一種煤礦開(kāi)采裂縫地區(qū)土壤水分運(yùn)移模擬方法，其特征在于：采用計(jì)算機(jī)利用hydrus軟件對(duì)煤礦開(kāi)采裂縫地區(qū)網(wǎng)格模型進(jìn)行模擬運(yùn)行的過(guò)程中，對(duì)獲取c點(diǎn)的含水量記作d點(diǎn)的含水量記作如果大于等于θs，且小于θs，則采用計(jì)算機(jī)利用hydrus軟件設(shè)定邊界bc和邊界fg為可變通量邊界，邊界de為可變通量邊界，邊界cd和邊界ef為可變通量邊界，邊界ao為自由排水邊界，邊界ob和邊界ga為零通量邊界；

30、如果大于等于θs，且大于等于θs，則采用計(jì)算機(jī)利用hydrus軟件設(shè)定邊界bc、邊界de和邊界fg為定水頭邊界，邊界cd和邊界ef為無(wú)壓滲流自由面邊界，邊界ao為自由排水邊界，邊界ob和邊界ga為零通量邊界。

31、上述的一種煤礦開(kāi)采裂縫地區(qū)土壤水分運(yùn)移模擬方法，其特征在于：步驟201，具體過(guò)程如下：

32、步驟a、在設(shè)定的總監(jiān)測(cè)時(shí)間中，按照設(shè)定采樣時(shí)間測(cè)量到的第1類測(cè)量數(shù)據(jù)至第9類測(cè)量數(shù)據(jù)發(fā)送至計(jì)算機(jī)；其中，第j個(gè)采樣時(shí)刻第i類測(cè)量數(shù)據(jù)記作aij；i為正整數(shù)，且1≤i≤i，i表示測(cè)量數(shù)據(jù)類型總數(shù)，j為正整數(shù)，且1≤j≤j，j表示采樣總數(shù)；

33、步驟b、采用計(jì)算機(jī)對(duì)j個(gè)采樣時(shí)刻第i類測(cè)量數(shù)據(jù)處理，獲取第i類測(cè)量數(shù)據(jù)的方差σi和均值μi，并對(duì)第i類測(cè)量數(shù)據(jù)進(jìn)行判斷，如果0＜αij-μi|＜3σi，則判定為正常值；反之，則aij判定為異常值；

34、步驟c、多次重復(fù)步驟b，得到第i類測(cè)量數(shù)據(jù)的正常值和第i類測(cè)量數(shù)據(jù)的異常值；其中，第i類測(cè)量數(shù)據(jù)的異常值的數(shù)量記作

35、步驟d、采用計(jì)算機(jī)判斷如果則將第i類測(cè)量數(shù)據(jù)的正常值的均值賦值給異常值，得到預(yù)處理后的第i類測(cè)量數(shù)據(jù)；

36、如果則選擇線性插值法、指數(shù)平滑法和三次樣條插值法分別進(jìn)行擬合替換異常值，得到三種擬合替換異常值后的第i類測(cè)量數(shù)據(jù)，并選擇三種擬合替換異常值后的第i類測(cè)量數(shù)據(jù)中均方根誤差最小的擬合替換異常值后的第i類測(cè)量數(shù)據(jù)記作預(yù)處理后的第i類測(cè)量數(shù)據(jù)；其中，表示向上取整。

37、上述的一種煤礦開(kāi)采裂縫地區(qū)土壤水分運(yùn)移模擬方法，其特征在于：步驟302，具體過(guò)程如下：

38、步驟3021、采用計(jì)算機(jī)設(shè)定土壤水吸力和土壤含水量的各個(gè)點(diǎn)擬合滿足三段式分段擬合函數(shù)，如下：

39、其中，θa(s)表示第一段擬合曲線的土壤含水量函數(shù)，θr為殘余含水量，θs為飽和含水量，s為土壤水吸力自變量，α為土壤孔隙第一相關(guān)參數(shù)，n為土壤孔隙第二相關(guān)參數(shù)，m為土壤孔隙第三相關(guān)參數(shù)，且θb(s)表示第二段擬合曲線的土壤含水量函數(shù)，s′表示土壤空氣進(jìn)氣時(shí)土壤水吸力，θc(s)表示第三段擬合曲線的土壤含水量函數(shù)，sa表示擬合曲線上第一個(gè)拐點(diǎn)的橫坐標(biāo)，sb表示擬合曲線上第二個(gè)拐點(diǎn)的橫坐標(biāo)；

40、步驟3022、采用計(jì)算機(jī)按照步驟3021中設(shè)定的三段式分段擬合函數(shù)對(duì)步驟301中的各個(gè)點(diǎn)進(jìn)行曲線擬合，并得到三段式分段擬合函數(shù)的參數(shù)；其中，所述三段式分段擬合函數(shù)的參數(shù)包括殘余含水量θr、飽和含水量θs、土壤孔隙第一相關(guān)參數(shù)α、土壤孔隙第二相關(guān)參數(shù)n和土壤孔隙第三相關(guān)參數(shù)m。

41、本發(fā)明與現(xiàn)有技術(shù)相比具有以下優(yōu)點(diǎn)：

42、1、本發(fā)明方法步驟簡(jiǎn)單，設(shè)計(jì)合理，解決目前土壤水分運(yùn)移模擬方法無(wú)法直接適用于開(kāi)采裂縫地區(qū)的問(wèn)題。

43、2、本發(fā)明煤礦開(kāi)采裂縫地區(qū)模型分割為裂縫域和基質(zhì)域，并針對(duì)裂縫域和基質(zhì)域分別建立了基質(zhì)域水分運(yùn)移方程和裂縫域水分運(yùn)移方程，從而有效地適應(yīng)裂縫域和基質(zhì)域的不同特性，有效地適應(yīng)裂縫地區(qū)土壤水分運(yùn)移情況，提高了水分運(yùn)移模擬準(zhǔn)確性。

44、3、本發(fā)明建立的基質(zhì)域水分運(yùn)移方程和裂縫域水分運(yùn)移方程對(duì)煤礦開(kāi)采裂縫地區(qū)網(wǎng)格模型進(jìn)行模擬運(yùn)行，能夠適應(yīng)不同需求裂縫下土壤水分運(yùn)移狀態(tài)模擬，適用范圍廣泛，模擬結(jié)果可靠。

45、4、本發(fā)明預(yù)先進(jìn)行測(cè)量數(shù)據(jù)的采集及預(yù)處理后輸入后續(xù)煤礦開(kāi)采裂縫地區(qū)網(wǎng)格模型進(jìn)行模擬運(yùn)行，提高了預(yù)測(cè)的準(zhǔn)確性。

46、5、本發(fā)明便于后續(xù)根據(jù)土壤含水量預(yù)測(cè)土壤水分運(yùn)移速率及水分運(yùn)動(dòng)路徑，為減少環(huán)境負(fù)面影響提供指導(dǎo)，也為后續(xù)礦區(qū)生態(tài)修復(fù)提供理論依據(jù)。

47、綜上所述，本發(fā)明方法步驟簡(jiǎn)單、設(shè)計(jì)合理，建立的煤礦開(kāi)采裂縫地區(qū)模型分割為裂縫域和基質(zhì)域，并對(duì)裂縫域和基質(zhì)域分別建立水分運(yùn)移方程進(jìn)行模擬，獲取煤礦開(kāi)采裂縫地區(qū)的土壤含水量，有效地適應(yīng)裂縫地區(qū)土壤水分運(yùn)移情況，提高了水分運(yùn)移模擬準(zhǔn)確性，便于后續(xù)根據(jù)壤含水量預(yù)測(cè)土壤水分運(yùn)移速率及水分運(yùn)動(dòng)路徑。

48、下面通過(guò)附圖和實(shí)施例，對(duì)本發(fā)明的技術(shù)方案做進(jìn)一步的詳細(xì)描述。