Моделирование виброакустических характеристик трубопровода с использованием метода конечных элементов

• ¹    Самарский государственный аэрокосмический университет ² Институт акустики машин

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

Вопрос изучения виброакустических характеристик трубопроводов заслуживает серьёзного внимания, так как он непосредственно связан с обспечением надёжности в условиях периодически действующих сил. На данный момент существует ряд теоретических работ посвящённых рассматриваемой тематике [1-6]. Одной из первых появилась работа [3]. В ней рассматривалась задача определения амплитуды параметрических колебаний прямолинейного участка трубопровода на податливых опорах с линейными упругими характеристиками с учётом демпфирования. Н.С. Кондрашовым было представлено аналитическое решение дифференциального уравнения, описывающего процесс установившихся параметрических колебаний. Полученное решение было ограничено случаем основного параметрического резонанса (®вын = 2®соб )• Кроме того, в указанной математической модели пренебрегалось распределением параметров пульсаций давления и скорости рабочей жидкости по длине исследуемого трубопровода. Последнее допущение существенно ограничивало область применения модели и было корректно лишь для низких частот, когда длина исследуемого трубопровода не превышала 1/8 длины волны акустических колебаний в рабочей жидкости, а механическая подсистема находилась в дорезонансной области. В работе [2] решалась задача определения вибропараметров колебаний криволинейного трубопровода, вызван- ных пульс ациями давления протекающей жидкости. В частности были получены результаты по определению возбуждающей силы. Однако, полученная методика справедлива только для низкочастотных колебаний, так как в указанной модели не учтено изменение параметров пульсаций давления и скорости жидкости по длине трубопровода. Кроме того не учитывалось влияние присоединённой гидравлической цепи, вида опор и некоторых геометрических параметров системы. В работе [6] рассмотрены методы снижения вибрации трубопроводов. Была разработана математическая модель, позволяющая производить расчёт влияния пульсаций давления на динамические характеристики трубопровода. За основу были приняты уравнения Н.Е. Жуковского, которые в дальнейшем были приведены к виду Матье-Хилла с периодическими коэффициентами. Однако, полученные соотношения позволили произвести только расчёт зон динамической неустойчивости. Таким образом, математическая модель работы [6] позволяет оценить динамику системы только на резонансных режимах. В работе [1] рассмотрены вопросы динамики и устойчивости трубопровода и его элементов. При этом учитывается большое число характеристик исследуемых физических объектов и внешних воздействий, таких как старение материалов, воздействие потока жидкости, вязкоупругого основания и сжимающего (растягивающего) продольного уси- лия. Однако, полученная автором аналитическая модель не позволяла производить расчёт параметров вибрации при вынужденных колебаниях от пульсаций рабочей жидкости, что представляет большой практический интерес. Сам автор отмечает, что решение подобных задач возможно только на базе численных методов расчёта. Развитием исследований Кондрашова Н.С. является работа [4], где представлено решение ранее полученного уравнения [3] для случая вынужденных колебаний трубопровода при возбуждении пульсациями рабочей среды произвольной частоты. При этом в учёт принималась и распределённость параметров пульсаций. Однако данное решение было справедливо только для однородного прямолинейного трубопровода. В работе [5] решена задача о поперечных колебаниях неоднородного прямолинейного трубопровода в условиях пульсаций давления протекающей жидкости. Предложена методика интегрирования полученной системы дифференциальных уравнений в частных производных с помощью приближённого метода Бубнова-Галеркина. Однако описанная модель была справедлива только для прямолинейного участка трубопровода. Кроме того, учёт сил трения производился с помощью подхода, изложенного в [3], который требует развития и уточнения.

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

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

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

Геометрические размеры исследуемого трубопровода следующие: длина прямолинейного входного участка: l 1 = 200 мм ; радиус колена R = 200 мм ; длина прямолинейного выходного участка l ₂ = 200 мм ; наружный и внутренний диаметры d н = 22 мм и d = 20 мм . вн

В качестве опор принимается жёсткая

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

В исследуемой системе моделируются следующие нагрузки: на входе задаются пульсации давления жидкости с частотой to и амплитудой А_п ; на выходе ‒ динамическая нагрузка гидравлической подсистемы. Моделирование производилось для двух крайних случаев динамической нагрузки: акустически закрытого конца (стоячая волна) и нагрузки на неотражающее сопротивление (бегущая волна) (рис.2).

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

Z

н

, pc И , pc ) .

x ■ — +1 У ■ — ■ j

,

s ( s где p — плотность рабочей жидкости;

c ‒ скорость звука в рабочей жидкости;

s ‒ площадь поперечного сечения трубопро- вода;

j = V—1 — мнимая единица;

,_       2 N x = (N2 +1)-(N2 - 1) cos 8 ’

‘=       ( N ² - 1 ) sin 8

y " ( N ² + 1 ) - ( N ² - 1 ) cos2 8 ’

N _ p max = 1 + V1 - а p min   1 - V1 - а

‒ коэффициент сто-

ячей волны;

X ) 2n --■------------

4 X

‒ фаза коэффициента отра- жения;

а — коэффициент звукопоглощения;

x0 ‒ координата ближайшего к концу трубы минимума давления.

Таким образом, по известному импедансу можно рассчитать коэффициент звукопоглощения. Для случая бегущей волны импеданс присоединённой гидравлической цепи равен волновому сопротивлению трубопровода Z=Zв. При этом амплитуда волны давления постоянна по длине трубопровода (при условии отсутствия трения), то есть xо=0, тогда а = 1. Для случая стоячей волны Z=^. Известно [8], что для рассматриваемой стоячей волны минимум пульсаций расположен в сечении, удалённом на четверть длины вол- ны от акустически закрытого конца, то есть x= Х/4, что по умолчанию реализовано в программном комплексе ANSYS.

Задача решается в трёхмерной постановке. Для моделирования жидкости используется восьмиузловой пространственный элемент FLUID30. Этот элемент предназначен для описания жидкой среды, а также её взаимодействия с твёрдым телом. Обычно применяется для исследования звуковых волн и различных динамических процессов. Узлы используемого элемента имеют 4 степени свободы: перемещения по координатам X , Y , Z и давление. Математическое описание элементов данного типа базируется на волновом уравнении:

J2 ■I^P-V 2 P = 0’ c о t где Р ‒ давление в жидкости; t ‒ время;

Рис.2. Краевые условия исследуемой системы

^ = { L } ' =

f d d d)

d x d y d z •

\ z /

В соответствии с методом конечных элементов волновое уравнение может быть записано в матричной форме:

[ m P H P } + k ][ P . } + P o [ R e 1 ^Г ( U e ) = { 0 } .

где [ м р ] = -2y J { N }{ N } T d ( vol ) — матрица c vol

масс жидкости;

[kP ]= J[B ]Г [Bd (vol) ‒ матрица жёсткости vol жидкости;

Po [Re ]= Po |{v Ы {V ' }d (S)

‒

матрица

S масс виброакустического взаимодействия;

{ P_e } - вектор узловых давлений;

{ v e } - вектор узловых перемещений;

{ V } - функция формы элемента давления;

{ v / } - функция формы элемента перемещения;

[ в ] = { ! № v } T ;

{ и } - единичный вектор нормали к поверхности S .

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

{ F e } - вектор внешних сил.

Нагрузка (коэффициент звукопоглощения) задаётся операцией IMPD. Математически это выражается добавлением слагаемого [ c eP ]{ P e } в левую часть волнового уравнения в матричной форме.

CP ]= - J{V}{V}d (S)‒ матрица дем-cS пфирования жидкости, определяющая условия диссипации на границе взаимодействия.

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

Для моделирования стенки трубопровода используется десятиузловой тетраэдральный (с промежуточными узлами на рёбрах) пространственный элемент с неплоскими гранями SOLID92. Этот элемент используется для трёхмерного моделирования твёрдых тел. Его узлы имеют 3 степени свободы: перемещения по осям X , Y , Z . Математическое описание данного типа элемента для рассматриваемого случая акустического анализа определяется уравнением:

[м ]V}+ [c ]{t7 }+[K Kv }_{f-}

Г^[ M e ^{]    [ 0 ]} || { e } 1

[M' ] [MP №Я

+

+

ГЕСе ]

e

[[0]

+

Г[Ke ]

e

[ [ о ]

[о]Ш }1 +

[CP Шр }/+

[ K f^- ] ]J { V e } | _| { F e } |

[KP jjliPe }J l{0}J ’

где [ M e ] - матрица масс твёрдого тела;

M ' ] = P o [ R e ] T ;

[ve ] - вектор узловых перемещений твёрдо го тела;

[ C e ] - матрица демпфирования твёрдого тела;

[ K e ] - матрица жёсткости твёрдого тела;

K f ] =- [ R e ] ;

[ м ] - матрица масс твёрдого тела;

[ С ] - матрица демпфирования твёрдого тела; [ K ] - матрица жёсткости твёрдого тела;

{ F a } - вектор внешней нагрузки.

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

c = I    c 0

D E

1 | вн ж

V 8 E ст тр

,

где c ₀ ‒ скорость распространения волн дав-

ления в жидкости, заключенной в трубопровод с абсолютно жесткими стенками;

D_BH - внутренний диаметр трубопровода;

• 5 ст - толщина стенки трубопровода;

E ж - модуль упругости рабочей жидкости;

E тр - модуль упругости материала трубопровода.

Величина скорости звука в жидкости при условии абсолютно жёстких стенок трубопровода принимается равной с= 1300 м/с . Аналитический расчёт для случая податливых стенок по вышеприведенной зависимости показал следующее значение скоро сти звука с= 1214 м/с. При численном моделировании скорость звука получилась равной с= 1202 м/с . Расхождение последних двух величин связана с погрешностью численных методов и не влияет на физическую адекватность получаемых результатов.

Геометрическое моделирование трубопровода и жидкости имеет ряд особенностей. Сначала производится построение поперечного сечения стенки трубопровода и заключённой в нём рабочей жидкости. Затем выполняется вытягивание построенных поверхностей вдоль заранее заданного контура при помощи операции VDRAG. Далее, при помощи операции VGLUE, производится склеивание объёмов жидкости по наружным поверхностям с соответствующими объёмами стенок трубопровода. При этом полученные объёмы стенки и рабочей жидкости имеют единую общую граничную поверхность. Это необходимо для дальнейшего задания взаимодействия пристеночных элементов со стенкой трубопровода. После этого на геометрические объекты наносится сетка конечных элементов. Следует отметить, что нанесение сетки конечных элементов на этапе построения сечения с дальнейшим вытягиванием вдоль оси трубопровода нецелесообразно, так как полученные при этом объекты не будут иметь общей единой поверхности, и моделирование взаимодействия стенки и жидкости станет невозможным.

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

А = —

А вх

.

Распределение пульсаций давления по длине (рис.3) является характерным для стоячей волны: на акустически закрытом конце реализуются максимальные амплитуды. Моделируемая частота колебаний на первом этапе была специально выбрана достаточно низкой ‒ ниже первой резонансной частоты как по гидравлической, так и по механической подсистеме. В этом случае амплитуды пульсаций давления монотонно уменьшаются по длине от выходного к входному сечению. Для оценки адекватности результатов численного моделирования пульсационного состояния были проведены аналитические расчёты по известным моделям [7, 8]. В соответствии с [8] отношение амплитуды пульсаций на акустически закрытом конце к амплитуде в текущем сечении определяется по соотношению:

А„ =         2,

\ 1 + cos ( 2 ■ £ ■ l )

Рис.3. Распределение пульсаций по длине трубопровода для случая акустически закрытого конца при частоте 150Гц где l ‒ расстояние от анализируемого сечения до акустически закрытого конца трубопровода;

£ - коэффициент фазы; е = ^— ;

с to - круговая частота колебаний.

Различие в результатах численного и аналитического моделирования не превышало 0,4% и объясняется погрешностями численных методов. Следует также отметить, что расчёт с использованием программного комплекса ANSYS показал софазный характер колебаний давления в трубопроводе с акустически закрытым концом на дорезонансных частотах, что также согласуется с известными аналитическими моделями [7, 8].

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

I М.) = 2nfx, c где x ‒ расстояние от рассматриваемого сечения до выхода исследуемого трубопровода.

В ходе моделирования также было получено распределение вибропараметров по длине трубопровода для обоих исследованных случаев (рис.4, 5).

Результаты представлены в виде областей равных значений для относительных виб- v м2 ■ с роскоростей v =---,-----.

А кг вх

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

Рис.4. Распределение виброскорости по длине трубопровода для случая акустически закрытого конца при частоте 150Гц носительной виброскорости для случая реализации бегущей волны составляет 7,46∙10-9, в то время как для нагрузки на акустически закрытый конец ‒ 9,09∙10-9.

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

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

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