Сравнение слабоИ жесткосвязных постановок при моделировании взаимодействия сверхзвукового потока и деформируемой панели Текст научной статьи по специальности «Математика»

УДК 519.63, 65 : 533, 539.3

СРАВНЕНИЕ СЛАБО- И ЖЕСТКОСВЯЗНЫХ ПОСТАНОВОК

ПРИ МОДЕЛИРОВАНИИ ВЗАИМОДЕЙСТВИЯ СВЕРХЗВУКОВОГО ПОТОКА

И ДЕФОРМИРУЕМОЙ ПАНЕЛИ

КОПЫСОВ СП., ТОНКОВ Л.Е., ЧЕРНОВА А.А.

Институт механики УрО РАН, 426067, г. Ижевск, ул. Т.Барамзиной, 34

АННОТАЦИЯ. Рассматриваются вопросы применения слабо- и жесткосвязных постановок задач численного моделирования взаимодействия деформируемого твердого тела и газа решаемых в инженерных комплексах программ. Численные результаты сопоставляются с известными экспериментальными и расчетными данными. Проводится анализ используемых алгоритмов связывания и интерполяции.

КЛЮЧЕВЫЕ СЛОВА: математическое моделирование, сопряжённая задача, сверхзвуковой поток, упругое деформирование, слабосвязная и жесткосвязная постановка.

Задачи численного моделирования взаимодействия деформируемого твердого тела и газа (жидкости) — Fluid Structure Interaction (FSI) последнее десятилетие привлекают к себе все большее внимание, что обусловлено как развитием методов, алгоритмов, вычислительных систем, так и растущей востребованностью результатов при решении широкого круга прикладных проблем. Численное моделирование междисциплинарных задач требует совместного решения уравнений динамики деформирования твердого тела и уравнений газодинамики. Существуют два основных подхода для связывания физических подзадач: монолитный и последовательный (разделенный), однако, в виду существующих ограничений монолитного подхода [1], для данного типа задач целесообразно применение последовательного, который позволяет очевидным образом учитывать специфику каждой задачи, но при этом требует реализации алгоритма сопряжения.

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

Создание и широкое распространение коммерческих программных продуктов численного моделирования (ANSYS, Solid WORKS, FlowVision и т.д.), обусловленное развитием численных методов и программных сред, приводит к необходимости проверки и сравнения алгоритмов связывания, схем интерполирования, применяемых в коммерческих продуктах. В работе, на примере задачи о взаимодействии сверхзвукового потока и деформируемой панели, выполнена валидация программного комплекса ANSYS (ПК ANSYS) с экспериментальными данными [2] и другими численными расчетами [3].

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

ПОСТАНОВКА ЗАДАЧИ СОПРЯЖЕНИЯ

В задачах взаимодействия газа с деформируемым телом на контактных поверхностях необходимо сформулировать дополнительные кинематические и динамические граничные условия. Пусть имеется физическая область Q ограниченная поверхностью 5Q. Область разделена на части занимаемые газом Q (t) и деформируемым твердым телом Qs (t), так

что Q(t) = Qs (t) U Qg (t). Граница раздела r(t) между Qs (t) и Q (t) существует и

определяется в каждый момент времени как r(t) = Qs (t) | Qg (t). Далее будем применять

обозначения rg и Г, если речь идет только о подобластях, занятых газом и деформируемым

телом соответственно.

Так как твердое тело предполагается непроницаемым для газа, то для произвольного момента времени скорости us (t) и ускорения us (t) на Г(?) для тела должны быть равны

скоростям Ug (t) и ускорениям U (t) газа на контактной границе:

us (t) = Ug (t), iis (t) = Ug (t). (1)

Эти равенства выражают условия кинематической совместности. Условия динамической совместности имеют вид

Psn s = Pg n g , (2)

где ps, pg - давление; ng, n s : ng =-n s - внешние нормали.

Следствием условий (1) и (2) является выполнение законов сохранения массы, импульса и энергии на границе Г, движущейся со скоростью U г. Для удовлетворения закона сохранения энергии на границе необходимо выполнение

JT Г+AtUrg ■ Pgngdtdrg = jr (us(t + At)-us(t))-Pns<. (3)

г g t г s

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

Для наглядности представления схем связывания разделённых задач введем обозначения: p — давление и иг — перемещение границы сопряжения. Тогда, на

временном шаге tn+1 перемещения иГ+1, полученные решением задачи деформирования с заданной на границе нагрузкой, обусловленной давлением газа ^Г+1, можно выразить через (неизвестную) функцию S: u Г+1 = S(^Г+1). Аналогично, решение газодинамической задачи для давления на интерфейсной границе, обозначим как ^Г+1 = F(uГ+1). Теперь, используя введенные обозначения, задачу сопряжения можно представить в виде системы нелинейных уравнений u^1 = S F (u^1) = S ■ F (u^1).

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

n+1 о/ n+1 n n+1 \ u s = S (u s , u s , Pg ) ,

n+1 77-/ n+1 n n+1 \

Pg = F (P , Pg , u s ) .

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

>r = F (pg+1, P, us),

/ \ (4)

n+1 al n+1 n „n+1 1 u. = SК , us , P ) .

Для решения системы уравнений (4) требуется единственный обмен данными на каждом временном шаге. В жесткосвязной постановке система решается итерационным

методом:

u

n+l,k+1 g

.n+1, k+1

s

-F (p S (u

n+1,k+1 g '

n+1,k s

n+1,k u n

n+1,k+1

), ).

(5)

В ПК ANSYS решение задачи деформирования (функция S) реализуется в конечно-элементном модуле ANSYS Transient Structural, а решение газодинамической задачи (функция F ) — в конечно-объемном модуле ANSYS CFX.

Последовательное решение каждой задачи на временном шаге tn+1 использует решение на предыдущем временном шаге и выполняется несколько итерационных приближений (обозначены индексом k) по схеме (5) до достижения сходимости по скоростям ||Ug+u+1 - Ug+u||, перемещениям ||uin+1-k+1 - un+u|| или давлению Цр^1*+1 - png+hk ||. Как правило,

в этом случае используются классические итерационные методы типа Гаусса—Зейделя или ньютоновского типа. Необходимо отметить, что постановка связывания в стандартных программных продуктах, как например в ANSYS, имеет свои характерные особенности.

РЕАЛИЗАЦИЯ АЛГОРИТМОВ СВЯЗЫВАНИЯ ДЛЯ ЗАДАЧ FSI В ANSYS

В ПК ANSYS имеются специализированные модули ANSYS MFS и ANSYS MFX. В ANSYS MFS связывания выполняются только для конечно-элементных модулей ANSYS Mechanical, а ANSYS MFX обеспечивает связывание конечно-элементных модулей ANSYS с конечно-объемным модулем CFX (Fluent). Модуль связывания ANSYS MFX содержит два основных алгоритма на несогласованных расчетных сетках: односторонний и двухсторонний [6]. Двухсторонний алгоритм передачи данных обеспечивает как слабо- так и жесткосвязную постановки задачи. Общая схема вычислительного процесса с использованием модуля MFX представлена на рис. 1.

Рис. 1. Схема вычислительного процесса при решении сопряженной задачи FSI в ANSYS

g

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

Алгоритмы интерполяции

Для передачи данных от одной физической подзадачи к другой в ПК ANSYS предусмотрено два типа интерполяции: «сохранение профиля» и консервативный [6, 7].

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

В данной работе использовалась интерполяция консервативного типа. Поверхностные элементы сетки на стороне получателя (газа) разделяются на полигоны, в терминологии ПК ANSYS -«пиксели» (рис. 2) таким образом, чтобы обеспечить требуемую подробность разбиения каждого элемента поверхности (по умолчанию 100*100 пикселей). Аналогичным образом разбиваются поверхностные элементы сетки на стороне отправителя [6]. Совокупность пересечений «пикселей» на сторонах отправителя и получателя (рис. 2) образует контрольную поверхность [8], на которой непосредственно осуществляется интерполяция значений.

Модели деформирования сетки

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

При деформировании сетки алгоритмом «регион с заданным движением» задаются перемещения узлов на подобласти. Движение всех узлов определено моделью движения сетки на основе упругой аналогии [8]. В этой модели смещения общей границы области или в отдельных подобластях вычисляются как: У-(лУд) = 0, где д - искомое смещение, Л - задаваемая жесткость сетки.

Рис. 2. Консервативная интерполяция

Итерационный процесс

Процесс согласования решения полученных систем в рассматриваемых модулях пакетных программ является итерационным на каждом временном шаге и в ANSYS имеет следующий вид:

фИ+1 = фП +а(фП+\ -фП),

Фп+1 I п+1,к+1 „п+1,к+11 I п+1,к „п+1,к I „ „ „ I

= {цж ,рг ], Ф = {цж ,рг ), а — параметр релаксации (по умолчанию а = 0,75).

Далее вычисляется а = |ФП+1 -Фп /\Фп 1 , здесь а — нормирующий коэффициент

1^2

изменения нагрузок на интерфейсе. При изменении нагрузок на интерфейсе «газ -конструкция» задается значение допустимой точности (К) и условие сходимости

определяется как: С = 1о§(а/ < 1. Задавая требуемую точность и число внутренних

1ое(10,0/К)

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

2

ЧИСЛЕННОЕ МОДЕЛИРОВАНИЕ

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

3750

ИГ

о

00

4-у

¡л

3405

250

Рис. 3. Схема ударной трубы

Под действием ударной волны, инициированной нестационарным расширением сжатого в камере высокого давления воздуха, возникает деформация упругой пластины, защемленной в рабочей области ударной трубы. Малая продолжительность испытания позволяет пренебречь влиянием турбулентности. В рассматриваемом диапазоне температур и давлений рабочее тело (воздух) можно считать совершенным газом с показателем адиабаты у = 1,41. Деформируемая панель, изготовленная из стали плотностью р = 7600 кг/м и модулем Юнга E = 220 ГПа, полагалась изотропной и однородной. Малые величины деформации и скорости деформации делают влияние демпфирования пренебрежимо малым, что позволяет решать задачу о деформировании, оставаясь в рамках физически и геометрически линейной теории упругости. Расчетная область представляет собой внутреннее пространство ударной трубы (рис. 2). Начальные данные соответствуют условиям проведения эксперимента [2, 3]: температура T = 293 K, скорость газа Ug= 0,

давление в камере (x < 0,75 м) р1 = 0,254 МПа, в канале р 2 = 0,1 МПа.

Используемые несогласованные неструктурированные расчётные сетки содержат 140000 гексаэдрических элементов и 530000 узлов для области газа и 4000 призматических элементов для деформируемой панели. Расчет проводился на кластере «Уран» Института математики и механики УрО РАН с использованием ANSYS Academic Research 14.0.

По результатам расчетов построены зависимости изменения давления в точке А (рис. 3) и перемещения правого края панели от времени при применении различных схем связывания и интерполирования (рис. 4 - 6). Как и в работах [2, 3], смещению кромки против потока соответствуют положительные значения на графиках.

Применение слабосвязанных постановок позволяет получить интересующие перемещения панели с требуемой точностью (рис. 5), однако приводит к эффектам, схожим с влиянием схемной вязкости — отмечено сглаживание второго пика давления: уменьшение амплитуды колебаний до 15 %, уменьшение частоты колебаний («вытягивание» кривых) на 23 % (рис. 4, кривая 4; рис. 6, кривые 4, 5).

Применение алгоритмов жесткого связывания позволяет повысить точность результатов моделирования (рис. 5) и получить физически верное распределение давлений во времени (рис. 4, кривая 3), однако, как было отмечено выше, требует значительных ресурсов.

Таким образом, из рис. 4, 6 видно, что для решения задач FSI в ПК ANSYS целесообразно, несмотря на значительное повышение требований к вычислительным ресурсам, использовать жесткосвязную постановку с предварительной оценкой числа требуемых итераций коррекции давления в модуле ANSYS CFX (kp ), а также количества итераций связывания ( k ).

|>-Па

1 -2 -3 4 /'♦ г /♦•у 1 *

I */ vAAj

[♦* ■ф VV

♦ ♦ / * •¿пф

9 ов?

JLZ

О 0,0002 0,0004 0,0006 0,0008 0,001 0,0012 0,0014 0,0016 0,0018 (, С

1 - эксперимент [2]; 2 - расчет [3]; 3 - расчет в ANSYS, жесткосвязанный алгоритм; 4 - расчет в ANSYS, слабосвязанный алгоритм Рис. 4. Графики зависимостей давления от времени

♦ 1 -3 О 4

♦ ♦ ♦ /

♦ч ♦ ■г

♦ "х

♦ ♦ ♦

О 0,0002 0,0004 0,0006 0,0008 0,001 0,0012 0,0014 0,0016 0,0018 t, с

1 - эксперимент [2]; 2 - расчет [3]; 3 - расчет в ANSYS, жесткосвязанный алгоритм; 4 - расчет в ANSYS, слабосвязанный алгоритм (kp = 1, k = 1); 5 - расчет в ANSYS,

слабосвязанный алгоритм с консервативным типом интерполяции (kp = 50, k = 1) Рис. 5. Графики зависимостей перемещения края пластины от времени

р. Па

с —2 — 3 4 -Б *

• Г

J °4 г VV

Кр* ♦ / i^'tf * * '£ у

ет:?

0 0,0002 0,0004 0,0006 0,0

0,001 0,0012 0,0014 0,0016 0,0018 t, с

1 - эксперимент [2]; 2 -расчет [3]; 3 - расчет в ANSYS, жесткосвязанный алгоритм; 4 - расчет в ANSYS, слабосвязанный алгоритм (кр = 1, к = 1); 5 - расчет в ANSYS,

слабосвязанный алгоритм с консервативным типом интерполяции (кр = 50, к = 1)

Рис. 6. Графики зависимостей давления от времени

Применение различных схем интерполяции при решении рассматриваемой задачи (рис. 6), как и изменение числа итераций на внутренних циклах (кр), не оказывает

значительного влияния на точность результатов моделирования, однако позволяет на 4 % уменьшить время расчета.

Проведенный анализ параллельной масштабируемости и эффективности организации вычислительных процессов в ПК ANSYS позволил выявить следующее — в сопряженной постановке существует возможность параллельного выполнения только ANSYS CFX. Так при слабосвязанной постановке один шаг по времени (А? = 10-6с) занимал 6 с, время процесса распределяется следующим образом: 2 с - на расчеты физических процессов течения газа и деформирования пластины, 3 с - интерполяция данных от одной физической подзадачи к другой и 1 с - передача данных между временными шагами. Увеличение количества итераций кр приводит к соответствующему увеличению временных затрат. При

жесткосвязной постановке один шаг по времени (кр = 50, к = 20) занимал 21 с, время

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

ЗАКЛЮЧЕНИЕ

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

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

Работа выполнена в рамках программы Президиума РАН № 18 при поддержке УрО РАН (проект 12-П-1-1034) и грантов РФФИ №№ 11-01-00275-а, 12-07-00080-а.

СПИСОК ЛИТЕРАТУРЫ

1. Петров И.В., Иванов А.И., Антонов A.H. и др. Взаимодействие сверхзвукового потока и упругой поверхности обтекаемого тела // Математическое моделирование. 1998. Т. 10, № 2. C. 15-24.

2. Giordano J., Jourdan G., Burtschell Y. et al. Shock wave impacts on deforming panel, an application of fluid-structure interaction // Shock Waves. 2005. V. 14. P. 103-110.

3. Копысов С.П., Кузьмин И.М., Тонков Л.Е. Моделирование взаимодействия сверхзвукового потока и деформируемой панели в ударной трубе // Вестник Удмуртского университета. Математика. Механика. Компьютерные науки. 2012. № 2. С. 156-165.

4. URL: http://www.mpcci.de/multiphysics-engineering.html (дата обращения: 19.09.2012).

5. Копысов С.П., Тонков Л.Е. Вычислительные алгоритмы и программное окружение для решения сопряженных задач FSI на гибридных вычислительных системах // Сб. трудов конф. «Теория управления и математическое моделирование». Ижевск : Изд-во ИжГТУ, 2012. С.31-33.

6. Galpin P.F., Broberg R.B, Hutchinson B.R Three-dimensional Navier Stokes predictions of steady state rotor / stator interaction with pitch chance // 3rd Annual conference of the CFD. Society of Canada, Banff, Alberta : Canada, Advanced Scientific Computing Ltd. 1995.

7. ANSYS® Academic Research, Release 12.1, Coupled-Field Analysis Guide. ANSYS Inc.

8. Аганин А.А., Кузнецов В.Б. Метод консервативной интерполяции интегральных параметров ячеек произвольных сеток // Труды семинара КФТИ «Динамика оболочек в потоке». Казань : Изд-во КФТИ КФ АН СССР. 1985. Вып. XVIII. C. 144-160.

9. ANSYS® Academic Research, Release 14, Help system. — ANSYS Inc.

COMPARISON WEAK AND STRONGLY COUPLING STATEMENTS FOR MODELING INTERACTION OF A SUPERSONIC STREAM AND THE DEFORMABLE PANEL

Kopysov S.P., Tonkov L.E., Chernova A.A.

Institute of Mechanics, Ural Branch of the Russian Academy of Sciences, Izhevsk, Russia

SUMMARY. Discusses the application of weakly and strongly coupling statements for numerical simulation of interaction between deformable solids and gas in standard programs are considered. Numerical results are compared with known experimental and estimates. Analyses the used algorithms of coupling and interpolation.

KEYWORDS: the mathematical modeling, the interfaced task, supersonic stream, elastic deformation, weak-coupling statement, strongly-coupling statement.

Копысов Сергей Петрович, доктор физико-математических наук, профессор заведующий лабораторией вычислительных и информационных технологий ИМ УрО РАН, тел. (3412) 21-45-83, e-mail: s.kopysov@gmail.com

Тонков Леонид Евгеньевич, кандидат физико-математических наук, старший научный сотрудник ИМ УрО РАН, e-mail: tnk@udman.ru

Чернова Алена Алексеевна, кандидат технических наук, научный сотрудник ИМ УрО РАН, e-mail: alicaaa@gmail.com