О новом подходе к расчету склонового стока Текст научной статьи по специальности «Математика»

ГЕОГРАФИЯ, ГЕОЛОГИЯ, ГИДРОЛОГИЯ

УДК 556.51.512

О НОВОМ ПОДХОДЕ К РАСЧЕТУ СКЛОНОВОГО СТОКА

© 2008 г. С.Л. Алита, В.С. Инюхин

The method suggested in this article is applicable to calculations of runoff rates for slopes of variable width and slope angle. One of the main advantages of our method is the fact that the values of integration steps over distance and time are independent of each other which makes it convenient for practical application.

Для осуществления оперативного прогноза паводков ливневого происхождения на горных реках на базе радиолокатора МРЛ-5 в ВГИ был создан метод дистанционного измерения количества и интенсивности осадков, выпадающих на водосборы горных рек. Полученная информация является исходной для математической модели, описывающей процессы дождевого стока.

При создании математической модели дождевого стока в реку определяющую роль играет выбор схемы интегрирования дифференциальных уравнений, описывающих процесс склонового стекания. В результате решения системы дифференциальных уравнений определяется пространственно-временное распределение интенсивности притока дождевой воды в реку. На основании этого распределения делаются выводы об изменениях уровней воды в притоках реки и самой реке.

При склоновом стекании роль инерционных сил и изменения глубин потока оказывается малой по сравнению с силой трения, и скорость потока можно представить в виде однозначной зависимости от глубины. Такое допущение приводит к модели течения воды называемой кинематической волной [1].

Наиболее известной из использующихся в гидрологии моделей стока с распределенными параметрами являются уравнения Сен-Венана [1-5]. Уравнения кинематической волны, описывающие неустановившееся течение воды по склону единичной ширины, можно представить в виде

дИ дд

— + — = г, (1)

3/ дх

где к - глубина потока, м; д - расход воды на единицу ширины склона, м2/с; t - время, с; х - расстояние от линии водораздела, м; г - интенсивность эффективных осадков, м/с.

Расход воды на единицу ширины склона можно определить по формуле

д = а-кт, (2)

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

а = С-у11 , т = 3, (3)

где С - коэффициент Шези; I - уклон склона.

После подстановки (2) и (3) в (1), уравнение преобразуется к виду

dh 3 dt 2

(4)

Полученное выражение (4) является уравнением в частных производных относительно к(хД) и не допускает интегрирования в общем виде. Численному интегрированию уравнения Сен-Венана посвящены десятки исследований, накоплен значительный опыт расчета неустановившегося течения с помощью этого уравнения. Однако многообразие решаемых задач и различия в исходной информации не позволяют выделить из существующих методов такой, который можно было бы считать универсальным для всех задач. Выбор того или иного метода зависит от свойств конкретной задачи, вычислительных мощностей, находящихся в распоряжении расчетчика, требуемой точности, а также объема исходной информации.

Существуют различные способы численного интегрирования уравнения (4), в основном это конечно-разностные методы [4, 5], которые основаны на аппроксимации производных, входящих в дифференциальное уравнение, разностными отношениями, причем в зависимости от способа этой аппроксимации эти методы могут быть разделены на две основные группы - явные и неявные схемы. В явных схемах определение искомых значений зависимых переменных производится для каждого из узлов разностной сетки по значениям, известным с предыдущего временного слоя. В отличие от этого в неявных схемах значения зависимых переменных находятся для всех узлов нового временного слоя одновременно путем решения системы алгебраических уравнений, включающих условия в граничных точках. Среди конечно-разностных методов особое место занимает метод характеристик, в котором исходное уравнение в частных производных преобразуется в систему двух обыкновенных дифференциальных уравнений. Метод характеристик обеспечивает высокую точность расчетов, его часто используют в качестве эталонного, однако практическое применение этого метода ограничено из-за того, что искомые величины находятся не в фиксированных точках фазовой плоскости х - а в тех, которые находятся в процессе вычисления. Шаги по времени при этом оказываются достаточно малыми.

Данная статья посвящена описанию метода расчета склонового стока, основанного на получении формул для численного интегрирования непосредственно

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

Рассмотрим склон постоянного уклона I, на который выпадают осадки интенсивностью г, м/с, и отсутствует инфильтрация воды в почву. Сначала, для простоты, будем считать ширину склона постоянной и равной Ь, м.

На рис. 1 изображен участок склона длиной йх, на который выпадает объем воды, равный Увып за время А.

Рис. 1. Схема для расчета склонового стока

Под воздействием выпавшего объема воды глубина на участке увеличится на йк. Сверху на участок за время й будет притекать объем воды, равный Упр, с верхнего участка, а на нижний участок за это же время стечет объем воды, равный Уот.

Изменение объема воды на участке вычислим по формуле, исходя из закона сохранения массы (объема) воды:

М = Уж +Уи -Уго ■ (5)

Объем выпавшей воды будет: Уш =с!Н-г-с11. (6)

где йБ - площадь элементарного участка; г - интенсивность эффективных осадков перпендикулярно склону.

Площадь элементарного участка склона определяем по формуле

с18 = Ь-с1х. (7)

С учетом (7) формула для расчета объема выпавшей на участок воды примет вид

Уаш =1' - Ъ ■ (Ье- сЬ . (8)

Объем воды, стекший с вышележащего участка: У~1д=Ъ-кй.11й.&, ' (9)

где кв, ив - соответственно глубина и скорость потока на входе в участок. Скорость потока определяем по формуле Шези [3]:

иа=Са-41а'ка . (10)

Подставляя в формулу (9) выражение (10), получим

(Н&)3/2 -Я. (11)

Аналогично определяем объем воды, стекший за время й на нижний участок склона:

у-1д =Ъ-С1 (к1 )3/2 а. (12)

где кн - глубина потока на выходе из участка.

Далее выразим суммарное изменение объема воды на участке через изменение глубины потока:

'ау^ъ-ах-аи . ' (13)

Подставляя выражения (8), (11), (12), (13) в формулу (5) и произведя преобразования, получим выражение для частной производной по времени глубины потока на участке:

3 3 _ _ 1

^=+Фа4П ■ 2 -щ 4П-■ (ы у2)—. (М)

от ах

а значение самой глубины будет иметь вид &г

h = Iiq + —— ■ dt

dt

(15)

Ah = rAt + (NäJTä (hä)2 - (к()2)

где к0 - глубина на момент начала текущего шага времени.

Переходя от бесконечно малых приращений к конечным приращениям, получаем формулу, пригодную для непосредственного численного интегрирования:

3 3

Ах

Для того чтобы распространить эту формулу на склон переменной ширины, следует задаться законом изменения ширины склона как функцией его длины; для определенности будем считать эту зависимость линейной: Ъ = к ■ х + йд - где Ь - текущая ширина склона; Ь0 - начальная ширина склона; х - расстояние до верхней границы водораздела; к - постоянный коэффициент.

Выражения для ширин верхней и нижней границ участка примут вид = к • х + Ь0 и = к -х + к ■ сЬс + Ь0 .

При вычислении площади участка будем исходить из средней ширины склона на участке:

= (к ■ х + ■ к ■ с1х + ¿о) • с!х.

С учетом всех вышесказанных замечаний формула (14) для склона переменной ширины примет вид

дН dt

= r+-

1

к -х + к ■1 • dx + br> 2 0

- ^-x + k-dx+bQ j^V/ ■■JTj'■

3~\

J_

dx

а формула (16) Ah = r • At + —

к -х + к — Ах+ br, 2 0

{■x + bo - fx + kAx + b0% Jl^iii ^ -

3

-s

Вычисления глубины потока как функции расстояния и времени сводятся к расчету значений глубины в узлах сетки, представленной на рис. 2.

Начальные и граничные условия задаются в виде А(0, х) = 0, h(t,0) = 0.

1

3

СО

X

m

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

Рис. 2. Схема численного интегрирования

На рисунке по горизонтальной оси откладываются расстояния от верхней границы водораздела, а по вертикальной - время от начала отсчета. Стрелки, проведенные между узлами, показывают очередность вычислений, т.е. из каких узлов берутся данные для расчета глубины в данном узле. Сначала рассчитываются значения к в горизонтальном ряду (для одного момента времени), затем происходит переход к следующему моменту и так далее (явная схема). Расчетные формулы для определения //(1,х) имеют вид

hi+1 j+1 ~hj+1"

■M,

Ч 7+1 •

(17)

где j+\ = r-At + -

k - x + k — Ax + bn 2 U

<(k^ + bo ^Cjjîjii j J- kx + kAx + bo ~Cj+1: '^ïij+lJ At

J Ах

В приведенных формулах i - индекс изменения глубины по времени; j - индекс изменения глубины по расстоянию.

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

r-säid = А-А -

где r - интенсивность осадков; Säiö - величина площади горизонтальной проекции склона; ык, Ик, Ък -скорость, глубина и ширина потока в нижнем створе склона.

Для иллюстрации работы метода был произведен сравнительный расчет склонового стока при помощи

уравнения Сен-Венана и предлагаемым методом прямого интегрирования. Для расчета был выбран склон длиной 200 м и шириной 1000 м, находящийся под уклоном 10о к горизонту, коэффициент шероховатости на протяжении всего склона был принят равным 0,03. На склон по всей его длине и ширине выпадают осадки интенсивностью 100 мм/ч. Рассматриваемый промежуток времени составил 260 с.

Интегрирование уравнения (1) производилось методом характеристик, который был выбран исходя из того, что он является одним из самых точных и часто используется в качестве эталонного для оценки точности других методов. В качестве эквивалентной характеристической системы обычных дифференциальных уравнений рассматривалась система:

dt 2 dh

dt '

В результате решения этой системы было получено уравнение характеристик, лежащих в плоскости x 0 t и выходящих каждая из своей начальной точки (x0,t0):

х = х0 +--) + r-t-r-tQj-

Г- Г

-—■h(t0)i,

Г

в каждой точке характеристики глубина потока находится как h(t) = Л(/0 ) + г ■ t - г ■ /0.

По результатам расчета был построен график значений глубины потока как функции времени и расстояния (рис. 3).

Рис. 3. Значения глубины как функции расстояния и времени

Пунктирными линиями показаны характеристики, сплошными горизонтальными линиями - интеграль-

1

ные линии метода прямого интегрирования. В местах пересечения характеристик и интегральных линий каждой точке соответствует два числа, верхнее является глубиной, рассчитанной методом характеристик, нижнее - глубиной, рассчитанной методом прямого интегрирования. Рассматриваемый интервал времени соответствует неустановившемуся процессу стекания.

Предлагаемый метод прямого интегрирования, хотя и относится к явным методам, но обладает хорошей устойчивостью к изменению исходных данных. В случае возникновения счетной неустойчивости можно поменять соотношение между величинами шагов по времени и расстоянию и тем самым устранить ее. Предлагаемый метод хорошо подходит для расчетов стока воды по склонам переменной ширины и уклона. Склон рассматривается как каскад последовательно соединенных склонов переменной ширины. Внутри каждого участка постоянными должны оставаться уклон и коэффициент шероховатости, при этом ниж-

ний створ верхнего участка рассматривается как верхний створ нижнего участка. Расчет стока по каждому из участков производится по формуле (17). К достоинствам предлагаемого метода можно отнести независимость друг от друга величин шагов интегрирования по расстоянию и времени, что делает данный метод удобным для практического применения.

Литература

1. Кучмент Л.С., Демидов В.Н., Мотовилов Ю.Г. Формирование речного стока. М., 1983.

2. Давыдов Л.К., Дмитриева А.А., Конкина Н.Г. Общая гидрология. Л., 1973.

3. Аполлов Б.А., Калинин Г.П., Комаров В.Д. Курс гидрологических прогнозов. Л., 1974.

4. Денисов Ю.М. Схема расчета гидрографа стока горных рек. Л., 1965.

5. Кучмент Л.С. Математическое моделирование речного стока. Л., 1972.

Высокогорный геофизический институт, г. Нальчик_17 августа 2007 г.