Создание вычислительных технологий для расчета ветровых и ударно-волновых воздействий на конструкции Вальгер Светлана Алексеевна

Содержание к диссертации

Введение

Глава 1. Математические модели и методы приближенного решения задач аэродинамики зданий 33

1.1. Математические модели аэродинамики 33

1.1.1. Модель сплошной среды и подходы к описанию движения сплошной среды 33

1.1.2. Уравнения движения сплошной среды 34

1.1.3. Турбулентные течения. Численные модели турбулентности 38

1.2. Граничные условия для задач моделировании атмосферных течений 51

1.2.1. Структура атмосферного пограничного слоя 51

1.2.2. Моделирование турбулентных характеристик в ПСА 58

1.2.3. Формирование профилей входных граничных условий 60

1.2.4. Характерные размеры расчетной области 61

1.3. Приближенные методы решения уравнений навье-стокса и уравнений эйлера 64

1.3.1. Метод конечных разностей 64

1.3.2. Метод конечных объемов 66

1.3.3. Численные методы решения уравнений Навье-Стокса несжимаемой среды. 71

1.3.4. Методы моделирования сжимаемых течений с высокими градиентами 73

1.4. Методы пространственной дискретизации расчетных областей 80

1.4.1. Классификация расчетных сеток 80

1.4.2. Методы построения расчетных сеток 83

1.4.3. Адаптация расчетных сеток 84

1.5. Программные средства реализации вычислительных алгоритмов 94

1.5.1. Технология обработки геодезических данных для построения крупномасштабной модели сложного рельефа местности 99

1.5.2. Набор подпрограмм адаптации расчетных сеток 102

Глава 2. Валидация и верификация методик численного моделирования 103

2.1. Тестирование алгоритма адаптации расчетной сетки 103

2.1.1. Решение 1D линейного и нелинейного нестационарных уравнений переноса 103

2.1.2. Решение 1D нестационарных уравнений Эйлера 106

2.1.3. Решение 2D нестационарных уравнений Эйлера 109

2.2. Расчет ударно-волнового воздействия на призму прямоугольного сечения 111

.2.1. Постановка задачи 112

2.2.2.Обсуждение результатов 116

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

2.3.1.Моделирование 2D течений в окрестности цилиндра 121

2.3.2.Моделирование 2D турбулентных течений в окрестности двух квадратных призм 126

2.4. 3D расчет обтекания потоком воздуха призмы с квадратным сечением 128

Глава 3. Моделирование течений воздуха в окрестности комплексов зданий и над местностью со сложной орографией 136

3.1. 3D расчет ударно-волнового воздействия на комплекс призм, имитирующих городскую застройку 136

3.2. Расчет ветрового воздействия на здание сложной формы 143

3.3. Расчет течения воздуха в приземном слое атмосферы над поверхностью со сложным рельефом 154

3.4. Расчет атмосферного течения воздуха в окрестности городского микрорайона 163

Заключение 171

Литература

• Турбулентные течения. Численные модели турбулентности
• Приближенные методы решения уравнений навье-стокса и уравнений эйлера
• Решение 2D нестационарных уравнений Эйлера
• Расчет ветрового воздействия на здание сложной формы

Введение к работе

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

Степень разработанности темы исследования. На основе выполненного обзора литературы можно выделить три основных подхода к оценке ветровых и ударно-волновых воздействий на здания и сооружения:

1. Проведение физических экспериментов и натурных испытаний. Экспериментальным исследованиям ветровых нагрузок на здания с использованием аэродинамических установок в рамках теории подобия посвящены работы Н. А. Попова, С. Д. Саленко, С. М. Горлина, С. В. Николаева, В. И. Терехова, А. И. Гныри, М. А. Березина, I. P. Castro, K. Mehta, T. Ishihara, K. Hibi и др. Экспериментальные исследования воздействия взрыва на строительные конструкции проведены в работах С. С. Набатова, В. А. Могилева, М. И. Ганопольского, А. С. Еременко, P. Smith, T. Rose и др. Эксперименты по изучению ветровых и ударно-волновых течений позволяют получить обоснованные данные о физике явлений, однако они являются дорогостоящими и существенно ограничены характеристиками экспериментальных установок (например, по способности воспроизводить толщину пограничного слоя, превышающего высоту обтекаемых зданий при оценке ветровых воздействий).

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

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

3. Моделирование, основанное на полных физико-математических моделях.
Математическое моделирование задач ветрового и ударно-волнового
воздействия на здания позволяет существенно дополнить физические
эксперименты, получить более детальную информацию о течении и учесть
факторы, которые не могут быть оценены с помощью нормативных мето
дик. Решению модельных задач обтекания плохообтекаемых конфигураций
в режимах, соответствующих нормальным атмосферным условиям, и
расчетам ветровой нагрузки на реальные здания и сооружения посвящены
работы А. М. Белостоцкого, С. И. Дубинского, С. А. Исаева, Е. Ю. Филатова,
С. В. Гувернюка, В. А. Гутникова, T. Stathopoulos, A. Baskaran, B. Blocken,
Y. Tominaga, A. Mochida, A. Kashef и др. Следует отметить, что в большин
стве работ для моделирования турбулентных течений используется прибли
жение логарифмического подслоя, что не позволяет описать отрывные зоны
и, следовательно, предсказать распределение коэффициента давления.
Отсутствуют результаты оценки влияния сложной орографии местности,
интерференционных эффектов, структуры отрывных течений на характер
ветровой нагрузки в плотной городской застройке.

Численному моделированию взрывных процессов посвящены работы

4. Мэйдера, Я. Б. Зельдовича, В. Н. Зубарева, Г. Т. Володина, В. Н. Охитина,
А. М. Ременникова и др. Сложная ударно-волновая структура течений,
инициированных детонацией зарядов конденсированных взрывчатых
веществ (далее - ВВ), требует качественного разрешения течения на разры
вах и в областях высоких градиентов, для чего используется адаптация
расчетной сетки. Разработке методов динамической адаптации расчетных
сеток посвящены работы В. Д. Лисейкина, Г. С. Хакимзянова, M. J. Berger,
J. F. Daunenhofer, W. Cao, H. D. Ceniceros и др.

Учитывая актуальность проблемы, а также опираясь на анализ имеющихся подходов к ее решению, целью данной работы является разработка комплексных вычислительных технологий для расчета внешних статических и динамических нагрузок на строительные объекты, обусловленных воздействиями воздушного потока при режимах обтекания в широком диапазоне амплитуд и временных характеристик нагрузок с учетом орографического фактора и интерференционных эффектов. Исходя из поставленной цели исследования, в работе решены следующие задачи: 1. Разработана и реализована вычислительная технология, основанная на полных физико-математических моделях механики сплошной среды, позволяющая применить комплекс инженерного анализа ANSYS Fluent для

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

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

3. Разработан и программно реализован алгоритм адаптации расчетной сетки в задачах с ударными волнами.

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

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

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

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

Методы исследования. В работе используются надежные методы численного моделирования задач механики сплошной среды, основанные на решении полных трехмерных нестационарных уравнений Навье-Стокса для расчета атмосферных течений воздуха. Для моделирования атмосферной турбулентности применяется подход осреднения уравнений Навье-Стокса по Рей-нольдсу, которые замыкаются с помощью двухпараметрических моделей турбулентности. Для численной аппроксимации уравнений используется метод конечных объемов повышенного порядка точности. Для расчета ударно-волновых течений, инициированных детонацией заряда конденсированного ВВ, используются трехмерные нестационарные уравнения Эйлера, для аппроксимации которых применяется метод конечных разностей повышенного порядка точности. Для повышения точности численного решения и увеличения ресурсной эффективности вычислений разработан и реализован эффективный алгоритм динамической адаптации расчетных сеток, основанный на методе эквираспределения. Моделирование проведено с использованием инструментов программных модулей ANSYS Fluent, ANSYS Autodyn, Civil 3D, дополненных набором пользовательских программ для создания физико-5

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

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

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

3. Разработан и использован для расчета течений с ударными волнами оригинальный программный комплекс решения уравнений Эйлера с динамической адаптацией расчетной сетки.

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

Положения, выносимые на защиту, в соответствии с пунктами паспорта специальности: (Пункт 3)

1.Численный алгоритм расчета одно- и двумерных уравнений невязкого сжимаемого газа в преобразованных криволинейных координатах с учетом динамической адаптации расчетной сетки на основе эффективного метода эквираспределения.

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

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

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

(Пункт 4)

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

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

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

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

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

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

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

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

Практическая ценность результатов исследования. Разработанные автором численные методики и алгоритмы могут быть использованы в проектной деятельности для уточнения данных о внешних аэродинамических нагрузках на строительные объекты; при разработке защитных мероприятий, нацеленных на повышение устойчивости сооружений при воздействии взрывов; в учебном процессе для подготовки специалистов в области проектирования зданий и сооружений, МЧС, планирования городской среды. Исследования были поддержаны грантами ФЦП «Научные и научно-педагогические кадры инновационной России» (2009-2013) № 01201275627 и № 01201275629; РФФИ № 115020610110; проектом в рамках задания Минобрнауки РФ № 2014/140 на выполнение государственных работ в сфере научной деятельности № 114091140003; проектом Аналитической ведомственной целевой программы «Развитие научного потенциала высшей школы» (2009–2011) № 2.1.1/11316; совместным грантом Минобрнауки РФ и DAAD (Германия) на 2013/2014 г № 11.9152.2014.

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

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

Апробация работы. Основные результаты работы представлены на 7 научных семинарах: кафедры информатики и прикладной математики НИУ МГСУ (рук. д.т.н., чл.-корр. РААСН А. М. Белостоцкий), «Информационно-вычислительные технологии» ИВТ СО РАН (рук. акад. РАН Ю. И. Шокин и д.ф.-м.н. В. М. Ковеня), «Механика вязкой жидкости и турбулентность» ИТПМ СО РАН (рук. д.ф.-м.н. В. В. Козлов), отдела термогазодинамики ИТФ СО РАН (рук. д.т.н. В. И. Терехов), лаборатории волновых процессов в ультрадисперсных средах ИТПМ СО РАН (рук. д.ф.-м.н. А. В. Федоров), научно-техническом Совете по физике атмосферы, океана и охраны окружающей среды ИВМиМГ СО РАН (рук. д.ф.-м.н. В. В. Пененко и д.ф.-м.н. В. И. Кузин), а также на всероссийских и международных научных конференциях: XI Всероссийский съезд по фундаментальный проблемам теоретической и прикладной механики (Казань, 2015), международная конференция «Потоки и структуры в жидкостях» (Калининград, 2015), Всероссийский семинар «Динамика Многофазных Сред» (Новосибирск, 2013, 2015), 6^th European Conference on Computational Fluid Dynamics (Барселона, 2014), 17-я Международная конференция по методам аэрофизических исследований ICMAR’14(Новосибирск, 2014), XII международная конференция «Забабахинские научные чтения» (Снежинск, 2014), Всероссийская научно-техническая конференция «Актуальные вопросы строительства» (Новосибирск, 2012–2015), «Вычислительные и информационные технологии в науке, технике и образовании –2013» (Усть-Каменогорск, 2013), Taiwan–Russia Bilateral Symposium on Civil Engineering (Тайпей, 2012), Всероссийский семинар с международным участием по струйным, отрывным и нестационарным течениям (Томск, 2012; Новосибирск, 2015) и др.

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

Публикации. Основные результаты по теме диссертации изложены в 19 печатных изданиях, из которых 4 статьи опубликованы в журналах, рекомендованных ВАК [1 - 4], 1 программа – в свидетельствах государственной регистрации программ для ЭВМ [5].

Структура и объём диссертации. Диссертация состоит из введения, трех глав, заключения и семи приложений. Общий объем диссертации составляет 220 печатных листов, включая 81 рисунок и 10 таблиц. Список цитированной литературы включает 251 наименование.

Турбулентные течения. Численные модели турбулентности

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

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

При описании движения сплошной среды в общем случае термодинамические параметры среды являются функциями координат х, у, z и времени t. Для однозначного описания движения сплошной среды может быть использованы Лагранжев или Эйлеров подходы. Подход Лагранжа, рассматривает термодинамические величины и скорость газа U для каждой частицы среды как функции времени t [83]. Независимыми переменными будут также признаки частицы, позволяющие единственным образом отличить её от остальных частиц (переменные Лагранжа). В качестве переменных Лагранжа можно выбрать координаты начального положения частицы в пространстве. Тогда движение частицы может быть описано: x g xlxlxltH = 1,2,3 (1.1) где xt - текущие координаты частицы, а х - координаты частицы в момент времени t = 0.

Альтернативным подходом к описанию движения сплошной среды является подход Эйлера, согласно которому параметры состояния сплошной среды в каждый момент времени t являются функциями от эйлеровых координат точки х1,х2,х3 в неподвижной системе координат. В подходе Эйлера U = U(x1,x2,x3,t) соответствует скорости частицы, находящейся в момент времени t в точке с координатами (х1, х2, х3).

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

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

Согласно закону сохранения массы скорость изменения массы для объема вещества V равна потоку массы, проходящему через поверхность S, ограничивающую объем V: где п - единичный вектор внешней к поверхности нормали, аv(u,v,w), где u,v,w -компоненты скорости по направлениям х, у, z соответственно. В соответствии с теоремой Гаусса-Остроградского [143], поток векторного поля через замкнутую поверхность S может быть выражен интегралом от дивергенции этого поля по объему V, ограниченному поверхностью S и выражение (1.2), следовательно, может быть преобразовано к следующему виду:

В случае, когда скорости движения среды существенно меньше скорости звука, эффектами сжимаемости среды можно пренебречь и принять p=const. В этом случае (1.5) приобретет вид:

Для расчетов течений газа, в которых силы вязкого трения вносят существенный вклад в структуру течения, поверхностные напряжения могут возникать в отличном от нормального направлении. В этом случае тензор напряжений можно представить как т = -pI+т, где -pI - нормальное давление на объем V, а т - тензор вязких напряжений. Тогда из (1.7) можно получить уравнения Навье-Стокса, описывающие движение вязкой несжимаемой жидкости:

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

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

Уравнения газовой динамики (1.5), (1.20-1.22) и (1.28), описывающие течение сжимаемого газа с изменением теплового состояния среды, должны быть дополнены уравнением состояния среды, которое определяется методами статистической физики [35] или по экспериментальным данным.

Для условий, далеких от фазовых переходов, среда описывается уравнением состояния идеального газа: P = P-R, (1.30) где R - универсальная газовая постоянная, для воздуха R = 287 дж/(кг-град). Используя соотношение R = Cp-Cv, (Ср, Су - удельные теплоемкости при постоянном давлении и при постоянном объеме), уравнение состояния можно записать в виде/?=(у–1)/ е, где у - показатель адиабаты Пуассона, у =C/CV =1.4.

Характеристики турбулентности. Одной из главных особенностей атмосферного движения воздуха является его неупорядоченный турбулентный характер, вследствие чего возникают случайные флуктуации скорости и других величин, описывающих атмосферную среду. Турбулентность атмосферных течений носит стохастический характер, что затрудняет описание её структуры аналитическими, численными и экспериментальными методами. Общие закономерности турбулентных течений описаны в фундаментальных работах [77,78, 80,144-146,147].

Под турбулентными течениями в общем смысле подразумевают нестационарные течения, гидродинамические характеристики которых испытывают хаотические изменения в пространстве и во времени. Турбулентность может быть описана в терминах вихревых структур как совокупность разномасштабных вихрей (рис. 1.1.,а-б). Крупные вихри с размерами, сопоставимыми с характерным линейным масштабом задачи L, в силу своей упорядоченности носят название когерентных структур. Когерентные вихри поглощают энергию основного потока. Их структура существенно зависит от типа течения. В результате существенного преобладания инерционных сил над вязкими происходит перемешивание объемов жидкости, движущихся с разными скоростями. Крупномасштабные вихревые структуры разрушаются, передавая энергию вихрям меньшего масштаба. Кинетическая энергия мелких вихрей диссипирует в тепло. Механизм передачи энергии от более крупных вихревых структур к мелким с последующей их диссипацией в тепло, впервые предложенный Л. Ричардсоном [148] и получивший развитие в работах Колмогорова [77, 144], носит название каскадного переноса энергии. Мелкомасштабная турбулентность слабо зависит от индивидуальных особенностей течения жидкости, т.е. является изотропной [8,9]. Линейный (tjk), скоростной (Uk) и временной (tk) масштабы наиболее мелких вихревых структур в соответствии с теорией Колмогорова [77, 144] определяются двумя параметрами: средней скоростью диссипации энергии є и коэффициентом вязкости//:

Приближенные методы решения уравнений навье-стокса и уравнений эйлера

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

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

Основной характеристикой турбулентного атмосферного течения являются флуктуации скорости ветра, т.е. мгновенные хаотичные отклонения скоростей ветра от средней величины. Использование RANS/URANS подхода, в рамках которого турбулентные пульсации оцениваются через осредненные параметры, требует адекватного задания характеристик турбулентности на границах расчетной области, имитирующей ПСА.

Один из первых подходов к описанию вертикального распределения турбулентных параметров в атмосферном пограничном слое в условиях нейтральной термической стратификации [207] основан на допущении постоянства сдвиговых напряжений по высоте:

Подход (1.85) - (1.86) позволяет явно задать уровень ТКЕ и скорость диссипации турбулентных пульсаций є в качестве параметров входного граничного условия. Однако натурные измерения интенсивности турбулентных пульсаций показывают, что допущение о постоянной по высоте ТКЕ не всегда применимо [205].

Альтернативный подход [103] основан на определении градиента интенсивности турбулентных пульсаций в ПСА из соотношения:

В (1.88) сделано предположение о равенстве флуктуационных компонент скорости ol (z) al(z) + a2W (z), что, исходя из данных [208], применимо для однородной поверхности с постоянной шероховатостью. В случае подстилающей поверхности сложного рельефа, можно предположить, что характеристики турбулентности в ПСА достаточно быстро изменяются с расстоянием от входной границы расчетной области. Это позволяет использовать изложенный выше подход в условиях, когда внешние границы расчетной области находятся на существенном расстоянии от исследуемого объекта. Как видно из (1.85-1.89) характерный масштаб турбулентных вихрей в обоих подходах предполагается равным расстоянию до поверхности z и, в конечном счете, определяется толщиной пограничного слоя 5. Это хорошо согласуется с данными [209] о широком диапазоне размеров вихревых структур в турбулентном ПСА, включающем как мелкие вихри с масштабом от 1-2 см вблизи поверхности земли, так и крупные структуры с масштабом до 10-1000 м в слое смешения.

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

Опишем методику оценки входных граничных условий для скорости, ТКЭ и диссипации турбулентных пульсаций с учетом изменения этих величин в пределах пограничного слоя на примере течения около плоской пластины. Моделирование выполнено на основе решения 2D стационарных осредненных по Рейнольдсу уравнений Навье-Стокса, дополненных к-со SST моделью турбулентности. 2D расчетная область представляет собой прямоугольник с размерами 1,6 м1,5 м и сеткой 200120 ячеек. Расстояние от первой расчетной ячейки до стенки пластины составляет 10-4 м.

Профиль скорости набегающего потока получен из экспериментальных данных [52] и с помощью метода наименьших квадратов аппроксимирован функцией u(z) =ucaz027, где z - вертикальная координата по пространству. Профиль ТКЭ набегающего потока k(z) получен из экспериментальных данных [197] с помощью линейной интерполяции. Значения диссипации турбулентных пульсаций вычислены по (1.89). На выходе из расчетной области для относительного давления использовано условие P = Pст -P0=0 атм. На верхней границе расчетной области использовано условие симметрии. На пластине задано условие прилипания для гладкой стенки.

На рисунке 1.8 представлены профили продольной компоненты скорости (а) и кинетической энергии турбулентности (б) во входном сечении (Inlet), и в сечениях x1 и x2, расположенных на расстоянии 5b (линия 1) и 10b (линия 2) от входного сечения. Хорошее соответствие профилей с разных сечениях по максимальным значениям ТКЭ и толщине пограничного слоя позволяет заключить, что граничные условия на входе заданы верно. Однако наполненность профиля средней скорости несколько растет при удалении от входного сечения, что может быть связано с шероховатостью подстилающей поверхности, а также может свидетельствовать о неравновесности пограничного слоя из-за наличия внешних возмущений в эксперименте.

В дальнейших 3D расчетах атмосферных течений в качестве входных граничных условий задавались профили газодинамических и турбулентных параметров с учетом толщины пограничного слоя.

Численное моделирование течений атмосферного воздуха в окрестности зданий и их комплексов осложняется тем, что масштабы расчетных областей, используемых в таких задачах, должны существенно превышать характерный размер обтекаемых препятствий. Размеры расчетной области не должны влиять на структуру течения в исследуемых тел. Практические рекомендации по выбору характерных размеров расчетной области приведены в работах [103,198-201, 202-203], согласно которым есть два возможных алгоритма: 1) определение лимитирующего фактора, характеризующего минимально допустимые расстояния от стенок исследуемого плохообтекаемого тела до внешних границ расчетной области; 2) определение лимитирующего отношения BR площади проекции наветренной стороны зданий на поперечное сечение расчетной области к площади этого поперечного сечения. В [103] приведены общие рекомендации по размерам расчетных областей для моделирования атмосферных течений в окрестности зданий, согласно которым: 1) расстояния от входной и боковых границ расчетной области до обтекаемого тела с характерным размером Нmax должны быть не менее 5Нmax; 2) расстояние от подветренной стороны обтекаемого тела до выходной границы расчетной области должно составлять как минимум 15Нmax; 3) высота расчетной области предполагается равной или превышающей 6Нmax. 4) отношение BR рекомендуется выбирать не превышающим 3%.

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

Решение 2D нестационарных уравнений Эйлера

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

Построение расчетной геометрии. Расчетная геометрия представляет собой некоторое математическое описание топологии и пространственного расположения узлов, ребер и граней геометрической модели. В настоящем исследовании расчетные области характеризуются большими размерами и сложной топологией. В некоторых случаях они формируются на основе данных лазерного сканирования и топографической съемки. Учитывая данные особенности, для построения расчетных областей требуется применение точных эффективных и относительно универсальных методов описания граней и ребер расчетной модели. В настоящей работе в качестве геометрического препроцессора использованы технологии CAD-системы Civil 3D [243], удовлетворяющие перечисленным выше требованиям. Civil 3D является системой автоматизированного 3D проектирования и поддерживает большинство форматов, используемых в работе проектировщика. Это позволяет формировать геометрическую модель для дальнейших расчетов на основе проектной документации. Следует также отметить, что численные расчеты задач внешней аэродинамики требуют использования технологий работы с большими объемами данных, например данными о рельефе местности. В настоящем исследовании разработана методика обработки данных геодезической съемки сложных рельефов местности для построения 3D расчетных моделей на основе инструментария Civil 3D, которая будет описана в п. 1.5.2.

Пространственная дискретизация геометрической модели. Построение расчетной сетки конечных объемов для объема несжимаемого вязкого газа проводится с использованием инструментария ПК ANSYS Meshing. Модуль Meshing позволяет автоматизировать процесс построения расчетной сетки структурированного и неструктурированного типов на основе гекса- и тетраэлементов с помощью методов, описанных в п. 1.4.

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

Построение физико-математической модели и методы решения. В качестве основного средства программной реализации исследуемых алгоритмов и методов в настоящей работе использован программный комплекс ANSYS [235]. ПК ANSYS – это многоцелевой пакет программ для численного моделирования физических процессов и явлений в области прочности, теплофизики, динамики жидкостей и газов, электромагнетизма. Кроме того, программный комплекс имеет широкие возможности для выполнения междисциплинарного анализа. Преимуществом программного комплекса является также и предоставление опций пользовательского программирования, что делает программный продукт актуальным для проведения исследовательских проектов. Для численного анализа течений жидкости и газа ANSYS включает CFD решатели, которые позволяют исследовать стационарные и нестационарные течения, сжимаемые и несжимаемые течения, невязкие, ламинарные и турбулентные течения, многокомпонентные и многофазные течения, течения с химическими реакциями, течения через пористые среды и т.д. Каждый решатель включает набор физико-математических моделей течения жидкости и газа, а также обширный набор численных методов, построенных на основе конечных объемов или конечных разностей. Также решатель управляет процессами выделения памяти для переменных и алгоритмом распараллеливания процесса решения.

В настоящей работе для расчета вязких несжимаемых течений воздуха в окрестности плохообтекаемых тел используется модуль ANSYS Fluent [236], программный код которого реализован на языке С [237]. В качестве основной математической модели используется описанная в п. 1.1 3D математическая модель, основанная на системе дифференциальных уравнений в частных производных и выражающая законы сохранения массы и импульсов. Управляющие уравнения в осредненной по Рейнольдсу постановке, дополненные полуэмпирическими моделями турбулентности, приведенными в п. 1.1.4, решаются с помощью численных методов, описанных в п. 1.3.3.

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

ANSYS Fluent включает возможности расширения функционала за счет пользовательских подпрограмм и функций на языке С, использующих стандартные библиотеки решателя. Загрузка подпрограмм в настоящей работе выполнена с помощью внутреннего интерпретатора решателя. Пример подпрограммы для формирования профилей скорости, ТКЭ и скорости диссипации ТКЭ на основе интерполяции экспериментальных данных кубическим сплайном приведен в Приложении A.

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

Численное решение 3D нестационарных уравнений Эйлера, описывающих течения невязкого сжимаемого газа, в настоящей работе выполнено с использованием ПК ANSYS AUTODYN. AUTODYN включает набор решателей для моделирования нелинейной динамики твердых тел, жидкостей, газов, а так же их взаимодействия. AUTODYN позволяет моделировать крайне сложные быстропротекающие процессы механики сплошной среды, такие как распространение ударных волн в среде. Для решения гидродинамических задач рассмотрен Эйлеров решатель, основанный на методе конечных разностей, хорошо зарекомендовавший себя в задачах с высокими градиентами течения и разрывами. Эйлеров решатель AUTODYN включает два метода повышенных порядков точности для расчета уравнений Эйлера: метод коррекции потоков и метод Годунова, описанные в п. 1.3.4. Существенным недостатком решателя является отсутствие функционала для адаптации расчетной сетки к высоким градиентам и разрывам решения, вследствие чего в настоящей работе предложены описанные в п. 1.5.2 программные реализации собственных алгоритмов адаптации применительно к решению 1D и 2D уравнений Эйлера.

Расчет ветрового воздействия на здание сложной формы

В главе 2 проведено численное моделирование ударно-волнового течения в окрестности отдельно стоящей призмы. По результатам расчета с помощью схем Годунова второго порядка точности и коррекции потоков FCT получена качественная структура течения: нерегулярное отражение УВ от подложки, формирование трехволновой конфигурации течения, отражение УВ от стенки призмы и дифракция УВ при обтекании призмы. Также проведено сравнение результатов расчетов с данными эмпирической функции CONWEP и приведенными в [245] экспериментальными данными по времени прихода УВ и пиковым давлениям на наветренную и подветренную стенки призмы. Показано, что для случая отдельно стоящей призмы численные расчеты, а также функция CONWEP дают хорошее соответствие с экспериментальными данными как по времени прихода УВ на наветренную стенку тела, так и по максимальному давлению на фронте УВ, приходящей на наветренную сторону тела. При этом получено, что схема коррекции потоков позволяет более точно разрешить фронты УВ по сравнению со схемой Годунова. Схема коррекции потоков также позволяет более точно разрешить вторичные пики волновой нагрузки, связанные с приходом на поверхность призмы вторичных отраженных УВ, в то время как схема Годунова, обладающая существенной схемной вязкостью, слабо разрешает такие особенности течения. Эмпирическая функция CONWEP не учитывает эффекты затенения и многократного переотражения УВ, вследствие чего сделан вывод о её неприменимости для сложных конфигураций и комплексов тел.

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

На рис. 3.1, а показана схема расчетной области, включающая 7 призм различного сечения и высоты. Максимальная высота призм (призмы А и G) составляет Hmax=0.45 м, в то время как наиболее низкие объекты D и F имеют высоту Hmin=0.3 м. Расстояние между соседними объектами меньше или сопоставимо по линейному масштабу с характерными поперечными размерами обтекаемых объектов, что характеризует данную конфигурацию застройки как плотную. Расчетная область (рис. 3.1, б) представляет собой объем воздуха с начальным распределением параметров среды, соответствующих нормальным атмосферным условиям. В отличие от тестовой задачи, описанной в п. 2.2 главы 2, моделирование ударно-волнового течения в окрестности данной конфигурации тел не может быть проведено в симметричной постановке и требует построения полной 3D расчетной области.

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

В начальный момент времени t0=0 сек в точке С (0.478 м;0.35 м;0.04 м) происходит детонация в ВВ массой MTNT = 0.016 кг тротилового эквивалента. Начальный процесс детонации ВВ в открытом пространстве описан в предположении осевой симметрии, затем полученные данные о параметрах течения интерполируются на 3D расчетную область. Движение воздушной среды описывается системой уравнений Эйлера, дополненной уравнением состояния идеального газа в соответствии с п. 1.1.2 и п. 1.3.4 настоящей работы. Течение расширяющихся продуктов детонации ВВ описывается с помощью гидродинамической модели материала [34] с использованием уравнения состояния Джона-Вилкинса-Ли (JWL).

Решение уравнений Эйлера выполнено с помощью метода конечных разностей. Интегрирование уравнений по времени осуществлялось с помощью явной схемы с соблюдением условия устойчивости схемы по числу Куранта. Для аппроксимации уравнений в работе использованы метод Годунова второго порядка точности по пространству и метод коррекции потоков (FCT), описанные в п. 1.3.4. Пространственная дискретизация 3D расчетной области включает 6 млн. расчетных гекса-элементов, характерный размер сеточного элемента составляет 1см1см1см.

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

Рассмотрим структуру ударно-волнового течения в окрестности комплекса призм на основе мгновенных полей статического давления на характерные моменты времени (рис. 3.3), а также зависимостей статического давления от времени в характерных точках - точках замеров Т1, Т11 и Т21 (3.4). На рис. 3.3 (аз) показаны мгновенные поля статического давления на характерные моменты времени t в горизонтальном сечении h=0.105 м. Сферическая УВ (1), сформированная в результате высвобождения большого количества энергии при детонации ВВ, распространяется в свободном пространстве и в момент времени t0.1 мсек. приходит на угол призмы F и отражается от него. Максимальное давление во фронте УВ 1 в момент отражения составляет 6 МПа. Через некоторое время УВ (1) приходит на стенки еще двух призм: стенку призмы C и острый угол призмы E. Отражение УВ (1) от стенок призм E и F носит нерегулярный характер, что приводит к формированию ножек Маха и тройных точек 3,4 и 5 (рис. 3.3, б). При отражении УВ 1 от стенки призмы С формируется УВ (2), которая движется по направлению к эпицентру взрыва. Отраженная от стенок призмы F УВ 6 (рис. 3.3, в) в момент времени 0.4 мсек начинает взаимодействовать с отраженной УВ 2 и образует ударно-волновую структуру 7 (рис. 3.3, г). При этом УВ 2 огибает острый угол здания C, взаимодействуя с первичной УВ 1 (рис. 3.3, в, г). В относительно узком "канале" между призмами E и F на момент времени 0.4 мсек также начинают интерферировать структуры 3 и 4, УВ 5 движется от эпицентра взрыва и огибает острую кромку призмы F (рис. 3.3, г). В момент времени t0.5 мсек УВ (1,2), образованная в результате взаимодействия УВ 1 и УВ 2, приходит на стенку призмы D и через некоторое время на стенку призмы B, в результате чего происходит формирование более мелких отраженных УВ, фронт которых движется по направлению к эпицентру взрыва (рис. 3.3, д). УВ 3 и УВ 4 формируют единый фронт (8) и движутся по направлению к стенке призмы G. С другой стороны, УВ 2, отраженная от поверхности призмы C взаимодействует с УВ 5. В эпицентре взрыва возникает сложная ударно-волновая структура, сформированная в результате взаимодействия УВ 7 и УВ, сформировавшихся из тройных точек 3 и 4. В момент времени t=0.9 мсек ударно-волновая структура 1,2 приходит в точку замера Т1 (рис. 3.3, е), максимальное статическое давление, приходящееся на призму B в точке Т1 составляет 260 КПа для расчета по методу коррекции потоков и 200 кПа при расчете со схемой Годунова второго порядка точности (рис.3.4, а). Значения давления, наблюдаемые в эксперименте для точки Т1 в момент времени прихода УВ на стенку, составляют порядка 300 кПа. УВ 8 приходит на стенку призмы G и отражается от нее регулярным образом (рис. 3.3,е). В это же время фронтальная УВ конфигурации 9, сформировавшаяся в результате нерегулярного отражения УВ 2 от стенки призмы B движется по направлению к стенке призмы A. На рис. 3.3., ж показано мгновенное распределение относительного статического давления на момент времени t1.3 мсек. Как видно из рисунка, отражение УВ 8 от стенки призмы G (10) взаимодействует с УВ 8. Часть УВ (10) распространяется в свободном пространстве к эпицентру взрыва, где взаимодействует с более слабой вторичной волной 7, а часть (8) в условиях стесненного пространства между призмами E и G продвигается вдоль их стенок, многократно отражаясь от них. Многократное отражение УВ 8 от стенок G и E обуславливает пики статического давления в точке Т21, расположенной на стенке призмы G в моменты времени t=2 мсек и 2,5 мсек (рис. 3.4, д). При этом за ударной волной 8 в этой зоне формируется зона разрежения с падением давления ниже одной атмосферы (до 43 кПа). В момент времени t=1.7 мсек УВ 9 приходит в точку Т11 на стенке призмы A, что сопровождается увеличением пиковой нагрузки в этой точке до 165 кПа (рис. 3.3,з). В точку Т1 на стенке здания B в этот же момент времени приходит вторичная ударная волна (11). Серия небольших локальных пиков ударно-волновой нагрузки, приходящихся на точку Т1 (рис. 3.4., а) связана с приходом слабых вторичных волн (11), возникших в результате интерференции УВ 7 и отраженных от стенок призм B и D более слабых ударных волн (рис. 3.3, з). На более поздних этапах времени структура течения характеризуется многократными переотражениями, дифракцией и интерференцией ударных волн, что приводит к возникновению локальных пиков нагрузки на стенки призм и ряда локальных зон разрежения за фронтами ударных волн (рис. 3.3, и, к).