Научная статья на тему 'Численное исследование динамики вращающейся изогнутой пружины с использованием конечного элемента в форме одного витка'

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

CC BY
165
23
i Надоели баннеры? Вы всегда можете отключить рекламу.
Ключевые слова
ВИНТОВАЯ ЦИЛИНДРИЧЕСКАЯ ПРУЖИНА / КОНЕЧНЫЙ ЭЛЕМЕНТ–ВИТОК / БОЛЬШИЕ ПОВОРОТЫ / ВЕКТОР ЭЙЛЕРА / МЕТОД ТЕНЕВОГО ЭЛЕМЕНТА / ГЕОМЕТРИЧЕСКАЯ НЕЛИНЕЙНОСТЬ

Аннотация научной статьи по физике, автор научной работы — Иванников В.В., Сорокин Ф.Д., Чжоу Су

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

i Надоели баннеры? Вы всегда можете отключить рекламу.
iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.
i Надоели баннеры? Вы всегда можете отключить рекламу.

A Numerical Study of the Dynamics of a Rotating Curved Spring Using a Single Coil as a Finite Element

A method is developed for numerical calculation of the motion of a curved coiled spring excited by the rotation of one of the grips. For the simulation purposes, the spring is replaced by a discrete set of coils, each of which is considered as a finite element. To construct a geometrically non-linear finite element, that is a coil that takes account of large displacements and rotations, small displacements are counted from an intermediate «shadow» position, in which the element is not deformed. To exclude special points, full rotations are represented as a combination of a tensor and a vector, where the tensor of the large rotation is constant at the integration step, while the variable small part of the rotation is described by the Eluer vector. The method of obtaining the tangent of the rigidity matrix, the generalized mass matrix and the hydroscopic matrix of the finite element (coil) are presented. Non-linear equations are obtained describing the motion of the spring model constructed of such finite elements and numerically integrated by the Newmark method. The calculations results are confirmed by the data obtained through full-scale experiments with a rotating curved spring.

Текст научной работы на тему «Численное исследование динамики вращающейся изогнутой пружины с использованием конечного элемента в форме одного витка»

УДК 534.11: 62-272.22 DOI 10.18698/0536-1044-2018-2-11-22

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

В.В. Иванников1, Ф.Д. Сорокин2, Су Чжоу2

1 ИКЦ ООО «Альфа-транзит», 141400, Московская область, Химки, Ленинградская ул., д. 1

2 МГТУ им. Н.Э. Баумана, 105005, Москва, Российская Федерация, 2-я Бауманская ул., д. 5, стр. 1

A Numerical Study of the Dynamics

of a Rotating Curved Spring Using a Single Coil

as a Finite Element

V.V. Ivannikov1, F.D. Sorokin2, Su Zhou2

1 Engineering Consulting Centre — OOO Alfa-Tranzit, 141400, Russian Federation, Khimky, Leningradskaya St., Bldg. 1

2 BMSTU, 105005, Moscow, Russian Federation, 2nd Baumanskaya St., Bldg. 5, Block 1

e-mail: vvivannikov@gmail.com, sorokin_fd@mail.ru, szleon602@gmail.com

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

Ключевые слова: винтовая цилиндрическая пружина, конечный элемент-виток, большие повороты, вектор Эйлера, метод теневого элемента, геометрическая нелинейность

A method is developed for numerical calculation of the motion of a curved coiled spring excited by the rotation of one of the grips. For the simulation purposes, the spring is replaced by a discrete set of coils, each of which is considered as a finite element. To construct a geometrically non-linear finite element, that is a coil that takes account of large displacements and rotations, small displacements are counted from an intermediate «shadow» position, in which the element is not deformed. To exclude special points, full rotations are represented as a combination of a tensor and a vector, where the tensor of the large rotation is constant at the integration step, while the variable small part of the rotation is described by the Eluer vector. The method of obtaining the tangent of the rigidity matrix, the generalized mass matrix and the hydroscopic matrix of the finite element (coil) are

presented. Non-linear equations are obtained describing the motion of the spring model constructed of such finite elements and numerically integrated by the Newmark method. The calculations results are confirmed by the data obtained through full-scale experiments with a rotating curved spring.

Keywords: coil cylindrical spring, finite element (coil), large rotations, Euler vector, shadow element method, geometrical non-linearity

Винтовые цилиндрические пружины можно ис- вращается в подшипнике, нагруженный момен-

пользовать не только как упругие элементы, но и том сил трения последнего [1].

как инструменты пружинных мельниц, вибро- Интенсивные вибрации вращающейся изо-

просеивателей, шнеков и других механизмов, гнутой пружины должны быть исследованы

предназначенных для обработки или перемеще- еще на этапе проектирования для того, чтобы

ния материала [1]. Например, спиральные гро- выполнить отстройку от резонансных режимов

хоты оснащены вращающимися цилиндриче- вращения или настройку на них, если они яв-

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

материал сквозь зазоры между витками, а также причинам.

дробящими крупные фракции материала. Ци- Цель работы — создание методики исследо-

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

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

пружины приложен крутящий момент от двига- вании витка пружины в качестве конечного

теля, в то время как другой ее конец свободно элемента.

Экспериментальный стенд и параметры пружины. Численное исследование проведено для пружины экспериментального стенда (рис. 1), разработанного Р.Н. Бадиковым [1].

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

• модуль упругости материала E = 2 • 1011 Па;

• коэффициент Пуассона материала ^ = 0,3;

• плотность материала р = 8 000 кг/м3;

• радиус цилиндра, описывающего пружину, R = 13,25 •10"3м;

• диаметр проволоки d = 2,6 • 10_3 м;

• угол подъема винтовой линии у = 3,5°;

• количество витков N = 38. Оба конца пружины закреплены в специальных захватах с винтовыми выемками. Захваты размещены в цилиндрических шарнирах. Вращение от электродвигателя передавалось на один из захватов, второй являлся пассивным (рис. 1, а). Стенд позволяет задавать произвольное положение захватов пружины в плоскости основания и создавать ее разнообразные начальные конфигурации. В ходе динамических испытаний пружина, установленная на стенде, вращалась с постепенно увеличивающейся угловой скоростью до момента наступления резонанса.

Представление поворота комбинацией фиксированного тензора и варьируемого вектора.

В кинематике известно очень много способов

Рис. 1. Внешний вид (а) и схема (б) экспериментального стенда (ф(£) — заданный закон вращения ведущего захвата; I — время)

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

Один из способов обхода указанной проблемы (инкрементный способ) представлен в работе [3], где предложено хранить поворот в виде двух частей — большой и малой. Большая часть должна сохраняться в виде тензора (или матрицы) и быть неизменной на шаге интегрирования, а малая часть представлена вектором Эйлера (вектор ориентации). Вектор Эйлера систематически используется во многих работах, в частности в [3-10]. Так как длина вектора Эйлера, имеющая смысл угла поворота, всегда мала, особая точка 2к, характерная для этого вектора, никогда не достигается.

Тензор поворота связан с вектором Эйлера Ф соотношением, которое можно рассматривать как тензорную функцию векторного аргумента [5-7]:

, . 1 - cos |Ф| sin |ф| L(d) = Ecos Ф +-^ФФ + ^-^Фх E, (1)

|Ф|2 б

где L(6) — функция, вычисляющая тензор поворота по заданному вектору Эйлера; E — единичный тензор; ФФ — диадное произведение векторов; Фх E — кососимметричный тензор с сопутствующим вектором Ф .

В начале шага итерационного процесса или шага интегрирования по времени положение и

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

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

Для перехода от малых перемещений и поворотов к большим в данной работе использован метод теневого элемента (corotational finite element) [3, 9, 12, 13], основанный на том, что деформации в материале пружины являются малыми.

Таблица 1

Алгоритм припасовывания малого поворота к большому

Этап Накопленный поворот Вектор Эйлера

Начало шага интегрирования Ro 0

Конец шага интегрирования Ro Ф

Припасовывание малого поворота к большому и сброс малого Ro ^ LW • R0 Ф^0

поворота

тензор Rm. Сопоставление радиусов-векторов и поворотов узлов КЭВ в состояниях (II) и (III) позволяет найти малые перемещения и повороты, вызванные деформациями

ua =

ub =■

ug - ub 2

ub - ua

■ + (R OT - E )• e zl2;

/ \ l -(Rm -E)• ezI^;

Рис. 2. Схема перемещения и поворотов КЭВ в базовом (0), исходном (I), теневом (II) и актуальном (III) состояниях

Для выделения малых перемещений и поворотов, вызванных незначительными деформациями, полное (большое) перемещение КЭВ рассмотрено как состоящее из нескольких этапов (рис. 2).

В исходном состоянии пружина и виток не деформированы, это положение (I). В актуальном состоянии (III), которого достигает виток в процессе движения, он деформирован, но так как деформации в материале проволоки малы, искажение конфигурации витка невелико. В связи с этим можно рассматривать переход из исходного состояния в актуальное как комбинацию жесткого перемещения из положения (I) в промежуточное теневое положение (II) с последующим наложением малых перемещений, обусловленных незначительными деформациями в материале (переход из (II) в (III)). Термин «теневое» объясняется тем, что состояние (II) реально не существует, оно конструируется искусственно.

Очевидно, что выбор теневого положения неоднозначен, но состояния (II) и (III) должны быть максимально близкими. При нахождении положения (II) в данной статье использован подход, предложенный в работе [3]: из тензоров фиксированных больших узловых поворотов Roя и Rofo с помощью специального алгоритма [3] формируется тензор среднего поворота Rm, а из векторов линейных перемещений узлов щ и щ — вектор среднего перемещения ит = = (щ + щ)/2. Тогда переход из (I) в (II) можно рассматривать как параллельный перенос КЭВ на вектор ит с последующим поворотом вокруг центра (средняя точка между узлами КЭВ) на

Щ -Щ-й; Щ -Щ + 0ъ

где — орт оси КЭВ в состоянии (I)

(см. рис. 2); I — расстояние между узлами КЭВ в исходном состоянии; Ща, Щ — векторы малых дополнительных узловых поворотов; Щ — вектор малого относительного поворота узлов

^0а = Ц-Щ^, R 0Ь = ЦЩ^ ).

Здесь и далее звездочка в верхнем индексе параметров указывает на малые перемещения и повороты (см. рис. 2), вызванные упругими деформациями.

По найденным малым упругим перемещениям (как и в работе [3]) получена потенциальная энергия деформации КЭВ

и=1 2

©а

ub ©b

[KI

ua

©а

ub ©b

[Kii] = [Rau ][Ko][Rau ] Г R II

[R all ] =

RII

R

RII

RII - Rт *R Ь

где [Кп] и [К0] — локальная матрица жесткости КЭВ из работы [11], после поворота КЭВ в положение (II) и в базовом положении (0) (см. рис. 2); ^ац] — полная матрица поворота из положения (I) в положение (II); R и — тензор поворота из положения (I) в положение (II); RI — тензор поворота из положения (0) в положение (I).

Здесь и далее квадратными скобками обозначены матрицы размером 12x12, при этом тензорам соответствуют матрицы 3x3.

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

выми силами, взятыми с обратным знаком, а вторые производные — касательными жесткос-тями. Таким образом, вектор неуравновешенных упругих сил Р и касательную матрицу жесткости можно представить в следую-

щем виде:

(„ >

P = U [K ]=0U.

r ~ > LAvtang J _ 2 >

oy oy2

Ua

Фа Ub

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

Фь

(2)

где y — вектор узловых перемещений КЭВ (тензоры поворота Ro<¡ и Rob фиксированы на шаге интегрирования, поэтому не входят в состав y).

Описанный алгоритм получения P и [K tang] практически полностью идентичен алгоритму, приведенному в работе [3]. Отличие заключается в матрице [K0], которая в работе [3] является общеизвестной матрицей жесткости прямолинейного стержневого конечного элемента. Здесь же [K0] — матрица жесткости КЭВ при малых перемещениях. Так как других различий нет, результат вычисления производных в формулах (2) будет тем же, что и в работе [3]:

где

P = -[H]T[Kn](y 0 +[H]y); [K tang] = [H]T[K n][H],

( E

(3)

[H]=2

0

-E 0

\ (

0

2E 0 0

y 0

-E о 1

0 0

E 0

0 2E /

E )• ezi -1 2

-Ф1

, , l -(Rm -E)• ez¡2

+Ф1

Следует отметить, что в работе [3] выполнено уточнение матрицы [К 1апё], которое сводится к добавлению слагаемых, учитывающих влияние внутренних сил на жесткости. Для КЭВ, который в отличие от обычного стержневого конечного элемента является «мягким», т. е. не содержит жесткостей, различающихся на порядки, необходимости в дополнительных слагаемых не возникает.

Построение матрицы обобщенных масс, вектора инерционных нагрузок и гироскопической матрицы КЭВ. С матрицей масс связана кинетическая энергия

( Va ^

т=1

2

«а

Vb

«b

v ь /

[R а ll ][Mi][R *all ]T

( V 1

«а

Vb

«Ь

v ь /

(4)

[Rail ] =

f

L(Фa ) • R

0a

L^a ) • R

0a

L^b) • R

0b

L^b) • R

0b

где Vа, vь — линейные скорости узлов; юа, щ — угловые скорости сечений; [И *ац ] — матрица поворота, учитывающая малые повороты (в отличие от [Кай ]); [М^ — локальная матрица масс КЭВ из работы [11], соответствующая положению (I) на рис. 2.

Линейные скорости Vа, Vь равны производным линейных перемещений, а угловые скорости юа, щ выражаются через производные векторов Эйлера с помощью тензора П.А. Жилина [4-6]:

( V 1

«a

V b

«b

V b /

= [B

all

i Él. i , >

dt

(5)

[Ball] =

В(ФЯ )

В(Фь )

sin Ф Ф- sin Ф 1 - cos Ф В(Ф) = E^-i-+ ' 1 , J 'ФФ + ——т^ФхE,

Ф

ФГ

Ф

где В(Ф) — функция, вычисляющая тензор П.А. Жилина по заданному вектору Эйлера.

Отметим, что название «тензор Жилина» не относится к общепринятым, например, в работе [10] использована матрица с названием «tangent operator», соответствующая тензору ВТ. Хотя П.А. Жилин не является автором этого тензора, но после публикации его работ тензор В превратился в повседневный инструмент очень многих исследователей. Поэтому исполь-

зование фамилии Жилина в названии тензора B вполне оправданно.

Подстановка выражения (5) в формулу (4) позволяет представить кинетическую энергию через обобщенные скорости

1 T T = - У T[M]y;

[M] = [Вя„ ] T[R *all ][M i][R *all ] T[BaH ],

где [M] — матрица обобщенных масс.

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

г 1 • =2 *<

(6)

где Шу — элементы матрицы [М]; , уу — компоненты вектора у; точкой обозначено дифференцирование по времени; индексы г и у пробегают значения от 1 до 12 (суммирование по повторяющимся индексам).

Уравнения Лагранжа второго рода имеют вид

( дТ >

dt

ду i

= Qi + А,

dyi dyi

(7)

где 0,1 — элементы вектора внешних нагрузок; Рг — элементы вектора внутренних упругих сил (вектор Р из соотношения (3)); Ф — диссипа-тивная функция Релея, Ф-0,5угСуУу (Су — обобщенные вязкости).

После подстановки выражения (6) в формулу (7) получим

+ у к

dMij

-у j

1 дм,

дук = Q> + P -

2 ду.

kj

у J =

(8)

где к — индекс суммирования (немой индекс).

Второе и третье слагаемые в левой части соотношения (8) после их переноса в правую часть можно трактовать как компоненты вектора динамических нагрузок Б:

dMj

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

с 1 . dMkj . .

F =- ук-г~ у j- у к —

2 дуг' дук

у j.

(9)

С учетом формулы (9) уравнения движения КЭВ принимают следующий вид:

Mi

d2 у d2

J + С,^ = Qt + Pt + Ft. (10)

dt

При использовании неявных методов численного интегрирования [10] требуется линеаризация нелинейных уравнений движения (10). Рассматривая уравнения (10) в два момента времени I и t + At и отбрасывая малые второго и последующих порядков, приходим к линеаризованным уравнениям движения КЭВ

Mij Ду J + (Gij + Cij )Ду J + Kij Дyj =

-с„у j +q, +Pt + Ft,

(11)

где величины без А соответствуют моменту времени £ Оу — элементы искомой гироскопической матрицы [С]; Ку — элементы матрицы [Ктд].

При завершении итерационного процесса метода Ньюмарка левая часть выражения (11) становится пренебрежимо малой и решение системы (11) на одном шаге сходится к решению исходных уравнений (10) (более подробное описание такого алгоритма приведено в работе [10]).

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

дЕ

Gj =-

дMiJ дMi,

= у к+ -

дук ду

■уп

ду j

1Г дMkJ дMJ

у к'

дуг' ду;

уп

С учетом симметрии матрицы обобщенных масс окончательное выражение для элементов гироскопической матрицы принимает вид

дШу

Gij =

1J ;

дук

у к + (Wji - Wj),

(12)

где Wij = (дM Jk %) у к; Wj, = /ду j) у к.

Вычисление производных по обобщенным перемещениям в формулах (9) и (12) выполнено с помощью компьютерной аналитики в пакете Wolfram Mathematica [14].

Составление ансамбля конечных элементов и учет граничных условий. Для составления полной модели пружины найденные матрицы касательных жесткостей [Kmg]> матрицы масс [M], гироскопические матрицы [G], а также векторы упругих сил P и инерционных нагрузок F отдельных КЭВ объединены в ансамбль конечных элементов с помощью традиционного для метода конечных элементов

(МКЭ) способа, при котором матрицы накладываются уголком друг на друга (аналогично работе [11]). В результате получены линеаризованные уравнения движения

[М£ ]ДУ + ([С £ ] + ]) ДУ + [К £ ]ДУ =

= Р£ + Р£"[М£ ]У - [С £ ]У, (13)

где нижний индекс «X» указывает на матрицы всего ансамбля КЭВ; У — полный вектор узловых перемещений пружины; ДУ — приращение вектора У; [Сх] — матрица демпфирования.

Вектор внешних нагрузок отсутствует в выражении (13), так как движение пружины задается кинематическим способом (см. рис. 1).

Размерности всех матриц и векторов в выражении (13) соответствуют набору из 38 КЭВ, в котором имеется 39 узлов по 6 степеней свободы в каждом. Таким образом, после составления ансамбля КЭВ модель имеет 234 степени свободы.

Учет закреплений и кинематического нагружения выполняли исключением фиксированных степеней свободы. В ведущем узле, вращение которого определяет двигатель (рис. 1, б), поворот вокруг оси шарнира задавали как функцию от времени, а остальные пять обобщенных перемещений — как нулевые значения. В свободном узле аналогичные пять обобщенных перемещений принимались нулевыми, а поворот вокруг оси шарнира остался обычной степенью свободы. Таким образом, после учета закреплений и кинематического возбуждения общее количество степеней свободы в модели пружины уменьшилось на 11 и стало равным 223. Строки и столбцы матриц в выражении (13), соответствующие фиксированным степеням свободы, были удалены.

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

[С £ ] = а[М£ ] + Р[К £ ],

где а = 0,78 рад/с; Р = 0,000925 с/рад.

Значения коэффициентов а и Р определены по методике работы [15], в которой коэффициент относительного демпфирования С, = 0,03.

Расчет частот собственных колебаний изогнутой пружины в случае отсутствия вращения.

При отсутствии вращения (оба захвата на рис. 1

Таблица 2

Значения частот собственных колебаний изогнутой пружины

Номер Расчет, Гц Эксперимент, Гц Погрешность, %

1 14,23 14,16 0,49

2 24,95 24,50 1,84

3 28,12 28,83 2,46

4 38,88 38,33 1,43

5 48,03 50,00 3,94

жестко зафиксированы) система (13) существенно упрощается и принимает вид

[Ms ]ДУ + [С£ ]ДУ + [К£ ]ДУ = 0. (14)

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

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

(-p2Ms + К s )ДУ = 0, (15)

где p — искомая круговая частота собственных колебаний; ДУ) — амплитудные значения приращений узловых перемещений относительно статического положения равновесия (форма колебаний).

Матрица демпфирования не учитывалась в выражении (15), так как ее влияние на собственные частоты весьма мало.

Для нахождения частот и форм к матрице [Ms]-1[Ks] применена процедура Eigensystem из пакета Wolfram Mathematica [14], предназначенная для вычисления собственных значений и собственных векторов квадратных матриц. Те же частоты определены экспериментально по методике работы [1]. Сопоставление результатов расчета и эксперимента для первых пяти частот собственных колебаний изогнутой пружины приведено в табл. 2.

Из табл. 2 следует, что результаты расчета собственных частот вполне удовлетворительно подтверждаются данными эксперимента.

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

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

[МЕ ]У + [С Е ]У = РЕ + .

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

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

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

На рис. 3 показаны экспериментальные и расчетные конфигурации пружины при первом и втором резонансных режимах.

Резонансная частота при наличии нелинейных эффектов определяется не очень четко, так как пружина как бы «вьется мелким бесом» и нахождение максимума амплитуды затруднено. В связи с этим для более надежного выделения резонансной частоты выполнен спектральный анализ зависимости перемещения от времени. Рассмотрена зависимость перемещения щ вдоль оси %2 (см. рис. 1) от времени для 18-го узла при первом резонансном режиме. Для наблюдений выбран узел с номером 18, так как он расположен наиболее далеко от захватов и поэтому имеет наибольшие перемещения. Угловая скорость вращения ведущего захвата

г -> 0,049 с

а

Рис. 3. Экспериментальная (слева) и расчетная (справа) конфигурации изогнутой пружины при первом (а) и втором (б) резонансных режимах

равнялась 89 рад/с, что соответствует первой собственной частоте / = 14,16 Гц, визуально (с использованием стробоскопа) наблюдаемой в эксперименте (см. табл. 2). Зависимость перемещения и2 от времени I приведена на рис. 4, а.

При спектральном исследовании непрерывную функцию и2(0 заменили массивом ее значений в равноотстоящих по времени точках, для чего последовательно с дискретным шагом по времени фиксировали значения = и(^). Полученный набор значений обработан дискретным преобразованием Фурье [16]. Результат преобразования Фурье представлен на рис. 4, б.

Как видно из рис. 4, б, максимальная амплитуда соответствует частоте /* = 14,01 Гц, которая несколько меньше /1. Расхождение /1 и /* можно объяснить нелинейными эффектами и влиянием демпфирования, а также неучтенным в данном расчете влиянием динамических эффектов на элементы матрицы жесткости [10].

Однако отмеченное расхождение небольшое и ожидаемое, поэтому можно утверждать, что спектральный анализ также подтверждает надежность разработанной методики.

Выводы

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

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

3. Результаты численного интегрирования нелинейных уравнений движения модели пру-

жины, построенной из разработанных элементов, практически не отличаются от данных натурных экспериментов. Совпадение наблю-

Литература

дается как по резонансным частотам, так и по конфигурациям пружины в резонансных режимах вращения.

[1] Бадиков Р.Н. Расчетно-экспериментальное исследование напряженно-деформирован-

ного состояния и резонансных режимов вращения винтовых пружин в пружинных механизмах. Дис. ... канд. техн. наук. Москва, МГТУ им. Н.Э. Баумана, 2009. 166 с.

[2] Журавлев В.Ф. Основы теоретической механики. Москва, Физматлит, 2001. 320 с.

[3] Попов В.В., Сорокин Ф.Д., Иванников В.В. Разработка конечного элемента гибкого

стержня с раздельным хранением накопленных и дополнительных поворотов для моделирования больших перемещений элементов конструкций летательных аппаратов. Тр. МАИ, 2017, № 92. URL: http://trudymai.ru/published.php?ID=76832 (дата обращения 6 ноября 2017).

[4] Сорокин Ф.Д. Особенности рационального способа описания больших поворотов в за-

дачах нелинейной динамики роторных машин. Нелинейная динамика машин. Матер. III Междунар. школы-конференции, Москва, 12-15 апреля 2016 г., Москва, ИМАШ РАН, 2016, с. 89-93.

[5] Жилин П.А. Векторы и тензоры второго ранга в трехмерном пространстве. Санкт-

Петербург, Нестор, 2001. 276 с.

[6] Жилин П.А. Рациональная механика сплошных сред. Санкт-Петербург, Изд-во Поли-

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

технического ун-та, 2012. 584 с.

[7] Елисеев В.В., Зиновьева Т.В. Механика тонкостенных конструкций. Теория стержней.

Санкт-Петербург, Изд-во Политехнического ун-та, 2008. 95 с.

[8] Сорокин Ф.Д. Прямое тензорное представление уравнений больших перемещений

гибкого стержня с использованием вектора конечного поворота. Известия Российской академии наук. Механика твердого тела, 1994, № 1, с. 164-168.

[9] Сорокин Ф.Д. Векторно-матричное описание больших поворотов при разработке ко-

нечных элементов механических систем. Компьютерные технологии анализа инженерных задач механики. Сб. тр. Междунар. науч. школы для молодежи, Москва, 913 ноября 2009 г., Москва, ИМАШ РАН, 2009, с. 106-113.

[10] Geradin M., Cardona A. Flexible multibody dynamics - A finite element approach. New York, Wiley, 2001. 327 p.

[11] Сорокин Ф.Д., Чжоу Су. Разработка конечного элемента в виде витка цилиндрической пружины. III Международная Школа-конференция молодых ученых Нелинейная динамика машин, Москва, 12-15 апреля 2016 г., Москва, ИМАШ РАН, 2016, с. 261265.

[12] Felippa C.A. A Systematic Approach to the Element-Independent Corotational Dynamics of Finite Elements. Department of Aerospace Engineering Sciences and Center for Aerospace Structures, USA, University of Colorado, 2000. 42 p.

[13] Felippa C.A., Haugen B. Unified Formulation of Small-Strain Corotational Finite Elements: I. Theory. Computer Methods in Applied Mechanics and Engineering, 2005, vol. 194, pp. 22852335.

[14] Дьяконов В.П. Mathematica 5.1/5.2/6. Программирование и математические вычисления. Москва, ДМК-Пресс, 2008. 574 с.

[15] Thomson W.T. Theory of Vibration with Applications. Cheltenham, Nelson Thornes, 2003. 546 p.

[16] Сергиенко А.Б. Цифровая обработка сигналов. Санкт-Петербург, Питер, 2002. 608 с.

References

[1] Badikov R.N. Raschetno-eksperimental'noe issledovanie napriazhenno-deformirovannogo sostoianiia i rezonansnykh rezhimov vrashcheniia vintovykh pruzhin v pruzhinnykh mekhanizmakh. Diss. kand. tekhn. nauk [Numerical and experimental investigation of

stress-strain state and the resonance modes of rotation of the helical spring in a spring. Cand. tech. sci. diss.]. Moscow, 2009. 166 p.

[2] Zhuravlev V.F. Osnovy teoreticheskoi mekhaniki [Foundations of theoretical mechanics].

Moscow, Fizmatlit publ., 2001. 320 p.

[3] Popov V.V., Sorokin F.D., Ivannikov V.V. Razrabotka konechnogo elementa gibkogo ster-

zhnia s razdel'nym khraneniem nakoplennykh i dopolnitel'nykh povorotov dlia modeliro-vaniia bol'shikh peremeshchenii elementov konstruktsii letatel'nykh apparatov [A flexible rod finite element with separate storage of cumulated and extra rotations for large displacements of aircraft structural parts modeling]. Trudy MAI [Trudy MAI]. 2017, no. 92. Available at: http://trudymai.ru/published.php?ID=76832 (accessed 6 November 2017).

[4] Sorokin F.D. Osobennosti ratsional'nogo sposoba opisaniia bol'shikh povorotov v zadachakh

nelineinoi dinamiki rotornykh mashin [Particularly rational way of describing large rotations in problems of nonlinear dynamics of rotor machines]. 3 Mezhdunarodnaia Shkola-konferentsiia Nelineinaia dinamika mashin [3 international School-conference on nonlinear dynamics of machinery]. Moscow, 12-15 April 2016, Moscow, IMASh RAN publ., 2016, pp. 89-93.

[5] Zhilin P.A. Vektory i tenzory vtorogo ranga v trekhmernom prostranstve [Vectors and second-

rank tensors in three-dimensional space]. Sankt-Petersburg, Nestor publ., 2001. 276 p.

[6] Zhilin P.A. Ratsional'naia mekhanika sploshnykh sred [Rational continuum mechanics].

Sankt-Petersburg, Politekhnicheskii un-t publ., 2012. 584 p.

[7] Eliseev V.V., Zinov'eva T.V. Mekhanika tonkostennykh konstruktsii. Teoriia sterzhnei [Me-

chanics of thin-walled structures. Theory of rods]. Sankt-Petersburg, Politekhnicheskii un-t publ., 2008. 95 p.

[8] Sorokin F.D. Priamoe tenzornoe predstavlenie uravnenii bol'shikh peremeshchenii gibkogo

sterzhnia s ispol'zovaniem vektora konechnogo povorota [Direct tensor representation of the equations of large displacement of the flexible rod by means of the nite rotation vector]. Izvestiia Rossiiskoi akademii nauk. Mekhanika tverdogo tela [Mechanics of Solids]. 1994, no. 1, pp. 164-168.

[9] Sorokin F.D. Vektorno-matrichnoe opisanie bol'shikh povorotov pri razrabotke konechnykh

elementov mekhanicheskikh system [Vector-matrix description of the big twists in the development of finite elements for mechanical systems]. Komp'iuternye tekhnologii analiza inzhenernykh zadach mehaniki. Sb. tr. Mezhdunar. nauch. shkoly dlia molodezhi [Computer technology of task analysis mechanical engineer. Proceedings of the International scientific school for young people]. Moscow, 9-13 November 2009, Moscow, IMASh RAN publ., pp. 106-113.

[10] Geradin M., Cardona A. Flexible multibody dynamics - A finite element approach. New York, Wiley, 2001. 327 p.

[11] Sorokin F.D., Su Chzhou. Razrabotka konechnogo elementa v vide vitka tsilindricheskoi pruzhiny [Development of a finite element in the form cylindrical coil springs]. 3 Mezhdunarodnaia Shkola-konferentsiia Nelineinaia dinamika mashin [3 International School-conference of young scientists Nonlinear dynamics of machinery]. Moscow, IMASh RAN publ., 2016, pp. 261-265.

[12] Felippa C.A. A Systematic Approach to the Element-Independent Corotational Dynamics of Finite Elements. Department of Aerospace Engineering Sciences and Center for Aerospace Structures, USA, University of Colorado, 2000. 42 p.

[13] Felippa C.A., Haugen B. Unified Formulation of Small-Strain Corotational Finite Elements: I. Theory. Computer Methods in Applied Mechanics and Engineering, 2005, vol. 194, pp. 22852335.

[14] D'iakonov V.P. Mathematica 5.1/5.2/6. Programmirovanie i matematicheskie vychisleniia

[Mathematica 5.1/5.2/6. Programming and mathematical calculations]. Moscow, DMK-Press, 2008. 574 p.

[15] Thomson W.T. Theory of Vibration with Applications. Cheltenham, Nelson Thornes, 2003. 546 p.

[16] Sergienko A.B. Tsifrovaia obrabotka signalov [Digital signal processing]. Sankt-Petersburg, Piter publ., 2002. 608 p.

Статья поступила в редакцию 30.11.2017

Информация об авторах

ИВАННИКОВ Владимир Вячеславович (Москва) — инженер отдела разработки математических алгоритмов и программного обеспечения. ИКЦ ООО «Альфа-транзит» (141400, Московская область, Химки, Ленинградская ул., д. 1, е-шаИ: vvivannikov@gmail.com).

СОРОКИН Федор Дмитриевич (Москва) — профессор кафедры «Прикладная механика». МГТУ им. Н.Э. Баумана (105005, Москва, Российская Федерация, 2-я Бауманская ул., д. 5, стр. 1, е-шаП: sorokin_fd@mail.ru).

ЧЖОУ Су (Москва) — аспирант кафедры «Прикладная механика». МГТУ им. Н.Э. Баумана (105005, Москва, Российская Федерация, 2-я Бауманская ул., д. 5, стр. 1, е-шаП: szleon602@gmail.com).

Information about the authors

IVANNIKOV Vladimir Vyacheslavovich (Moscow) — Engineer, Department of Mathematical Algorithms and Software Development. Engineering Consulting Centre — OOO Alfa-Tranzit (141400, Russian Federation, Khimky, Leningradskaya St., Bldg. 1, e-mail: vvivannikov@gmail.com).

SOROKIN Fedor Dmitrievich (Moscow) — Professor, Department of Applied Mechanics. Bauman Moscow State Technical University (105005, Moscow, Russian Federation, 2nd Baumanskaya St., Bldg. 5, Block 1, e-mail: sorokin_fd@mail.ru).

ZHOU Su (Moscow) — Postgraduate, Department of Applied Mechanics. Bauman Moscow State Technical University (105005, Moscow, Russian Federation, 2nd Baumanskaya St., Bldg. 5, Block 1, e-mail: szleon602@gmail.com).

Ю.в Никифоров, A.A. Казакова, М.6. Алехина

ПРОЦЕССЫ ДИФФУЗИИ И АДСОРБЦИИ В ИНЖЕНЕРНЫХ ЗАДАЧАХ

ПРИМЕРЫ РАСЧЁТА

В Издательстве МГТУ им. Н.Э. Баумана вышло в свет учебное пособие Ю.В. Никифорова, A.A. Казаковой, М.Б. Алехиной

«Процессы диффузии и адсорбции в инженерных задачах. Примеры расчета»

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

Для студентов МГТУ им. Н.Э. Баумана, обучающихся по направлениям подготовки 16.03.03 и 16.04.03 «Холодильная, криогенная техника и системы жизнеобеспечения», 16.05.01 «Специальные системы жизнеобеспечения».

По вопросам приобретения обращайтесь:

105005, Москва, 2-я Бауманская ул., д. 5, стр. 1. Тел.: +7 499 263-60-45, факс: +7 499 261-45-97; press@bmstu.ru; www.baumanpress.ru

i Надоели баннеры? Вы всегда можете отключить рекламу.