Модель – wastr 6



жүктеу 120.63 Kb.
Дата12.04.2019
өлшемі120.63 Kb.

Модель – WASTR - 6.

Основные уравнения


Алгоритм основан на решении системы уравнений, описывающих движение влаги и солей в насыщенной и ненасыщенной зонах почвогрунта. Предполагается наличие следующих ионов:

1) в почвенном растворе – Ca2+, Mg2+, Na+, Cl-, SO42-, HCO3-;

2) в почвенно-поглощающем комплексе – Ca2+, Mg2+, Na+;

З) в кристаллической форме возможно существование гипса – СaSO42H2O

Модель может, быть использована для прогнозов влагосолепереноса в почвогрунтах.

В одномерной постановке уравнения модели имеют следующий вид:



; (19)

; (20)

; (21)

; (22)

; (23)

; (24)

; (25)

; (26)

; (27)

; (28)

; (29)

; (30)

Уравнение (19) описывает одномерное движение влаги в зовах неполного и полного водонасыщения почвогрунта. Уравнение (20) определяет скорость влагопереноса. Уравнения (21) – (25) списывают конвективно-диффузионный перенос рассматриваемых ионов с учетов ионообменной сорбции и растворения - кристаллизации солей. Уравнения (26), (27) — изотермы катионного обмена Б. П. Никольского. Уравнение (28) задает постоянство емкости поглощения во времени. Уравнение (29) описывает кинетику растворения-кристаллизации гипса. Здесь приняты следующие обозначения:

H - обобщенный потенциал почвенной влаги, м; H = P + x,

P - капиллярный потенциал, м;

X - вертикальная координата (x = 0 на поверхности почвы, положительное направление вниз);

W - объемная влажность – отношение объема воды в дм3 к общему объему в дм3, % или доли единицы;

- капиллярная влагоемкость, м-1;

k(W) - коэффициент влагопроводности, м/сут;

e - функция отбора, влаги корнями растений, сут;

V - скорость влагопереноса, м/сут;

С1, С2, С3, С4, С5, С6 – концентрации Ca2+, Mg2+, Na+, Cl-, SO42-, HCO3- в почвенном растворе, мг-экв/л;

N1, N2, N3 – содержание Ca2+, Mg2+, Na+ в ППК, мг.экв/100 г. почвы;

N15 – концентрация кристаллического гипса в почве, мг.экв/100г почвы;

С15 – концентрация гипотетической соли CaSO4 в почвенном растворе, мг-экв/л (в программе предусмотрена связка ионов в гипотетические соли);

СН15 - концентрация предельного насыщения раствора гипсом, мг-экв/л. Значение СН15 вычисляется в программе по эмпирической формуле Ю.Н. Чиркина в зависимости от содержания в растворе хлоридов и сульфатов (С24, С25, С34, С35 – содержание в растворе гипотетических солей MgSO4, Na2SO4, MgCl2, NaCl)

; (31)

D* – коэффициент конвективной диффузии ионов, м2/сут;

 – объемная масса скелета почвы, г/см3;

k31, k12 – константы изотерм ионного обмена Na - Ca и Ca Mg

N0 – сумма обменных оснований, мг.экв/100 г почвы;

15 – коэффициенты скорости растворения гипса, сут;



t – время, сут.

Начальные и краевые условия


Начальные условия задают исходное распределение потенциалов влаги (или влажности), концентрации ионов Ca2+, Mg2+, Na+, Cl-, SO42 и гипса по координате

; или ;

; ; ; (32)

; ; ;

где: Н0(х), W0(х) - начальные значения потенциала (м) или влажности (дм3/дм3);



Сn0(х) - начальное содержание растворенных ионов Ca2+, Mg2+, Na+, Cl-, SO42 (мг-экв/л или %)

N150(x) - исходное содержание гипса в почвогрунте (мг-экв/100 г или %).

Краевые условия задают значения потоков влаги и солей на верхней и нижней границах области расчета.

Для уравнения влагопереноса на верхней границе предусмотрено задание двух видов краевых условий:

1) условие 1го рода (поверхностный полив)



, ; (33)

где: H1(t) – слой воды на поверхности почвы (м) при поливе затоплением и по полосам; при поливе по бороздам следует принимать H1 = - 0,2 м.

2) условие 2го рода (осадки, испарение, полив дождеванием)

, ; (34)

где: Q1 – величина потока через поверхность (м/сут)

В программе WAS052 величина Q1 на заданный интервал времени T0 вычисляется по формуле:

; (35)

где: Oc, Nn, Ec, Tp - соответственно осадки, поливная норма, суммарное испарение и транспирация, (м3/га) за период Т0 (сут) Если Q1 <0 , то испарение превышает осадки, при Q1 > 0 осадки и поливная норма превышают испарение.

Для уравнения солепереноса на верхней: границе условно 3го рода Данквертса-Бреннера

, , i=1, , 5; (36)

Cni – концентрация соответствующих ионов в поливной воде (мг-экв) при поливе; в межполивной период Cni= 0.

На нижней границе области расчета возможно задание следующих краевых условий:

1) условие 1го рода – потенциал

, ; (37)

где: H3(t) – заданный потенциал, соответствующий определенней влажности, или уровень грунтовых вод (от поверхности почвы - со знаком минус) Это условие для расчетов используется редко, как правило, для обсчета лизиметрических данных.

2) условие 2го рода - поток

, ; (38)

где: Q2(t) – заданный поток через нижнюю границу;



Q2(t) < 0 – отток за пределы области расчета (переток через нижележащие слои, боковое растекание);

Q2(t) > 0 – приток в области расчета (напорное питание, приток с соседней территории)

Q2(t) = 0 – отсутствие потока (водоупор или равенство притока и оттока).


  1. условие 3го рода – поток, зависящий от H; в выражении (10) :

; (39)

В программе WAS052 данное условие моделирует вертикальный отток в центре междренья при работе горизонтального дренажа в зависимости от его параметров.


4) условие 2го рода - свободное стекание влаги

, ; (40)

Данное условие используется при глубоких грунтовых водах, когда их подъем не прогнозируется, т.е. предполагается, что вода проходит через нижнюю границу под действием гравитационных сил. Для уравнения солепереноса на нижней границе принято:



, ; (41)

т.е. полагается, что диффузионный поток ионов мал.


Параметры массопереноса

Параметры уравнения влагопереноса


Для уравнения влагопереноса необходимо задавать два основных параметра:

  • зависимость капиллярного потенциала от влажности (основная гидрофизическая характеристика – ОГХ)

  • зависимость влагопроводности от влажности.

При неполном влагонасыщении капиллярный потенциал является функцией влажности. Он отрицателен в зоне неполного насыщения. На границах зон неполного и полного насыщения он равен нулю. Ниже уровня грунтовых вод он положителен, величина его меняется по гидростатическому закону. Вид кривых P(W) приведен на рис.5.2. Определение ОГХ производится для каждого литологического слоя экспериментально на образцах почвогрунта. Для этого могут быть использованы методы мембранного и пластинчатого пресса, капилляриметрические и другие. Затем экспериментальная зависимость P(W) должна быть аппроксимирована функцией:

Р(W ); (42)

где: Wп – полная влагоемкость (меньше значения пористости на величину содержания защемленного воздуха);

W0 – нижний предел иссушения почвогрунта (максимальная гигроскопичность);

hk* – условная высота капиллярного поднятия.

Параметры W0, hk* должны находиться на основе аппроксимации экспериментальных кривых P(W). Если таковых не имеется, то априори можно принять W0 – максимальная гигроскопичность.



; (43)

где: hk* – высота капиллярного поднятия;  = 1,2 (чем тяжелее грунт, тем больше )

Экспериментальная зависимость влагопроводности K(W) должна быть аппроксимирована формулой С.Ф. Аверьянова

; (44)

где: K0 - коэффициент фильтрации, m - пористость. Значение W0 принимается тем же, что и в зависимоcти P(W). В результате аппроксимации находится параметр Nkw, который зависит от механического состава почв (безразмерный, т. е. число) (рис. 4.3).

Для расчета в программу закладываются значения m, Wп, W0, k0, hk*, nkw. Экспериментально определенные гидрофизические характеристики при изыскательских работах.
Для расчета отбора влаги корнями растений в программе используется формула, учитывающая зависимость транспирации от влажности и распределения корневой системы

ek – интенсивность транспирации из всего корнеобитаемого слоя, м/сут;

e(W) – функция, учитывающая влияние влажности на отбор

; (45)

где: – влажность завядания (Wз~1,4Wмг),



Wппв – предельно-полевая влагоемкость;

f(x) – функция распределения корневой системы (рис. 4.3)

; (46)

x, x1, лпараметры распределения, определяемые в программе по коду культуры;

fкор – мощность корнеобитаемой зоны.

В качестве исходных данных вводятся значения Wз, Wппв, fтр (fтр – доля транспирации от суммарного испарения).


Параметры солепереноса


Для расчета задаются следующие параметры солепереноса: доля «сквозных» пор – ; приведенный параметр гидродисперсии - *; константы изотерм катионного обмена k31, k12, коэффициент скорости растворения 12.

Параметр гидродисперсии (приведенный) c* и "активная" пористость mak определяются по экспериментальным данным промывки монолитов для каждого литологического слоя (по иону Cl-) Коэффициент конвективной диффузии в программе вычисляется по формуле:



; (47)

где: V - скорость фильтрации. Для ориентировочных расчетов можно определять c по графику приведенному на рис. 4.2. Значение c*, используемое в программе, связано с c соотношением:



; (48)

mak - активная пористость.

Koэфициенты изотерм катионного обмена рассчитывают на основе экспериментальных данных совместного определения катионов Ca2+, Mg2+, Na+ в почвенном растворе и ППК. Навески почвы заливаются растворами солей (обычно CaCl2 или NaCl) различной концентрации (1, 5, 10 г/л) и выдерживаются 1-2 сут до установления катионного равновесия. Затем раствор отфильтровывается и выполняется анализ на содержание в нем и ППК почвы ионов Ca2+, Mg2+, Na+. По результатам эксперимента для каждого горизонта составляется таблица, в которую записывают значения С1, С2, С3 (мг-экв/л), N1, N2, N3 (мг-экв/100г), С1/С2. Полученные значения наносятся на графики: на одном по оси абсцисс откладывают значения С3/С1, а по оси ординат N3/N1 , на другом – соответственно С1/С2 и N1/N2. Подученные области точек, аппроксимируются прямыми линиями. По тангенсу угла наклона прямой на первом графике определяется k31, а на втором k12.

Коэффициенты изотерм катионного обмена составляют: k31~0,04–0,08; k12~1–3.

Величина коэффициента интенсивности растворения гипса 15 варьирует в довольно широких пределах (0,01 – 10 сут). Подбираться по данным промывки монолитов, содержащих гипс, путем сопоставления расчетных и экспериментальных выходных кривых ионов SO42-.


Численная реализация модели на ЭВМ


Система уравнений 19 – 25 решается численным методом конечных разностей. Область расчета [0,L] разбивается на узлы равномерной или неравномерной сеткой . Задается число узлов сетки kx и их координаты i. Первый узел располагается на поверхности почвы 1 = 0, последний – на нижней границе kx= L. Внутри области расчета строение почвогрунта может быть неоднородным. Задается число литологических слоев ks и верхняя граница каждого слоя Hsi, тогда Hs1=1 (рис. 4.4).

Для каждого литологического слоя задаются свои гидрофизические и гидрохимические характеристики. Исходные значения потенциалов или влажности и концентрации солей задаются в узлах, число которых kj может не совпадать с kz. По числу узлов задаются координаты Hj, в которых в начальный момент времени заданы исходные эпюры. Программа их аппроксимирует по координате x - ломаной (потенциалы, влажности и концентрации Cj) и ступенчатой (N0, N15) функциями.

Для прогнозирования влагосолепереноса каждый расчетный год разбивается esl на отдельные периоды (задается число периодов в году).

Nпер и их продолжительность Т0i, i = 1, 2, Nпер (с определенными значениями испарения Eci, осадков Oci, поливных норм Nni и нижнего краевого условия). Каждый период программа разбивает неравномерной сеткой по времени и производит расчет. За начальное распределение компонент для каждого периода принимается распределение, полученное на конец предыдущего периода. В программе предусмотрено несколько режимов расчета, определяемых переменными MR, MH, MC. Уравнение влагопереноса может решаться без уравнения солепереноса (при MH = 1, MC = 0), и наоборот (при MH = 0, MC = 1), а также совместное их решение (при MH = 1, MC = 1).

Вводимые параметры приведены ниже - построчно.

Форматы ввода:

1) целые переменные KX, KS, KJ, MR, MHW, MPF, G2, MD, MC, MGP, MIZ, MAK, MCP, ND, NM, NPEP, KE вводятся в формате 12 (целое число, заснимающее две позиции), NP, KP, KCP - в формате 11 (целое число, занимающее одну позицию), NG - в формате 14 (целое число, занимающее четыре позиции). Сокращения будут пояснены ниже при пояснении схемы (Рис. 4.6) ввода исходной информации для программы WAS061 отражающей Модель –2.

2) остальные переменные вводятся в формате E6.2 - вещественное, число со знаком или без него, с десятичной точкой (обязательно).

ПРИМЕР ВВОДА ИСХОДНОЙ ИНФОРМАЦИИ К ПРОГРАММЕ WASO61

(исходная информация формируется в 8-м групп,

описание которых приводится ниже)

1) Исходные данные,формирующие расчет

kx ks kj mr mhw mpf g2 md mh mc mgp miz mak mcp

11 4 7 1 1 0 3 2 1 1 1 1 1 2

kx -число узлов разностной сетки

ks -число литологических слоев

kj -число узлов в исходной эпюре

mr -режим счета:1-водно-солевой,2-промывка

mhw-признак исходной эпюры:h-потенциал,w-влажность

mpf=0

g2 -признак нижнего краевого условия:2-поток,3-горизонтальный



дренаж,4-свободное стекание

md - признак горизонтального дренажа:0-нет,1-однослойная схема, 2-двухслойная схема

mh -признак уравнения влагопереноса:0-нет,1-есть

mgp-признак размерности исходных концентраций:0-мг.экв/л жидкая

фаза, мг.экв/100г-твердая фаза, 1-%-жидкая и твердая фазы

miz=1


mak-признак учета активностей:0-не учитывается,1-учитывается

mcp-число вариантов качества оросительной воды

2)Координаты точек: 50 02.50 03.50

3.Параметры слоев (число слоев ks= 4)

а)водно-физические параметры

hs пор пв ппв вз мг кф нк G nkw nu

00.00 0.566 0.509 0.400 0.181 0.129 0.790 02.00 1.160 6.000 1.000

00.99 0.561 0.505 0.444 0.155 0.111 0.110 02.50 1.190 6.500 1.000

01.60 0.563 0.533 0.400 0.169 0.121 0.380 02.25 1.310 7.000 1.000

02.60 0.558 0.530 0.400 0.122 0.087 0.380 02.00 1.310 7.000 1.000 hs -верхняя граница слоя,м

пор-пористость

пв -полная влагоемкость

ппв-предельная полевая влагоемкость

вз -влажность завядания

мг -максимальная гигроскопичность

кф -коэффициент фильтрации,м/сут

нк -условная высота капилярного поднятия,м

G -об емная масса скелета почвы

nkw-показатель степени в формуле влагопроводности Аверьянова

б)химические параметры

каппа LMD KNA-CA KCA-MG ППК В15 В16

1.000 1.640 0.030 2.180 28.51 0.020 0.015

1.000 1.550 0.020 1.970 37.64 0.020 0.015

1.000 0.620 0.017 1.230 35.93 0.020 0.015

1.000 0.440 0.017 1.230 28.43 0.020 0.015 каппа -доля сквозных пор

LMD -параметр гидродисперсии,м

KNA-CA -коэффициент изотермы обмена,Na-Ca

KCA-MG -коэфициент изотермы обмена, Ca-Mg

ППК -сумма обменных оснований,мг.экв/100г

B15 -коэффициент интенсивности растворения гипса

В16 -то же для кальцита

4) Параметры исходной эпюры:

H(j) j=1,...,kj -значения потенциала -1.500-1.500-1.500-1.500-1.500-1.500-1.500

C1(j) J=1,...,KJ -концентрация Ca в поровом растворе

0.005 0.005 0.048 0.008 0.010 0.010 0.010

C2(j) j=1,...,kj -концентрация Mg в поровом растворе

0.00090.00290.00350.00300.00470.00470.0047

C3(j) j=1,...,KJ -концентрация Na в поровом растворе

0.025 0.023 0.090 0.034 0.027 0.022 0.022

C4(J) j=1,...,kj - концентрация Cl в поровом растворе

0.013 0.013 0.026 0.014 0.013 0.013 0.013

C5(j) j=1,...,kj -концентрация SO4 поровом растворе

0.048 0.048 0.290 0.072 0.054 0.054 0.054

C6(j) j=1,...,kj -концентрация гипса в твердой фазе

0.00

C7(j) j=1,...,kj- концентрация кальцита в твердой фазе



02.00 02.00 02.00 03.00 03.50 04.00 04.00

C8(j) j=1,...,kj - парциальное давление CO2 в почвенном воздухе

0.00750.00750.00750.00750.00750.00750.0075

5) Данные для расчетов дренажа

B DDP HDP TB KФ PH

600.0 1.000 2.500 11.00 0.380 56.0

В- междренное расстояние,м

DDP-эффективный диаметр дрены,м

HDP-глубина заложения дренажа,м

TB- мощность покровного слоя от поверхности земли до водоупора

КФ -коэффициент фильтрации покровного слоя,м/сут

PH -проводимость подстилающего слоя,м2/сут

6) Качество поливной воды ( количество вариантов mcp=2)

Сп1 Сп2 Сп3 Сп4 Сп5

16.00 10.40 15.18 7.200 24.94

16.00 10.40 15.18 7.200 24.94

Cп1- концентрация ,мг.экв/л, в поливной воде Ca

Сп2- то же Mg

Сп3- то же Na

Сп4 -то же Cl

Сп5 -то же SO4

6) Задание календарного времени на начало засчета

11 03 1961

7) Исходные данные для сельскохозяйственной культуры

Nпер KE

07 2


Nпер-число расчетных периодов в году

KE -код сельскохозяйственно культуры:

1-кукуруза

2-подсолнечник

3-озимая пшеница

4-яровая пшеница

5-свекла

6-люцерна

7-картофель

8-томаты


9-яблоня

10-виноград

11-хлопчатник

8) Статьи водного баланса и данные о корневой системы.

(число расчетных периодов Nпер=7)

KP KCP TO ESM FTR OS PNOR HKOR

KP-признак печати результата (0-печати нет,1-печатаются эпюры в узлах Х)

KCP- номер варианта качества воды

TO - расчетный интервал времени

ESM- суммарное испарение за время ТО, м3/га

FTR- доля транспирации от ESM

OS - осадки за период ТО минус поверхностный сток, м3/сут



PNOR-поливная норма за период ТО, м3/сут

HKOR-глубина корневой системы, м

Достарыңызбен бөлісу:


©kzref.org 2019
әкімшілігінің қараңыз

    Басты бет