Применение метода граничных интегральных уравнений для моделирования аэродинамической устойчивости мостовых конструкций Текст научной статьи по специальности «Физика»

ПРИМЕНЕНИЕ МЕТОДА ГРАНИЧНЫХ ИНТЕГРАЛЬНЫХ

УРАВНЕНИЙ ДЛЯ МОДЕЛИРОВАНИЯ АЭРОДИНАМИЧЕСКОЙ УСТОЙЧИВОСТИ МОСТОВЫХ КОНСТРУКЦИЙ

А.С. ДОРОГАН, канд. техн. наук

Дальневосточный государственный университет путей сообщения, Хабаровск

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

КЛЮЧЕВЫЕ СЛОВА: МГИУ, МГЭ, аэроупругость, аэродинамическая устойчивость.

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

В работе [1] изложен альтернативный подход решения связанных систем, основанный на многозональном МГЭ. Предполагается квази - динамическое решение для кратковременных воздействий шквала At = 1-3 мин. Приведенный там расчетный аппарат позволяет находить полные значения перемещений и усилий: U(x) = u(x) • At (м), P(x) = p(t) • At (Н).

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

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

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

Наиболее общее решение получается, если применять уравнение движения 2М и 3М упругих тел с учетом демпфирования. Для их решения будем применять численный метод прямого интегрирования по времени.

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

- пружин с коэфф. жесткости K;

- поглотителей колебаний с коэфф. вязкого (линейного) затухания D. Уравнение равновесия Навье под действием динамических и демпфирующих сил в перемещениях в тензорной форме [3]:

dV ,„2 „2ч F д2ы . D ды. K . ч (1)

Cp--^ + (C2 -C2)--J— + ^ = +—^--^ + ^ut(x,y,z,t) i,i = 1,2,3, (1)

дх.дх. dxidxj p dt p dt p

где p - плотность материала; F - объемные силы; D - коэффициент вязкого затухания (т/с); К - коэффициент жесткости пружин (т/с2); Cs, Cp - скорости распространения продольных и поперечных волн; Cs = -Jo7p , Cp = + 2G)/p ,

C /C > 1; G - упругие константы.

p s

Для нестационарного движения возможны следующие формулировки. Применение гармонического возбуждения f (X, t) = P(®) • emt. Реакция конструкции зависит от частоты ш внешнего воздействия. Используется преобразование Фурье для уравнения движения (T.A.Cruse , F.J. Rizzo, R.E. Bellman, А.А. Ходжибоев). Гармоническое воздействие характерно для механизмов.

Применение импульсного возбуждения f (X,t) = а • P(t) •5(X). Происходит воздействие импульса P(t), действующего из точки в направлении а . Такой подход основан на зависящих от времени функций влияния (W.J. Mansur) и будет применен ниже. Уравнение (1) необходимо интегрировать, чтобы получить перемещения. Для этого дополнительно определяются:

- начальные условия при t = to в области S: ы(x, t) = ы°(x); ыi(x,t) = v°(x);

- граничные условия: ut(x,t) = ui(x,t) наГ1; pi(x,t) = pi(x, t) на Г2.

Известно решение (1) без включения 2-го и 3-го членов в правой части, когда учитываются только силы инерции. Тогда ГИУ для нестационарной задачи имеет вид [4]:

X

Cy (%) ы. (%, t) + JJ Pj(%, х, t) • ы. (x, t) • dr (x) • dt

to г

X X (2)

= J J U*. , х, t) • py (x, t) • dr(x) • dt + J J U* , z, t) • bj (z, t) • dS (z) • dt +

tor to S

P • J[U*-(d,z,t) • Uj(Z,z,t) + Uj(%,z,t) • ыj (%,z,t)]|io dS(z), % e Г, z e S, t e T

S

где Cy = 0.5 Sj- - для гладкой границы; U*.(%, x, t),P*.(%, x, t) - зависящие от времени функции флияния Ахенбаха перемещений и усилий [4]. Последняя строка здесь учитывает влияние начального условия R1 при to. ГИУ (2) линейное по своей структуре и требует только интегрирования по времени.

Предлагается учет всех членов ДУ (1). В этом случае с учетом сил демпфирования и сил жесткости после интегрирования по частям окончательно получается следующее ГИУ:

C^&ufat) + f \p;{^xj)- Ujixj) ■ сПхусе = ЙЛ zeS, isT

U^xj)-Pj(xj)-dr(x) m-\ |byez,f)-b,(zj)- dS(z)- dt+

t,T 15 w

-p ■ J[E/*i£, zj) - z,f) - ul(f, Z-. 0 - zj)\ dS(z) -

■z

-fz,i) Dj{z) uj{z,t\ -dSizHzA)-DXz)-utzJ)-dSiz)-^

5 4 5

где Kj - параметр жесткости; Dj - параметр демпфирования. Видно, что это ГИУ становится дважды нелинейным относительно функции u на S в правой части. Ниже дается специальная процедура его численного решения.

Пространственная и временная дискретизация. Пространственная дискретизация задана L граничными элементами и Mвнутренними ячейками области.

Обычно считают функции u и p постоянными вдоль выбранных временных шагов. Общий интервал времени T делится на временные шаги At = T/N ; N -число шагов. Временная сеть имеет равные шаги: tn = n - At, n = 1, 2,...,N

Тогда дискретизация ГИУ (2) в тензорной форме будет:

N L

Cj (I) u j (|, At) + J P^I, X, At) - Фк - uk (x, At) - dr (x) = (4)

n=1 e=1 Ге

N L N M

=ZZ J U (I,X, At) - Фк - pk(x, At) - dr(x) J U* (I,z, At) - bj (z, At) - dS(z) +

n=1 e=1 Ге n=1 m=1 Sm

NM

P-TTJ [U* (I, z, At) - Uj (|, z, At) + U* (|, z, At) - Uj (|, z, At)]|At dS (z), | e Г, z e S

n=1 m=1 Sm

где Ге - граничный элемент; Sm ячейка области; Ф - базисные функции пор-ка к.

Видно, что сначала выполняется суммирование по граничным элементам, а потом - по временным шагам. Это определяет схему численного решения.

В матричной форме система уравнений МГЭ, соответствующая ГИУ (2), имеет вид:

H - U = G - P + B + R1, (5)

где H и G - глобальные матрицы влияния; U, P - глобальные вектора узловых перемещений и усилий; B - вектор объемных сил; Ri - вектор начальных условий; находится из предыдущего временного шага.

Аналогично выглядит матричная система для ГИУ (3):

H - U = G - P + B + R1 + (K + D) - U + R 2, (6)

что К - вектор жесткости, D - вектор демпфирования; R2 - вектор начальных условий. Здесь нелинейные члены перенесены в правую часть и считаются известной функцией, т.е. это возбуждающие силы.

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

На каждом шаге временного цикла требуется решение нелинейной задачи. Используется алгоритм метода простой итерации вида:

H - U,+1(z, tn) = G - P + B + R1 + (K + D) - U (z, tn) + R 2, (7) где задано начальное условие:

H - U0( z,0) = G - P + B + R1. (8)

При больших значениях коэфф. демпфирования (для ж.б. конструкций) возникает система с большим затуханием колебаний, что приводит к большей нелинейности. Требуются более точные методы решения нелинейной задачи: метод Эйткена - Стеффенсона или метод предиктор-корректора Ньюарка (Newark predictor-corrector [5]).

Схема 1 повременного численного решения ГИУ (2) состоит из шагов: Шаг 1. Задание начальных полей перемещений uo (z) = 0, uio (z) = 0; z e S; to = 0. Шаг 2. Начало временного цикла с шагом At и дискретом tn = n • At, n = 1, 2,...,N Шаг 3. Определение граничных скоростей усилий p(S, At) из ГИУ (см. (2) [1]) для At; p (S, At) = p (S, At) •At.

Шаг 4. Решение граничной задачи (5) и определение перемещений на гран-це Г.

Шаг 5. Решение внутренней задачи в точках коллокаций области S: перемещений

u (z, tn) и их скоростей u (z, tn) = u (z,tn)/ At.

Шаг 6. Накопление значений поля перемещений и его скоростей:

U (z) = U (z) + u (z, tn), U (z) = U (z) + U (z, tn)

Шаг 7. Идти к шагу 2, если tn < T.

Схема 2 повременного численного решения ГИУ (3): Шаг 1. Задание начальных полей перемещений uo (z) = 0, uio (z) = 0; z e S; to = 0. Шаг 2. Начало временного цикла с шагом At и дискретом tn = n • At, n = 1, 2,...,N Шаг 3. Определение граничных скоростей усилий p (S, At) из ГИУ (см. (2) [1]) для At; p (S, At) = p (S, At) •At.

Шаг 4. Решение линейной граничной задачи и определение перемещений на границе Г по матричной формуле (8).

Шаг 5. Решение линейной внутренней задачи в точках коллокаций области S: перемещений u (z, tn) и их скоростей u (z, tn) = u (z,tn)/ At.

Шаг 6. Сохранение значений поля скоростей перемещений uk(z) итерации k - 1.

Шаг 7. Решение нелинейной граничной задачи и определение скоростей перемещений и

усилий на границе Г по матричной формуле (7).

Шаг 8. Решение нелинейной внутренней задачи итерации k в точках области S. Шаг 9. Контроль условия сходимости итерацион. процесса:

max | U'+l - Uj | < е';

j

если условие выполняется - решение закончено, иначе идти к шагу 6. Шаг 10. Накопление значений поля перемещений и его скоростей:

U (z) = U (z) + u (z, tn), U (z) = U (z) + u (z, tn)

Шаг 11. Идти к шагу 2, если tn < T.

Программа Aero Zon (Azon) - это программа для решения 2М и 3М задач аэроупругости. Основана на методе граничных элементов (МГЭ). Имеет возможность учета флаттера.

Укрупненный алгоритм её решения:

1. Начинается цикл по скорости ветра V = 10, 20, ..., 100 м/с; задается ускорение в порыве и общая его длительность T.

2. Составляется квазидинамическая многозональная МГЭ модель А (рис. 1, а), которая включает ветровой и структурный макеты (пролетное строение и пилон).

3. Для модели А находится её граничное значение p(At) для всех временных и граничных узлов на внешнем контуре.

4. Составляется динамическая многозональная МГЭ модель Б (рис. 1, б), которая включает только структурный макет, где заданное p(At) учитывает воздействие ветра.

5. Для модели Б решается задача упругой динамики по временным и граничным узлам с накоплением.

6. Строится график вертикальных перемещений граничных узлов для времени Т и анализируется его поведение. Если происходит нарастание амплитуды колебаний, то скорость ветра характеризует начало изгибно-крутильных колебаний и Vcr = V.

По такой схеме выполняется аналогичное МКЭ моделирование связанного флаттера. В Azon не требуется значительного сгущения ГЭ сетки вокруг обтекаемых конструкций, как в МКЭ. При угле атаки а = 0° анализируется плоский (классический) флаттер, а при а = 10° - срывной флаттер, когда критическая скорость Vcr(a = 10°) < Vcr(a = 0°) в 2-3 раза.

Демпферы в вантовых мостах активно применяются с 90-х годов. Их демпфирующая способность определяется требуемым логарифмическим декрементом [6]. Внутренние демпферы применялись для вантовых мостов во Владивостоке. Демпферы установили также для подвесок 3-го висячего моста через п. Босфор.

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

Особенности выбора моделей для расчета динамических колебаний 2М и 3М тел с применением многозального МГЭ (рис. 1, б):

- для учета сосредоточенных подвесок и вант создаются 2М и 3М ячейки ограниченного размера;

- соответствующие вантам / подвескам пружины и демпферы задаются своими паспортными параметрами демпфирования Dx, Dy, Dz, и жесткости Кх, Ку, К;

- сосредоточенные воздействия корректируются (растягиваются) по ячейке:

поверхности d = D /А, к = К/А, А - площадь;

объема d = D / V, к = К/ V, V - объем.

Заключение. Предложено решение ГИУ колебаний упругого тела (конструкция моста) с учетом сил демпфирования и жесткости, характерных для подвесок и вант. Такая ГИУ дополняется ГИУ многозонального тела (газ + конструкция) в скоростях для решения связанной задачи аэроупругости. Предложен алгоритм численного моделирования флаттера с учетом демпфирования в виде нелинейного временного интегрирования.

S2

-F-

-ANY^

Si

а)

I I I I I I I

Dy D=0; К=0

]Ру Dy

Кх

Рл

б)

51—WV-1

Кх

Рис. 1. Многозональные модели расчета: а) ветровой и структурный макет; б) структур-

ный макет

.....,,r,(rrmrmffff;

'«¡¡тГтщ^

шшве......""""'

liliiijj

„.."Г1"

¡Uli!

Предложенное ГИУ (3), основанное на зависящих от времени функциях влияния, имеет более удобную и универсальную форму численной реализации.

На рис. 2 показаны вертикальные колебания для различных скоростей ветра висячего моста Lo = 3.3 км через Мессинский пролив [5].

Рис. 2. История вертикальных перемещений висячего моста при скорости ветра

30, 40 и 50 м/с.

Л и т е р а т у р а

1. Дороган А.С. Применение МГИУ теории упругости для моделирования аэроупругости мостовых конструкций // Строительная механика инженерных конструкций и сооружений. - 2014. - № 6. - С. 30-35.

90 100 TIME ls«)

2. Казакевич М.И. Аэродинамика мостов. - М.: Транспорт, 1987. - 240 с.

3. Кин Н. Тонг. Теория механических колебаний. - М.: Машгиз, 1963. - 350 с.

4. Dominquez J., Alarcon E. Elastodynamics // Progress in BEM. Berlin: Springer Verlag, 1983. - Pp. 213 - 257.

5. Brotton D.M. Time integration in suspension bridge dynamics // Int. J. Numer. Meth. Eng. - Vol. 20, № 4, 1984, pp. 715 - 732.

6. Constantinou M.C., Tsopelas P. Toggle-brace-damper seismic energy dissipation systems using fluid viscous dampers // J. Struct. Eng., Vol. 127, №2, 2011, pp. 105 - 112.

1. Dorogan, A.S. (2014). The application of BEM elasticity for modeling of the wind flow - bridge structures interaction, Struct. Mech. of Eng. Construct. and Buildings, № 6, pp. 30-35.

2. Kazakevich, V.I. (1987). AerodinamikaMostov, Moscow: Transport, 240 p.

3. Kin N. Tong (1963). Theory of Mechanical Vibration, Moscow: Mashgiz, 350 p.

4. Dominquez, J., Alarcon, E. (1983). Elastodynamics, Progress in BEM. Berlin: Springer Verlag, pp. 213 - 257.

5. Brotton, D.M. (1984). Time integration in suspension bridge dynamics, Int. J. Numer. Meth. Eng., Vol. 20, № 4, pp. 715 - 732.

6. Constantinou, M.C., Tsopelas, P. (2011). Toggle-brace-damper seismic energy dissipation systems using fluid viscous dampers, J. Struct. Eng., Vol. 127, №2, pp. 105 - 112.

THE APPLICATION OF A METHOD OF THE BOUNDARY INTEGRATED EQUATIONS FOR MODELING OF AERODYNAMIC STABILITY OF BRIDGE

CONSTRUCTION

The work is devoted the modeling of oscillations of bridge construction taking into account the forces of damping and rigidity. The boundary integral equation (BIE) in velocities of a multizonal body (gas + construction) is supplemented BIE oscillations of an elastic body (construction) that allows to find the increase of amplitude of their oscillations.

KEY WORDS: MBIE, BEM, aeroelasticity, aerodynamic stability.

References

A.S. Dorogan, Khabarovsk