| Результаты этапа: По данным МЧС наводнения — самое распространенное стихийное бедствие, причем как по площади распространения и по повторяемости, так и по причиняемому ущербу. По количеству человеческих жертв наводнения находятся на втором месте после землетрясений [Государственный доклад…, 2018]. Вскрытие северных рек Европейской территории России часто сопровождается заторами льда. Причиной этого является более позднее вскрытие рек в нижнем течении. Крупные реки этого района — Северная Двина, Печора, Онега — текут с юга на север, и вскрытие их, следовательно, начинается с верховьев. Перемещаясь, лёд попадает в районы, где река ещё не вскрылась и ледяной покров достаточно прочен, что приводит к розливу реки на пойму, где зачастую находятся населенные пункты.
Полностью избежать наводнений невозможно, возможно лишь ослабить или предотвратить негативное воздействие. Чтобы оценить возможный ущерб, выявить возможные зоны затопления, влияние отдельных факторов, воздействующих на реку во время прохождения паводка, всё чаще и активнее применяются методы математического моделирования. Эта отрасль в настоящее время стремительно развивается благодаря появлению достаточно мощных вычислительных машин, ГИС-технологий. Во многих зарубежных странах особое внимание уделяется разработке технологий, направленных на прогнозирование и оповещения населения о наводнениях [Buizza et al., 1999; Demeritt et al., 2007; Bogner, Kalas, 2008; Alfieri et al., 2011; Cloke et al., 2013].
На данный момент в России система предупреждения населения в случае наводнений развита слабо, несмотря на очевидную ее необходимость, особенно в районах, где наводнения случаются регулярно. Только в XXI веке произошло множество наводнений, повлекших за собой человеческие жертвы.
Одним из наиболее подверженных затоплению районов северной части ЕТР является бассейн Северной Двины, где регулярно формируются ледовые заторы, приводящие к подъему уровня воды на 9 и более метров. В данной работе на примере слияния рек Сухона и Юг, где расположен город Великий Устюг рассмотрены особенности прохождения наводнений и возможности прогнозирования характеристик наводнений на основе синтеза различных типов гидрологических моделей.
В рамках данного проекта за два года были решены следующие задачи:
1) Произведена адаптация всех используемых в работе моделей (STREAM_2D, ECOMAG и MSFR_WG) для исследуемого участка реки Северная Двина.
2) Разработана методика верификации двумерной гидродинамической модели STREAM_2D при условии нехватки гидрологических данных на основе автоматической классификации пикселей космического снимка.
3) Разработана методика ансамблевого краткосрочного прогноза характеристик осенне-летних паводков и долгосрочного прогноза весеннего половодья.
4) Получены значения калибровочных параметров модели STREAM_2D (коэффициента шероховатости и величину повышения отметки дна) в зависимости от мощности ледового затора на основе моделирования затопления за многолетний период.
5) На основе совместного использования динамико-стохастического моделирования, физико-математического моделирования стока и гидродинамического моделирования рассчитаны наиболее неблагоприятное соотношение метеоусловий, приводящее к прохождению 0,1 % расхода воды, с целью дальнейшего определения характеристик подобного наводнения.
6) Была создана программа на языке Python, объединяющая модель формирования стока ECOMAG и двумерную гидродинамическую STREAM_2D, с возможностью визуализации полученных результатов.
7) Были разработаны и созданы подпрограммы модели STREAM_2D, позволяющие получать осредненные характеристики затопления по руслу и пойме исследуемого участка, а также позволяющие рассчитать время стояния воды.
8) Были смоделированы и проанализированы катастрофические наводнения 1998 и 2016 года, в том числе с применением детальной расчетной сетки.
9) Результаты расчетов были визуализированы в виде графиков и тематических карт.
Материалы и методы.
В данной исследовании были использованы различные типы гидрологических моделей, каждая из которых нуждается в определенных исходных данных.
Гидродинамическая модель STREAM_2D.
В основу создания гидродинамической модели исследуемого участка вошли материалы экспедиционных исследований, проводимые сотрудниками научно-исследовательской лаборатории эрозии почв и русловых процессов им. Н.И. Маккавеева и кафедры гидрологии суши географического факультета МГУ имени М.В.Ломоносова [Отчет о научно-исследовательской работе, 2013]. В ходе данных полевых работ были получены:
1. Промеры глубин на реках Сухона и Юг, а также от города Великий Устюг и до города Котлас;
2. Уклоны водной поверхности рек Сухона, Юг и Малая Северная Двина
3. Измеренные расходы воды и скорости течения рек Сухона, Юг и Малая Северная Двина
4. Морфометрические створы через долины рек Сухона, Юг и Малая Северная Двина
5. Положение и высотные отметки дорожных насыпей, наличие водопропускных отверстий.
6. Для создания цифровой модели рельефа использовались данные о рельефе с крупномасштабных топографических карт (1:25000) и планов города Великий Устюг масштаба 1:10000.
В качестве входных данных для модели STREAM_2D использовалась архивная информация среднесуточных режимных наблюдений (за уровнем воды и расходом воды) на гидрологических постах реки Сухона (гидрологические посты в деревне Каликино и городе Великий Устюг), реки Юг (гидрологические посты в деревне Гаврино, селе Кичменский Городок и поселке городского типа Подосиновец) и реки Малая Северная Двина (гидрологические посты в деревне Медведки и городе Котлас), полученная по материалам гидрологических ежегодников за период с 1980 — 2007 год. Данные с 2008 по 2018 год были получены с помощью информационного ресурса АИС ГМВО Федерального Агентства Водных Ресурсов [https://gmvo.skniivh.ru]. Данные об уровнях воды по гидрологическим постам за 2019 год были предоставлены Вологодским ЦГМС. Для определения типа подстилающей поверхности, с целью задания коэффициента шероховатости, использовались разновременные оптические и радиолокационные космические снимки, предоставленные информационным порталом Геологической Службы США (снимки спутника Landsat) и инженерным-технологическим центром СканЭкс.
В основе модели STREAM_2D лежит система уравнений Сен-Венана, известные как уравнения мелкой воды. Ситуация, когда глубина акватории намного меньше ее горизонтальных размеров, достаточно распространена в случае русел и пойм равнинных рек.
В качестве начальных условий для моделирования берутся начальная поверхность дна Z (x, y, 0), соответствующие ей мгновенные поля скорости v (x, y, 0), глубины h (x, y, 0). На жидких границах задаются расходы воды и наносов и/или уровни водной поверхности. Для учета ледовых явлений необходимо дополнительно задать значение плотности льда; и подготовить файлы со значениями толщины льда и коэффициента шероховатости льда в каждой ячейке в различные моменты времени. Это позволяет учесть в расчете нестационарность покрытия льдом водной поверхности, толщины и шероховатости льда. Коэффициенты шероховатости льда необходимо подбирать в процессе калибровки модели на основе наблюдавшихся случаев заторообразования путем сопоставления рассчитанных и наблюденных уровней в период затора.
Для расчетов применялись гибридные треугольно-четырехугольные сетки нерегулярной структуры. Для русла использовалась четырехугольная, прямоугольная сетка, а для поймы – треугольная. Такие сетки хорошо адаптируются под плановые очертания расчетной области и особенности течения. Таким образом, была построена расчетная сетка, содержащая более 50 тыс. ячеек, с длинами сторон от 100 до 150 м вдоль русла и 40–50 м – поперек, на пойме длины сторон треугольных ячеек составили от 100 до 250 м.
Помимо описанной выше расчётной сетки отдельно была построена детальная сетка на город Великий Устюг с учетом всех жилых и хозяйственных построек. В России применение детальной расчетной сетки для затопления городских территорий в настоящее время не распространено, первый опыт подобного моделирования был получен в ИВП РАН под руководством В.В. Беликова для городов Ярославль и Крымск [Норин, 2016]. В качестве исходной информации о жилых и хозяйственных постройках, расположенных в черте города, были использованы данные Open Street Map, выбранные из ГИС-проекта, приобретенного в компании «Next GIS». Согласно выборке, сделанной в программе ArcGIS, в черте располагается 3276 жилых и хозяйственных построек различного назначения (рисунок 1).
Рисунок 1. Полигональный SHP-файл жилых и хозяйственных строений города Великий Устюг
Программа для построения сетки AutoMesh-2D была разработана в Шаньдунском университете (г. Цзинань, Китай). Изначально данная программа была написана для построения сетки, которая используется при формовании металлических деталей [Ma et al, 2011].
Процесс построения расчетной сетки в программе AutoMesh-2D состоит из пяти этапов:
1) Разделение области построения на отдельные участки, то есть на два класса: «дом» и «область затопления».
2) Разбиение линий на интервалы.
3) Непосредственно построение четырехугольной нерегулярной сетки по алгоритму Looping-Quard [Lee et al, 2003].
4) Ручная корректировка неправильно построенных ячеек.
5) Адаптация созданной детальной расчетной сетки для модели STREAM_2D.
Конечный, «чистовой» вариант сетки содержит 122 404 ячейки и 129 535 узлов (рисунок 2 и 3). Остальные параметры модели были взяты из предыдущей версии, таким образом, отдельная калибровка не выполнялась.
Рисунок 2. Фрагмент детальной модели города Великий Устюг в редакторе STREAM_2D
Рисунок 3. Узлы детальной сетки города Великий Устюг, построенной в программе Auto Mesh-2D и открытой в редакторе STREAM_2D.
Физико-математическая модель формирования стока ECOMAG.
Для моделирования стока рек Сухона и Юг использовался ИМК ECOMAG, включающий в себя различные базы данных. Они делятся на три 4 типа:
1. Картографическая основа. Для моделирования речной сети бассейна реки Северная Двина, который включает в себя бассейны рек Сухона и Юг и выделения элементарных водосборов была выбрана наиболее популярная модель, находящаяся в открытом доступе, SRTM (Shuttle radar topographic mission) (рисунок 4). SRTM является совместной миссией "Шаттл" (международный проект) между НАСА, DoD / NGA (Национальное агентство геопространственной разведки - ранее NIMA), DLR (Немецкий аэрокосмический центр - Deutsches Zentrum für Luft- und Raumfahrt) и ASI (Agenzia Spaziale Italiana). Разрешение полученной цифровой модели рельефа составило 30х30 м.
Рисунок 4. Цифровая модель рельефа бассейна реки Северная Двина
2. Тематические карты. В ИМК ECOMAG используется ряд тематических карт: почвенная карта (масштаба 1:2500000 разработанная коллективом авторов под руководством М.А. Глазовской и представлена в электронном виде Почвенным институтом РАН) (рисунок 5), ландшафтная карта стран бывшего СССР в масштабе 1:2500000 составлена коллективом авторов под руководством В.А. Николаева и представлена в электронном виде Почвенным институтом РАН (рисунок 6).
Рисунок 5. Карта почв бассейна реки Северная Двина
Рисунок 6. Карта ландшафтов бассейна реки Северная Двина
3. Метеорологические данные. В качестве входной информации для ИМК ECOMAG использовались среднесуточные данные о температуре воздуха, дефиците влажности воздуха и сумме осадков по 68 метеостанциям. Калибровка модели проводилась на основе данных за 1969-1984 гг., верификация – за период 1985-2019 гг. Сведения о метеорологических характеристиках были предоставлены лабораторией региональной гидрологии ИВП РАН.
4. Гидрологические данные. Для калибровки ИМК ECOMAG использовались данные среднесуточных расходов воды реки Сухона (гидрологический пост в деревне Каликино), реки Юг (гидрологические посты в деревне Гаврино и поселке городского типа Подосиновец) и на реке Малая Северная Двина (гидрологический пост в деревне Абрамково). Большинство из них ведут наблюдение за стоком с 1966 по настоящее время, за исключением поста в деревне Гаврино, который был закрыт в 1989 году. Использованные материалы гидрологических наблюдений были получены из архива кафедры гидрологии суши географического факультета МГУ имени М.В. Ломоносова и из архива лаборатории региональной гидрологии ИВП РАН
Стохастический генератор погоды MSFR_WG (Multi-Site FRagment-based stochastic Weather Generator).
Основой для получения продолжительных выборок данных о метеорологических элементах являлась стохастическая модель, «генератор погоды». Впервые термин «генератор погоды» упомянул в своей работе К.Ричардсон [Richardson, 1981], впоследствии он стал нарицательным, им стали называть класс стохастических моделей, основанных на моделировании псевдовременных рядов. Ричардсон описал модель, которая воспроизводила совместный ход максимальной и минимальной температуры, суточных сумм осадков и суммарной солнечной радиации. По мере развития вычислительных ресурсов в 1990-ые годы появились различные стохастические модели, основанные на принципах работы, заложенные К.Ричардсоном, например, модель WXGEN [Sharpley, Williams, 1990], В 1991 году в Великобритании была создана модель LARS-WG [Racsko et al., 1991], которая воспроизводила внутригодовой ход тех же величин, что и в модели К.Ричардсона: суточные осадки, суммарная солнечная радиация, суточные максимумы и минимумы температур воздуха. Количество осадков в ней задавалось на основе эмпирических распределений, для описания температуры воздуха использовалась множественная модель авторегрессии 1-ого порядка, однако в отличие от модели К.Ричардсона для задания сезонного хода параметров применялись ряды Фурье. Модель CLIGEN [Arnold, Elliot, 1996] была разработана для нужд Министерства сельского хозяйства США, позднее она была дополнена [Stockle et al., 1999]. Помимо вышеописанных четырех параметров в модели CLINGEN воспроизводится среднесуточная скорость ветра, подчиняющаяся распределению вероятностей Вейбулла, а также минимальные и максимальные значения дефицита влажности, оцениваемые на основе данных о температуре с помощью уравнения линейной регрессии. Для калибровки генератора погоды MSFR_WG, предназначенного для моделирования псевдовременных рядов среднесуточных значений температуры воздуха, дефицита влажности и суммы осадков, используемых в качестве входных данных в ИМК ECOMAG, использовалась архивная информация по 35 метеостанциям, расположенных в пределах бассейна р. Северная Двины и в его буферной зоне за период наблюдений с 1971 по 2012 год.
В монографии «Динамико-стохастическое моделирование талого стока» А.Н. Гельфан [Гельфан, 2007] отмечает, что подобные генераторы погоды хорошо воспроизводят основные климатические особенности изучаемого региона за счет использования большого числа параметров, задаваемых по фактическим данным, однако также отмечает два недостатка:
Проблема обеспечения генератора погоды исходной метеорологической информацией, которая позволяла бы получать устойчивые оценки статистических параметров.
Вторая проблема использования генератора погоды заключается в существовании локальных метеоусловий, которые не характерны для обширной площади, но тем не менее будут фиксироваться метеостанциями.
В данной работе для моделирования выборок метеорологическх характеристик на территории бассейна реки Северная Двина использовалась стохастическая модель MSFR_WG (Multi-Site FRagment-based stochastic Weather Generator), представляющая собой «генератор погоды», содержащий в себе вложенные генераторы, которые описывают последовательности метеорологических переменных с временным осреднением год, месяц и сутки [Гельфан, Морейдо, 2014]. Этот генератор погоды предназначен для моделирования методом Монте-Карло многолетних временных рядов метеорологических характеристик - «входных» величин гидрологической модели (среднесуточных значений температуры и дефицита влажности воздуха, интенсивности осадков) с учетом взаимной временной и пространственной статистической связности между ними. Пространственная и временная изменчивость метеорологических характеристик учитывается в генераторе погоды использованием модифицированного метода фрагментов Г.Г. Сванидзе [Сванидзе, 1964, 1977].
Ниже коротко изложен алгоритм разработанного метода.
Рассчитываются нормализованные поля (пространственные фрагменты) среднесуточных значений метеорологических переменных по данным наблюдений на 35 станциях за 43 года (1971-2012 годы). Для каждого года наблюдений, но нормированным относительно среднего по пространству в рассматриваемом году значению рассматриваемой метеорологической характеристики строятся пространственные фрагменты.
На следующем этапе выполняется моделирование методом Монте-Карло многолетних искусственных последовательностей средних по пространству среднегодовых значений искомых метеорологических переменных с учетом их взаимной корреляции.
Заключительный шаг алгоритма – расчет полей среднесуточных значений метеорологических переменных по соответствующим пространственным фрагментам и с использованием смоделированных среднегодовых, средних по площади значений каждой метеорологической переменной. Фрагменты выбираются случайным образом методом латинского гиперкуба из созданного ранее набора фрагментов.
Используемый в работе метод пространственных фрагментов – малопараметрический по сравнению с методами, обычно применяемыми в подобных генераторах погоды и основанными на использовании аналитических пространственных корреляционных функций, часто плохо аргументированными. Этот генератор погоды имеет всего 9 параметров, перечисленных в таблице 1. Здесь же даны выборочные оценки параметров, вычисленные по данным многолетних метеорологических наблюдений в рассматриваемом бассейне за 1971-2012 гг.
Таблица 1 - Список параметров модели MSFR_WG и их значения для водосбора р. Северной Двины.
Методика верификации двумерных гидродинамических моделей на основе сравнения площади затопления, полученной в ходе моделирования и в ходе автоматической классификации космических снимков.
Верификация – проверка модели на независимом материале (в англоязычной литературе этот этап называется также валидация). Она проводится на основе сопоставления расчетных и фактических данных об уровнях воды, уклонах и скоростях течения, не использовавшихся при калибровке.
Одним из источников проверки правильности работы модели являются данные дистанционного зондирования Земли. Сравнение площадей затопления позволяет оценить насколько точно откалибрована гидродинамическая модель.
Определение площади затопления выполнялось при помощи программы ERDASIMAGINE 2014. Методика обработки космических снимков зависит главным образом от типа снимка. Наиболее просто обрабатывать снимки, полученные в видимом диапазоне. В данном случае хорошо подходит как автоматическая классификация, так и классификация с предварительным обучением.
Алгоритм, использующийся для классификации без обучения (безэталонной классификации), базируется на кластерном анализе. Для формирования кластеров используется формула минимального спектрального расстояния. Кластеризация начинается с произвольно заданных значений (средних) или средних значения взятых из существующих сигнатур. После отнесения всех возможных пикселов к одному из классов, центры классов сдвигаются и процесс повторяется полностью сначала (следующая итерация) [Книжников и др., 2011].
Процесс продолжается до тех пор, пока не будет достигнуто максимальное количество итераций или достигнут максимальный процент пикселов, не изменивших свой класс (предел сходимости (ПС) - convergencethreshold). Например, если ПС=0,95, это значит, что процесс кластеризации закончится, как только количество пикселов, не поменявших свой класс между итерациями, достигнет 95%. Другими словами, если только 5% или меньше пикселов поменяют свой класс, процесс классификации закончится (центры классов будут установлены равным тем, которые участвовали в кластеризации на последней итерации). При работе со снимками был использовал процент сходимости равный 98. Инициирующие (начальные) центры классов могут рассчитываться несколькими путями:
1. Путем вычисления заданного количества средних значений (количество классов) по всему изображению вдоль диагонали, распределение центров классов равномерное;
2. Путем вычисления заданного количества средних значений (количество классов) по растру, рассчитанному по методу главных компонент.
Распределение центров классов по главным осям также может рассчитываться по-разному:
1. По методу стандартных отклонений. Чем большее количество стандартных отклонений задается, тем большее количество классов может попасть в "хвост" распределения (при ПС=1 большее количество классов попадет в основную массу пикселов);
2. Автоматически, в зависимости от количества классов, на основе предположения, что данные распределены нормально (Гауссово-распределение).
Данный метод был опробован на снимке участка Северной Двины, сделанный в 1998 году отечественным спутником Ресурс-Ф на пике половодья (рисунок 7).
Рисунок 7. Космический снимок моделируемого участка реки Северная Двина спутника Ресурс-Ф
С помощью описанного выше метода ISODATA была получена маска воды (рисунок 8). В качестве параметров было указано:
Количество классов — от 10 до 20.
Количество итераций — 15.
Минимальный размер кластера (в процентах) от общего числа пикселей в изображении, который удалялся в конце каждой итерации – 0,5.
Порог сходимости — 0,98.
Минимальное расстояние между классами, которое допускает существование двух кластеров как самостоятельные
Рисунок 8. Результат автоматической классификации воды методом ISODATA (на примере участка р. Северная Двина)
После проведения классификации, как и в предыдущем случае, были отобраны классы, попадающие на воду, им также было присвоено значение 1 (остальным классам 0). Далее полученные таким образом границы воды были автоматически векторизованы и подгружены в ArcGIS 10.1
Определение калибровочных параметров модели STREAM_2D (коэффициента шероховатости и величину повышения отметки дна) в зависимости от мощности ледового затора на основе моделирования затопления за многолетний период.
Как показывает практика, для отдельных гидрологических ситуаций можно с большой точностью подобрать калибровочные параметры гидродинамической модели. Однако, для воспроизведения заторных уровней воды в течение всего периода наблюдений необходима разработка дополнительных методов уточнения параметров модели.
Для этой цели на первом этапе были проанализированы все заторные уровни воды за весь период наблюдений с 1980 до 2016 года. Было использовано сценарное моделирование без ледового модуля. Принимая во внимание, что уровни открытой воды воспроизводятся моделью с высокой точностью [Agafonova et al.,2017; Zelentsov et al.,2018], максимальные наблюдаемые уровни наводнений, связанных с ледовыми заторами, сравнивались с моделируемыми стоковыми уровнями за те же даты, и заторные уровни рассчитывались как их разность.
Выяснилось, что в разные годы величина заторного подъема уровня воды колеблется в широких пределах – от 1 до 3,8 м. Из 18 проанализированных наводнений в семи случаях заторные уровни воды были выше 2,5 м, из них в трех случаях – выше 3 м (1992, 2013, 2016 годы). В семи случаях заторные уровни составили от 1,5 до 2,5 м, а в остальных – менее 1,5 м. Несколько случаев ледовых заторов с незначительными заторными уровнями, такие как 1981, 1986, 2000, 2009 годы, не были включены в анализируемый набор данных о ледовых заторах.
Описанные результаты позволяют обосновать использование различных групп параметров модели, отражающих дополнительное увеличение коэффициента шероховатости и уменьшение поперечного сечения потока в зависимости от уровня воды, вызванной ледовым затором.
Существующие классификации мощности заторов [Бузин, 2014; Донченко, 1987] были адаптированы для ключевого участка с учетом величины заторного повышения уровней воды за различные годы, что позволило выделить 4 градации величины заторных подъемов уровня воды: небольшой затор (<1,5 м), средний затор (1,5–2,5 м), мощный затор (2,5–3,5 м) и катастрофически мощный затор (>3,5 м).
Рассматриваемый временной отрезок был разделен на две равные части: период калибровки (1980–1998 гг.) и период верификации модели (1999–2016 гг.). Расчетный период задавался с 1 апреля по 30 июня, чтобы начало и конец любого половодья попадали в расчет. За это время было выявлено 18 лет, когда происходило образование ледового затора, приводившего к повышению уровня воды. Для лет, когда заторное повышение составило менее 1 м, например, 1981, 2000, 2009 гг., ледовый блок модели не подключался.
Для периода, когда русло свободно ото льда, задавался коэффициент шероховатости русла 0,025. Данное значение было получено в процессе калибровки модели на основе экспедиционных исследований 2013 года. Для периода, предшествующего началу половодья и вскрытию реки, устанавливался коэффициент шероховатости, равный 0,035. и среднемноголетнее значение толщины льда, равное 0,7 м.
Для периода калибровки параметры ледового модуля, обеспечивающие лучшее совпадение моделируемых и наблюдаемых уровней воды для поста Великий Устюг, определяли по заторным уровням (таблица 2). Заторные уровни менее 1,5 м (тип 1) лучше воспроизводятся при дополнительном коэффициенте шероховатости льда 0,045 и параметре толщины льда 1 м; для средних (тип 2) ледовых заторов параметр толщины льда следует увеличить до 1,5 м. Мощные заторы (тип 3) с уровнями воды более 2,5 м и менее 3,5 м хорошо воспроизводились при коэффициенте шероховатости льда 0,075 и параметре толщины льда 2 м. Последние наводнения наблюдались в период калибровки в 1984, 1992 и 1998 годах.
С использованием описанного набора параметров было получено удовлетворительное совпадение моделируемых и наблюдаемых суточных уровней воды для поста Великий Устюг. (рисунок 9).
Рисунок 9. График сопоставления отметки уровня воды по гидрологическому посту г. Великий Устюг за период калибровки модели (1980 - 1999 гг.): 1 — фактические данные, 2 — смоделированные данные
Сводная информация по типовым параметрам приведена в таблице 2.
Таблица 2. Итоговые параметры ледового блока модели STREAM_2D для разных градаций заторов
Для оценки качества моделирования уровней воды был использован критерий Нэша-Сатклиффа (2).
Сравнение производилось как для каждого из двух периодов, так и отдельно для дней с ледовыми заторами и дней, когда русло было свободно ото льда. Большую роль в корректной работе модели играет правильное установление времени формирования и разрушения ледового затора в параметрах модели.
Для калибровочного периода, для дней, когда наблюдался ледовый затор, было получено значение критерия Нэша-Сатклиффа, равное 0,72, что можно считать удовлетворительным результатом, однако для всего калибровочного периода оно составило 0,94, что говорит о том, что в целом модель откалибрована хорошо.
Для периода верификации (рисунок 10), когда наблюдалось заторное повышение уровня более чем на 3,5 м во время двух выдающихся наводнений 2013 и 2016 г., критерий Нэша-Сатклифа для дней, когда наблюдался ледовый затор, составил 0,79, что является также удовлетворительным результатом, а в целом за период верификации он составил 0,91.
Рисунок 10. График сопоставления отметки уровня воды по гидрологическому посту г. Великий Устюг за период верификации модели (1999 - 2016 гг.): 1 — фактические данные, 2 — смоделированные данные
Некоторое снижение значения критерия для всего периода 1999-2016 в целом связано с ухудшением качества данных о стоке для реки Юг в связи с закрытием гидрологического поста Гаврино.
Результаты проверки работы модели сведены в таблице 3.
Таблица 3. Результаты проверки качества работы модели STREAM_2D при воспроизведении уровней воды по г.п. Великий Устюг с помощью статистического критерия Нэша-Сатклиффа
Максимальные уровни воды для каждого года в течение периода наблюдений были проанализированы отдельно, поскольку они являются наиболее важной характеристикой для определения опасности наводнения. Модель воспроизводит их достаточно хорошо; значение эффективных критериев Нэша-Сатклиффа для максимальных уровней воды за период 1980–2016 гг. равно 0,92, среднее значение абсолютных отклонений моделируемых максимальных уровней воды от наблюдаемых составляет 0,18 м.
Анализ временного ряда максимальных уровней воды за 1980–2016 гг. показал, что заторные уровни определяются максимальными годовыми уровнями воды в 15 случаях из 18 рассмотренных ледовых заторов. Максимальные уровни воды в результате особо крупных наводнений (например, 1998 г.) обычно соответствуют высоким заторным уровням воды и (более 2,5 м) и стоковому уровню, превышающему средние значения. Недавнее наводнение 2016 года характеризуется как повышенным стоковым уровнем, так и катастрофическим значением дополнительного заторного подъема уровня воды.
Предложенные значения коэффициентов шероховатости для каждого из трех возможных периодов половодья (период перед вскрытием, установление ледового затора и освобождение русла ото льда) позволяют с удовлетворительной точностью воспроизводить характеристики наводнения на протяжении всего рассматриваемого периода с 1980 по 2016 год. Полученные значения параметров могут быть использованы в оперативной практике прогнозирования зоны затопления, основанной на методе математического моделирования.
Гидродинамическое моделирование катастрофических наводнений 1998 и 2016 года.
1998 год.
Наводнение 1998 года стало одним из самых крупных за период инструментальных наблюдений в городе Великий Устюг. 3 мая 1998 года около 7 утра начался резкий подъем уровня воды на реке Сухона, а в 12 часов дня произошла подвижка льда, и начался интенсивный ледоход. Спустя четыре часа в районе деревни Аристово образовался ледяной затор. К исходу 3 мая образовался еще один, но уже более мощный затор южнее деревни Бобровниково на реке Малая Северная Двина, в полночь уровень воды в реке Сухона составил 804 см, начался потоп. Тогда водой залило площадь более 12 тысяч гектар. В зоне затопления оказался Великий Устюг и 24 населенных пункта с общей численностью населения 38960 человек. Всего от водной стихии пострадало 3 686 домов и зданий и около восьми тысяч человек. Паводком было смыто шесть автодорожных и восемь трособлочных подвесных мостов, размыто и разрушено большое количество межрайонных дорог. Из-за воздействия воды вышли из строя 18 электроподстанций, 4 насосных и одна газораспределительная станции. Повреждены линии электропередачи, упали опоры ЛЭП, повреждены линии и средства связи, телефонные кабеля. Прекратили работать пять школ, шесть детских садов, один детский дом и четыре производственные объекта. Уровень воды в Сухоне поднялся почти на 10 м. Общий ущерб от наводнения тогда оценили более чем в 300 миллионов рублей.
Чтобы подробней изучить динамику затопления в период этого страшного наводнения была использована модель участка реки Северная Двина от города Великий Устюг до города Котлас.
В качестве исходных данных использовались гидрографы половодья на период, 37 дней (с 24.04 по 30.05). По реке Сухона были использованы данные о расходах воды, полученные с гидрологического поста Каликино. На реке Юг гидрологический пост в деревне Гаврино был закрыт в 1989 году, поэтому расходы реки Юг были получены с помощью модели формирования стока ECOMAG. Максимальный расход по реке Сухона составил 5560 м3/с, а по реке Юг 3183 м3/с. На нижней границе модели, ниже по течению города Котлас для каждого момента времени был задан уровень, для пика половодья он составил 46,76 метров. Для воспроизведения моделью заторных уровней 1998 г. на всем протяжении заторных участков коэффициент шероховатости повышался до 0,15, глубина уменьшалась на 3,5 м. Условие затора вводилось с 235 по 300 час (10 – 12 сутки от начала расчета).
На основе проведенного моделирования была построена карта глубин максимального затопления рассматриваемого участка (рисунок 11). На реке Сухона согласно построенной карте преобладают глубины от 9 до 12,5 метров. Максимальные отметки на реке Сухона наблюдаются перед границей начала ледового затора.
Согласно модельному расчету г. Великий Устюг подвергается затоплению с восточной стороны, это же подтверждает космический снимок. Глубина затопления в городской черте составляет около 1,5 метров, в некоторых местах эта отметка достигает 2 метров. При такой глубине происходит частичное затопление первых этажей зданий, порча мебели, бой стекол. Жители покидают первые этажи домов. Полная невозможность передвижения без плавсредств. Массовые короткие замыкания электропроводки, травмы электротоком. Затопление автотранспорта [Методика…, 2003].
На реке Юг глубина потока существенно меньше, до узла слияния на реке преобладают от 7 до 8,5 метров, а в локальных участках до 9,5 метров.
Глубина затопления правобережной поймы ниже по течению узла слияния составляет от 3 до 5 метров (рисунок 11).
Рисунок 11. Карта глубины затопления моделируемого участка реки Северная Двина во время выдающегося половодья 1998 года
В районе поселка Приводино наблюдается подтопление окраины поселка, глубина затопления составляет 0,7-1 метр. Согласно таблице, составленной МЧС России для оценки тяжести разрушения, нанесенного негативным воздействием вод [Методика…, 2003], при такой глубине происходит подтопление цоколей и первых этажей зданий, затопление сельхозугодий. Город Котлас остается нетронутым затоплением, в то время как происходит обширный разлив на левый берег реки. На карте глубин отчетливо видны затопленные протоки, глубина воды в которых на момент пика половодья достигает 6,5 метров.
Помимо карты глубины затопления была построена карта распределения скоростей течения (рисунок 12).
Рисунок 12 - Карта распределения скоростей течения на моделируемом участка реки Северная Двина во время выдающегося половодья 1998 года
На реке Сухона скорости равны от 0,1 м/с до 2,5 м/с. Максимальная скорость течения наблюдается в наиболее узком участке русла на верхней границе модели. По мере удаления от стрежня реки скорость течения снижается за счет уменьшения глубины и возникновения трения о берег. По мере приближения к верхней границе затора скорость течения падает до 1,6-1,8 м/с.
На затопленной территории города Великий Устюг скорость течения составляет 0,1 – 0,5 м/с. За счет резкого увеличения шероховатости подстилающей поверхности и появления препятствий в виде жилых и хозяйственных построек в виде домов скорость течения резко падает. На территории затопленной поймы скорость течения также, как и в городе составляет 0,1-0,5 м/с.
На участке между двумя ледовыми заторами скорость течения достигает 2,5 м/с. На месте, где находится основное русло скорость течения меняется от 1,6 м/с до 2,5 м/с, в то время как на месте протоки скорость существенно ниже, там она не превышает 1,6 м/с.
Ниже по течению деревни Приводино, поток, преодолев ледовый затор, вновь начинает набирать скорость, достигая своего максимального значения в наиболее узкой части русла. По окончанию сужения, поток разливается на левобережную пойму, при этом скорость течения заметно падает. На затопленной пойменной территории скорость течения составляет 0,1-0,5 м/с, на месте затопленных во время половодья рукавов, скорость течения возрастает до 1 м/с.
Кроме анализа глубин и распределения скоростей течения был произведен анализ продолжительности затопления при помощи созданной ранее подпрограммы «Flood time».
При проведении расчётов с использованием подпрограммы "Flood time" основная модель STREAM_2D настраивается таким образом, чтобы выдача результатов производилась с необходимым нам шагом, в данном случае шаг составлял 1 час. То есть за каждый расчётный час пользователь получает txt-файл, в котором прописаны отметки водной поверхности. Эти отметки в подпрограмме «Flood time» сравниваются с высотными отметками незатопленного рельефа. Если их разность будет больше нуля, это означает, что земля находится под водой, то есть затоплена. Затем ведется общий подсчет, сколько раз за время расчёта ячейка была затоплена, это и будет равняться продолжительности затопления. В итоге формируется txt-файл, в котором указывается количество часов затопления для каждой ячейки.
Карта продолжительности затопления отображает продолжительность стояния воды в конкретной ячейке модели (рисунок 13). Так как рассматривался весьма продолжительный период, более месяца, то для построения карты был выбран временной интервал 24 часа, то есть 1 сутки.
Рисунок 13. Карта продолжительности затопления на моделируемом участке реки Северная Двина во время выдающегося половодья 1998 года
Наибольшая продолжительность наблюдается в русле реки, по мере удаления от русла реки продолжительность стояния воды снижается.
В черте города Великий Устюг вода во время наводнения 1998 года стояла от 2 суток на южной окраине города, имеющей достаточно высокие абсолютные высотные отметки до 25 дней. Учитывая тип застройки города, а именно преобладание малоэтажных кирпичных и деревянных домов, то согласно классификации МЧС, продолжительность затопления более 10 суток приводит к сильному разрушению дома.
Продолжительность затопления левобережной поймы реки Юг составляет от 24 до 28 дней в зависимости от рельефа.
Продолжительность затопления деревни Красавино составляет около 3 часов, что не столь критично, за это время возможна порча сельхозугодий.
Продолжительность затопления деревни Приводино составляет около 2 дней, за это время согласно классификации МЧС, учитывая глубину затопления, происходит сильное разрушение деревянных жилых и хозяйственных построек.
На нижней границе модели, на левом берегу, напротив города Котлас, вода стоит около 8 дней. Так характерный тип застройки для данного места — это малоэтажные деревенские дома, то такая продолжительность стояния воды достаточна для нанесения существенного ущерба и приведение в негодность жилых построек.
2016 год.
Наводнение 2016 года также было одним из самых масштабных за период наблюдений в городе Великий Устюг. Причиной произошедшего наводнения 2016 года стал комплекс сложившихся метеоусловий в осенне-зимний период 2015–2016 гг (рисунок 14). С ноября по декабрь температура воздуха на территории бассейна р. Северная Двина превышала среднемноголетнее значение. Это привело к тому, что в начале декабря вследствие теплых дождей река вскрылась и начался ледоход, который продолжался более двух недель. Все это привело к образованию протяженных наторошенных участков и наличию большого количества шуги подо льдом на этих реках. Одной из особенностей зимнего ледохода в нижнем течении р. Сухона осенью 2015 г. было неоднократное образование заторов льда. 23 декабря 2015 г. у г. Великий Устюг был зафиксирован самый высокий уровень воды для декабря за весь период наблюдений с 1881 по 2015 гг. — 674 см над нулем поста. Ледостав установился при высоких уровнях воды, которые сохранялись и на начало февраля и были сопоставимы со средними уровнями воды в период весеннего ледохода. В первой декаде апреля (6–7 числа) уровни воды р. Сухона в среднем и нижнем течении начали расти. К 9–10 апреля лед потемнел, наблюдались закраины, и 13 апреля начался весенний ледоход у г/п Березовая Слободка. 14 апреля лед остановился у льдоулавливающей плотины в заторе, начался резкий рост уровней воды.
Рисунок 14. Изменение температуры воздуха (t) и осадков (x) по м/с Великий Устюг и уровня воды (H) р. Сухона в осенне-зимний период 2015–2016 гг.
Задержка льда дамбой была кратковременной и неполной. Начало ледохода в нижнем течении проходило на высоких уровнях воды и повсеместно сопровождалось навалами льда. Днем 14 апреля начался стремительный рост уровней воды у г/п Каликино. К моменту остановки льда в заторе ниже п. Подсосенье уровень достиг отметки 1099 см над нулем поста, превысив исторический максимум 1953 г. За 6 часов уровни воды выросли на 4 м. Резкий рост уровней воды обеспечил прорыв затора к 11 утра 15 апреля. Волна прорыва затора привела к подвижкам и мощному ледоходу в районе г. Великий Устюг. Высокие уровни воды привели к выходу льда на набережную города, гидрологический пост сразу оказался завален льдом и был расчищен уже после образования затора (рисунок 15). Ледоход в городе был непродолжителен, к 19 часам 15 апреля лед уперся в осенний затор на р. Малая Северная Двина и остановился в 40 км ниже города. Вода пошла в обход заторных скоплений, затопив правый берег, в том числе д. Дымково. Этот сценарий во многом обеспечил устойчивость образовавшегося затора, который простоял 90 ч. Уровни за этот период практически не менялись (рисунок 16).
Рисунок 15. Затор льда на р. Сухона у г. Великий Устюг в апреле 2016 г. (фото Н.Л. Фроловой)
Рисунок 16. Ход уровней воды и ледовой обстановки на р. Сухона, г/п Каликино и г/п Великий Устюг
С 14 апреля в Великоустюжском районе Вологодской области был введен режим ЧС. Проводились попытки разрушения затора взрывными работами, а также бомбометанием силами ВКС, с помощью СУ-34. После разрушения затора уровни в районе города упали на 0,5 м за 1 час. После 20 апреля в г. Великий Устюг начались аварийно-восстановительные работы, воду с затопленных территорий удалось откачать только к 20 м числам мая.
Для моделирования наводнения 2016 г. была использована как обычная модель (от слияния рек Сухона и Юг до города Котлас), так и модель с детальной сеткой на территорию города Великий Устюг.
Согласно восстановленным с помощью ECOMAG данным о стоке, максимальный суммарный расход рек Сухона и Юг в 2016 г. за период половодья составил более 9000 м3/с (27 апреля). Данные о расходе были использованы в качестве начальных условий в модели STREAM_2D. Совместное использование модели ECOMAG и STREAM_2D позволили определить ход уровней воды в пределах исследуемого участка, оценить максимальную глубину и площадь затопления. Расчеты проводились за период с 1 апреля по 30 июня. Качество работы модели оценивалось путем сопоставления графиков хода уровня воды (рисунок 17).
Рисунок 17. Ход фактического уровня и уровня, полученного по данным моделирования ECOMAG — STREAM_2D (на примере половодья 2016 года, гидрологический пост г. Великий Устюг)
На рисунке 17 видно, что было достигнуто хорошее соответствие даже на период стояния ледового затора. Несмотря на то, что по посту Медведки смоделированный уровень превышает фактический, данная ситуация допустима, так как при составлении прогноза уровня воды для планирования спасательных операций лучше завысить отметку воды, нежели занизить. Помимо сопоставления смоделированного уровня с фактическим по гидрологическому посту города Великий Устюг, сравнение уровня было проведено по гидрологическим постам деревни Медведки и города Котлас (рисунок 18).
Рисунок 18. Ход фактических уровней и уровней, полученных в результате моделирования ECOMAG — STREAM_2D по контрольным точкам (половодье 2016 г.)
По результатам расчетов по модели STREAM_2D были построены карты на момент пика половодья, а именно карта глубин и карта скоростей течения (рисунок 19).
На карте глубин видно, что затоплению подверглась юго-восточная часть города, в отдельных местах глубина потока достигала 7–7,5 м. Анализ скоростей течения показал, что на затопленной территории скорость течения не превышает 0,4 м/с, а в некоторых районах и вовсе близка к нулю. Наибольшие скорости на исследуемой территории наблюдаются в районе д. Красавино, в месте сужения русла, то есть в месте сжатия потока. Площадь затопления соответствует имеющимся данным, полученным с помощью космических снимков.
Рисунок 19. Карта глубины затопления (м) (слева) и скорости течения (м/с) (справа) на момент пика половодья 2016 г.
Для детального изучения динамики потока в черте города Великий Устюг был проведен расчет пика половодья на детальной расчетной сетке.
В качестве временного отрезка для моделирования наводнения было выбрано три дня, когда наблюдались наибольшие уровни воды, а именно 16, 17 и 18 апреля.
На момент пика половодья, то есть 18 апреля, по результатам проведенного моделирования была построена карта глубин (рисунок 20) и скоростей течения (рисунок 21).
Рисунок 20. Карта глубин в районе города Великий Устюг на момент пика половодья 2016 года (18 апреля).
На карте глубин видно, что затоплению подвергся пригород Великого Устюга, а именно деревня Юдино, расположенная на северо-восточной окраине города, причем в центральной части деревни, например, на улицах Парковая и Тихая, глубина достигала 6 м, то есть были полностью затоплены двухэтажные здания. В самом городе Великий Устюг затоплению подверглась относительно низкая центральная часть, от ул. Красной, которая вблизи набережной регулярно подтапливается, до ул. Цветочной. В наиболее низко расположенной части города глубина затопления, согласно данным моделирования, составляла 1,5–2 м. Однако на большей части затопленной территории глубина затопления составила 0,5 м. Следует отметить, что во время наводнения 2016 года юго-восточная часть города оказалась отрезанной от остальной части города.
Из-за установления ледового затора затоплению подверглись правобережные населенные пункты, причем в отдельных местах ширина фронта подтопления составляла около 1 км. В зону подтопления попали деревни Коромыслово, Барсуково, Ивашево, Дымково и Добрынино. Наибольшая глубина затопления отмечена в районе деревень Барсуково и Ивашево, где она достигала 5 м. Глубина в районе деревни Коромыслово, находящейся в месте слияния рек Сухона и Юг, составила около 0,5 м.
Для сравнения границ затопления, полученных с помощью модели STREAM_2D, с фактическими данными был использован космический снимок Landsat-8, сделанный на момент пика половодья. С его помощью была получена маска затопленной территории с помощью автоматической классификации по методу K-means. Согласно данным космической съемки, площадь затопления составила 2 км2, а по данным моделирования —1,85 км2. Таким образом, относительная ошибка моделирования составила 7,5%. Данный результат можно оценить, как хороший.
По результатам моделирования была построена карта скоростей течения (рисунок 21).
Рисунок 21. Карта скоростей течения в районе города Великий Устюг на момент пика половодья 2016 года (18 апреля).
Максимальная скорость течения на стрежне реки достигает 1 м/с, однако по мере направления к берегу она угасает и в урбанизированных районах уменьшается до нуля. Отдельного внимания заслуживает центральная часть города, где ложбина выступает в роли русла. В этом месте наблюдается движение воды, причем скорость течения достигает 0,4 м/с.
Чтобы понять, как движется поток воды, и оценить, с какой стороны вода будет проникать в жилые районы, была построена карта скоростей течения в совокупности с векторами направления течения. Фрагмент такой карты представлен на рисунке 22.
Рисунок 22. Карта скоростей течения в районе города Великий Устюг на момент пика половодья 2016 года (18 апреля) с наложением векторов направления течения.
На карте видно, что вода в центральную часть города проникает со стороны реки Сухоны, причем по мере продвижения в город вода постепенно распространяется на северо-восток и на юго-запад, вглубь жилого квартала.
Расчет наиболее неблагоприятных сценариев прохождения половодья на основе динамико-стохастического моделирования совокупности метеоусловий, модели формирования стока ECOMAG и двумерной гидродинамической модели STREAM_2D.
В данной работе предпринята оценка максимально возможных расходов воды для реки Северная Двина в районе слияния у города Великий Устюг для реализации сценария наихудшего прохождения половодья.
Применение комплекса моделей для изучения наводнений стокового генезиса обусловлено двумя факторами:
Наличие реальных случаев прохождение больших расходов воды по реке Северная Двина в районе города Великий Устюг без непосредственного образования ледовых заторов. К таким примерам можно отнести половодье 1984 года, при котором суммарный расход воды рек Сухона и Юг достиг 4090 м3/с.
При прохождении расходов воды 1% обеспеченности и менее ледовые заторы не образуются, так как вода выходит на пойму, площадь поперечного сечения потока увеличивается и как следствие увеличивается пропускная способность русла.
Проблема экстраполяции эмпирической кривой в область, находящейся за областью, освещенными фактическими данными возрастает по мере уменьшения процента обеспеченности. А.В. Христофоров в 1993 году теоретически обосновал, что даже при идеальном соответствии аналитической и эмпирической кривой, погрешности в области редких повторяемостей могут достигать порядка нескольких десятков процентов [Христофоров, 1993].
Проблема расчета расходов редкой обеспеченности до сих пор полностью не решена в прикладной гидрологии, однако предложенный подход делает упор на статистическое обоснование (формирование временных метеорологических рядов большой продолжительности), а не на перестраховку в намеренном завышении параметров, например, установление Cs/Cv равное 3-4.
С помощью модели MSFR_WG были смоделированы соответствующие современному климату (так как использовались параметры, определенные за 1971-2012 гг.) искусственные ряды среднесуточных значений метеорологических характеристик: интенсивности осадков, температуры и влажности воздуха длиной 1000 лет для каждой станции на рассматриваемой территории. Для проверки качества разработанного генератора погоды сравнивались статистические характеристики, определенные по фактическим и искусственным рядам метеорологических переменных: средние многолетние величины за год, месяц и сутки; стандартные отклонения среднегодовых и среднемесячных величин; коэффициенты корреляции между среднегодовыми (среднемесячными) значениями на отдельной станции и соответствующими значениями на других станциях; пространственные корреляционные функции значений осадков, температуры и влажности воздуха, осредненных за разные временные интервалы (год, месяц, сутки) (рисунок 23).
Рисунок 23. Норма и среднеквадратичное отклонение среднегодовых величин суммы осадков, температуры и влажности воздуха, оцененные по фактическим и смоделированным данным
Расчеты показали, что статистики искусственных рядов искомых метеорологических переменных оказались в пределах стандартного отклонения оценок соответствующих характеристик по рядам наблюдений за этими метеорологическими переменными на сети станций. Этот вывод оказался справедлив для разных временных осреднений искомых переменных – от суток до года.
Удовлетворительные результаты проверки стохастического генератора погоды MSFR_WG позволяют использовать смоделированные временные ряды метеорологических характеристик для задания пространственно распределенных случайных входов в описанную выше детерминистическую гидрологическую модель.
Разработанная динамико-стохастическая модель формирования речного стока р. Северной Двины проверялась по данным наблюдений в замыкающих створах рек Сухона и Юг, в узле слияния которых находится город Великий Устюг, подверженный частым затоплениям в период весеннего половодья. Повышение уровня воды в реке во время наводнений обусловлено, наряду с притоком талых вод, частыми заторами льда. В работе [Agafonova et al., 2017] показано, что вероятность катастрофического подъема воды в период весеннего снеготаяния возрастает на рассматриваемом участке реки вследствие формирования ледового затора на 30-35%. Численные эксперименты были организованы следующим образом:
С помощью генератора погоды были смоделированы искусственные последовательности среднесуточных метеорологических характеристик длиной 1000 лет, которые являются входной информацией в модель ECOMAG. Затем по гидрологической модели рассчитан 1000-летний временной ряд ежедневных расходов воды в замыкающих створах рек Сухона и Юг, по которому получены 1000 значений максимальных за год расходов воды. По полученным выборкам оценены статистические характеристики максимальных расходов и построены кривые вероятности превышения, которые сравнивались с кривыми, построенными по наблюденным рядам максимальных расходов в рассматриваемых створах рек.
Из смоделированных искусственных гидрографов стока в замыкающих створах рек Сухона и Юг выбирались такие сочетания, которые формировали максимальный, со средней повторяемостью 1 раз в 1000 лет, гидрограф стока в узле слияния указанных рек.
С помощью двумерной гидродинамической модели STREAM_2D [Belikov et al., 2015], в качестве граничного условия в которой задавался рассчитанный по динамико-стохастической модели гидрограф максимального стока, рассчитывался ход уровня воды и динамика затопления поймы в створе г. Великий Устюг. Ниже изложены результаты численных экспериментов
На рисунке 24 сравниваются кривые распределения максимальных расходов воды в замыкающих створах рек Юг и Сухона, построенные по выборкам фактических величин максимального стока и величин, рассчитанных с помощью динамико-стохастической модели. Использовались фактические данные о максимальных расходах воды р. Юг за 77 лет (с 1931 по 2007 годы) и по р. Сухоне за 127 лет (с 1881 по 2007 годы). В таблице 4 приведены статистические характеристики максимального стока, оцененные по фактическим и рассчитанным рядам.
Рисунок 24. Вероятности превышения максимальных за год расходов воды, оцененные по ряду фактических значений (желтые точки), по 1000-летнему ряду, рассчитанному с помощью динамико-стохастической модели (красная линия). Черным пунктиром показана кривая трехпараметрического гамма-распределения, построенная по фактическим данным.
Таблица 4. Статистические характеристики максимальных расходов стока рек Юг и Сухона, оцененные по фактическим и рассчитанным с помощью динамико-стохастической модели рядам стока
Расчеты показали, что динамико-стохастическая модель позволила с удовлетворительной точностью воспроизвести средние значения стока на обеих реках, а также стандартное отклонение максимального стока р. Сухона: для перечисленных статистик отклонения фактических и рассчитанных значений не выходят за пределы 95%-го интервала выборочной погрешности оценок по фактическому ряду стока. Однако стандартные отклонения рассчитанных величин максимального стока р. Юг статистически значимо завышают стандартные отклонения соответствующих фактических величин (см. таблицу 4). Анализ результатов показал, что для р. Юг оказались занижены максимальные расходы высокой обеспеченности (т.е. максимальные расходы низких половодий; см. также верхний рисунок 24), в то время как максимальные расходы малой обеспеченности (наиболее высоких половодий) получились близкими к наблюденным, что позволяет использовать результаты расчетов в области малых обеспеченностей максимальных расходов для гидродинамического моделирования максимального затопления в створе г. Великий Устюг (в узле слияния рек Сухона и Юг). Анализ имеющихся данных наблюдений показал, что наиболее критический сценарий прохождения половодья – это одновременное прохождение пиков по реке Сухона и Юг. Такой сценарий наблюдался, например, во время формирования последнего выдающегося наводнения 2016 года (рисунок 25).
Рисунок 25. Гидрограф реки Сухона — пост д. Каликино (1) и восстановленный по данным ECOMAG гидрограф реки Юг — пост Гаврино (2) во время половодья 2016 года
Аналогичный сценарий был выбран по результатам динамико-стохастического моделирования из 1000 лет: для реки Юг максимальное значение составило 5994 м3/с, а для реки Сухона — 7860 м3/с. Соответствующие гидрографы стока задавались в качестве граничных условий в гидродинамической модели, с помощью которой рассчитывались поля глубин затопления (рисунок 26) и скоростей потока. Показано, что при рассмотренном сценарии прохождения максимальных расходов воды на реках Сухона и Юг глубина затопления в черте города будет составлять 0.5 м. среднее значение скорости течения в пределах моделируемой области составит 0.8 м/с, средняя глубина затопления составит 4.5 м, а площадь затопления будет 300 км2. Максимальная отметка поверхности уровня воды по посту Великий Устюг будет составлять 57,36 метров; для сравнения уровень воды при прохождении 1% расхода воды при наличии ледового затора на 2 метра больше и составляет 59,39 м.
Рисунок 26. Карта глубины затопления г. Великий Устюг, рассчитанных по граничным условиям, задаваемым по результатам динамико-стохастического моделирования максимального стока в узле слияния рек Сухона и Юг
В качестве наиболее экстремального сценария было рассмотрено прохождение 0,1% расхода воды при условии образования ледового затора воды в узле слияния рек Сухона и Юг. В действительности такой вариант прохождения половодья маловероятен, так как, анализируя ряд расходов воды с 1980 года, можно отметить, что установление ледового затора происходит на подъеме половодья, то есть до момента достижения его пика. Учитывая этот факт, вероятность такого наводнения была определена как 0,001%.
Для того, чтобы выявить влияние ледового затора на характеристики наводнения в черте города Великий Устюг были построены графики изменения площади затопления (км2), средней глубины затопления (м), средней скорости течения (м/с) в зависимости от расхода воды. Для этого, помимо вышеописанного расчета, с помощью модели STREAM_2D были рассчитаны значения этих характеристик при прохождении наводнения стокового генезиса, сформированным максимальным расходом воды 5% и 1% обеспеченности и при прохождении наводнения заторного генезиса, сформированным максимальным расходом воды 5% обеспеченности. Стоит отметить, что для вычисления расхода воды различной обеспеченности использовались различные кривые зависимости расхода воды от уровня, а именно для стокового и заторного генезиса. Результаты моделирования сведены в таблице 5.
Таблица 5. Значения характеристик наводнений заторного и стокового генезиса при различной обеспеченности
Максимальная отметка уровня воды при установление ледового затора и прохождение расхода воды 0,1 % согласно проведенному расчету составила 61,07 м. Стоит отметить, западная и центральная, то есть историческая часть города расположена на высоте 60-62 м. Таким образом, оставаясь не затопленной во время прохождения самых крупных наводнений за инструментальный период наблюдения, она будет затоплена при образовании ледового затора и прохождение 0,1% расхода воды.
Сравнивая отметки уровня воды стокового и заторного генезиса, образованные расходом воды одной обеспеченности можно отметить, что заторные уровни превышают стоковые на 2 и более метра (рисунок 27).
Рисунок 27. Смоделированные на основе двумерной гидродинамической модели STREAM_2D уровни воды стокового (1) и заторного (2) генезиса у г. Великий Устюг
Наличие ледового затора приводит к увеличению не только уровня воды, но и других перечисленных характеристик. Так, например, площадь и глубина затопления превышают аналогичные по повторяемости стокового генезиса в полтора раза и более, а скорость течения возрастает в 1,2 раза.
Программная оболочка на языке Python, позволяющая на основе задания климатических параметров, получать характеристики затопления.
Для хранения данных средствами СУБД PostgreSQL была спроектирована и реализована схема хранения, учитывающая потребность в отслеживании результатов отдельных циклов расчетов и хранении промежуточных данных с необходимыми пространственными привязками (рисунок 28).
В результате получены следующие таблицы:
launch – таблица учета запуска циклов расчетов (таблица 6);
eco_input – таблица идентификаторов исходных данных модели ECOMAG;
meteo_station – таблица, содержащая метеостанции, на которых выполнялись наблюдения (таблица 7);
value_type – таблица возможных типов наблюдаемых значений (таблица 8),
meteo – таблица, содержащая результаты наблюдений на метеостанциях (таблица 9);
eco_meteo – развязочная таблица для отношения «многие ко многим» для таблиц eco_input и meteo;
hydro_station – таблица, содержащая гидропосты, для которых рассчитываются расходы (таблица 10);
hydrograph – таблица для хранения гидрографов, расчитанных моделью ECOMAG (таблица 11);
eco_out – таблица идентификаторов рассчитанных значений модели ECOMAG;
eco_hydro – развязочная таблица для отношения «многие ко многим» для таблиц hydrograph и eco_out;
eco_out_stat – таблица идентификаторов статистических значений результатов расчетов модели ECOMAG;
eco_stat_value – таблица статистических значений результатов расчетов модели ECOMAG для отдельных гидропостов (таблица 12);
stream_input – таблица идентификаторов исходных данных модели STREAM_2D;
stream_hydro – развязочная таблица для отношения «многие ко многим» для таблиц hydrograph и stream_input;
stream_out – таблица идентификаторов рассчитанных значений модели
stream_point – таблица точек, для которых рассчитываются значения в модели STREAM_2D (таблица 13);
stream_value – таблица значений, рассчитанных моделью STREAM_2D (таблица 14);
stream_out_value - развязочная таблица для отношения «многие ко многим» для таблиц stream_out и stream_value;
stream_out_stat - таблица идентификаторов статистических значений результатов расчетов модели STREAM_2D;
stream_stat_point – таблица статистических значений результатов расчетов модели STREAM_2D.
Рисунок 28. Схема хранения, учитывающая потребность в отслеживании результатов отдельных циклов расчетов и хранении промежуточных данных с необходимыми пространственными привязками
Таблица 6. Столбцы таблицы launch
Таблица 7. Столбцы таблицы meteo_station
Таблица 8. Столбцы таблицы value_type
Таблица 9. Столбцы таблицы meteo
Таблица 10. Столбцы таблицы hydro_station
Таблица 11. Столбцы таблицы hydrograph
Таблица 12. Столбцы таблицы eco_stat_value
Таблица 13. Столбцы таблицы stream_point
Таблица 14. Столбцы таблицы stream_value
Таблица 15. Столбцы таблицы stream_stat_point
Реализованный на языке Python модуль объединения цикла расчетов в автоматическом режиме выполняет следующую последовательность обработки данных (рисунок 29):
Преобразование исходных данных и сохранение в БД.
Подготовка исходных данных модели ECOMAG.
Запуск расчета модели ECOMAG.
Преобразование результатов расчетов, выделение статистических характеристик, сохранение в БД.
Подготовка исходных данных модели STREAM_2D на основе расчетов ECOMAG.
Запуск расчета модели STREAM_2D.
Преобразование результатов расчетов, выделение статистических характеристик, сохранение в БД.
С учетом современных тенденций, при проектировании данного программного модуля была применена архитектура, позволяющая запускать отдельные его компоненты в различных средах, наиболее подходящих, например, для конкретной конфигурации используемой модели.
Рисунок 29. Последовательность обработки данных
С помощью языка программирования Python (с использованием модуля Dash) было разработано аналитическое приложение, визуализирующее результаты работы расчетного модуля.
Данная программа имеет очень простой интерфейс, что позволяет пользователю легко разобраться со всеми ее функциями.
Основной функционал программы:
В панели “Настройка” задаются параметры моделей ECOMAG и STREAM_2D.
Выбирается исследуемый период наблюдений
Происходит работа расчетного модуля
Отображение результатов расчетов и основных характеристик, сохраненных в БД, пользователю:
Максимальный смоделированный уровень
Максимальный расход
Максимальное значение скорости
График смоделированного и фактического уровней для анализа работы модели
Карты смоделированных скоростей и глубин затопления, полученных в результате работы STREAM_2D
Итоговый вид программной оболочки представлен на рисунке 30.
Рисунок 30. Интерфейс программного расчетного комплекса.
Методика ансамблевого краткосрочного прогноза характеристик осеннее-летних паводков и долгосрочного прогноза весеннего половодья.
Ансамблевое прогнозирование в работе основывается на воспроизведении с помощью модели формирования стока (ECOMAG) возможных условий стока на водосборе р. Северная Двина – г. Великий Устюг путем подстановки за период заблаговременности прогноза реализованных в предыдущие годы метеоусловий. Эти условия накладываются на сформированное в конкретный год состояние водосбора к дате выпуска прогноза, учтенное с помощью внесения в модель фактических метеоданных (рисунок 31). В качестве параметров, задающих метеорологические условия за период заблаговременности, были использованы среднесуточная температура воздуха, суточная сумма осадков и суточный дефицит влажности воздуха за 1966-2014 гг. (49 сценариев). Выходными данными модели является ансамбль гидрографов данного года.
Рисунок 31. Схема ансамблевого моделирования гидрографов речного стока
На основе анализа даты начала и конца половодья по ансамблевым сценариям за различные годы дата выпуска прогноза была назначена на 30 марта, а период заблаговременности выбран в течение 30 марта – 30 июня (31 сутки).
Ансамблевый прогноз характеристик весеннего половодья в работе выполнен для 1993 г. и 2012 г (рисунок 32). Такой выбор связан с наблюдением высоких значений максимальных расходов воды, отсутствием заторов льда на реках в рассматриваемые годы, их соответствием современным условиям формирования стока, и низким значением относительной ошибки моделирования максимального расхода (таблица 16). Это позволяет сконцентрировать внимание исключительно на весеннем стоке. В дальнейшем данные могут быть привлечены к прогнозу затопления территории вблизи г. Великий Устюг.
Рисунок 32. Ансамбль гидрографов половодья в створе р. Сухона – д. Каликино:
а) 1993 г.; б) 2012 г.
Таблица 16. Характеристики исследуемых половодий по данным створа р. Сухона – д. Каликино
Использование ансамбля равновероятных гидрографов даёт возможность оценки разнообразных характеристик половодья, среди которых в работе для прогноза были выбраны: максимальный расход, объем, продолжительность, дата начала, продолжительность подъема, дата наступления максимального расхода. Прогнозирование каждой из них имеет важное значение в решении вопросов управления водными ресурсами и защиты населения и объектов хозяйственного использования от наводнений.
Особенностью ансамблевого моделирования является возможность получения прогноза в вероятностной и детерминистической формах. Первая может быть получена путем оценки одной характеристики половодья по разным ансамблевым сценариям, вторая – с помощью осреднения характеристики по ансамблю или же моделирования на основе климатических норм метеорологических величин (Кучмент, Гельфан, 2009). Однако ансамбль гидрографов также позволяет получить средний из них, т.е. осредненный вариант развития половодья. В случае простого вычисления среднего арифметического значения расхода воды за каждую дату половодья могут наблюдаться большие погрешности в виду того, что при различных метеорологических условиях (по различным сценариям) на одну и ту же дату могут приходиться подъем, максимум или спад половодья. Из-за этого средний гидрограф, полученный таким способом, может кардинально не соответствовать рассматриваемой ситуации. В работе использован подход к осреднению ансамблевых гидрографов, позволяющих проводить более корректное осреднение расхода воды на каждую дату.
На основе выделенного по каждому ансамблевому сценарию периода половодья выполнялся пересчет абсолютных временных координат гидрографов половодья в относительные (от 0 до 1) таким образом, чтобы максимальный расход приходился на абсциссу 0,5 (рисунок 33).
Рисунок 33. Ансамбль гидрографов половодья 2012 г. приведенный к относительным координатам по времени
Далее предстоит вычислить среднее арифметическое значение расхода воды при различных Kt. Это представляет сложность, т.к. непосредственно снять с графика значения по всему ансамблю в каждой абсциссе невозможно. Был рассмотрен вариант осреднения максимального и минимального значения расхода воды в данной абсциссе (рисунок 34). Значение расходов воды снималось с графика. Шаг по оси Kt был выбран равным 0,05. Полученные таким образом средние максимальные значения расхода воды отличались от фактических на 579 и 178 м3/с в 1993 и 2012 гг. соответственно.
Рисунок 34. Средний гидрограф половодья 2012 г. (красным), полученный путем осреднения ансамбля
Для перехода от относительных координат времени к абсолютным необходимо задать продолжительность подъема и спада половодья. В этом заключается одно из преимуществ метода, т.к. можно оперировать продолжительностью интересующей обеспеченности. В данном случае для расчетов были использованы значения обеспеченностью 50%.
Получение вероятностного прогноза на основе ансамблевого моделирования заключается в оценке функции распределения вероятности исследуемой характеристики, полученной по разным сценариям. Форма представления результатов прогноза зависит от нужд потребителя, однако наиболее распространенной является определение доверительных интервалов в 90 и 85%. Наглядные иллюстрации могут быть приведены на графиках интегральных функций распределения вероятности рассматриваемых характеристик. В работе для прогнозируемых величин также был оценен 75% и 98% доверительный интервал. Первый определяет диапазон значений характеристик половодья соответствующий (в вероятностном отношении) прогнозному интервалу детерминистической формы прогноза, ограниченному допустимой погрешностью. Однако при ансамблевом подходе 75%-ый доверительный интервал обусловлен различием возможных сценариев стока, а не погрешностью методики, что наблюдается при прогнозе по математическому ожиданию. Второй интервал показывает общий диапазон изменения исследуемых величин, его оправдываемость связана с разнообразием реализуемых сценариев метеоусловий периода заблаговременности прогноза.
Ниже представлены результаты прогноза характеристик весеннего половодья в 1993 и 2012 гг. (таблицы 17 и 18).
Таблица 17. Результаты ансамблевого прогноза характеристик половодья р. Северная Двина – г. Великий Устюг в 1993 г.
В 1993 г. обеспеченность фактического значения максимального расхода воды по ансамблю гидрографов составила 15%. Оно входит в 90% доверительный интервал (3150-4790 м3/с) (рисунок 35). Значения объема, продолжительности половодья, даты наступления максимального расхода воды в 1993 г. отвечают доверительному интервалу 75%.
Значение характеристик половодья, осредненное по ансамблю, показывает относительную ошибку определения максимального расхода воды, равную 21%. В соответствии практике гидрологических расчетов подобную ошибку можно считать допустимой [Евстигнеев и др., 2016]. Наибольшая точность наблюдается при определении объема, даты начала половодья и даты наступления максимального расхода воды. Эту же закономерность можно проследить в 2012 г.
Рисунок 35. Интегральная функция распределения долгосрочных прогнозов максимального расхода воды по ансамблю фактических метеоусловий
для половодья 1993 г.
В 2012 г. наблюдалось крайне высокое значение максимального расхода воды половодья. Фактическое значение составило 5490 м3/с, что соответствует обеспеченности 5% (по ансамблю). Таким образом, эта величина не попадает в доверительный интервал 90% (2800-5100 м3/с), однако описывается 98% доверительным интервалом (рисунок 36).
Таблица 18. Результаты ансамблевого прогноза характеристик половодья р. Северная Двина – г. Великий Устюг в 2012 г.
Объем половодья и дата наступления максимального расхода воды в 2012 г. отвечают доверительному интервалу 75%. Среднее значение максимального расхода по ансамблю гидрографов имеет относительную погрешность 30%, что связано с выдающимся фактическим значением. Большая погрешность также характерна для продолжительности подъема половодья.
Рисунок 36. Интегральная функция распределения долгосрочных прогнозов максимального расхода воды по ансамблю фактических метеоусловий
для половодья 2012 г.
На основе анализа аномалии температуры воздуха исследуемого бассейна за 1-25 марта и 1-10 апреля была предпринята попытка уточнения вероятностного прогноза путем исключения некоторых сценариев (около 15) метеоусловий периода заблаговременности. При этом не учитывались годы, имеющие наибольшие различия в аномалии температуры воздуха за рассматриваемый период с исследуемым годом.
После сокращения количества сценариев в 1993 и 2012 гг. для максимального расхода воды, объема, продолжительности половодья и продолжительности его подъема не наблюдаются какие-либо значимые изменения (динамика в пределах 5% вероятности распределения как в сторону высоких, так и низких значений). Так, например, для максимального расхода воды в 1993 г. обеспеченность фактического значения по данным с уточнением составила 17% (14% без уточнения), для 2012 г. – 3% (по уточнению за 1-25 марта) и 8% (по уточнению за 1-10 апреля) вместо 5% (без уточнения).
Значительное улучшение результатов прогноза на основе анализа данных за 1-25 марта (т.е. до выпуска прогноза) было отмечено для даты наступления половодья. В 2012 г. это обусловило переход фактического значения из доверительного интервала 80% в 75% (увеличение обеспеченности от 21% до 27%). В 1993 г. обеспеченность фактического значения изменилась от 65% до 50%.
Обеспеченность фактического значения даты наступления максимального расхода воды после исключения некоторых сценариев изменяется в 1993 г. с 60% на 67 и 72%, в 2012 г. с 60% на 68 и 55% по данным 1-25 марта и 1-10 апреля соответственно.
Таким образом, на основе оценки исследуемых лет можно отметить, что уточнение прогноза даёт положительные результаты только для даты начала половодья. Возможно, к этому процессу следует привлекать также анализ рядов осадков.
Оценка эффективности методики выполнялась с помощью расчёта показателя RPSS (ranked probability skill score) (3), который показывает преимущество одной методики вероятностного прогноза относительно другой (в данном случае относительно прогноза по климатической норме) [Морейдо, 2015].
При положительном значении RPSS сравниваемая прогностическая методика даёт новую информацию относительно прогноза по климатической норме. Идеальным прогноз является при значении RPSS равным 1.
Для полученной методики при анализе по максимальному расходу воды половодья значение RPSS составило 0,12. В соответствии с этим можно заключить, что методика даёт новую информацию в сравнении с климатическим прогнозом.
Относительная характеристика RPS показала наибольшие погрешности вероятностного прогноза в годы наблюдения максимальных расходов низкой обеспеченности (1974, 1979, 1993, 2000, 2012 гг.).
Оценка эффективности детерминистической формы прогноза по характеристике максимального расхода, полученной путем осреднения ансамбля, дала неудовлетворительный результат (S/б = 0,96).
Таким образом, была получена и исследована методика составления прогноза характеристик весеннего половодья р. Северная Двина в районе г. Великий Устюг на основе 1993 и 2012 г. Полученный по ней вероятностной прогноз отвечает практическим требованиям качества. Можно отметить актуальность дальнейшего рассмотрения способов уточнения с помощью привлечения других метеорологических характеристик и варианта получения детерминистической формы прогноза путем моделирования на основе климатических данных. |
