«ЧИСЛЕННОЕ МОДЕЛИРОВАНИЕ ГИДРОФИЗИЧЕСКИХ ПРОЦЕССОВ В СТРАТИФИЦИРОВАННЫХ ОЗЁРАХ ...»

РОССИЙСКАЯ АКАДЕМИЯ НАУК

СИБИРСКОЕ ОТДЕЛЕНИЕ

ИНСТИТУТ ВЫЧИСЛИТЕЛЬНОГО МОДЕЛИРОВАНИЯ

На правах рукописи

Белолипецкий Павел Викторович

ЧИСЛЕННОЕ МОДЕЛИРОВАНИЕ ГИДРОФИЗИЧЕСКИХ

ПРОЦЕССОВ В СТРАТИФИЦИРОВАННЫХ ОЗЁРАХ

05.13.18 – математическое моделирование, численные методы и комплексы программ Диссертация на соискание учёной степени кандидата физико-математических наук

Научный руководитель:

член-корреспондент РАН, доктор физико-математических наук В.В. Шайдуров Красноярск -

СОДЕРЖАНИЕ

Введение

Глава 1. Озёра. Обзор математических моделей гидротермодинамики озер

1.1. Стратифицированные озёра. Озеро Шира

1.2. Перемешивание и течения

1.3. Тепломассообмен

1.4. Методы решения задач гидротермики

1.5. Трёхмерные модели

1.6. Двумерные в вертикальной плоскости модели

1.7. Плановые модели

1.8. Одномерные модели

Выводы

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

2.1. Одномерная модель

2.2. Двумерная в вертикальной плоскости модель

2.3. Гидростатическое приближение

Выводы

Глава 3. Комплекс компьютерных программ. Численные эксперименты

3.1. Комплекс программ ShiraModel

3.2. Результаты численных расчетов

Выводы

Заключение

Литература

Введение Антропогенное эвтрофирование, т. е. увеличение первичной продуктивности озер в результате обогащения их вод биогенными элементами, главным образом фосфором, стало во второй половине XX в. повсеместным явлением, охватившим около 90 % всех озер мира, включая крупнейшие из них. Большая часть озер мира расположена в зоне умеренных широт северного полушария. Общей особенностью природных вод этой зоны в естественных условиях является слабая насыщенность их фосфором, который, как правило, и является одним из основных факторов, лимитирующих развитие автотрофных организмов - первичных продуцентов. Рост поступления биогенов, а также увеличение сброса загрязнений, как следствие промышленного роста и интенсификации сельского хозяйства, привели к развитию процесса антропогенного эвтрофирования озер. Необходимость борьбы с антропогенным эвтрофированием водоемов и их загрязнением, принявшими глобальный характер, стимулировала проведение широкого круга исследований в области лимнологии, математического моделирования, экономики, связанных с проблемой сохранения, восстановления и эффективного использования природных ресурсов озер. Экологическое состояние водных объектов зависит от большого числа разнообразных факторов и процессов: гидрофизических, гидробиологических, гидрохимических, метеорологических и антропогенных. Гидрофизические процессы в значительной мере формируют среду обитания гидробионтов, определяют перенос и седиментацию веществ, интенсивность загрязнения и самоочищения водоёмов.

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

В настоящее время в прикладной математике при решении задач на ЭВМ сложилась технологическая цепочка, описывающая основные этапы вычислительного эксперимента: объект исследования – физическая модель – математическая модель – численный алгоритм – программа – расчёт на ЭВМ, сравнение с экспериментальными и другими данным [19, 79, 91].

Решение задачи начинается, как правило, с формализации объекта и построения соответствующей математической модели, основным требованием к которой является адекватное описание физических процессов, протекающих в исследуемых системах. Однако охватить всё многообразие физических явлений чрезвычайно трудно. Необходимо упростить проблему и рассмотреть только основные процессы, первостепенность которых определяется, в первую очередь, свойствами изучаемой системы и кругом моделируемых физических явлений. Известен ряд общих требований, предъявляемых к каждой математической модели: соответствующая система уравнений должна быть замкнутой и непротиворечивой; модель должна описывать широкий класс физических явлений и допускать конструирование реализуемого численного алгоритма.

На втором этапе вычислительного эксперимента для решения входящих в математическую модель уравнений при различных краевых условиях используется основной аппарат вычислительной математики – численные методы. Разрабатывается численный алгоритм и проводится его исследование.

На третьем этапе составляется вычислительная программа, реализующая построенный алгоритм. Далее проводятся вычисления на ЭВМ, выполняется анализ результатов, сопоставление их с теоретическими прогнозами, натурными данными и результатами физических экспериментов.

Одним из критериев выбора алгоритма, используемого при численном моделировании той или иной физической задачи, является объем вычислительной работы, необходимый для его реализации. Существует правило, что этот объем должен быть пропорционален реальным физическим изменениям, происходящим в моделируемой системе. Если алгоритм требует большого количества тяжелой вычислительной работы для расчета слабого эффекта или очень медленного физического процесса, то от такого «затратного» алгоритма следует, отказаться, выбрав более эффективный. Примером могут служить решения нестационарных задач, с шагом по времени (выбор которого диктуется условиями устойчивости) много меньшим масштаба реального изменения решения. Из вышесказанного можно сделать вывод, что создание адекватных математических моделей на основе «незатратных» алгоритмов является важной научной задачей.

Математические модели для расчета течений и температурного режима в озерах стали создаваться уже достаточно давно вслед за моделями циркуляции океана в 60-х годах. В 70-е годы были разработаны первые математические модели экосистем для американских Великих озер, озер Байкал и Ладожского. Но каждый гидрологический объект – река, озеро и т.д. – является частью географического ландшафта и одновременно комплексом взаимосвязанных физических явлений, происходящих внутри данного объекта. Взаимодействие же этого объекта с окружающей средой определяется метеорологическими условиями, разностью свойств и характеристик данного скопления воды и окружающей среды, физическими процессами, происходящими в данном объекте. Тесная связь окружающей среды и гидрологического объекта, их взаимное влияние делают задачи гидрофизических расчётов достаточно разнообразными и сложными, а методы их решения часто зависят не только от необходимой точности решения этих задач, но и от индивидуальных особенностей рассматриваемого объекта, не говоря уже об ограничениях, накладываемых имеющимися в распоряжении исследователя материалами наблюдений и трудоёмкостью расчётов.

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

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

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

Кроме того, высокая температура сама оказывает существенное воздействие на физические и химические свойства воды и связанные с ними биологические процессы. В частности, при повышении температуры снижается растворимость кислорода. Для озера Шира имеются данные измерений сотрудников Института биофизики СО РАН по распределениям температуры, солёности и других компонентов водной системы в глубоководной части водоема. Из этих данных следует, что концентрации исследуемых величин зависят в основном от глубины, а не от горизонтальных координат. Хотя к настоящему времени созданы трёхмерные модели термогидродинамических процессов в водоёмах, они являются весьма ресурсоёмкими и не всегда удобными в использовании. Применение таких моделей значительно ограничивает скорость проведения вычислительных экспериментов.

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

Целью работы является разработка одномерной вертикальной и двумерной в вертикальной плоскости математических моделей гидрофизических процессов в озёрах, создание компьютерных программ для моделирования гидрофизических и гидробиологических процессов в озёрах, проведение численных расчётов по озеру Шира.

Для достижения данной цели были поставлены и решены следующие задачи:

1. Разработать одномерную и двумерную модели гидрофизических процессов в озёрах с использованием «незатратных» алгоритмов.

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

3. Провести численное моделирование гидрофизических процессов в озере Шира.

Научная новизна диссертационного исследования заключается в разработке двумерной в вертикальной плоскости математической модели гидрофизических процессов в озёрах на основе «незатратного алгоритма», проведение численного моделирования гидрофизических процессов в озере Шира.

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

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

Работы выполнены совместно с Институтом биофизики СО РАН и при поддержке грантов РФФИ–ККФН (проект № 05-01-97700_ р_енисей), РФФИ-НВО (проект № 05-05НВО), Министерства образования и науки Российской Федерации и Американского фонда гражданских исследований и развития (грант RUX0-002-KR-06, программа «Фундаментальные исследования и высшее образование»), интеграционного проекта СО РАН № 24 «Роль микроорганизмов в функционировании живых систем», РФФИ (проект № 07а).

На защиту выносятся следующие результаты:

1. Разработана двумерная в вертикальной плоскости математическая модель гидрофизических процессов в озёрах.

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

3. Разработан комплекс компьютерных программ, реализующий одномерные и двумерные гидрофизические модели. Компьютерная модель позволяет исследовать динамику вертикальных распределений гидрохимических и гидробиологических компонентов водной экосистемы.

4. Проведено численное моделирование гидрофизических процессов в озере Шира.

Представление материалов диссертационной работы на конференциях. Основные результаты работы были доложены на IV Всероссийской конференции молодых ученых по математическому моделированию и информационным технологиям (Красноярск, 2003), на Международной научной студенческой конференции «Студент и научнотехнический прогресс» (Новосибирск, 2004), на конференции-конкурсе молодых учёных ИВМ СО РАН (Красноярск, 2006), на VII Всероссийской конференция молодых ученых по математическому моделированию и информационным технологиям (Красноярск, 2006), на конференции-конкурсе молодых учёных ИВМ СО РАН (Красноярск, 2007), на VIII Всероссийской конференция молодых ученых по математическому моделированию и информационным технологиям (Новосибирск, 2007), на конференции-конкурсе молодых учёных ИВМ СО РАН (Красноярск, 2008).

Публикации. По теме диссертации опубликовано 7 работ, в том числе (в скобках в числителе указан общий объём этого типа публикации, в знаменателе – объём, принадлежащий лично автору) 2 статьи в издании, рекомендованном ВАК для представления результатов докторских диссертаций [15, 10] (21/15 стр.), 3 – в трудах международных и всероссийских конференций [12, 13, 14] (12/8 стр.), 2 – в материалах конференций молодых учёных [16, 18] (11/9 стр.) Личный вклад автора. В работе [10] автор участвовал в постановке задач, конструировании численных алгоритмов решения двумерных задач ветровых течений в стратифицированных водоёмах, выполнил программную реализацию и провёл численные эксперименты. В работе [14] участвовал в конструировании численного алгоритма для одномерной вертикальной модели и его программной реализации. В работе [16] автор участвовал в постановке задачи и реализации численного алгоритма решения задачи протекания.

Структура и объём диссертации.

Работа состоит из введения, трёх глав, заключения и списка литературы из 101 наименований. Полный объём диссертации составляет 112 страниц, включая 39 рисунков, 5 таблиц.

Обзор содержания работы.

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

Первая глава содержит описание объекта исследования - озера Шира, и краткий обзор работ по математическому моделированию гидрофизических процессов в водоёмах.

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

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

В третьей главе рассматривается конкретный объект исследования озеро Шира и комплекс программ, в котором были реализованы описанные во второй главе численные алгоритмы. Приводятся сравнения результатов расчетов по одномерной и двумерной моделям с натурными измерениями.

В заключении приведены основные результаты, полученные в диссертационной работе.

Автор выражает благодарность Шайдурову В.В. за поставленные задачи и научное руководство, а также Белолипецкому В.М., Дегерменджи А.Г., Прокопкину И.Г. за ценные советы и плодотворное обсуждение результатов работы.

Глава 1. Озёра. Обзор математических моделей гидротермодинамики озер 1.1. Стратифицированные озёра. Озеро Шира Озеро - это водоем, окруженный сушей. По размерам озера варьируют от очень крупных, таких как Каспийское море и Великие озера в Северной Америке, до крошечных водоемов площадью несколько сотен квадратных метров и даже меньше. Вода в них может быть пресной, как в оз. Верхнем, или соленой, как в Мертвом море. Озера встречаются на любых высотах, от самой низкой на Земле абсолютной отметки на поверхности суши -408 м (Мертвое море) и до самой высокой (в Гималаях). Некоторые озера не замерзают круглый год, тогда как другие, например оз. Ванда в Антарктиде, большую часть года скованы льдом. Многие озера существуют постоянно, а другие (например, оз.

Эйр в Австралии) — лишь изредка заполняются водой. Несмотря на разнообразие, озера всех типов имеют ряд общих физических, химических и биологических характеристик и подчиняются многим общим законам. Изучением озер во всем их многообразии и во всех аспектах занимается одна научная дисциплина — озероведение, или лимнология.

Вода является превосходным растворителем и поэтому в озерных водах содержится много растворенных веществ. Примечательно, однако, что подавляющая масса этих веществ в большинстве озер представлена ограниченным числом соединений, а именно, положительно заряженными ионами (катионами) кальция, магния, натрия и калия и отрицательно заряженными ионами (анионами), состоящими из углерода и кислорода (бикарбонаты), серы и кислорода (сульфаты) и хлора (хлориды) (обе группы ионов перечислены в порядке убывания их содержания). Эти семь ионов составляют от 90 до 95% общего количества растворенных веществ в водах большинства озер, а их суммарная концентрация, обычно измеряющаяся в миллиграммах на литр (мг/л), характеризует соленость (минерализацию) воды. Другие вещества, например элементы питания растений (азот и фосфор) и металлы (железо и марганец), присутствуют в существенно меньших количествах, так что их концентрации измеряются в микрограммах на литр (мкг/л). В бессточных озерах испарение приводит к изменению состава солей. Озера называются хлоридными, сульфатными или карбонатными в зависимости от того, какие анионы накопились в них в наибольшем количестве под воздействием испарения или атмосферных осадков.

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

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

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

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

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

Движение воды в озерах значительно отличается от высокоамплитудных приливоотливных и мощных океанических течений. Только в таких крупнейших озерах, как Верхнее и Мичиган, существуют постоянные течения, но даже в них практически отсутствуют приливо-отливные колебания (их амплитуда в оз. Верхнем составляет лишь 3 см). Под воздействием температурного градиента, впадающих водотоков и ветров в озерах совершается движение воды. Когда в стратифицированное озеро впадает река, либо в поверхностном слое, либо на средних глубинах возникает стоковое течение. Поверхностные течения формируются, когда воды притока имеют меньшую плотность, чем воды самого озера, как, например, летом при впадении р. Иордан в Тивериадское озеро. Среднеглубинные течения образуются, если водоток устремляется вниз к слоям, соответствующим его собственной плотности. Если плотность водотока выше плотности любого слоя озерной воды, он опустится на дно и образует придонное течение. При этом возможно даже формирование подводного русла, как, например, при впадении р. Роны в Женевское озеро. Под влиянием ветра возникает несколько типов движений озерных вод. Во-первых, происходит турбулентное перемешивание верхних слоёв жидкости. Другой тип движения воды происходит, когда ветер постоянно дует над поверхностью озера. Поскольку вода перемещается по ветру, уровень воды в дальнем конце озера несколько поднимается, что приводит к формированию компенсационного течения — либо вдоль берегового, если озеро мелкое, либо, в более глубоких озерах, противоположно направленного и проходящего на некоторой глубине от поверхности.

Хотя химический состав озера важен для всех организмов, о чем свидетельствуют, например, специализированные виды растений и животных, обитающие в соленых озерах, именно растения, осуществляющие фотосинтез, имеют наибольшее влияние. В процессе фотосинтеза солнечная энергия используется для превращения углекислоты и воды в углеводороды и кислород. При этом помимо диоксида углерода и воды в фотосинтезе участвуют еще 18-20 химических элементов, и уменьшение содержания любого из них ниже оптимальной потребности существенно замедляет процесс фотосинтеза. Эта так называемая гипотеза лимитирующей роли питательных элементов, выдвинутая в середине 19 в. Юстусом Либихом, до сих пор используется при характеристике водных экосистем. В пресных водоемах большинство питательных элементов присутствует в количествах, превышающих потребность в них, однако два из них — азот и фосфор — относительно редки. Именно эти элементы, порознь или совместно, лимитируют процесс фотосинтеза, или первичную продукцию. Более того, поскольку некоторые синезеленые водоросли способны связывать атмосферный азот, превращая его в аммоний и используя в процессе фотосинтеза, а фосфор не имеет такого источника, то последний становится наиболее важным лимитирующим элементом. В результате многие существенные характеристики озер, как, например, суммарный прирост первичной продукции или обилие водорослей, находятся в прямой зависимости от содержания фосфора в озерах. Поэтому озера классифицируют по этому показателю. Выделяются олиготрофные озера (с низким содержанием питательных веществ), мезотрофные (со средним содержанием) и эвтрофные озера (с высоким содержанием питательных веществ).

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

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

Однако озера редко пребывают в равновесном состоянии. Напротив, они часто используются как источники воды для орошения, питьевой воды, для сельскохозяйственных нужд или же для сброса таких продуктов современной цивилизации, как сточные воды предприятий, ливневые и сельскохозяйственные стоки. Озера загрязняются все возрастающим количеством пестицидов, гербицидов и попадающих в воду из воздуха органических соединений, таких, как полихлорированные бифенилы, а также кислотными дождями, образующимися в результате выбросов загрязняющих веществ двигателями автомобилей и тепловыми электростанциями. В них проникают чуждые им виды растений и животных, заносимые рыбаками на днищах судов и иными случайными способами. Угрожающие размеры принимает эвтрофикация, или избыточное обогащение озер питательными веществами из антропогенных источников, которое наносит значительный экологический ущерб. В некоторых случаях большие, имеющие хозяйственное значение озера, находятся даже под угрозой полного исчезновения. Так, например, объем воды в Аральском море (крупном соленом озере) сократился в настоящее время вдвое вследствие разбора на орошение вод впадающих в него Амударьи и Сырдарьи. В результате его соленость возросла почти в три раза (с 9,6-10,3 г/л до 27-30 г/л). Обнажившиеся участки морского дна развеваются пыльными бурями, что приводит к выносу солей и пестицидов и осаждению их в пределах близлежащих заселенных территорий. Загрязнение озер очень серьезная проблема. Например, чтобы снизить эвтрофикацию водоемов, во многих странах приняты законы по ограничению концентрации фосфора в водах, прошедших через очистные станции и которые могут попасть в озера. Сформировалась целая наука о восстановлении озер, базирующаяся главным образом на эмпирических соотношениях, связывающих такие показатели, как обилие водорослей и прозрачность воды, с концентрациями фосфора в озерных водах. В некоторых регионах регулируется забор воды из озер. Тщательно изучается применение пестицидов.

Уникальным объектом для математического моделирования является озеро Шира вследствие его малого биоразнообразия. Это позволяет использовать математическое моделирование для описания динамики водной экосистемы. Озеро Шира относится к солёным, стратифицированным водоёмам с высоким содержанием сульфата. Круглогодично в его воды поступает органическое вещество аллохтонного и автохтонного происхождения.

Органика интенсифицирует процессы окисления и восстановления серы. Сложившийся комплекс физико-химических условий функционирования озера Шира не имеет аналогов.

Озеро Шира представляет собой бессточное озеро без островов, в которое впадает одна речка Сон. В силу малости притока все влияние реки сосредоточено в приустьевой зоне.

На рис.1.1 представлена батиметрия озера.

Основные характеристики озера: озеро имеет овальную форму, длина – 9.4 км, ширина – 5 км, площадь водного зеркала 34.7 кв. км, средняя глубина – 11.2 м, максимальная глубина – 24 м, естественный среднемноголетний водосток – 1.6 млн. куб. м/год, подземный водообмен составляет 9% от общего водопоступления. Минерализация воды – 16-19 г/л. В летний период в озере формируется плотностная стратификация. В эпилемнионе более высокая температура воды (15 – 250С) и низкая солёность (15-17 г/л), в гиполемнионе – более высокая соленость (18-19 г/л) и холодная вода (2-30С). Наблюдаются заметные изменения во времени толщины верхнего квазиоднородного слоя и заглубления термоклина. Температурный режим озера зависит от метеоданных: температуры и влажности воздуха, силы и направления ветра, облачности. Характерные распределения температуры и солёности зимой, весной, летом и осенью представлены на рис. 1.2.

Рис. 1.2. Характерные распределения температуры и солёности в озере Шира 1.2. Перемешивание и течения Под перемешиванием понимают обмен между отдельными слоями или объемами водоема. Внешние воздействия, изменение объема воды под влиянием разных причин, формирующие течения и перемещающие значительные массы воды, перемещают вместе с ними и заключенные в них взвеси, химически растворенные вещества, газы, запасы тепла. При передвижении происходит смешение различных объёмов вод в месте их соприкосновения и выравнивание физико-химических и других характеристик в слоях, вовлеченных в перемешивание.

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

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

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

Другой фактор - динамическое (или вынужденное) перемешивание, возникающее при наличии вертикальных и горизонтальных градиентов скорости, создающихся за счет трения на поверхностях раздела: внутри водной толщи, у дна и поверхности воды. Основной причиной течений и интенсивного перемешивания в озерах чаще всего бывает действие ветра. С началом ветра возникает перемешивание поверхностных слоев за счет образования волн и опрокидывания их гребней. Постепенно перемешивание захватывает все большую глубину, так как трение воздуха о воду и давление ветра на тыловую сторону волн создает ветровое течение. На поверхности раздела между верхним движущимся слоем и подстилающим его неподвижным возникают вихри с горизонтальными осями, перпендикулярными направлению движения, способствующие распространению перемешивания в глубину. Чем сильнее и продолжительнее ветер, тем на большую глубину распространяется перемешивание. Максимальные градиенты скорости, обусловленные ветром и волнением, наблюдаются в поверхностных слоях, где динамическое перемешивание наиболее сильно. На мелководьях и у берегов наибольшие градиенты скорости наблюдаются у дна, и интенсивность перемешивания увеличивается ко дну. Глубина проникновения ветрового перемешивания зависит от размеров и глубины озера, от его термического состояния, а для соленых озер – от вертикального распределения солености.

Мелководные озера перемешиваются до дна. На глубоководных озерах перемешивание обычно захватывает только поверхностные слои до глубины 30 м, но за счет возникновения глубинных противотечений, связанных с поверхностными ветровыми течениями и уплотнения при смешении, перемешиванием может быть охвачена и более значительная толща. Так, на Байкале глубина динамического перемешивания при длительных штормовых периодах осенью захватывает толщу до 200-300 м [27]. Осень и весна наиболее благоприятны для развития ветрового перемешивания в пресноводных водоемах благодаря малым градиентам температуры (а, следовательно, и плотности) в водной толще. Летом, когда вертикальное распределение температуры характеризуется прямой стратификацией за счет интенсивного прогрева поверхностных слоев, глубина проникновения и скорость распространения ветрового перемешивания ограничиваются слоями с большими градиентами плотности (слоем температурного скачка).

Важное значение имеют ветровые течения. Теория ветровых течений, возникающих при постоянном ветре над беспредельным и однородным по плотности водоемом, разработана В. Экманом (1905 г.) и развита в работах В. В. Шулейкина и В. Б. Штокмана.

Дальнейший вклад был сделан В. М. Маккавеевым [57] и А. В. Караушевым [43], рассматривавшими конкретные условия формирования течений в водоёмах при нагоне. Однако связь между скоростью течения и скоростью ветра, создавшего эти течения, на основе теоретических расчетов еще не может быть получена для любого водоема. Трудность решения этого вопроса связана с разнообразием условий передвижения воды в ограниченной и сложной по форме котловине озера, поэтому для реальных озер эту зависимость пока еще приходится устанавливать по данным непосредственных измерений.

По данным В. К. Давыдова, скорость ветрового течения на поверхности Онежского озера при скорости ветра от 1,0 до 4,2 м/с изменялась от 3 до 18 см/с. По данным Г. Р. Юнусова, скорость нагонных течений, наблюдавшихся в проливе на оз. Балхаш, достигала 0,87м/с, доходя при отдельных нагонах до 1,07- 0,95 м/с в начале и уменьшаясь к концу нагона до 0,21 м/с. А. Н. Охлопкова [70] для Ладожского озера получила зависимость скорости поверхностного течения от скорости установившегося ветра для участков с различными глубинами на открытой части озера, причем зависимость имеет линейный характер. Для заливов, где сказывается влияние окружающих берегов, зависимость имеет параболический характер. Неодинаковый вид зависимости между скоростями ветра и течения даже в пределах одного озера объясняется тем, что скорость поверхностного течения в озерах зависит не только от скорости ветра, но и от продолжительности действия, его направления, распределения глубин в озере, наличия островов, изрезанности берегов и их высоты, создающей «ветровую тень», характера погоды в предшествующий период (направления и скорости предшествовавших ветров, наличия постоянных течений и т. п.).

Для определения напряжения трения ветра используются различные формулы [83].

Наиболее очевидным эффектом воздействия ветра на поверхность водоёма является генерация волн. Помимо колебательного волнового движения ветер вызывает также стационарное (в определённом смысле) горизонтальное движение поверхностного слоя примерно в том же направлении, что и у ветра, причём за счёт внутреннего трения это движение передаётся глубоким слоям. Поверхностное течение определяет скорость и направление дрейфа плавающего материала, такого как нефтяная плёнка или твёрдый мусор. Ниже поверхности скорость ветрового течения изменяется с глубиной: как правило, вектор скорости уменьшается по глубине и поворачивается. Дрейф плавучего объекта, имеющего осадку несколько метров, например, буя или судна, зависит от среднего течения в этом слое, которое может значительно отличаться от поверхностного течения. Берег ограничивает движение воды, поэтому вода может здесь накапливаться и на большой площади у берега возникает наклон уровня. Наклон уровня порождает горизонтальные градиенты давления в воде, а это в свою очередь приводит к возникновению течений, которые накладываются на первоначальный ветровой дрейф. Градиентные течения могут быть направлены против течения, вызванного непосредственно воздействием ветра, и в этом случае на некоторой глубине направление течения изменяется на обратное. Таким образом, возникает двухслойный поток: в верхнем слое вода течёт в направлении ветра, в нижнем, компенсационном, потоке – в противоположном направлении [22].

При перемене направления и скорости ветра глубинное течение изменяется медленнее, чем поверхностное, что объясняется большой инерцией процесса; направление и скорость течения по вертикали могут изменяться в достаточно широких пределах [69].

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

Например, на водохранилищах формирование сгоннонагонных течений усложняется наличием постоянного движения воды к плотине. Ветровые напряжения могут менять характер распределения течений, и даже их направление (рис. 1.3). Случаи I и II отвечают относительно большому транзитному расходу через водоём или слабому ветру. При этом в случае, когда ветер дует по направлению к плотине, в водоёме создается сточно-дрейфовое или дрейфовое течение (I). В случае (II) ветер, дующий от плотины, не может вызвать обратное течение, хотя в верхней части водоёма уровни все же повышаются, но течение в водоёме сохраняет градиентный характер. В случае III вследствие малого расхода либо сильного ветра течение на поверхности у плотины принимает характер дрейфового (ветрового), а у дна и верхней части водоёмах сохраняется градиентное стоковое течение. В случае IV при малом или очень сильном ветре нагон к плотине будет настолько велик, что приведет к возникновению компенсационного течения в приплотинной части водоёма [43].

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

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

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

где U i, p - средние значения; U i', p ' - пульсационные составляющие. После подстановки этих выражений в уравнения Навье-Стокса, записанные относительно мгновенных значений, получается система уравнений осредненного турбулентного потока, уравнения Рейнольдса [56]:

Здесь x j - координаты в прямоугольной системе (в дальнейшем будем считать, что x3 вертикальная координата, направленная вниз); U i - компоненты вектора скорости; плотность; p - давление; t - время; Fi - компоненты вектора объемных сил; - кинематическая вязкость, здесь и далее используется тензорное правило суммирования по повторяющемуся индексу.

Для описания течений используются различные модели на основе этой системы.

Система уравнений Рейнольдса является незамкнутой, так как содержит кроме четырех искомых функций U i, p девять неизвестных величин турбулентных напряжений U 'jU i'.

Проблемам замыкания системы уравнений Рейнольдса посвящена обширная литература (см., например, [65]). Можно, используя только уравнения Навье-Стокса и неразрывности (не привлекая дополнительной теоретической и экспериментальной информации), построить последовательность систем уравнений для осредненных значений гидродинамических характеристик потока, однако каждая из этих систем оказывается незамкнутой (цепочка уравнений Келлера-Фридмана). Поэтому системы уравнений для осредненных характеристик потока приходится замыкать, привлекая дополнительные зависимости, основанные на экспериментальных данных, аналогиях и т.п. (полуэмпирические теории турбулентности). Наиболее разработаны методы замыкания уравнений Рейнольдса. Так, согласно известной концепции Буссинеска о турбулентной вязкости принимается, что турбулентные напряжения определяются по формуле где Tij - составляющие коэффициента турбулентной вязкости, характеризующие статистические свойства поля пульсационных скоростей. В общем случае коэффициент турбулентной вязкости является тензорной величиной, для определения которой нужны дополнительные соотношения. Поэтому формула по существу вводит вместо неизвестных U 'jU i' новые неизвестные Tij. Далее считают, что коэффициенты Tij включают молекулярные и турбулентные составляющие. Предполагается, что режим турбулентных флюктуаций изотропен в горизонтальной плоскости и анизотропен в вертикальной. Коэффициенты горизонтального турбулентного обмена полагают постоянными. Так как изменение плотности воды в водоёмах существенно меньше самой плотности, то стратифицированное течение может быть описано в приближении Буссинеска, т.е. в уравнениях движения переменная плотность заменяется некоторым ее средним постоянным значением везде, кроме членов, описывающих архимедову силу.

С учетом указанных допущений уравнения Рейнольдса в приближении Буссинеска имеют вид [86]:

Здесь Кориолиса, К1, К2 – коэффициенты горизонтального турбулентного обмена, Кz – коэффициент вертикальной турбулентной вязкости, здесь и далее знак осреднения у величин Ui, p опущен.

Для замыкания системы уравнений (1.2.3) используются различные полуэмпирические модели турбулентности, например, уравнение энергии турбулентности, уравнение скорости диссипации турбулентной энергии и т.д. [46, 67, 76]. Коэффициенты горизонтального турбулентного обмена полагают постоянными.

К основным процессам, формирующим турбулентную структуру вод глубоких водоемов, относятся ветровые течения, обрушение поверхностных и внутренних волн, конвективные процессы, трение о дно, речной приток и др. [81]. Многообразие механизмов приводит к определенным трудностям в создании методов расчета характеристик турбулентности.

Методы, основанные на алгебраических соотношениях, позволяют определить локальные характеристики турбулентного движения независимо от его предыстории. Модель для описания коэффициента турбулентной вязкости K z, предложенная Прандтлем, основана на понятии пути смешения. Согласно теории Прандтля коэффициент турбулентной вязкости определяется соотношением где U - горизонтальная составляющая скорости, z - вертикальная координата, l - длина пути смешения, характеризующая масштаб турбулентности в данной точке. Установлено, что для пристеночного пограничного слоя однородной жидкости l = z, где = 0.4 постоянная Кармана, а в остальной области величина l пропорциональна толщине пограничного слоя. Для расчётов турбулентного переноса тепла или массы коэффициент турKz ное число Шмидта).

Плотностная стратификация может существенно изменить распределение масштаба турбулентности. Этот эффект обычно учитывается эмпирическими зависимостями длины пути смешения от числа Ричардсона Манком и Андерсеном [100] на основе данных натурных измерений для стратифицированных течений предложена связь между коэффициентом турбулентного обмена импульса и температуропроводности в зависимости от числа Ричардсона:

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

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

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

Наиболее распространенной полуэмпирической моделью является k модель турбулентности, основанная на использовании уравнений переноса кинетической энергии ( k ) и скорости диссипации этой энергии ( ), уравнения для которых используются в виде [76]:

где С помощью k модели можно рассчитывать сдвиговые слои и ограниченные рециркуляционные течения.

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

Анализ всех членов в уравнении баланса турбулентной энергии показывает [59], что в верхней области турбулентного слоя преобладающими по отношению к диффузии турK z k булентной энергии ] являются диссипация турбулентной энергии и сдвиговый источник турбулентности Если пренебречь вкладом диффузионного члена в уравнение энергии турбулентноcµ k здесь K z - коэффициент вертикального турбулентного обмена, с = 0.05, h - глубина верхнего квазиоднородного слоя, которая определяется по первой от поверхности расчетной точке, в которой K z ( z ) min K z.

Численные решения уравнений экмановского пограничного слоя (с применением формулы (1.2.4)), удовлетворяющие соответствующим граничным условиям, хорошо согласуются с решениями, полученными по динамическим уравнениям турбулентности [59]. С учетом этого, для вычислительной экономичности математических моделей, коэффициент вертикального турбулентного обмена можно определять по формуле (1.2.4).

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

1.3. Тепломассообмен Распространение теплоты в движущейся сплошной среде описывается уравнением [56]:

где Т – температура, U i - компоненты скорости, aij - компоненты тензора молекулярной температуропроводности, Fтепл - источники (стоки) теплоты.

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

Обычно предполагается, что главные направления составляющих aij совпадают с направлением координатных осей - это дает ортотропную молекулярную температуропроводность ai. Для несжимаемой жидкости с учетом уравнения неразрывности уравнение (1.3.1) может быть записано в виде Турбулентность потока учитывается путем представления переменных в уравнениях (1.3.1), (1.3.2) в виде суммы соответствующих осредненных (по ансамблю реализаций) величин и пульсационных составляющих Подставляя эти выражения для мгновенных величин U i и T в (1.5.2), получаем Пульсационный член в (1.3.3) обычно представляется в виде где bi - компоненты турбулентной температуропроводности, т.е. в виде, аналогичном закону Фурье в теории молекулярной теплопроводности.

В результате после подстановки (1.3.4) в (1.3.3) получается полуэмпирическое уравнение турбулентного теплопереноса Здесь K1T = a1 + b1, K 2T = a2 + b2, K zT = a3 + b3 - эффективные коэффициенты температуропроводности, включающие молекулярные и турбулентные составляющие, знак осреднения у величин U i, T опущен.

При решении (1.3.5) могут быть использованы следующие граничные условия:

- на части S1 границы области S задается температура T = T * ;

- на части S 2 границы области S задается нормальный тепловой поток Qтепл = Q *.

Применительно к водоёмам необходимо задать граничные условия на свободной поверхности воды, "твердой" границе (дно, боковая поверхность) и "жидкой" границе (водозабор, водовыпуск, границы с другими водоемами и водотоками), а также внутренние источники и стоки тепла.

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

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

В тепловом балансе в соответствии с направлением поступления тепла в водоем могут быть выделены три группы потоков [69, 56]:

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

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

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

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

Здесь F - изменение количества тепла в водоеме.

При расчетах теплового баланса удобно объединять в один член ряд составляющих, общих по физическому смыслу или характеру теплообмена. Общепринятым является уравнение теплового баланса, записываемое в сокращенном виде как где R = Fр Fи + Fа - радиационный баланс водоема; B = Fпр Fст ± Fдно ± F - теплоаккумуляция водоема.

В этом случае пренебрегают рядом малозначащих составляющих теплового баланса, такими, как Fос, Fгр, Fб. Величины, связанные с особенностями гидрологического режима водоема Fдин или сезона Fл, учитываются при необходимости.

Соотношение величин, входящих в уравнение теплового баланса, непостоянно и зависит от времени года, климатических и метеорологических условий, размеров и высоты расположения водоема над уровнем моря и др. Составляющие Fр, Fа, Fдин характеризуют приток тепла к воде; Fи, Fст являются потерями тепла. Члены P, Fис, Fос, Fдно, Fгр могут менять знак в зависимости от условий. При более высокой температуре воздуха по сравнению с температурой воды P будет характеризовать приток тепла, при обратном соотношении - потери. Аналогично может быть определен знак других составляющих в соответствии с соотношением температуры воды водоема и осадков, поверхности дна, грунтового питания, притоков.

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

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

Таким образом, в первую половину сезона, свободного ото льда, водоем накапливает тепло, во вторую - расходует его. Амплитуда этих величин зависит от конкретных метеорологических условий года.

Проблеме теплообмена через свободную поверхность водоемов посвящена обширная литература (см., например, [77, 8, 75, 69]). Основными факторами, определяющими плотность теплового потока на свободной поверхности водоема, являются: суммарная (прямая и рассеянная) солнечная радиация с учетом отражения, эффективное излучение водной поверхности, испарение, конвективный теплообмен между водой и воздухом.

Выражение для плотности теплового потока на свободной поверхности водоема можно записать в виде где FR - составляющая, характеризующая суммарный радиационный теплообмен (включая излучение); FE - плотность теплового потока за счет испарения; FC - плотность конвективного теплового потока.

Имеются различные зависимости для определения составляющих радиационного теплообмена, тепловых потоков, связанных с испарением и конвекцией. Например, в [101] предлагается следующая формула для расчета суммарного радиационного теплообмена]:

FR = 0.94 FSC (1 0.65n 2 ) + 5.17 10 13 (Ta + 273.15) 6 (1 + 0.17n 2 ) 5.5 10 8 (TS + 273.15) где FSC - коротковолновая солнечная радиация при ясном небе; n - облачность в долях единицы; Ta и TS - температура воздуха на высоте 2 м и водной поверхности соответственно, °С.

В ряде работ (см., например, [75]) приводятся зависимости для определения FR, учитывающие ряд дополнительных факторов, влияющих на радиационный теплообмен (местоположение и размеры водоема, альбедо поверхности, ярусность расположения облаков и т.д.). Как правило, эти факторы учитываются введением различных поправочных коэффициентов. В качестве примера приведем формулу сравнительно простого вида где F0 - суммарная солнечная радиация на уровне моря при альбедо, равном нулю; a p альбедо поверхности воды; I - эффективное излучение при безоблачном небе, определяемое как функция влажности и температуры воздуха; k и c - коэффициенты, зависящие от широты места расположения водоема.

Для расчета теплоотдачи испарением используют зависимость где f (W z ) - так называемая ветровая функция; e M - давление насыщенного пара при температуре водной поверхности; W z и e z - скорость ветра и давление пара на высоте z над поверхностью (обычно принимают z=2 м).

В литературе встречаются различные формулы для определения ветровой функции (см., например, [56, 8]); некоторые авторы учитывают зависимость f (W z ) не только от скорости ветра, но и от соотношения между температурами воздуха и водной поверхности [51].

Расчет теплообмена конвекцией обычно выполняют исходя из соотношения Боуэна, полученного в предположении подобия процессов тепло- и массопереноса:

где e - давление пара на высоте 2 м над водоемом.

Из (2.25) и (2.26) следует Таким образом, составляющие плотности теплового потока FE и FC часто записывают в виде где ae и aC - коэффициенты теплоотдачи испарением и конвекцией соответственно.

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

Тепловые потоки через дно и боковые границы водоема обычно пренебрежимо малы по сравнению с потоком через свободную поверхность и поэтому на этих границах полагается ( n - нормаль к поверхности) На водовыпуске может быть задана температура сбрасываемой в водоем воды Более сложным является вопрос о задании условия для температуры на водозаборе.

Часто предполагают, что конвективный перенос теплоты существенно больше диффузионного, что приводит к условию Другой путь для формулировки граничного условия на водозаборе - идентификация функции, описывающей тепловой поток на границе по результатам измерений температуры в водоеме.

Уравнение переноса субстанции аналогичное уравнению для температуры:

или в дивергентном виде Здесь С - переносимая субстанция, K1С, K 2С, K zС - коэффициенты диффузии, FС - внутренние источники (стоки). Для постановки конкретных задач к уравнению (1.3.6) добавляются начальные и граничные условия.

1.4. Методы решения задач гидротермики Для изучения гидротермического режима (гидродинамических характеристик и полей температур) в водоёмах в зависимости от целей исследования, особенностей водоема, полноты исходной информации, возможностей вычислительных средств, традиций и т.п.

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

В общем случае уравнения гидродинамики и теплопереноса следует рассматривать совместно, так как плотность воды зависит от температуры.

Возможны два способа решения задачи:

- уравнения гидродинамики и теплопереноса решаются совместно, в результате чего определяются все характеристики гидротермического процесса;

- задача решается в два этапа: сначала определяется поле скорости при известной температуре, а затем решается уравнение теплопереноса.

Если в первом случае исходить из трехмерных полуэмпирических уравнений турбулентного потока и теплопереноса, не используя более сложных моделей турбулентности, то будем иметь систему уравнений гидротермики относительно шести неизвестных U i, p,, T [56] которая решается при соответствующих граничных и начальных условиях. Во втором случае сначала решается система уравнений гидродинамики, а затем с использованием полученного поля скоростей и уравнение теплопереноса.

Если коэффициенты турбулентного обмена, содержащиеся в уравнениях (1.4.1), определяются с привлечением какой-либо модели турбулентности (однопараметрической, двухпараметрической и т.п.), то количество неизвестных и уравнений в системе (1.4.1) увеличивается.

Некоторое упрощение трехмерных уравнений (1.4.1) может быть достигнуто также за счет каких-либо приближений – медленных течений, гидростатики и т.д.

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

Иногда используется еще один способ видоизменения исходной системы типа (1.4.1) - выделение баротропной и бароклинной составляющих [60].

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

Уравнения гидродинамики и теплопереноса, служащие основой для изучения гидротермических процессов в водоёмах, могут быть решены аналитическими методами лишь в отдельных частных случаях при существенной схематизации физической и геометрической сторон явления (последовательное рассмотрение задач гидродинамики и теплопереноса, линеаризованные уравнения, постоянные коэффициенты, стандартные области и т.п.). Если удается свести исходные уравнения к некоторым типам классических уравнений математической физики (уравнениям Лапласа, Пуассона, и т.п.), то тогда известными методами (разделения переменных, интегральных преобразований и др.) можно получить аналитическое решение. Однако при этом нередко форма аналитического решения оказывается достаточно громоздкой и затруднительной для получения численных результатов.

Все это приводит к тому, что при решении задач гидродинамики и теплопереноса (в частности, при изучении гидротермических процессов в водоёмах) широко используются различные численные методы с использованием ЭВМ.

Для численного решения задач гидродинамики и теплопереноса в той или иной постановке в настоящее время применяются различные методы: метод конечных разностей (МКР), метод конечных элементов (МКЭ), метод граничных элементов (МГЭ), различные их комбинации и т.п. В результате применения численных методов исходное дифференциальное уравнение (или система уравнений) с граничными условиями сводится к системе алгебраических уравнений.

Начиная с конца 60-х - начала 70-х годов 20 века все большее распространение при решении самых разнообразных задач механики сплошных сред (в том числе гидродинамики и тепломассопереноса) получает метод конечных элементов [63]. По сравнению с широко распространенным МКР преимущества МКЭ обусловлены его вариационноразностной сущностью, большей гибкостью при задании положения узлов и типов элементов для областей сложной геометрической конфигурации, простотой исследования неоднородных сред, широкими возможностями использования стандартных программ на всех этапах решения и т.п.

Вопросы применения МКЭ к исследованию гидротермических процессов в водоемах-охладителях изложены в [81, 101].

В последние годы для решения аналогичных задач начал применяться также метод граничных интегральных уравнений (МГИУ) или МГЭ, позволяющий понизить размерность задачи на единицу и во многих случаях, как показывают сопоставительные расчеты, оказавшийся весьма эффективным. В частности, этот метод не обладает такими недостатками классических МКР и МКЭ, как трудность описания бесконечных и полубесконечных областей, необходимость решения систем уравнения высокого порядка и определения сложных характеристик объекта (матриц "жесткости", масс и т.п.).

Для исследования гидрофизических процессов в стратифицированных водоемах рассматриваются математические модели различного уровня. При выборе математических моделей для исследования гидрофизических процессов необходимо учитывать основные особенности изучаемого объекта. Так для расчета конвекции глубоких озер типа Байкал используется модель, в которой учитывается сжимаемость воды и два параметра Кориолиса [89]. В работах О.Ф. Васильева, В.И. Квона и их сотрудников были предложены и реализованы двумерные и трёхмерные модели [24, 25, 29, 44]. Для морей и глубоких озер разработаны модели на основе уравнений гидродинамики турбулентного движения и распространения солености и температуры с учетом влияния сил Кориолиса, ветрового напряжения на поверхности и трения на дне [4, 5]. Для вытянутых в продольном направлении бассейнов используется двумерная в вертикальной плоскости (x,z) модель с двухслойной температурной стратификацией [36, 37]. Двумерные xz-модели часто используются для моделирования стратифицированных течений в водоемах [95, 98]. Двумерные (в горизонтальной плоскости) так называемые уравнения мелкой воды используются в случае, если горизонтальные размеры водоемов гораздо больше глубины (Онежское озеро, Каспийское море, Чудское озеро). Получаемые при этом осредненные по глубине значения скоростей позволяют определить горизонтальную циркуляцию вод [1, 7, 37, 71]. В тоже время в ряде работ [5, 35, 37, 52, 29] рассматриваются методы, позволяющие находить скорости в зависимости от глубины, используя предварительно вычисленные величины осредненных по глубине скоростей и возвышения свободной поверхности. Для оценки вертикального распределения температуры в слабопроточных водоемах применяются одномерные модели [41, 93].

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

В работах О.Ф. Васильева, В.И. Квона и их сотрудников [29, 44, 24, 25] были предложены и реализованы некоторые трехмерные модели, учитывающие ряд характерных особенностей водоёмов и позволяющие получить подробную информацию о характере гидротермических процессов с учетом большого числа различных факторов. Так как изменение плотности воды в водоёмах существенно меньше самой плотности, то стратифицированное течение может быть описано в приближении Буссинеска, т.е. в уравнениях движения переменная плотность заменяется некоторым ее средним постоянным значением везде, кроме членов, описывающих архимедову силу. Глубина водоёмов обычно бывает значительно меньше плановых размеров, поэтому течение в водоёме можно рассматривать как течение в мелководном бассейне, когда вертикальная составляющая скорости много меньше горизонтальных составляющих.

В результате из уравнений (1.4.1) была получена система, аналогичная уравнениям, описывающим течение в бароклинном океане [61]:

Для решения конкретных задач использовались различные упрощенные варианты системы (1.5.1).

Если скорость течения или ее пространственные изменения невелики, то в уравнениях движения можно пренебречь нелинейными инерционными членами; можно также считать, что отклонения свободной поверхности воды от ее невозмущенного горизонтального положения малы по сравнению с глубиной, можно опустить члены, характеризующие процессы горизонтального обмена импульсом и теплотой. С учетом указанных допущений система (1.5.1) упрощается:

Для системы (1.5.2) в [5] принимались следующие граничные условия:

на твердой части боковой поверхности U n = 0 ;

на жидкой части боковой поверхности, где вода втекает в водоем (водовыпуск), на жидкой части боковой поверхности, где вода вытекает из водоёма (водозабор), на свободной поверхности воды при x3 = Уравнения (1.5.2) при граничных условиях (1.5.3) могут быть численно решены для заданных (из тех или иных соображений) коэффициентов турбулентного обмена или с использованием различных моделей турбулентности.

Коэффициенты вертикального турбулентного обмена в (1.5.2) определяются на основе двухпараметрической модели турбулентности [5]:

где Здесь - скорость диссипации энергии турбулентности; T = E 2 / - коэффициент турбулентной вязкости; K = T + M ; = T + M - суммарный коэффициент вязкости;

c1, c 2,, - эмпирические коэффициенты.

Для системы уравнений (1.5.4) принимались следующие граничные условия:

на свободной поверхности воды на дне где y 0, y 0 - шероховатости водной поверхности и дна соответственно; k, c - некоторые постоянные.

Первое условие (1.5.5) определяет поток энергии, вызванный обрушением поверхностных волн. Связь между диссипацией и энергией турбулентности в условиях (1.5.5) и (1.5.6) вытекает из известного соотношения для [65] при замене в нем масштаба турбулентности соответственно шероховатостью водной поверхности y 0 и шероховатостью дна y 0.

Задача (1.5.2) – (1.5.6) решалась численно с помощью метода дробных шагов, неявной конечно-разностной схемы, методов простой и матричной прогонок [44].

Описанная математическая модель применялась, в частности, для изучения гидротермического режима (полей скоростей и температуры) водохранилища Экибастузской ГРЭС-1 (площадь водного зеркала 18,5 км2, объем воды 860 м3, средняя глубина 4,6 м).

Расчет выполнялся для следующих условий: циркуляционный расход 140 м3/с; температура воды на водовыпуске 36° С; равновесная температура 25,2 °С; суммарный коэффициент теплообмена на водной поверхности 52 Вт/(м2°С); скорость ветра 3 м/с. Стационарная задача решалась методом установления. Полученные результаты показали, что в данном случае распределение температуры по глубине почти однородно; это может быть объяснено относительно небольшой глубиной водоёма и воздействием ветра. Течение в поверхностном слое воды следует в основном направлению ветра, осредненное по глубине течение имеет одну большую циркуляцию с центром в средней части водоема. Придонное течение качественно похоже на среднее по глубине, но в нем наблюдается компенсационное течение, противоположное по направлению поверхностному.

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

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

Система уравнений (1.5.1) записывается в виде Коэффициенты вертикального турбулентного обмена в (1.5.7) определяются на основе двухпараметрической модели турбулентности (1.5.4) – (1.5.6).

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

В работе [50] задача для уравнений (1.5.7) решалась без учёта теплообмена между водоёмом и атмосферой.

1.6. Двумерные в вертикальной плоскости модели При изучении течения воды и процессов теплообмена в относительно узких водоемах задача может быть упрощена путем перехода к двумерной (в вертикальной плоскости) модели с применением операции осреднения по ширине водоема. Такая двумерная модель температурно-стратифицированного течения в водоеме вытянутой формы была построена в [24, 25] на основе осреднения трехмерных уравнений турбулентного течения в мелководном водоеме. При выводе предполагалось также, что кривизна водоема в плане невелика, его поперечное сечение плавно изменяется по длине, коэффициенты турбулентного обмена - скалярные величины и что можно пренебречь продольными диффузионными потоками импульса и теплоты по сравнению с вертикальными. В упрощенном варианте такой модели дополнительно предполагается (как и в случае трехмерной задачи), что нелинейные инерционные члены пренебрежимо малы, а свободная поверхность водоема близка к горизонтальной ( 0 ). Поперечное сечение водоема схематизируется в предположении, что периметр сечения включает донный горизонтальный участок и криволинейные боковые поверхности x 2 = b1 ( x1, x3 ) и x2 = b2 ( x1, x3 ).

При этих предположениях после осреднения соответствующих трехмерных уравнений по координате x 2 и исключения давления задача сводится к решению системы где Здесь U 1, U 3, T и - средние по ширине водоема значения скорости, температуры и плотности.

Напряжение трения боковой поверхности принимается в виде где - коэффициент сопротивления трения.

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

Если коэффициенты турбулентного обмена определяются с привлечением какойлибо модели турбулентности, то система уравнений (1.6.1) соответственно расширяется, добавляются граничные и начальные условия, аналогичные по смыслу условиям, приведенным выше при формулировке трехмерной задачи. В данном случае на входном и выходном сечениях водоема достаточно задать расход (или уровень водной поверхности) и температуру в тех точках, где вода втекает в рассматриваемую область; в качестве динамического условия на дне задается "условие прилипания", которое может рассматриваться как предельный случай приведенного выше более общего "условия скольжения". В [24] приводятся результаты численных расчетов некоторых нестационарных стратифицированных течений, выполненных методом конечных разностей в описанной (а также в несколько уточненной) постановке. В частности, был выполнен расчет неустановившегося движения, возникающего при внезапном удалении вертикальной перегородки, отделяющей холодную воду от теплой; при этом оказалось, что описанная математическая модель хорошо воспроизводит наблюдаемые в физическом эксперименте особенности явления.

В ряде работ переход от трехмерной задачи к плоской (в вертикальной плоскости) осуществляется без выполнения операции осреднения по координате x 2, т.е. предполагается, что все характеристики потока несущественно зависят от x 2 и, следовательно, можно пренебречь их производными по этой координате по сравнению с другими слагаемыми исходных уравнений. Так, в [49] для описания течения стационарного стратифицированного потока использовалась следующая двумерная модель:

Для замыкания уравнений движения используется уравнение баланса турбулентной энергии, принимаемое в форме где кроме ранее принятых обозначений k ET - коэффициент турбулентной диффузии энергии Е, t = cµ E LT, LT - линейный масштаб, - температурный коэффициент объемного расширения; c и c1 - эмпирические постоянные ( c = 0.2, c1 = 3.93 ).

На основе анализа экспериментальных данных были приняты следующие функциональные зависимости для коэффициентов обмена:

где Re T, Re лок - турбулентное и локальное числа Рейнольдса соответственно, Pr = M число Прандтля.

Система уравнений (1.6.2), (1.6.3) решалась методом конечных разностей при следующих граничных условиях:

В ряде работ плоские задачи теплопереноса (или аналогичные им задачи массопереноса) в турбулентном потоке решаются методом конечных элементов (см., например, [2, 3] и др.). Так, в [3] рассматривается двумерное нестационарное уравнение тепломассопереноса вида при соответствующих граничных и начальных условиях.

1.7. Плановые модели Во многих случаях при изучении гидротермических процессов в мелководных водомах используется плановая (двумерная в горизонтальной плоскости) модель гидродинамики и теплопереноса, т.е. рассматриваются осредненные по глубине водоема характеристики [34, 53, 55, 62, 81]. Выводу уравнений плановой задачи, их анализу, методам решения посвящена обширная литература (см. например, [30]).

Уравнения плановой задачи гидродинамики (так называемые уравнения мелкой воды), получаемые осреднением по вертикали соответствующих трехмерных уравнений, могут быть записаны в виде Здесь где si и bi - составляющие напряжения трения на поверхности воды и на дне p a - давление у поверхности; ( x1, x 2 ) - заданное распределение плотности воды; C - коэффициент Шези; - коэффициент ветрового напряжения; W - скорость ветра; - угол между осью x1 и направлением ветра, ось x3 направлена вверх.

Первое уравнение в (1.7.1) представляет собой уравнение сохранения импульса, второе - уравнение неразрывности. В выражениях (1.7.3) первые слагаемые представляют силы трения на свободной поверхности (обусловленные действием ветра), вторые - касательные напряжения на дне водоема (слагаемые, учитывающие силы Кориолиса, здесь и далее опущены).

Другая распространенная форма записи уравнений мелкой воды (силы вязкости здесь не учитываются) где U 1 ( x1, x 2, t ), U 2 ( x1, x 2, t ) - горизонтальные компоненты скорости вдоль оси x1 и x 2 соответственно, осредненные по глубине водоема.

В уравнениях (1.7.1) требуется найти неизвестные функции q1, q 2 и. В начальный Пренебрегая соответствующими слагаемыми в уравнениях (1.7.1) и (1.7.4), получаем различные частные случаи плановой задачи (стационарную, без учета нелинейных инерционных членов и т.п.).

Вопрос о постановке граничных условий для нелинейной системы двумерных уравнений гидродинамики и различных ее упрощенных вариантов неоднократно обсуждался в литературе. Было, в частности, показано, что для линеаризованной системы (1.7.4) при сохранении конвективных членов количество условий, которые нужно задавать на границе, зависит от знака нормальной к границе составляющей U n. На участках границы, где жидкость вытекает из области ( U n > 0 ) следует задавать одно условие, а на участках, где втекает ( U n < 0 ) - два. Если конвективные члены отбросить, то на всех участках границы следует задавать по одному условию.

При решении системы (1.7.1) обычно используют следующие граничные условия [56]. На береговой границе Г Б полагают q n = 0. Если на границе Г Б задан расход жидкости, например, на водовыпуске или водозаборе ТЭС, то q n = q n. Кроме того, на водовыпуске можно задать условие q s = 0, где q s - касательная составляющая расхода к граничному контуру. На жидкой границе Г Ж (например, условной границе, отделяющей море от океана) в общем случае (т.е. при учете вязкости) необходимо задать нормальные и касательные силы. Если слагаемыми, учитывающими вязкость, пренебречь, то граничные условия сводятся к следующим:

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

Уравнение плановой задачи теплопереноса относительно осредненной по глубине водоема температуры T (t, x1, x 2 ), получаемое из трехмерного уравнения (вывод приводится, например, в [90]), может быть записано в виде где Fs - плотность теплового потока на свободной поверхности; a - суммарный коэффициент температуропроводности, учитывающий здесь помимо процессов молекулярного и турбулентного теплообмена также гидродинамическую дисперсию (конвективную диффузию), обусловленную вертикальными градиентами осредненной скорости и температуры.

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

При этом следует иметь в виду, что Fs зависит от температуры поверхности, в то время как в другие слагаемые уравнения (1.7.5) входит средняя по глубине температура.

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

Для уравнения (1.7.5) задается начальное условие и граничные условия, например, следующего вида:

где T Г - температура на внешней стороне границы; k Г - коэффициент теплообмена через границу.

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

Уравнение (1.7.5) было использовано для моделирования гидротермического режима ряда водоемов-охладителей [81]. Рассматривался стационарный случай. Численная модель строилась на основе метода конечных элементов в формулировке Галеркина с использованием двумерных треугольных симплекс-элементов.

Для проверки эффективности численной реализации и программы была рассмотрена следующая модельная задача [56]. Для прямоугольной области численным и аналитическим путем решалось уравнение с постоянными коэффициентами при условиях первого рода на границе. Среднеквадратичное отклонение при 45 узлах и 64 элементах составило 1%. Для более редкой сетки несколько выше. Например, при использовании 15 узлов и 16 элементов - 3,4%. Проведенные численные эксперименты показали, что для рассмотренной задачи влияние изменения в довольно широких пределах густоты конечно-элементной сетки на результаты численного решения.

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

1.8. Одномерные модели Для оценки различных параметров гидротермического режима в водоемах и водотоках часто применяются (ввиду их относительной простоты) одномерные математические модели гидродинамики и теплопереноса.

Одномерные модели могут быть получены из трехмерных уравнений путем осреднения соответствующих величин по двум координатам (глубине и ширине или длине и ширине водоема) или непосредственно из рассмотрения одномерной задачи (т.е. в предположении независимости процессов от двух других координат).

Так, для описания движения жидкости в вытянутом водоеме используется хорошо известная система одномерных уравнений сохранения импульса и неразрывности (одномерные уравнения Сен-Венана), записанная в виде [33]:

где ( x, t ) - возвышение свободной поверхности; b( x ), S ( x) - ширина и площадь поперечного сечения водоема соответственно; Q ( x, t ) - расход; K = C R S - модуль расхода;

R - гидравлический радиус; C - коэффициент Шези.

Одномерные модели теплопереноса широко используются для изучения распределения температуры по длине или глубине водоёмов. Так, для вытянутого мелководного водоёма может быть использовано одномерное уравнение теплопереноса:

где T ( x1, t ) - осреднённая по площади поперечного сечения водоема температура; FV и FS - соответствующим образом осредненные объемные и поверхностные тепловые потоки.

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

Одномерные модели гидродинамики и теплопереноса (при соответствующих граничных и начальных условиях) могут рассматриваться совместно (система одномерных уравнений гидротермики) или последовательно, без привлечения или с привлечением различных моделей турбулентности и т.д. Важное методическое и практическое значение имеют различные одномерные модели термоклина (т.е. модели, получаемые в предположении горизонтальной однородности теплофизических процессов в водоемах). Существует много различных одномерных моделей термоклина. В одних моделях используются постоянные коэффициенты турбулентного обмена [54, 77, 93]. В других моделях коэффициенты турбулентного обмена связываются со стратификацией с помощью числа Ричардсона [100, 58]; при этом, например, в [41] коэффициенты турбулентного обмена определяются на основе двухпараметрической модели турбулентности. Одномерная модель такого типа описывается системой уравнений где К уравнениям (1.8.2) добавляются граничные условия:

на свободной поверхности на дне Для замыкания системы (1.8.2) необходимо определить дополнительно коэффициент турбулентной температуропроводности K T или отношение KT K z.

Одномерная модель для вертикально стратифицированного водоема переменной площади сечения (без применения моделей турбулентности) использовалась во многих работах. Так, в [101] нестационарное уравнение теплопереноса с учетом конвективной составляющей записывается в форме где S, B - площадь сечения, нормального к вертикальной оси, и ширина водоема соответственно; U вх, U вых и Tвх, Tвых - горизонтальные скорости и температуры втекающего в водоем и вытекающего из него потоков; FR - плотность потока поглощенной солнечной радиации; a - суммарный коэффициент температуропроводности (принимаемый постоянным по глубине). Второе слагаемое в правой части учитывает горизонтальный обмен, третье - радиацию.

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

на дне водоема на свободной поверхности где a - коэффициент теплообмена с атмосферой; T p - равновесная температура.

Выводы Рассмотрены основные особенности стратифицированных озёр. Подробно описан объект исследования – солёное озеро Шира. Отмечены основные процессы, происходящие в озёрах. Для моделирования турбулентных течений применяются уравнения Рейнольдса.

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

Глава 2. Комплекс математических моделей и численных алгоритмов для исследования гидрофизических процессов в стратифицированных водоемах 2.1. Одномерная модель Согласно натурным данным в глубоководной части оз. Шира, параметры водной экосистемы в основном изменяются только в вертикальном направлении. Поэтому на первом этапе предлагается для математического моделирования вертикальной структуры водоёма использовать одномерную модель.

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

Формирование температурного режима в непроточном стратифицированном водоеме осуществляется вследствие перемешивания верхнего слоя и теплообмена с атмосферой.

Предполагается, что характеристики воды достаточно однородны в горизонтальном плане водоема. Основные процессы, влияющие на распределение температуры – это турбуT лентное перемешивание, моделируемое членом ( КT ) ; нагрев воды за счёт солнечz z ной радиации, описываемой членом Fтепл ; изменение температуры воды, связанное с теплообменом на водной поверхности. Остальные процессы и теплообмен с дном являются не существенными.

Задача для температуры в одномерном приближении формулируется в виде [92]:

с граничными условиями Здесь T – температура воды, KT(z) – коэффициент вертикального турбулентного обмена, Fn - полный тепловой поток через свободную поверхность, FI – приходящая коротковолновая радиация, – коэффициент поглощения излучения, - параметр, определяющий часть коротковолновой радиации, проникающей на глубину ( 0 1 ), cp – удельная теплоемкость воды, 0 – характерное значение плотности воды, H – глубина водоема.

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

где Fэф – эффективное длинноволновое излучение, Fис – теплоотдача испарением, Fт – конвективный теплообмен, FI – приходящая коротковолновая радиация.

Известны различные формулы для вычисления составляющих теплового потока через свободную поверхность. В данной работе использовались следующие с соотношения. Коротковолновая радиация вычисляется по формуле где n =1.11-1.23 в зависимости от влагосодержания атмосферы, N0 – балл общей облачности в долях единицы, hc – высота солнца в градусах, – плотность воздуха, =0.94, k – широта местности в градусах, t=0,1,…,23 – местное астрономическое время; tn=12 – полуденное местное время; 1 – склонение солнца, d – порядковый номер суток с начала года).

Длинноволновое излучение турбулентный обмен между водной поверхностью и атмосферой поток тепла, обусловленный испарением здесь TS – температура поверхности воды (0С), Tа – температура воздуха (0С), N0 – облачность (доли единицы), eS – давление насыщенного пара при данной температуре водной поверхности (Мб), eа – давление насыщенного пара в атмосфере, измеренное на одной высоте с температурой (Мб), W2 – скорость ветра на высоте 2 м (м/сек), - относительная влажность воздуха.

Напряжение ветра определяется по формуле Давтян [83]:

где a – плотность воздуха, W 2 = (wx, wy ) – вектор скорости ветра на высоте 2 м (м/с).

Аналогично задача ставится для определения вертикального распределения солености:

Здесь S – соленость воды, KS(z) – коэффициент вертикального турбулентного обмена для солености, FSH – массообмен с дном, FS – поток через свободную поверхность. В озере Шира мы принимаем для солёности массообмен с дном и атмосферой равными нулю.

Для решения нестационарных задач необходимо задать начальные распределения температуры и солености - T (0, z ) = T 0 ( z ), S (0, z ) = S 0 ( z ). При проведении численных экспериментов для этого использовались натурные данные.