Разработка математической модели и расчет напряженного состояния вращающегося диска Текст научной статьи по специальности «Механика и машиностроение»

УДК 539.4:678.07

А. Х. Валиуллин

РАЗРАБОТКА МАТЕМАТИЧЕСКОЙ МОДЕЛИ И РАСЧЕТ НАПРЯЖЕННОГО СОСТОЯНИЯ ВРАЩАЮЩЕГОСЯ ДИСКА

Ключевые слова: вращающийся диск, напряженное состояние, метод конечных элементов, составной диск, основной диск,

покрывной диск.

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

Key words: rotating disk, the state of stress, finite element method, composite disk, the primary disk, disk coating.

A variant of the calculation of a rotating disk finite element method Taken triangular element with linear stress distribution. Built stiffness matrix element and form a global stiffness matrix of the disc. A program and calculated the impeller as a composite disk assuming rigid connection of the ground and coating disks.

Вращающийся диск является основным элементом различных центробежных машин и силовых установок: турбина двигателя, рабочее колесо центробежного насоса или компрессора и т.п. Важное место в проектировании этих машин занимают расчеты на прочность, что отмечалось в работах [1], [2]. Скорости вращения обычно очень велики (на периферии диска линейная скорость может достигать звуковой), поэтому возникают большие центробежные силы и, соответственно, большие напряжения и деформации в материале диска. В данной работе приводится расчет напряженно-деформированного состояния ротора, состоящего из основного и покрывающего дисков, соединенных между собой лопатками. Задача решается методом конечных элементов в упругой постановке, влияние температуры не учитывается. Так как сам диск и приложенная к нему нагрузка симметричны относительно оси вращения, то и картина напряженно-деформированного состояния также осесимметрична, поэтому объемная задача сводится к плоской задаче. Задача решается в перемещениях на основании вариационного принципа Лагранжа - условия экстремума дополнительной работы, или - что при справедливости закона Гука равносильно -экстремума потенциальной энергии системы.

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

Для определения узловых перемещений составляется система уравнений

И-{}={}, (1)

где [а] - матрица жесткости размерами п х п , п -число узловых перемещений, {и}- матрица-столбец

искомых перемещений, {{}- матрица-столбец внешних сил, приложенных в узлах.

Входящая в это уравнение глобальная матрица жесткости [а] формируется сложением элементарных матриц в определенном порядке, который зависит от взаимного расположения элементов. Для этого выполняется сквозная нумерация всех элементов и узлов. Порядок нумерации ясен из рисунка. 1, где приведены глобальные номера узлов, а в кружочках показаны номера элементов. Узлы нумеруются и внутри каждого элемента: 1, 2, 3 - наименьший глобальный номер входящих в элемент узлов будет первым локальным, следующие (против часовой стрелки) два узла получают второй и третий локальные номера. Глобальные номера узлов образуют массив / (к,}) здесь / - глобальный номер узла, } -

номер элемента, к - локальный номер данного узла; например, узлы, входящие в элемент номер 7 (рис. 1), имеют следующие глобальные номера: / (1,7)= 5, / (2,7) = 6, / (3,7) = 8.

89

т Ж

0/ /©

4

(3)

4

1

5 (2)

1 (1)

2 (2

2 (1)

Рис. 1

Перемещениям узлов также присваиваются сквозные глобальные номера. В каждом узле возникает два перемещения: горизонтальное и/ и вертикальное У . В одномерном массиве перемещений и (п) горизонтальное перемещение и/ получает нечетный номер 2/ -1 и обозначается

7

8

1

2

и2,_1, а вертикальное перемещение получает четный номер и обозначается и2,.

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

1 <•

и = - ¡еЕ еСУ, (2)

2 V

получим соответствующее выражение в матричном виде.

Выразив деформации через перемещения

ди

дv

ди дv

ех =-г-, ег = ~г-, е =—, Ухг = — + — , (3) дх дг у дг дх

а напряжения через деформации из обобщенного закона Гука

сгх = ^[("1 _ у)ех + уег + уе( ], сгг = ^[("1 _ у)ег + уе1 + уех ],

где

И = '

= /и\( _ у)( + уех + уег], ^хг = М2ухг ,

Е „ 1 _ у _ -у2

(4)

2 =

, Е _ модуль

1 _ у _ 2у2' " 2(1 + у) упругости, у _ коэффициент Пуассона материала диска, получим матричную форму обобщенного закона Гука:

М=[с]-{е}. (5)

Здесь {сг} _ вектор напряжений, [о]_ матрица упругости, {е}_ вектор деформаций:

С} =

с/

сг ,, [ф

с

тхг

1 _ у у у 0

у 1 _ у у 0

у у 1 _ у 0

0 0 0 2

{} =

(6)

с-г

[Ухг

При линейном законе распределения перемещений внутренних точек элемента их можно выразить через узловые перемещения: и = Ы1и1 + Ы2и2 + Ызиз,

V = + N2v2 + Ызуз. (7)

Входящие в эти формулы коэффициенты N

называют функциями формы элемента, они

совпадают с L _ координатами точки [3]. Они определяются так:

1

N = Ц = 2А (2 + Ь х+СУ).

(8)

где А площадь элемента, а,- = Х+1У+2 _ X ,+2У+1,

Ь = у+1 _ у+2 , С/ = Х\+2 _ Х\+1, Х \, у _ координаты узлов, правило круговой перестановки: 1+1=2, 2+1=3, 3+1=1 и т.д.

Выполнив операции (3) над функциями (7) с учетом формул (8), выразим деформации ех, еу, уху через перемещения узлов:

1

А (Ь1и1 + Ь2и2 + Ь3и3 ^

= ег = (c1V1 + ^ 2 + 3 ^

уху = ухг =

ех = х 2А

'у 2А

1

= 2А (си + + Сзиз + + Ь^2 + bзVз).

Найдем выражение тангенциальной деформации через узловые перемещения.

Если считать эту деформацию постоянной в пределах элемента, то

е = ^ = ^ + ^ + у8 = У1 + У2 + Уз.

уср уз уз уз При представлении этой формулы в виде е( = + к2у2 + к^з все коэффициенты одинаковы

к = к2 = кз = 1/уз. (9)

Если функцию 1/ у разложить в ряд Тейлора и сохранить два первых члена, для коэффициентов формулы (9) получим следующие выражения:

к, = ^

уз

2 _^2 + У)

, = 1, 2, з, (10)

при сохранении трех членов разложения 5 з

к''=^" (2у+У+1+У+2)+ +(зу2+у+1+у+ 2+2уу+1+

+ У+1У+2 + 2У+2 У).

(11)

В результате получим следующее матричное выражение деформаций через перемещения узлов:

2А

Ь1 0 Ь2 0

0 С1 0 С2

0 2Ак1 0 2Ак2

С1 Ь1 С2 Ь2

и1

Ьз 0 " v1

0 сз и2

0 2Акз v 2

сз Ьз . из

V з

= [в].{и}.

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

и = 1 ¡ 1жуо{и}Т -[в]Т ■[о]-[в]-{и}-СА =

А (12) = 1 {и}Чв]Тф].[вИи}.2,у о А.

Входящая в эту формулу квадратная матрица размерами 6 х 6 является матрицей жесткости элемента:

[А ] = В]Т-[й]-[в]. уоА, (13)

здесь у0 = уз / з _ ордината (радиус) центра тяжести элемента. Заметим, что в формулу матрицы жесткости элемента не входит множитель 2ж . При

£

х

составлении системы уравнений (1) этот множитель исключается и из матрицы-столбца правой части.

Определим часть матрицы-столбца внешних сил, приходящуюся на данный элемент. На единицу объема диска действует центробежная сила F = рсо2г = ра>2у, которая на перемещении V совершает работу

М = | FvdV = | Fv • 2лу • dA.

v а

Опуская, как и в левой части системы (1), множитель 2я , запишем формулу работы внешних сил виде

М = {} •{Т},

/1 = Гз = 4 = 0,

(14)

где

f2 =

f. =

fa =

Apw2

30 Apw2

30 Apw2 30

(2У,2 + ys2 - V2Y3), ( + ys2 - Y/,), (2У32 + ys2 - YY).

Рис. 2 - Схема рабочего колеса: 1 - основной диск, 2 - покрывной диск, 3 - лопатка; напряжения (МПа): а) в ступице, б) в основном диске, в) в тарелке покрывного диска, г) в цилиндрической части покрывного диска

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

Формирование глобальной матрицы ведется одновременно с вычислением матрицы жесткости элемента, по мере вычисления коэффициентов очередной элементарной матрицы они рассылаются

по соответствующим местам глобальной матрицы жесткости. Параллельно формируется и вектор правой части.

В результате получается система уравнений (1), которая решается методом Гаусса. После определения узловых перемещений поэлементно выполняется определение напряжений:

М=М-[в]-{и}. (15)

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

В тестовой задаче были опробованы различные схемы разбиения диска на элементы, получено, что измельчение сетки повышает точность расчета, на основании этого приняты оптимальные параметры сетки: для каждого диска около 150 приблизительно равновеликих треугольников (без тупого угла), при этом число узлов получается около 100 . Здесь также решалось. как влияет форма представления тангенциальной деформации или - по какой из формул (9), (10), (11) нужно определять коэффициент к/. Расчеты показали, что для принятых здесь геометрических размеров, когда размеры элементов малы по сравнению с радиусом элемента, различия не существенны, поэтому следует применять самую простую из упомянутых формул, то есть формулу (9).

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

ar2

При расчете дисков по отдельности были приняты следующие параметры сетки: в основном диске число элементов 152 , число узлов 111, в покрывном диске - 116 и 90 .

Далее выполнен расчет колеса в сборе. Здесь важно правильно учесть жесткость лопаток, ведь именно они передают вращательное движение от основного диска к покрывному диску. Здесь рассмотрен самый простой случай, когда жесткость лопаток бесконечно велика, то есть покрывной диск как бы насажен неподвижно на основной диск и они составляют одно целое. В этом случае сечение составного диска разбивается на 301 элемент, число узлов составляет 200 . Эпюры эквивалентных

экв

напряжений приведены на рис. 3. Максимальное напряжение возникает в массивной ступичной части колеса в крайнем правом примыкающем к валу элементе (445 МПа), минимальное - в крайнем левом элементе, тоже примыкающем к валу, только с тонкой стороны ступицы (61 МПа). По сравнению с отдельно вращающимися дисками в составном диске напряжения изменились так: в тарельчатой части покрывного диска они заметно уменьшились, а в ступичной части основного диска увеличились.

Расчет выполнен при следующих данных:

- материал - сталь 34ХН3ВЛ, плотность р = 7850 кг / мз, модуль упругости Е = 200 ГПа, коэффициент Пуассона у = 0,з, предел текучести сТ = 700 МПа.

- наружный диаметр колеса 700 мм,

- диаметр вала 1з0 мм,

- внутренний диаметр покрывного диска (входной диаметр) з8з мм,

- остальные размеры на рис. 2 изображены пропорционально,

- угловая скорость вращения со = 1000 с 1.

Литература

1. Максимов, В.А. Компрессорное и холодильное машиностроение на современном этапе / В.А.Максимов, А.А.Мифтахов, И.Г. Хисамеев/ Вестник Казан. технол. ун-та. - 1998. № 1. - С. 104 - 113.

2. Демьянушко, И.В. Расчет на прочность вращающихся дисков/ И.В. Демьянушко, И.А. Биргер. - М.: Машиностроение. 1978. - 247 с.

3. Сегерлинд, Л. Применение метода конечных элементов/ Л.Сегерлинд. - М.: Мир, 1979. - 392 с.

Рис. 3 - Расчетная схема колеса в сборе и эпюры напряжений

© A. Kh. Valiullin, канд. техн. наук, доц., проф. каф. «Теоретическая механика и сопротивление материалов» КНИТУ, tmsm@kstu.ru.

© A. Kh. Valiullin, candidate of PhD, associate professor, department of theoretical mechanics and strength of materials, Kazan National Research Technological University, tmsm@kstu.ru.