ПРОВЕДЕНИЕ НЕЛИНЕЙНЫХ РАСЧЕТОВ ВОЗДУХООПОРНЫХ ОБОЛОЧЕК

Полный текст

Воздухоопорными оболочками (сооружениями) называются мягкие безмоментные оболочки, форма которых поддерживается за счет поддержания избы- точного давления воздуха, непрерывно нагнетаемого под оболочку [5, 6]. В связи с малой толщиной материала, его изгибная жесткость пренебрежимо мала и в оболочке возникают только мембранные усилия. «Пик популярности» пневматических конструкций пришелся на 70-е года XXвека, после всемирной выставки «Экспо-70» в г. Осака, Япония. Последующий энергетический кризис [с. 84, 5], в результате которого резко возросли эксплуатационные расходы воздухоопорных сооружений, на время замедлил их распространение. В последние годы в связи с развитием химической промышленности интерес к воздухоопорным оболочкам, и к текстильной архитектуре вообще, вновь возрастает во всем мире [13]. В нашей стране подобные сооружения возводятся ежегодно, поэтому отечественные исследователи вновь обращаются к этой теме [6, 11]. Подобные конструкции, как правило, применяют- ся в качестве покрытий большепролетных сооружений различного назначения. Конструкции из тканей и тросов обладают намного большей податливостью, чем сооружения из традиционных материалов, поэтому при действии рас- четных нагрузок происходит сложное взаимодействие конструкции с окружающей средой, которое приводит как к изменению формы конструкции, так и к изменению самой действующей нагрузки (иногда и по величине, и по направлению). До недавних пор, основным источником знаний о подобных явлениях служили экспериментальные исследования, которые, как правило, проводились для нужд авиационной и космической промышленности. Разработанные аналитические и полуэмпирические методы решения подобных задач преимущественно применимы для узкого класса задач с простейшей геометрией и рядом ограничений, накладываемых на постановку задачи и граничные условия [2]. С развитием вычислительной техники и расчетных программ появилась возможность численно смоделировать взаимодействие сооружения с окружающей средой (англ. fluid-structuralinteraction) с минимальным количеством упрощений. Проведение подобных связанных расчетов в строительной отрасли является достаточно новым и актуальным направлением [2]. Например, расчет воздухоопорных оболочек с тросовым усилением в связанной постановке пока не был встречен в литературе. Можно выделить два типа расчетов взаимодействия тела со средой (жидкостью или газом): 1) Односторонний расчет (1-wayFSI) - производится расчет абсолютно жесткого тела (первоначальной неизменной формы) в потоке жидкости или газа, а затем полученное на поверхности тела поле давлений прикладывается в виде внешней нагрузки, на действие которой тело рассчитывается как деформируемое. 2) Двусторонний связанный расчет (2-wayFSI)- представляет собой последовательность односторонних расчетов, выполняемых с малым шагом по времени, так чтобы приращения давлений и перемещений точек тела были достаточно малы. Его суть сводится к тому, что на каждой итерации последовательно выполняются два расчета - динамический расчет конструкции и аэродинамический. Результаты динамического расчета передаются в виде перемещений узлов расчетной сетки, после чего производится решение задачи аэродинамики для модели с измененной геометрией. Результаты этого расчета в виде поля давлений прикладываются к конструкции в виде нагрузок и производится новый динамический расчет с учетом этих дополнительных нагрузок, в результате которого уточняются значения перемещений узлов сетки. После этого опять выполняется уточнение поля давлений по результатам аэродинамического расчета и т.д. пока относительные невязки решений не достигнут заранее заданной точности, в нашем случае 0.01. Для автоматизации таких расчетов применяются специальные программ- ные средства. В программном комплексе Ansys имеются инструменты для решения всех поставленных задач: решение задач строительной механики - Ansys- Mechanical, задач гидрогазодинамики -Ansys Fluent и др., а также средство для- проведение связанных расчетов - SystemCoupling. Таким образом, полный дву- сторонний связанный расчет можно произвести в рамках одного программного комплекса, что достаточно просто осуществить с помощью платформы An- sysWorkbench, которая и была использована в данной работе. Рис. 1.Схемапроекта: а) односторонний расчет;б) двусторонний сопряженный расчет Для корректного обмена данными, на двух расчетных областях должны быть назначены совпадающие интерфейсные поверхности. Передача данных производится после автоматической интерполяции, следовательно, узлы расчетных сеток на интерфейсных поверхностях не обязаны совпадать, что позволяет использовать различную дискретизацию интерфейсных поверхностей. Для проведения одностороннего расчета не требуется использование надстройки SystemCoupling, импорт нагрузки выполняется внутри модуля Ansys- Mechanical. Поведение воздухоопорных оболочек в потоке воздуха интересовало мно- гих ученых как в нашей стране, так и за рубежом, поскольку именно ветровая нагрузка является для таких сооружений основной, и нередко приводит к их разрушению [с. 247, 4]. Эксперименты по продувке жестких и мягких моделей в аэродинамических трубах проводились многими отечественными и зарубежными исследователями [1, 3, 4, 8, 9, 14]. Однако, в мелкомасштабных моделях удается добиться только частичного подобия модели реальным условиям. Эксперимент, поставленный В.П. Поляковым, описанный в нескольких статьях [3, 8], «отличается наибольшей широтой» [4], поскольку продувалась крупномасштабная модель усеченной сферической оболочки диаметром 4.2 м и высотой 3.36 м в аэродинамической трубе с площадью поперечного сечения 264 м2.Размеры оболочки сопоставимы с размерами реальных сооружений, а исследуемые скорости ветра (12-45 м/с) - с расчетными скоростями ветра в различных ветровых районах. Помимо больших размеров модели следует отметить, что набор данных для сопоставления представлен достаточно полно: размеры оболочки, параметры материала, направление раскроя, скорость потока и внутреннее давление в раз- личных сочетаниях, перемещения после деформации, аэродинамические коэффициенты, а также наибольшие главные усилия в зоне наибольшего разряжения.Однако, для еще более детального и точного моделирования не хватало де- формационной кривой для материала оболочки(прорезиненной капроновой ткани арт. 51-060 толщиной 0.6 мм), а также геометрических параметров рас- считываемой области - размеров поддерживающей конструкции и положения оболочки в потоке. В связи с этим, в численной модели использованы следующие упрощения - рассмотрена область в виде прямоугольного параллелепипеда с «удаленным» подоболочечным пространством, а материал оболочки принимался ортотропным, но линейно-упругим. Задача аэродинамики решалась в Ansys Fluentс помощью метода конечных объемов (МКО), в котором расчетная область (рис. 2), включающая в себя пространство вокруг оболочки, разбивается на множество локальных элементов (расчетных ячеек). Для каждого локального элемента записывается система за- конов сохранения массы, импульса и энергии в интегральной форме, которая затем преобразуется к системе алгебраических уравнений относительно искомых величин - плотности, скорости, температуры и др. в центрах расчетных ячеек [12]. Интенсивность турбулентности принималась равной 1%, поскольку в аэро-динамических трубах она как правило мала. В качестве модели турбулентности принималась довольно широко распространенная двухпараметрическая модель k-?. Шероховатость поверхности оболочки принималась равной 1 мм. Несмотря на то, что рассматривалось действие установившегося потока, расчет производился как для переходного процесса, поскольку изменение фор- мы оболочки и поиск ею равновесного положения сопровождается затухающими колебаниями, то есть в начальной стадии является нестационарным. Шаги по времени проводились до тех пор, пока оболочка и поток вокруг нее не стабилизировались. Рис. 2. Расчетные модели для решения задачи:а) аэродинамики; б) строительной механики. Определение напряженно-деформированного состояния (НДС) оболочки (задача строительной механики) решалась в модуле Ansys Mechanical методом конечных элементов (МКЭ). Оболочка моделировалась плоскими треугольными конечными элементами (КЭ)типа Shell181 с включенной опцией KEYOPT(1)=1, учитывающей только мембранную группу усилий. Последний параметр играет ключевую роль в достижении сходимости решения задачи. Расчет производился в геометрически нелинейной постановке (активна опция «LargeDeflection») с помощью метода Ньютона-Рафсона. Матрица демпфирования задавалась в виде суммы матрицы масс и матрицы жесткости, умноженных на коэффициенты, равные 0.1. Поскольку производился квазистационарный расчет, величина демпфирования влияла только на скорость затухания колебаний в переходном режиме, поэтому величина коэффициента 0.1 была определена путем тестовых расчетов. Избыточное давление воздуха под оболочкой, необходимое для ее существования, моделировалось в виде постоянной «следящей» нагрузки, т.е. всегда направленной по нормали к оболочке. Корректная оценка влияния избыточного давления на оболочку воздухоопорных сооружений (далее ВОС) является обязательным условием при проведения любых инженерных расчетов данных сооружений.Особенности расчета воздухоопорных оболочек на действие избы- точного давления воздуха, те же, что описаны в докладе [7]. Следует отметить одну дополнительную особенность проводимых расчетов: поскольку толщина оболочки крайне мала по сравнению с ее размерами, то пространственное тело заменялось двухмерной сферической поверхностью, не имеющей толщины, которая также служила интерфейсом для передачи данных между решателями. Однако, в этом случае, в связи с ограничениями надстройки SystemCoupling, отсутствует возможность приложить внутреннее давление в виде нагрузки на внутреннюю поверхность оболочки. Моделирование оболочки объемными КЭ потребовало бы неоправданно сильного сгущения сетки (размер элемента не должен отличаться от толщины оболочки более чем на порядок). В результате, было принято решение моделировать оболочку двумя поверхностями:внешняя является интерфейсной и служит для передачи данных при выполнении связанного расчета, а внутренняя - используется для приложения избы- точного давления. Для совместной работы обеих оболочек между ними устанавливалось контактное взаимодействие, что дополнительно усложняло задачу. Кроме этого, близкое размещение (на расстоянии меньше 1 мм) двух поверхностей выходило за пределы точности системы геометрического моделирования. Поэтому, оболочки размещались на расстоянии 6 мм друг от друга, а толщина каждой из них принималась равной 3 мм. Суммарная толщина оболочки hу величилась в 10 раз, поэтому модули упругости материала Eбыли уменьшены в 10 раз, чтобы значение жесткости при растяжении, равное произведению Eh не изменилось. Как показали расчеты, такая замена привела к незначительным искажениям результатов. На рис. 2,а в статье [3] приведены относительные перемещения оболочки по нормали к поверхности при действии избыточного давления в 150 и 500 Па. Увеличение нагрузки более чем в три раза привело к росту перемещений лишь в полтора раза, что свидетельствует об увеличении модуля деформации с рос- том напряжений, то есть физической нелинейности. Поскольку вид деформационной кривой при дальнейшем росте напряжений в работах [3, 8] не приведен, то было принято решение использовать линейно-упругую ортотропную модель материала.Модули упругости подбирались так, чтобы деформации оболочки в потоке были сопоставимы с результатами эксперимента. Отметим, что перемещения оболочки связаны не столько с ее растяжением, сколько с изменением формы оболочки из-за действия воздушного потока. Соотношение жесткостей по направлению основы и утка было задано таким же, как и в источнике [3], т.е. 13/12. Модуль упругости в направлении основы (меридиональном) принимался равным E1 = 200 МПа, а в направлении утка (окружном) - E2 = 185 МПа.При двусторонних расчетах, эти модули были уменьшены в 10 раз, потому что в 10 раз возрастала толщина материала. Следовательно, полученные в результате расчета напряжения также следовало увеличить в 10 раз. При расчетах воздухоопорных оболочек большое значение имеет отношение избыточного давления воздуха под оболочкой p к динамическому давлению потока воздуха q, определяемому как: , (1) где - плотность воздуха при температуре 15?С; v - скорость потока, м/с. Отношение ветровой нагрузки к динамическому давлению (скоростному напору) называется аэродинамическим коэффициентом и обозначается СР.Отрицательные значения СР соответствуют отсосу. Линия пересечения вертикальной плоскости, проходящей через центр сферы параллельно потоку и поверхности сферы называется главным меридианом. На рис. 3 показано распределение величины СР по плоскости главного меридиана для одностороннего и двустороннего расчетов (без учета (а) и с учетом (б) деформаций оболочки). Рис. 3.Распределение аэродинамического коэффициента в плоскости главного меридиана (ветер слева): а) Односторонний; б) Двусторонний связанный расчет В эксперименте исследовалось распределение аэродинамического коэффициента по главному меридиану при внутреннем давлении p = 1000 Па и скоростях потока v = 12, 20, 30, 40 и 45 м/с. А также при скорости потока v = 40 м/с изменялись значения внутреннего давления: p = 500, 1000, 1500, 2000 Па. Широтный угол ? отсчитывается по часовой стрелке (если смотреть на рис. 3) от горизонтальной оси направленной навстречу потоку и проходящей через центр сферы. Некоторые графики распределения СР показаны на рис. 4-7. Можно заметить, что результаты одностороннего расчета слабо изменяют- ся в зависимости от скорости потока. В результате двусторонних расчетов минимальное значение с ростом скоростного напора СР уменьшается от -1.23 до - 1.38 (в эксперименте от -1.16 до -1.4), а различие между результатами одностороннего и связанного расчета увеличивается. Результаты двустороннего связанного расчета при больших скоростях ветра лучше соответствуют эксперименту, как количественно, так и качественно (заметен второй минимум в районе 100?, который отсутствует на результатах одностороннего расчета). В целом, соответствие эксперименту аэродинамических коэффициентов, полученных в результате двустороннего расчета можно считать удовлетвори- тельным. При испытании оболочки, перемещения точек ей поверхности измерялись методом стереофотограмметрии [3]. В зоне стереоскопического видения была расположена передняя часть оболочки, где наблюдаются наибольшие прогибы оболочки вовнутрь. Рис. 4. Распределение аэродинамического коэффициента вдоль главного меридиана при p = 1000 Па и v = 20 м/с. Рис. 5.Распределение аэродинамического коэффициента вдоль главного меридиана при p = 1000 Па и v = 40 м/с. Рис. 6.Распределение аэродинамического коэффициента вдоль главного меридиана при p = 1000 Па и v = 45 м/с Общий вид изменения формы оболочки в потоке показан на рис. 8. Сравнение деформаций показано на рис. 8, в), где они изображены в натуральную величину. Рис. 7.Распределение аэродинамического коэффициента вдоль главного меридиана при p = 2000 Па и v = 40 м/с. Рис. 8. Деформированная оболочка (поток направлен слева направо): а) Фотография эксперимента [3]; б) Результат двустороннего расчета; в) Сопоставление деформаций при p=1000 Па, v=40 м/с В статье [3] в табл. 2.3 приведены формулы для определения величин глав- ных растягивающих усилий, возникающих в области максимального разрежения (вблизи полюса оболочки) по безмоментной теории оболочек вращения для различных соотношений p/q. Для наглядности, продублируем табл. 2.3 (см. табл. 1). Табл. 1 p/q 0.5 0.5-1 1.0-2.0 2.0-5.0 2N1/R p + 3.65q p + 3.14q p + 3.29q p + 2.61q Для сравнения результатов с Ansys произведем расчет главных напряжений в табличной форме для p = 1000 Па и скоростях потока v = 45, 40, 30 и 20 м/с. Переход от главных натяжений к главным напряжениям осуществляется по формуле: ?1 = N1/h,, (2) где N1 - главное растягивающее усилие (натяжение), определяемое по табл. 1; h - толщина оболочки в расчетной модели, h = 6 мм = 0.006 м. Результаты расчета главных напряжений в зоне наибольшего разряжения (отсоса) в Ansys обозначим ?1-wayи ?2-way соответственно занесем в табл. 2 в отдельные строки.В последней строке сравниваются значения ?2-wayи ?1. Распределение напряжений показано на рис. 9. Нетрудно заметить, что разница расчета максимальных напряжений составляет менее 10%, что свидетельствует о получении достаточно качественных результатов. Напряжения, полученные в ходе одностороннего расчета, также достаточно хорошо согласуются с результатами эксперимента и двустороннего расчета, особенно при меньших скоростях потока, из чего можно сделать вывод, что при стационарных режимах и небольших изменениях формы оболочки односторонний расчет дает решение близкое к двустороннему. Таблица 2 p, Па 1000 1000 1000 1000 v, м/с 45 40 30 20 qпо (1), Па 1244.4 983.2 553.1 245.8 p/q 0.804 1.017 1.808 4.068 N1, Н/м 5152.66 4446.46 2960.51 1723.61 ?1, МПа 0.859 0.741 0.493 0.287 ?2-way, МПа 0.87 0.75 0.47 0.31 ?1-way, МПа 0.84 0.69 0.48 0.31 ?, % 1.31 1.20 4.98 7.91 Рис. 9. Распределение главных напряжений вблизи полюса сферы: а) p=1000 Па, v= 20 м/с, б)p=1000 Па, v=40 м/с Выводы:проведенные расчеты в программном комплексе Ansys показали достаточно хорошую сходимость результатов с данными, полученными в ходе экспериментального исследования. Расхождения результатов отчасти можно объяснить упрощением закона деформирования материала до линейного и не- достатком информации по конфигурации экспериментальной установки. В работе показано, что при больших скоростях потока учет деформаций оболочки приводит к изменению распределения аэродинамического коэффициента и, как следствие этого, изменению НДС оболочки. Кроме того, учет деформаций оболочки в потоке обязателен при расчете нестационарных процессов, выявлении режима колебаний оболочки в неустановившемся сильно-турбулентном потоке. В связи с этим намечены следующие направления исследования работы воздухоопорных оболочек: ? разработка и применение физически-нелинейной ортотропной модели материала на основе экспериментальных данных; ? уточнение характеристик демпфирования конструкции; ? анализ существующих и выбор оптимальной модели турбулентности; ? учет совместной работы оболочки с усилением в виде тросовой сети.* * Работа опубликована по рекомендации оргкомитета Молодежного международного научно-практического семинара молодых ученых и студентов «Современные исследования в области прикладных инженерных наук», 6-8 декабря 2016 г., грант РФФИ №16-38-10439мол_г.

Об авторах

НИКОЛАЙ АНДРЕЕВИЧ МОКИН

ФГБОУ ВО «Национальный исследовательский Московский государственный строительный университет» (НИУ МГСУ)

Автор, ответственный за переписку.
Email: mokiavelli@mail.ru

аспирант кафедрыcтроительной и теоретической механики

129337, г. Москва, Ярославское шоссе, д. 26

Список литературы