CC BY
58
УЧЕНЫЕ ЗАПИСКИ Ц А Г И Том VII 19 7 6
№ /
УДК 534.121.014.2:538.601.155
ПРИМЕНЕНИЕ МЕТОДА КОНЕЧНЫХ ЭЛЕМЕНТОВ ДЛЯ ИССЛЕДОВАНИЯ УСТОЙЧИВОСТИ ТРЕУГОЛЬНЫХ ПЛАСТИН В СВЕРХЗВУКОВОМ ПОТОКЕ
В. А. Вислоух, В. П. Кандидов, С. С. Чесноков
Для исследования динамической устойчивости пластин предлагается треугольный конечный элемент, который используется совместно с известным прямоугольным. Аэродинамические воздействия сверхзвукового потока вычисляются на основе постоянной сетки разбиения пластины на конечные элементы. Учитывается распределение скоса потока в обратном конусе Маха. Рассмотрены примеры, приведено сравнение с экспериментальными данными.
В настоящее время метод конечных элементов находит применение в задачах динамической устойчивости пластин, обтекаемых сверхзвуковым потоком [1, 2]. При исследовании пластин сложной геометрической формы обычно используются треугольные элементы с шестью и более обобщенными координатами в каждом узле [3]. Их практическое применение требует значительного объема вычислений.
В настоящей работе предлагаются треугольные элементы с четырьмя обобщенными координатами на узел, которые можно использовать совместно с известными совместными прямоугольными элементами [4]. Развивается алгоритм, согласно которому аэродинамические воздействия вычисляются с учетом распределения скоса потока в обратном конусе Маха. При этом используется та же постоянная сетка разбиения, что и при расчете конструкции. При вычислении обобщенных аэродинамических и инерционных матриц, для понижения порядка уравнений динамики вводятся укороченные системы базисных функций.
На фиг. 1 изображена схема пластины, состоящая из прямоугольных и треугольных элементов. Рассмотрим к-й треугольный элемент в локальной системе координат Юу. Введем обобщенные координаты, равные значениям функции поперечного перемещения г(Е, т), <) и ее производным в узлах:
дг дг д2 г
®в. ?л = ^ > *я= й * Хп== д%дт\ ' л =1,2..., (1)
где п — номер узла элемента.
Эти величины образуют вектор обобщенных координат элемента ?(/). имеющий 12 компонентов:
дГ= (да, <р, V, т }. (2)
Форму прогиба элемента возьмем в виде
г (?, 1), *) = 1Г ад (*), (3)
где хГ(Е. ч)— строка базисных функций; а — числовая матрица порядка 12 X 12.
В качестве системы базисных функций выберем следующий степенной ряд
ХГ={1, £> т). £2. £*1» Ч2. £3> £2,1. ^12. ",13> ^Ч3}- (4)
Ограничимся случаем однородной пластины с жесткостью на изгиб О.
Используя известное (см. [5]) выражение для плотности энергии, а также (3), представим потенциальную энергию элемента следующим образом:
Ур = ЧТ Кч> К = аг 11 в (I, т)) йЫф,
£
здесь /С — матрица жесткости элемента;
е(£, 1\) = В
д2 х д2 хг
-ар"
д2Х д2х
дтг]2
йт,2
^ д2 у2 б2 ХГ ,
'^"аі2 5Ї2" +
есть матрица плотности потенциальной энергии, « — площадь элемента.
Матрица К определяет вектор обобщенных сил упругости ф, сопряженный <?. Его компонентами являются поперечные силы упругости Р и моменты ЛГ, Т, Н в узлах;
О = КЧ.
Для вычисления узловых сил инерции воспользуемся принципом виртуальных перемещений:
(5)
где р — нагрузка на единицу площади. В
,и
случае сил инерции р = — рЛ —(л—толщина пластины, р — массовая плотность). Здесь перемещение возьмем в виде (3), используя укороченную систему базисных функций
х[(1 ■Ц) = (1’ 5. Ч/- <6)
Тогда, после преобразований получим
$т (<) = Щ (0, Л1 = а[ 11 7л х[ й&Ь)«!,
где Ж — укороченная матрица масс, — числовая матрица 3X12, соответствующая системе функций (6).
Прямоугольный элемент в каждом узле имеет координаты (1); его вектор обобщенных координат содержит 16 компонентов. Матрица жесткости такого элемента описана в работе [4]. Для получения укороченной матрицы масс используется система базисных функций вида [7]:
І2 (5. 1)):
{1, £, ч, &)}•
(7)
Числовую матрицу порядка 4x16, соответствующую этой системе, обозначим а2.
Вычисление аэродинамических воздействий основывается на интегральном соотношении [6] между потенциалом возмущений Ф (л:, у, *) и скосом потока V (£'* V, *)» которое справедливо для малых чисел Струхаля.
Ж, 1 ГГ У(6', V, ,8.
Е'
здесь 2 — площадь в обратном конусе Маха с вершиной в точке (х, у)', З2 = Мд — 1, М0 — число М.
Рассмотрим узел /-го элемента, в обратном конусе Маха которого находится к-Ь элемент. Стационарная часть скоса потока на к-м элементе с учетом (3) и (6) выразится через его обобщенные координаты:
У (Є, Чі і)= £/«аі?(0, (9)
ос
где и — скорость потока.
Аэродинамическое демпфирование в приближенных выражениях (8) и (9) не учитывается, так как при малых числах Струхаля оно, как показано в [6], лишь незначительно увеличивает критическую скорость флаттера.
Переходя в (9) к координатам £', ц', подставим выражение для скоса в (8). После интегрирования по площади 5* 2 получим значение потенциала в узле /-го элемента. Полагая в (9) (х, у) равными последовательно координатам всех узлов /-го элемента, определим вектор потенциалов в узлах /-го элемента за счет скоса на &-м.
Ф т = иа№д{і)і (10)
здесь (?;•* — числовая матрица, зависящая от взаимного расположения элементов / и Распределение потенциала внутри /-го элемента'представим в виде, подобном (3), используя укороченную систему базисных функций прямоугольного элемента (7).
Ф(Е, ч. *) = г1а2Ф((). (11)
Переходя от потенциала к давлению, с учетом (10), (11), получим зависимость плотности аэродинамической нагрузки на /-м элементе от вектора обобщенных координат к-го.
д'12 -
(6, 1), о = 2р0 т — а2 0;к Я ((). (12)
Вектор аэродинамических сил ($а определим из принципа виртуальных перемещений (5), учитывая (3) и (7). После интегрирования по площади /-го элемента получим
Я а (0 = 2р0 и2 аВ)к Ц ((),
где а — размер /-го элемента вдоль по потоку; р0—плотность воздуха, В^— матрица аэродинамического влияния, определяющая узловые силы на /-м элементе за счет возмущения, созданного £-м.
При вычислении матриц аэродинамического влияния прямоугольного элемента на прямоугольный или треугольный в (9) используются базисные функции (7), а в соотношении (11) — функции (7) или (6) соответственно.
Для получения уравнения движения пластины образуем вектор ее .обобщенных координат
- Яъ — та}>
где — вектор узловых смещений размерности Я. X 1; /? — число свободных узлов модели.
Суммируем силы упругости, инерции, аэродинамического воздействия в каждом узле схемы с учетом взаимного расположения элементов. В результате получаем матричную систему уравнений движения пластины в потоке газа. Вследствие использования укороченных систем ^базисных функций в матрицах и Въ отличны от нуля только блоки и В^, связывающие силы и смещения в узлах. Переменные ^ не являются независимыми и их можно
исключить [6]. Полагая №ъ(() = '№ех*, где X = 6 + ш— частотный параметр, получаем задачу на собственные значения:
= (»з)
здесь х=2ро£.3а0/£> — приведенное аэродинамическое давление; а2 = X2 рА/^/О — приведенная частота; £ — размер пластины по потоку; с0—скорость звука. Редуцированная матрица жесткости получается из исходной путем преобразований, подобных описанным в [7]. Таким образом, число степеней свободы пластины [порядок матричного уравнения (13)] понижается до числа свободных узлов.
При анализе колебаний пластины без потока в уравнении (8) следует положить х = 0. В качестве примера рассмотрим пластину в форме равнобедренного треугольника. Безразмерные частоты а вычислены при различных разбиениях пластины на элементы. В частности, на фиг. 1 изображено разбиение 4X4. Расчетные данные, а также частоты пластины, составленной из элементов с шестью степенями свободы на узел [1] и экспериментальные [8], приведены в таблице.
С укороченной матрицей масс Совместимые элементы [1] Эксперимент [7]
Разбиение 4X4 5X5 3X3 4X4
Порядок матриц 10 15 40 65
1 тон 6.039 6.065 6.292 6.288 5.93
2 тон 23.81 23.65 23.92 23.91 23.5
3 тон 33.35 33.31 33.34 33.31 32.6
4 тон 59.81 59.44 57.32 57.25 55.7
Из таблицы видно, что результаты, полученные с укороченными матрицами масс, хорошо согласуются с экспериментальными данными. По первым трем тонам предложенные элементы обеспечивают точность не хуже, чем совместные, при значительно меньшем порядке матричного уравнения движения (соответственно 15 и 40).
Для исследования динамической устойчивости пластины в потоке использовалась схема, имеющая 21 степень свободы, что соответствует разбиению 6X6-Аэродинамические матрицы были получены для нескольких значений числа М0= 1,7; 2; 2,5; 3. Собственные значения и векторы задачи (13) вычислялись на ЭЦВМ БЭСМ-4 при различных %. С ростом приведенного аэродинамического давления частоты первых двух тонов в1( о2 сближаются, оставаясь чисто мнимыми. При критическом значении чкр они сливаются и в дальнейшем становятся комплексными: мнимые части совпадают, действительные отличаются знаком.
Происходит потеря устойчивости типа флаттера. Вычисленные критические значения обозначены на фиг. 2 треугольниками. Там же приведены экспериментальные значения из работ [1] (отмеченные кружками) и результаты расчетов методом конечных элементов на основе поршневой теории [1] (сплошная линия). Видно, что полученные результаты хорошо согласуются с экспериментальными и несколько отличаются от вычисленных по поршневой теории без учета аэродинамического демпфирования. Развитый алгоритм определения аэродинамических сил расширяет область применения метода конечных элементов к задачам устойчивости на случай, когда М0<»£/£/<С 1.
ЛИТЕРАТУРА
1. Olson М. D. Some flutter solutions using finite elements. AIAA Journ., vol. 8, N 4, 1970.
2. Кандидов В. П., ЧесноковС. С. Расчет устойчивости прямоугольных пластин в потоке воздуха методом конечных элементов. Вестник МГУ, сер. „Физика, астрономия', № 5, 1972.
3. Каупер и др. Применение высокоточных треугольных элементов изгибаемых пластин в статических и динамических задачах. Ракетная техника и космонавтика, т. 7, № 10, 1969.
4. Carson W. G., Newton R. Е. Plate byckllng analysis using a fully compatible finite element. AIAA Journ., vol. 7, N 3, 1969.
5. О г и б а л о в П. М. Изгиб, устойчивость и колебания пластинок. М., изд. МГУ, 1958.
6. БисплингхоффР. Л., Эшли X. Аэроупругость. М., Изд. иностр. лит., 1958.
7. Vysloukh V. A., Kandidov V. P., Chesnokov S. S. Reduction of the degrees of freedom in solving dynamic problems by the finite element method. Intern. Journ. for Numerical Meth. in Engng. vol. 7 N 2, 1973.
8. Gustafson P. N. at all. An experimental study of natural vibrations of cantilevered triangular plates. Journ. of the Aeronaut. Sci vol. 20. N 4, 1953.
Рукопись поступила I0j VII 1974 г.
