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

УДК 62-233.3

БОТ: 10.18698/2308-6033-2023-9-2306

Модель шлицевого соединения для задач роторной динамики газотурбинных двигателей

© И В. Николаев1, М.К. Леонтьев1, В В. Попов2

Московский авиационный институт, Москва, 125993, Россия 2МГТУ им. Н.Э. Баумана, Москва, 105005, Россия

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

Ключевые слова: роторная динамика, шлицевые соединения, центрирование валов, аналитическая модель, верификация

Введение. Шлицевые соединения являются одним из наиболее распространенных способов передачи крутящего момента в роторах авиационных газотурбинных двигателей (ГТД). Под действием динамических сил, возникающих в роторных системах, в шлицевых соединениях часто появляются перекосы, которые, в свою очередь, меняют и упругие характеристики шлицевого соединения — радиальную и угловую жесткость. Эти изменения приводят к изменениям и динамических характеристик всей роторной системы. С учетом этого задача о влиянии шлицевых характеристик на динамику роторной системы является весьма актуальной.

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

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

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

Разработанная модель шлицевого соединения может быть включена в состав всей роторной модели, что позволяет провести ее анализ не только с точки зрения механических характеристик (распределение контакта, давление, зазоры между шлицами), но и судить о влиянии соединения на собственные частоты системы в условиях действующих динамических сил и перекоса [4-6].

Математическая модель шлицевого соединения. Движение любой механической системы подчиняется известному принципу наименьшего действия Гамильтона [7], вследствие чего это движение можно описать уравнением Лагранжа второго рода:

где Ь — лагранжиан системы; ^ — время; у, у — векторы обобщенных координат и их первых производных по времени; Жгх1 — вектор обобщенных внешних сил; ^с1атр — вектор диссипативных сил; Ж— вектор сил трения.

При раскрытии лагранжиана как разницы между кинетической и потенциальной энергией ( Ь = Т - и ) уравнение (1) примет вид

(1)

— —---1--= Ж - Ж + Ж

я \ду) су су 6x1 "атр "

й (СТ Л ст си

(2)

Уравнение (2) описывает систему в динамическом равновесии. С учетом этого левую часть можно представить как сумму инерционных и внутренних сил:

F + F = F - F + F^ '

inerc inner ext damp fric >

(3)

% J öy

F. = —

inerc

dt

öU

Fn

öy

inner

Здесь ГШегс — вектор инерционных сил; ГШег — вектор внутренних сил.

Если в начальный момент времени векторы обобщенных координат и их производные имеют вид у, у, у, а в следующий момент они получают приращения у + Лу, у + Лу, у + Лу, то можно провести линеаризацию левой части уравнения (3) в момент существования у, у, у. С учетом этого общее уравнение движения примет вид

[М] Лу + ([С] + [О]) Лу + [К] Лу = Гех1 - ГШег - ¥ШеГс - ¥аЬтр + Г^, (4)

—р

где [М] — матрица масс, [М] = ——у^; Лу, Лу, Лу — векторы приращений обобщенных координат и их первых и вторых производных по времени; [С] — матрица демпфирования, [С] = —-^ттР; [О ] — матрица инерционных и гироскопических эффектов

—Г -Г

[О ] = ; [к ] — матрица жесткости, [К ] = .

Матрицы и векторы, описывающие инерционные и диссипатив-ные свойства соединения, могут быть вычислены так же, как в работах [8, 9]. Уравнение (4) является исходным для применения неявных методов интегрирования, таких как метод Ньюмарка [10].

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

ортов — е^, е^, еЗ1, и ef, е^, е| соответственно. Орты е^ и е|

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

Движение системы происходит вследствие воздействия внешних нагрузок, задаваемых вектором сил Гех;. В процессе движения каждый

из центров приобретает векторы линейных перемещении иА и ив, тензоры больших накопленных поворотов Rл и Rв и малые векторы Эйлера дополнительного поворота &л и [8].

Таким образом, вектор состояния соединения содержит 12 степеней свободы — по три линейных перемещения и три компонента вектора малых дополнительных поворотов отдельно для вала и втулки. В процессе движения, вследствие упругости зацепления, образуются линейные смещения центра втулки относительно вала ( Ли ^ на переднем торце и Ли6 на зад-

\ А В

нем торце), а также угловое смещение между ортами е1 и е1 :

Рис. 1. Исходное положение шлицевого соединения

Ли

Г -

1 Ео ■ еВ

2

1 Ео ■ еВ

2

А

Лф - агссоБ

(Ео■еВ

Е

о 'С1

(5)

(6)

(7)

где Ео - Е - ез ® ез — плоский единичный тензор; Е — единичный тензор; I — длина шлицев; Лф — угол относительного поворота.

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

Рис. 2. Система координат для шлица

Рис. 3. Пересечение геометрических профилей шлицев на торцах

хе = г I

(8)

Поскольку профили шлицев — эвольвенты, их геометрия описывается следующими соотношениями:

■ (sin ф - ф cos а cos (а + ф));

ye = r (cos ф + ф cos а sin (а + ф)).

Здесь xe, ye — координаты профиля; r — радиус делительной окружности шлицев; а — угол давления; ф — параметр эвольвенты, принимающий значения < ф < фтах.

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

(

фт = arctan

ГХ0 cos а + y0yjх^ + y^ -r2 cos2 а

Л

2 Г 2 : 2 2 2 )Vx

+ а.

(9)

rx0y0cos a-x0^ x0 + y0 -r cos a

Здесь x0, y0 — координаты эвольвентного профиля в системе координат, начало которой находится в центре шлицевого соединения.

Величину взаимопроникновения в точке контакта находим как разницу:

Í

Ь0(Дф) • L0 | |. L0(v) •

с

2га

Lo| у +

fx (ф )1 ^

e (фтах )[+Дц

tye (фтах ) J Хе ( фт ) 1

ye (9m )Jy

(10)

Здесь Ь0 — матрица поворота; г — количество шлицев; I — номер шлица, для которого определяется взаимопроникновение; у — угол, соответствующий половине толщины зуба на делительной окружности; А — зазор по боковым поверхностям между шлицами; Аи — относительное смещение.

Значение параметров 5 и фт определялись из решения задачи о минимуме расстояния между шлицами с помощью метода Ньютона:

5 = min

Хмф (фмф) - Хв (фв )) - (Умф (фмф) - Ув (фв ))

(11)

где индекс «мф» указывает на внутреннюю часть шлицевого соединения, индекс «в» — на внешнюю часть.

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

5 =

параметр характеризует зазор между поверхностями шлицев внутренней и внешней части.

В результате взаимопроникновения каждая пара зуб—углубление приобретает упругую энергию:

^ш.п

0,

кш.п ( 5)

2

5 > 0;

5 < 0,

(12)

где кшп — жесткость шлицевой пары, определяемая как величина,

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

к -ш.п

1

с ■ с А . с Б ^конт °изг ^изг

(13)

Контактная податливость рассчитывается согласно теории линейного контакта двух упругих цилиндров:

с 4/ (14)

(1 +1 ^

л

Е

Е2

у

где |1, |2 — коэффициенты Пуассона внутренней и внешней частей шлицев соответственно; Еь Е2 — модуль упругости внутренней и внешней частей шлицев.

Для прямозубых зацеплений наиболее распространенным методом определения изгибной податливости являются балочная модель и взятие интеграла Мора [11, 12]. В связи с тем, что в эвольвентном шлицевом соединении линия зацепления в общем случае не сонаправ-лена с осью шлицев и действует не по всей ширине шлица, изгибную податливость рассчитывают более сложным способом.

Шлиц рассматривается как толстостенная пластина (рис. 4).

Рис. 4. Схема шлица как толстостенной пластины

Толщина зуба рассчитывается так:

h = 2 (ye sin y-xe cos у). (15)

В классической теории изгиба тонких пластин, построенной на гипотезах Кирхгофа, упругая энергия пластины записывается следующим образом [13]:

Í я2,,, .,2 -2 ^

1 г г , , О W о w о w

п = -ll Mx-т- + M -+ 2M -

п , J J x 2 y ~ 2 ХУ

СхСу, (16)

Л.2 У ^ 2 ху ^ л.

^ дх ду дхду)

где Мх, Му, Мху — изгибающие и крутящий моменты; н — поперечное перемещение; х, у — координаты в плоскости пластины.

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

2 \\кОЪу2Схйу, (17)

где G — модуль сдвига; у — деформация поперечного сдвига; к — геометрический коэффициент для сдвига.

В теории толстых пластин, построенной на гипотезе Тимошенко, упругая энергия является комбинацией (16) и (17) с учетом сдвиговых добавок в углах поворота:

1 утах 1

ип = - | |(к^у^ +кхуОку2ху + Дх2 + 2щхХу + 2(1 -ц)хХу))СхСу. (18)

2 о о

О-Эх о^ y 1 Í0-9x о^ y ^

Здесь X x = , X y , Xxy =

- + -

кривизны и круче-

^х ' Лу ^ ' Гчху _ ~ -

дх дУ 2 ^ ду дх

Он _ _ дн

~ ^ ^х, уху _ ~

дх ду

ции поперечного сдвига; кХ2 и ку — геометрические коэффициенты для сдвига (оба принимаются по 6/5); ц — коэффициент Пуассона;

ние пластины; yxz =--+dx, yxy =--l-$y — осредненные деформа-

D — цилиндрическая жесткость, D = -

Eh3

12 (1V У

Для расчета пластины воспользуемся методом Ритца. Перемещение и повороты разыскиваются в виде полиномов от координат:

н _ а0у + а1 ух + а2ух2 +а3у2 +а4у2х + а5у2х2;

2 2 2 2 2 _а6у + а7ух + а8ух +а9у +а10у х + а11у х ; (19)

^у _а12у + а1зух + а14ух2 +а15у2 у2х + а17у2х2, где а0,..., а17 — искомые коэффициенты полиномов.

После подстановки выражений (19) в (18) и вычисления интегралов с помощью компьютерной аналитики энергия деформаций пластины превращается в квадратичную форму коэффициентов а0,..., а17.

Для реализации метода Ритца кроме энергии деформаций ип необходим еще потенциал внешних сил Уп, который вычисляется следующим образом. Для каждого шлица зона контакта определяется путем линейного соединения точек взаимопроникновения профилей. В результате на боковой поверхности шлица образуется линия контакта. Зона контакта разделяется на п участков (было принято п = 10), для которых записывается потенциал внешних сил ¥}, действующих на шлиц (рис. 5) в точке контакта с номером у.

К (, ^ ). (2°)

]=1

Решение задачи об изгибе шлица как толстой пластины получается путем минимизации суммы (18) и (20) относительно коэффициентов а°,..., а17 с помощью метода Ньютона.

Линия контакта

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

п F■ 1

^изг = ^ Т У ^изг = ~, . (21)

>1 ^ (xj,Уj) ^изг

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

и = 1^,. (22)

г=1

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

К„

д 2и

дУг дУi

ди

дУ г ''

(23)

где у;, у| — компоненты вектора обобщенных координат у.

Учет сил трения в контакте проводится на основе оценки относительной скорости скольжения профилей относительно друг друга (рис. 6).

Точка контакта

Рис. 6. Кинематическая схема для расчета относительной скорости скольжения

Вектор сил трения имеет вид

F

тр I-

- { N ■

J тр г

(V 2? - V В)

V 2? - V В

(24)

Здесь / — коэффициент трения; Ыг — нормальная сила; ¥2 —

= ¥а + ша у2 =Ув + шВ хг[ — скорость скольжения в точке

контакта (V?, Vв — скорость в центре масс внешней и внутренней части, тх, Г2 — радиус-вектор от центра масс до точки определения силы трения, ш а , ш В — угловые скорости внешней и внутренней частей шлицев).

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

Нагружение системы происходит по внутренней поверхности шлицев через мастер-узел, соединенный с узлами поверхности нагружения ЯБЕ-3 элементами.

Таблица 1

Параметры модели шлицев

Параметр Значение Параметр Значение

Модуль упругости Е, МПа 210 000 Наружный диаметр кольца, мм 158.5

Коэффициент Пуассона д 0,3 Внутренний диаметр кольца, мм 135

Плотность р, кг/м3 7850 Длина шлицев 1, мм 55

Модуль шлицев т, мм 2,5 Зазор по боковым поверхностям, мкм 160

Число шлицев г 52

Рис. 7. Общий вид модели

Для верификации рассматривались два варианта нагружения — силой в 1 Н и 100 Н и изгибающим моментом, равным 1 Н-м и 100 Н-м, в радиальных направлениях. Жесткостные характеристики соединения, полученные при использовании двух подходов, представлены в табл. 2. По результатам можно заключить, что разница значений между аналитической моделью и КЭ-моделью не превышает 2 %.

Таблица 2

Результаты расчетов аналитической модели и КЭ-модели шлицевого соединения

Нагрузка Радиальная жесткость соединения, Н/м Угловая жесткость соединения, Н-м/рад

КЭ-модель Аналитическая модель КЭ-модель Аналитическая модель

1 Н, 1 Н-м 2,23-1010 2,261-1010 1,17-107 1,194-107

100 Н, 100 Н-м 3,182-1010 3,241-1010 5,47-107 5,58-107

В шлицевом соединении при появлении перекосов можно наблюдать эффект центрирования при приложении крутящего момента (рис. 8). Описание данного эффекта подробно представлено в работах [6, 14]. Возникновение центрирования также подтверждается результатами конечно-элементного анализа шлицевого соединения (см. рис. 8).

Рис. 8. Процесс центрирования шлицевого соединения с перекосом под действием крутящего момента

Поскольку данный эффект может оказывать влияние на динамическое поведение системы, была проведена верификация эффекта центрирования путем сравнения двух моделей с перекосом и дальнейшим нагружением крутящим моментом. Рассматривались радиальный и угловой перекосы величиной 120 мкм и 0,07° соответственно. Согласно результатам, представленным на рис. 9 и 10, можно заключить, что с возрастанием крутящего момента в соединении происходит уменьшение величины несоосности (эффект центрирования) и соединение возвращается в обычное состояние. Следует также отметить, что разница между аналитической и конечно-элементной моделями не превышает 7 %.

Рис. 9. Результаты моделирования процесса центрирования при радиальном перекосе аналитической и КЭ-модели

— мкэ

— Модель

30 Мкр, Н м

Рис. 10. Результаты моделирования процесса центрирования при угловом перекосе аналитической и КЭ моделей

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

Верификация аналитической модели показала хорошую сходимость с результатами КЭ-анализа. Разница в жесткостных характеристиках не превышала 2 %. Приведены также результаты по центрированию соединения вследствие действия крутящего момента при перекосе. Наибольшая разница между характеристиками центрирования составила 7 %.

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

ЛИТЕРАТУРА

[1] Dupertuis C., Ligier J.-L. Contact pressure in misalignment spline couplings. Mechanics & Industry, 2020, vol. 21, art. no. 510, pp. 1-12.

[2] Xiangzhen Xue, Qixin Huo, Jian Liu, Jipeng Jia. Nonlinear dynamic load analysis of aviation spline coupling with mass eccentricity and misalignment. Advances in Mechanical Engineering, 2021, vol. 13 (2), pp. 1-19.

[3] Hong J., Talbot D., Kahraman A. A stiffness formulation for spline joints. Journal of Mechanical Design, April 2016, vol. 138, art. no. 043301 (pp. 1-8).

[4] Yunbo Hu, Huibin Wang, Yuanqiang Tan, Likuan Jiang, Shengqiang Jiang. Study on the meshing force of misaligned gear coupling. In: International Forum on Mechanical, Control and Automation (IFMCA 2016). Advances in Engineering Research, vol. 113. Atlantis Press, 2017.

URL: https://www.atlantis-press.com/proceedings/ifmca-16/25871874

[5] Cura F., Mura A. Analysis of a load application point in spline coupling teeth. Journal of Zhejiang University SCIENCE A, 2014, vol. 15 (4), pp. 302-308.

[6] Cura F., Mura A. Experimental and theoretical investigation about reaction moments in misaligned splined couplings. Mechanical Systems and Signal Processing, 2014, vol. 45, pp. 504-512.

[7] Мусалимов В.М., Сергушин П.А. Аналитическая механика. Уравнение Ла-гранжа второго рода. Свободные колебания. Санкт-Петербург, СПбГУ ИТМО, 2007, с. 53.

[8] Попов В.В., Сорокин Ф.Д., Иванников В.В. Конечный элемент гибкого стержня с раздельным хранением накопленных и дополнительных поворотов для задач нелинейной динамики конструкций летательных аппаратов. Труды МАИ, 2018, № 99. URL: https://trudymai.ru/published.php?ID=91790

[9] Попов В.В., Сорокин Ф.Д., Иванников В.В., Дегтярёв С.А. Разработка пространственной модели зубчатой передачи с раздельным хранением накопленных и дополнительных поворотов для решения нелинейных задач динамики авиационных трансмиссий. Труды МАИ, 2020, № 112.

DOI: 10.34759/trd-2020-112-7

[10] Geradin M., Cardona A. Flexible Multibody Dynamics. John Wiley & Sons, Ltd, 2001, pp. 31-33.

[11] Wan Z., et al. An improved time-varying mesh stiffness algorithm and dynamic modeling of gear-rotor system with tooth root crack. Engineering Failure Analysis, 2014, vol. 42, pp. 157-177.

[12] Karpat F., et al. Effect of rim thickness on tooth root stress and mesh stiffness of internal gears. In: ASME 2014 International Mechanical Engineering Congress and Exposition. Montreal, Canada, 2014, 7 p.

[13] Образцов И.Ф. Строительная механика летательных аппаратов. Москва, Машиностроение, 1986, с. 114-117.

[14] Luciano D.I. Analisi F.E.M. su accoppiamenti scanalati: Tilting Moment e Friction Moment. Politecnico di Torino, 2020.

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

Ссылку на эту статью просим оформлять следующим образом:

Николаев И.В., Леонтьев М.К., Попов В.В. Модель шлицевого соединения для задач роторной динамики газотурбинных двигателей. Инженерный журнал: наука и инновации, 2023, вып. 9. http://dx.doi.org/10.18698/2308-6033-2023-9-2306

Николаев Илья Витальевич — аспирант кафедры «Конструкция и проектирование двигателей» Московского авиационного института. e-mail: nikolka831@mail.ru

Леонтьев Михаил Константинович — д-р техн. наук, профессор кафедры «Конструкция и проектирование двигателей» Московского авиационного института. e-mail: lemk@alfatran.com

Попов Валерий Васильевич — старший преподаватель кафедры «Прикладная механика» МГТУ им. Н.Э. Баумана. e-mail: vvppov.bmstu@gmail.com

Spline coupling in the rotor dynamics problems of a gas turbine engine

© I.V. Nikolaev1, M.K. Leontiev1, V.V. Popov2

1Moscow Aviation Institute (National Research University), Moscow, 125993, Russia 2Bauman Moscow State Technical University, Moscow, 105005, Russia

Involute spline couplings are widely used in the aircraft gas turbine engine rotors. In the rotor dynamics problems, they are usually simulated by a link connecting the shafts. Rigidity characteristics of such a connection depend on a set offactors: forces acting from the connected shafts section leading to the shafts distortion in the spline coupling; torque transmitted by the splines; gaps between the splines; friction forces in the contact pairs. Taking into account all these factors is quite labor-intensive even in the finite element formulation. In addition, spline coupling is a nonlinear element, which characteristics depend on the engine operating mode. The paper presents an analytical model making it possible to account for all operating factors of a spline coupling. The model was verified using the finite element calculations.

Keywords: rotor dynamics, spline coupling, shaft alignment, analytical model, verification

REFERENCES

[1] Dupertuis C., Ligier J.-L. Contact pressure in misalignment spline couplings. Mechanics & Industry, 2020, vol. 21, art. no. 510, pp. 1-12.

[2] Xiangzhen Xue, Qixin Huo, Jian Liu, Jipeng Jia. Nonlinear dynamic load analysis of aviation spline coupling with mass eccentricity and misalignment Advances in Mechanical Engineering, 2021, vol. 13 (2), pp. 1-19.

[3] Hong J., Talbot D., Kahraman A. A stiffness formulation for spline joints. Journal of Mechanical Design, April 2016, vol. 138, art. no. 043301 (1-8 p.)

[4] Yunbo Hu, Huibin Wang, Yuanqiang Tan, Likuan Jiang, Shengqiang Jiang. Study on the meshing force of misaligned gear coupling. International Forum on Mechanical, Control and Automation (IFMCA 2016). Advances in Engineering Research, vol. 113. Atlantis Press, 2017. Available at: https://www.atlantis-press.com/proceedings/ifmca-16/25871874

[5] Cura F., Mura A. Analysis of a load application point in spline coupling teeth. Journal of Zhejiang University SCIENCE A, 2014, no. 15 (4), pp. 302-308.

[6] Cura F., Mura A. Experimental and theoretical investigation about reaction moments in misaligned splined couplings. Mechanical Systems and Signal Processing, 2014, no. 45, pp. 504-512.

[7] Musalimov V.M., Sergushin P.A. Analiticheskaya mekhanika. Uravnenie Lagranzha vtorogo roda. Svobodnye kolebaniya [Analytical mechanics. Lagrange equation of the second kind. Free oscillations]. Saint Petersburg, 2007, SPbGU ITMO Publ., p. 53.

[8] Popov V.V., Sorokin F.D., Ivannikov V.V. Konechnyi element gibkogo ster-zhnya s razdelnym khraneniem nakoplennykh i dopolnitelnykh povorotov dlya zadach nelineynoy dinamiki konstruktsiy letatelnykh apparatov [Finite element of a flexible rod with separate storage of accumulated and additional rotations for the problems of nonlinear dynamics of aerial vehicles structures]. Trudy MAI, 2018, no. 99. Available at: https://trudymai.ru/published.php?ID=91790

[9] Popov V.V., Sorokin F.D., Ivannikov V.V., Degtyarev S.A. Razrabotka pros-transtvennoy modeli zubchatoy peredachi s razdelnym khraneniyem nakoplen-nykh i dopolnitelnykh povorotov dlya resheniya nelineynykh zadach dinamiki aviatsionnykh transmissiy [Developing a spatial model of a gear transmission with separate storage of the accumulated and extra turns for solving non-linear problems of aviation transmissions dynamics]. Trudy MAI, 2020, no. 112. DOI: 10.34759/trd-2020-112-7

[10] Geradin M., Cardona A. Flexible Multibody Dynamics. John Wiley & Sons, Ltd, 2001, pp.31-33.

[11] Wan Z., et al. An improved time-varying mesh stiffness algorithm and dynamic modeling of gear-rotor system with tooth root crack. Engineering Failure Analysis, 2014, vol. 42, pp. 157-177.

[12] Karpat F., et al. Effect of rim thickness on tooth root stress and mesh stiffness of internal gears. In: ASME 2014 International Mechanical Engineering Congress and Exposition. Montreal, Canada, 2014, 7 p.

[13] Obraztsov I.F. Stroitelnaya mekhanika letatelnykh apparatov [Construction mechanics of the air vehicles]. Moscow, Mashinostroenie Publ., 1986, pp. 114-117.

[14] Luciano D.I. Analisi F.E.M. su accoppiamenti scanalati: Tilting Moment e Friction Moment. Politecnico di Torino, 2020.

Nikolaev I.V., Postgraduate, Department of Engine Design and Engineering, Moscow Aviation Institute (National Research University). e-mail: nikolka831@mail.ru

Leontiev M.K., Dr. Sc. (Eng.), Professor, Department of Engine Design and Engineering, Moscow Aviation Institute (National Research University). e-mail: lemk@alfatran.com

Popov V.V., Senior Lecturer, Department of Applied Mechanics, Bauman Moscow State Technical University. e-mail: vvppov.bmstu@gmail.com