УДК 502/504:532.5
Б. Ф. НИКИТЕНКОВ, А. В. ЕВГРАФОВ
Федеральное государственное образовательное учреждение высшего профессионального образования «Московский государственный университет природообустройства»
ИСПОЛЬЗОВАНИЕ НЕЛИНЕЙНОГО УРАВНЕНИЯ
ПАРАБОЛИЧЕСКОГО ТИПА
ДЛЯ РАСЧЕТОВ СКЛОНОВОГО СТОКА
Представлено преобразование системы уравнений Сен-Венана к нелинейному уравнению параболического типа и его конечно-разностная аппроксимация. Научная новизна состоит в специальной неявной схеме его решения, обеспечивающей устойчивость при малых слоях воды. Данная математическая модель может быть использована для практических расчетов прогнозирования поверхностного притока воды к руслам рек и оценки влияния хозяйственной деятельности при водосборе на водные ресурсы.
Склоновый сток, математическое моделирование, уравнения Сен-Венана, приток воды к руслу реки, геоинформационная система (ГИС).
There is given a transformation of the Saint-Venant equations system into the nonlinear parabolic equation and its finite difference approximation. The scientific novelty is a special implicit scheme of its solution which provides stability at small water layers. The given mathematical model can be used for practical forecasting estimations of the surface water inflow to river channels and assessment of the influence of the economic activity on water resources at the catchment basin.
Overland flow, mathematical simulation, Saint-Venant equations, water inflow to the river channel, geo-informative system (GIS).
Целью данной работы является совершенствование модели склонового стока в виде сплошного слоя, начатое в [1].
С помощью геоинформационной системы «Учебная», созданной на кафедре общей и инженерной экологии ФГОУ ВПО МГУП, производится:
оцифровка контуров земельных угодий, горизонталей, высотных отметок, линии водотока и линии водораздела. Высотные отметки берегов вводить не нужно, поскольку алгоритмы, заложенные в ГИС, позволяют автоматически наметить на линии реки необходимое количество створов и определить высотные отметки берегов интерполированием;
оцифровывка вручную с занесением высотных отметок точек, если масштаб рабочих карт 1:100 000 не позволяет качественно найти контур водораздела аналитически. Геоинформационная система, используя сплайн-интерполяцию, доводит количество точек до 300 и, после того как отметки x, y и z контура водосбора введены, определяет их высотные отметки;
с помощью дополнительных горизонталей соединение точек на берегу реки с точками на границе водосбора;
на завершающих этапах программа разбивки водосбора на полосы по линии максимального градиента. С каждой такой полосы вода стекает на участок, заключенный между соседними створами (рис. 1).
Система дифференциальных уравнений Сен-Венана для фрагмента водосборной площади записывается так:
дх дх ю
дИ
+ gra — + дх
+ga
V
Q\Q\ , 1 BqQ _ С1)
ю
K2
д< + дQ _ ^
дх дх
где X - коэффициент Буссинеска (количества движения); ф - средний по сечению расход воды, м3/сут; т — время (независимая переменная); П — коэффициент Кориолиса (кинетической энергии); х — координата (независимая переменная); ю - площадь живого сечения; g -ускорение свободного падения; Н — глубина; К — коэффициент расхода, К2 = w2C2R (С —
дг БдЯ дюЯ2
=--+ —г
дх ю К
(3)
дЯ + Б — _ Бд.
дх дх
(4)
Рис. 1. Модель водосбора с полосами склонового стекания
коэффициент Шези; Я — гидравлический радиус); ¿д - уклон склона; В - ширина полосы; д - интенсивность поверхностного стокообра-зования (разность между интенсивностью осадков д0 или снеготаяния и интенсивностью впитывания I, перехватом листвой и замкнутыми понижениями, а также другими потерами, м/сут.
Выражение в квадратных скобках представляет инерционную часть уравнения, однако этими членами нередко пренебрегают. Это связано во многом с оценкой их величины и существенным превышением остальных членов над инерционными для руслового течения. При движении воды по водосборной площади вообще трудно представить значительные величины ускорений, так как талые или дождевые воды движутся, с одной стороны, по весьма шероховатой поверхности, с другой - незначительным средним слоем, величина которого вообще сравнима с колебаниями микрорельефа, гасящими ускорение. В литературе встречаются модели, в которых движение на водосборах принимается ламинарным.
В основу предлагаемой модели склонового стока положено уравнение Сен-Венана без инерционных членов и уравнение баланса массы [2, 3].
Условие, отражающее верхнюю границу водосбора, на которой отсутствует приток воды, можно записать так:
е!,=о =
На нижней границе фрагмента водосбора, примыкающего к реке, из допустимых типов граничных условий для уравнений параболического типа можно задать следующие условия:
« = о = °; Н = о = 0 и (5).
На начало расчета (т.е. снеготаяния или выпадения осадков) т = 0, транзитные расходы, как и слой воды, равны нулю: = ° = 0; Н\; = ° = 0.
Уравнение (3) приводим к такому виду:
Я _-
К 2д
дг
(5)
БдС2Я + дЯ дх Подставив (5) в (4), получим нелинейное уравнение параболического типа:
Б дг _ д дх дх
К2 д
БдС Я + дЯ
дг дх
(6)
_ Б (до -1 )•
Расход « в дивергентной части уравнения (6) можно определить в итерациях по формуле (5) или, разрешив (5) как квадратное уравнение, получить модифицированную формулу Шези с учетом изменения массы в виде дождя и инфильтрации:
Я _-С2Ш + 2д
С 2ЯБд 2д
дг
- К2—.
дх
(7)
для уравнения (6) была разработана неявная конечно-разностная схема с обеспечением устойчивости при малых слоях воды (рис. 2) - трехточечная по длине I + 1, I, I - 1 (нумерация узлов сетки по координате х) и двухточечная по времени.
Показатель степени / относится к предыдущему моменту времени, для которого все параметры системы уравнений известны; ] + 1 - к последующему моменту с интервалом А!;, который
№ 5' 2010
Е
X Хп - 1 Хп - 2 Хп - 3
^П-1
г
п — 2
Ьп_ 1 Ьп- 2 Ьп- 3 Ь п = 0
у у у
Х У 1 Х: Х !
+ 1
^ 1
г 1-1
Ь1 + 1 Ь Ь1
• //•
—•-
г.
Ь3 Ь2
2 Ь1 ь = о
го
ей
К М
ей &
¡4 «
М
X &
ф м
Рис. 2. Расчетная схема для построения конечно-разностной аппроксимации системы
необязательно постоянен при решении и определяется сходимостью и устойчивостью схемы.
Весь водосборный фрагмент состоит из п участков, для каждого из которых записывается конечно-разностная аппроксимация уравнения (6). Для обеспечения однородности алгоритма решения задачи верхняя и нижняя части фрагмента дополнены участками нулевого размера.
для конечно-разностной аппроксимации уравнение (6) выглядит так:
дю „ дг д ¡—\ дг
дх
„дг д ¡—\ дг „/ гч
1 -_Бдг _дх(К)дг+Б(до-1), <8)
где дЮ _ Б (г« - 3 ) / Дх;
дх
дх
(9)
А К )£+в(до - / >
дх дх
1
к
К>+У2 2(22/+ ) К 1,1 +1 , ,
к + к
1+1
-К,+12 2-(¿Г -¿Л)1 + (10)
к -1+к + в, (70+12 -1Г).
Учитывая, что расчетные значения находятся в центрах участков, для аппроксимации дивергентной части уравнения (8) требуются значения параметров потоков К на границах участков, интегральные по расчетному интервалу времени дт, т.е. отнесенные к середине интервала времени, если шаг по времени мал (рис. 3). В расчетах, приведенных ниже, шаг равен 10 с.
Для решения методом прогонки сумма (9) и (10) преобразуется к следующему виду:
а1 г1++1 -Р1 3+1 + у 1 ,
где а 1 _
Ь
КН1/2
К 1,1+1
Ь1 + Ь1+1
КШ-12
_ 1
1 Ь1 |_ ь1-1 + Ь1 Р1 _ а1 + у 1 + Б / Дх;
81 _-^ + Дд0+1/2 - БД Дх
При прямой прогонке определяются коэффициенты:
^+12
р _
а.
и и _
у и. . - 8.
11 1 -1_1
1 Р1 -У1 р -1 1 Р4-у 1Р-1 а при обратной - искомые уровни:
т1 _
+ и.
1 1+1 1
г,
г,
ь + 1
Ь.
г,
ь - 1
Хь-1
Хь
Хь + 1
Рис. 3. Пояснение к построению конечно-разностной схемы
Расчетная интенсивность впитывания определяется по следующей зависимости:
Ь
Ь
ь-1
I = (уо - V! )е-ат + у^
где V - начальная скорость впитывания; V 1 -предельная скорость впитывания; а - параметр; т - время от начала впитывания.
В ходе численных экспериментов было установлено, что предлагаемая методика расчета склонового стока правильно реагирует на изменение величин входных параметров (уклона, интенсивности дождей, шероховатости и др.), а результаты расчета удовлетворяют необходимым балансовым соотношениям.
Ординаты кривых на рис. 4 с некоторыми результатами расчетов нормированы соответственно на 75, 100 и 50 % от своих максимальных значений.
Элементы моделирования новой экологической ситуации на этапе расчета склонового стока состояли в задании на участках расчетного фрагмента различных шероховатостей и параметров инфильтрации, соответствующих различным типам земельных угодий. В табл. 1
100
80-
к я
¡60-
а>
3 а
я
«
§ 40-
а
з §
а о X
20
1 2 / Вариант [^лина с] Ширина Ширина расчета ;лона 10 по верх1 по низу 1 00 м 1 100 м 100 м
/ У Уклон и Число рг Величин Цлитель Гип оса/ 01 [Счетны; а осад ко ГОСТЬ ОС) ков 1,0 точек 1 в 20 мм щков 3,( 02 ) ч
/ Л \ Шерохо! Ларамет Ларамет атость С рвпиты рвпиты ,05 *ания у0 ¡ания V] 0,01 0,001
/ 3 / \ Время Объем Во ко в< окончат впиты в >й прито ия сток ания 20, Ь 1593,с 1 5,5 ч 3 % от ( м:< юадков
/ 4 Макси Норм* \ мальньп рующий [расход множип 0,18 м 8, 'ель 0,6( 'с
3 4 5 6 7
Время от начала выпадения осадков
10
Рис. 4. Результаты расчетов склонового стока: 1 - расход с полосы за время стока, м3/с; 2 - суммарный сток с полосы, м3; 2' - суммарный сток с полосы (аппроксимирующая кривая), м3; 3 - слой осадков от начала дождя
приведены результаты модельных расчетов, показывающие, как сказывается на стоке с фрагмента изменение уклона полосы и коэффициента шероховатости.
Результаты модельных расчетов, позволяющие проследить связь суммарного склонового стока и его продолжительности со слоем осадков и длительностью их выпадения при неизменных уклонах и шероховатостях, представлены в табл. 2.
Можно отметить, что сток с полосы прямо связан со слоем осадков (в гидрологии эту связь описывают
№ 5' 2010
коэффициентом поверхностного стока) и в конечном счете мало зависит от продолжительности выпадения дождя.
Продолжительность стока в представленном диапазоне дождей с увеличением слоя дождя или продолжительности дождя растет слабо. Последнее обстоятельство можно объяснить так: при более затяжных дождях с меньшей интенсивностью последние порции дождя на верхних участках водосбора успевают впитаться, и сток прекращается раньше. В итоге общее время стока
(83
Таблица 1
Влияние уклона полосы и коэффициента шероховатости на суммарный объем склонового стока W и его продолжительность T
ж м л л сум г " max
(слой осадков - 20 мм, время выпадения - 3 ч)
h Коэффициент шероховатости
0,02 0,05 0,08
Ж™,м3 ^тах» Ч WCVM, м3 ^тах? Ч Ж™, м3 ^тах? Ч
0,001 1583 5,6 1097 6Д 724 6,3
0,005 1755 5,2 1505 5,7 1277 6,0
0,01 1785 5,0 1594 5,5 1418 5,8
0,02 1789 4,9 1643 5,4 1505 5,6
Таблица 2
Влияние слоя дождя и продолжительности выпадения осадков на суммарный объем склонового стока W и его продолжительность T
сум max
(при уклоне 0,01 и коэффициенте шероховатости 0,05)
Слой дождя, мм Продолжительность дождя, ч
1 2 3
Ж™, м3 Тщах» Ч Ж™,м3 Тщах? Ч Ж™,м3 Тщах» Ч
10 778 5,1 719 5,3 633 5,4
20 1772 5,2 1705 5,3 1594 5,5
30 2770 5,2 2699 5,4 2577 5,6
увеличивается не на то время, на которое продолжительнее дождь, а на меньшее.
Выводы
Создана и отлажена программа расчета склонового стока для фрагмента водосбора. Краткая характеристика модели: 11 входных параметров, характеризующих условия водосбора и интенсивность поступления воды на него (см. рис. 4); специальный параметр (для расчета дождевого склонового стока с фрагмента водосбора); распределенные параметры, часть значений которых задается едиными для всего фрагмента и позволяет вести непрерывный расчет а) количества воды на водосборе, б) потерь на инфильтрацию, в) стока с фрагмента в виде расхода и г) суммарного стока от начала дождя.
В отличие от работ на эту тему [2, 3], моделирование осуществлено на принципиально новом научно-техническом уровне - на базе геоинформационной системы «Учебная», созданной на кафедре общей и инженерной экологии ФГОУ ВПО МГУП с использованием оригинальных алгоритмических ходов. Научно доказана применимость данной ГИС для решения подобного рода гидрологических задач.
Предлагаемая программа расчета поверхностного притока к русловой сети после включения ее в более общую модель, учитывающую русловую трансформацию, может стать хорошим инструментом прогнозирования стока вдоль всего водотока в зависимости от антропогенных изменений на водосборе.
1. Евграфов А. В. Моделирование интенсивности склонового стока с водосборных бассейнов малых рек с использование геоинформационной системы: дис. ... канд. техн. наук. - М.: ФГОУ ВПО МГУП, 2003. - 190 с.
2. Виноградов Ю. Б. Математическое моделирование процессов формирования стока: опыт критического анализа. -Л.: Гидрометеоиздат, 1986. - 312 с.
3. Кюнж Ж. А., Холи Ф. М., Вервей А. Численные методы в задачах речной гидравлики. - М.: Энергоатомиздат, 1985. -256 с.
Материал поступил в редакцию 17.03.10. Никитенков Борис Федорович, доктор технических наук, профессор, заведующий кафедрой «Общая и инженерная экология»
Евграфов Алексей Викторович, кандидат технических наук, доцент
Тел. 8 (499) 976-09-37