Численное определение температурного поля в системе подшипников на общем валу с учетом скорости его вращения

Математические заметки СВФУ Июль—сентябрь, 2014. Том 21, № 3

УДК 621.89:536.24

ЧИСЛЕННОЕ ОПРЕДЕЛЕНИЕ ТЕМПЕРАТУРНОГО ПОЛЯ В СИСТЕМЕ ПОДШИПНИКОВ НА ОБЩЕМ ВАЛУ С УЧЕТОМ СКОРОСТИ ЕГО ВРАЩЕНИЯ Р. С. Тихонов, Н. П. Старостин

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

Ключевые слова: подшипник скольжения, математическая модель, трение, теплопроводность, температура, уравнение теплопроводности, численное решение.

Введение

Математическое моделирование нестационарного теплового процесса в системе подшипников на общем валу используется для решения различных прикладных задач, в том числе в методе тепловой диагностики трения, позволяющего определять фрикционные тепловыделения и, соответственно, моменты трения по замерам температур путем решения граничной обратной задачи [1]. Температурное поле в системе подшипников ранее определялось при допущении однородного распределения температуры по сечению вала вследствие достаточно высоких скоростей вращения (более 5 рад/с) и рассмотрения вала как одномерного стержня [1—4]. Скорости вращения вала в подшипниках скольжения, при которых необходимо учитывать движение вала при расчете температурного поля, определялись в работе [5]. Влияние вращения вала на температурное поле учитывалось для одного подшипника скольжения в плоском случае в работе [6], в трехмерном случае — в работе [7]. В данной работе рассматривается задача определения нестационарного температурного поля в системе подшипников на общем валу, вращающемся с невысокой скоростью, при которой необходимо учитывать движение вала (рис. 1).

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

Работа выполнена при финансовой поддержке Министерства образования и науки Российской Федерации в рамках государственного задания на выполнение НИР на 2014—2016 гг. (проект №3047).

© 2014 Тихонов Р. С., Старостин Н. П.

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

Рис. 1. Расчетная схема системы подшипников скольжения: 1 — вал; 2 — вкладыш (втулка); 3 — обойма.

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

Нестационарное температурное поле в подшипниках описывается двумерным квазилинейными уравнениями теплопроводности для втулок с корпусами:

дТ 1 д ( дТк дЬ г дг V дг

1 д Л дТ\ + — — ( — I

о I 'чк I

г2 дф \ дф у

Л2к < г < Я4к, —п < ф < п, 0 < Ь < Ьг

к = г = 2,3,

а для вала — трехмерным уравнением с конвективным членом, учитывающим его вращение:

С1

ди 1 д

дЬ

ди_\ J_д_

г дг \ дг ) г2 дер

0 < г < Я1, —п < ф < п, 0 < Ь < Ьш.

ди д ди

+ +

В зонах трения втулок с валом задаются условия фрикционного тепловыделе-

А1

ди

дг

— Л2

дТк

г=К 1

дг

= Як(ф,Ь), |ф| < фо,

Г = В-2к

1

dk

Zfc

J U(Ri,p,z,t) dz = Tk(R2k,p,t).

На свободных поверхностях вала, втулок и обойм задаются условия конвективного теплообмена:

А

dU

A2k

dr

dTk

-ai(U(Ri,p,z,t) - To),

r=Ri

A^k

dr

dTk

dr

= a2(Tk(R2k,p,t) - To), |p| > tpo

r=R2

= аз(Tk(R4k,p,t) - To), -n < p < n.

r=R4k

На концах вала задаются условия I и III рода: dU

Ai

dz

= -ai

z=L

В центре вала задается условие ограниченности теплового потока: По угловой координате выполняются условия периодичности

dTk

dp

dU

dTk

dp

п dp

dU

Tk(r, n, t) = Tk(r,n,t),

dp

U

Начальные распределения температур в элементах узла трения примем равными и однородными:

Тк

1. Численное решение

Задача (1)—(11) решалась методом конечных разностей, позволяющим решать уравнения с переменными коэффициентами, нелинейные уравнения и широко используемым для решения различных прикладных задач [9-15]. Наличие конвективного члена в уравнении теплопроводности (2), учитывающего вращение вала, приводит к определенным трудностям при численном решении поставленной задачи. Использование монотонных локально-одномерных разностных схем суммарной аппроксимации позволяет обеспечить выполнение принципа максимума, т. е. при любых шагах т и ¡г^ по временной и угловой переменным может быть определено приближенное решение. Из множества приближенных решений выберем решение, удовлетворяющее условию «слипания» при стремлении нулю шага по времени. В дальнейшем разработанный алгоритм численного определения температурного поля намечается использовать для решения краевых задач при реализации алгоритма решения обратной задачи теплообмена по определению фрикционного тепловыделения. Поэтому затраты машинного

Zk-1

п

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

Пусть имеется достаточно детальная пространственная сетка. Для этой пространственной сетки вычислительными экспериментами определим шаг по времени, обеспечивающий сходимость решения. Варьировалось число Куранта (7 = ти/Нр, и = Ях характеризующее связь шага по угловой переменной со скоростью движения вала и шагом по времени. Поскольку для всех подшипников скольжения вращающийся вал является общим, достаточно рассмотреть случай одного подшипника. В последующем, найденный таким образом шаг по времени будет использован для расчетов температур в системе с нескольких подшипников.

т,°с

38

24 -1-р-г---г--1-1

О 150 300 450 600 750 t, сек

Рис. 2. Расчетные зависимости максимальных температур в зоне трения при различных числах Куранта 7: 1 — Y = 36; 2 — Y = 12; 3 — 7 = 2; 4 — 7 = 1; 5 — 7 = 1/8.

На рис. 2 представлены результаты расчетов температур в подшипнике скольжения в зависимости от числа Куранта при следующих геометрических размерах: Яхк = 12 мм, Я2к = 13 мм, Я3к = 16 мм, Я4к = 12 мм, к =1. Материалом для вала и обоймы служит сталь, а для втулки — наполненный фторопласт Ф4К20. Скорость вращения вала п рад/с, угол контакта — 30°. Удельная интенсивность тепловыделения постоянна и равна ^ = 67кВт/м2. Расчеты показывают, что при 7 < 1 решения сходятся. Для практических расчетов шаг по времени можно определить из условия 7 = 2, поскольку при 7 < 2 значения температур меняются в пределах одного градуса.

2. Сопоставление расчетных и экспериментальных температур

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

Т,°с -| 42 -40 -38 -36 -34 -32 -30 -28 -26 -

24 -I................-......т......................»........................................................................т.............................................,.......................,...........................г...........................,..................- ................т

0 300 600 900 1200 1500 1800 2100 2400 2700 3000 3300 1, сек

Рис. 3. Расчетная зависимость температуры от времени во внутренней точке втулки на расстоянии 0,5 мм от зоны трения.

(Значком I отмечены доверительные интервалы экспериментальных температурных данных)ю

Испытания проводились на машине трения СМТ-1 для одного подшипника. Температуры регистрировались с помощью медь-константановых термопар диаметром 0,1 мм с использованием многоканального устройства «ТЕРМОДАТ» в 5 точках втулки на расстоянии 0, 5 мм от зоны трения. Угол контакта составлял 60°. Трехкратное повторение эксперимента с регистрацией температур позволило определить доверительный интервал для зависимости температуры от времени. На рис. 3 показано, что расчетная зависимость температуры при р = 0 лежит в пределах разброса экспериментальных данных. Подобные результаты получены и для других точек замера температур, что свидетельствует об адекватном описании теплового процесса в рассматриваемом подшипнике скольжения.

В результате расчетов определены эффективные коэффициенты в уравнении теплопроводности, обеспечивающие близость расчетных и экспериментальных температурных данных. При этом эффективные коэффициенты отличаются от значений справочных теплофизических характеристик материалов на величину, не превышающем разброс свойств композиционного полимерного материала и стали. Эффективный коэффициент теплопроводности Л2 наполненного фторопласта Ф4К20 составил 0, 34 Вт/(м • °С), объемная теплоемкость С2 равна 2, 2 • 106 Дж/(м3 • °С), для стали Лх = 46 Вт/(м • °С), С =3,48 • 106 Дж/ (м3 • ° С). Коэффициент теплообмена вращающего вала вычислялся по фор-

муле [16]:

«1 Nu = 0, 95

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

3. Моделирование теплового процесса в системе подшипников

Разработанный алгоритм решения задачи был обобщен для систем подшипников скольжения. В качестве примера рассматривалась система из четырех одинаковых подшипников из наполненного фторопласта Ф4К20. Интенсивности тепловыделения в подшипниках задавались в виде функций времени:

Qi((p,t) = 3851,656(£ + l)i|cos<p|, Q2(<P,t) = 215'89"^оооо' | cos<p|,

л (t-600)2

Q3{<P,t) = 3851, 656тг| COS Q4{(p,t) = 124'3194" зеоооо I costal.

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

20 ............. ,-1-г—-1-1-1-1-1-1

0 200 400 600 800 1000 1200 1400 1600 t, сек

Рис. 4. Зависимости максимальных температур в зоне контакта подшипников скольжения от времени: 1 — температура 1-го подшипника, 2 — 2-го подшипника, 3 — 3-го подшипника, 4 — 4-го подшипника.

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

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

Т.'С

26,5 H-1-1-1-1-1-1-1-1-1-1-1-1

-180 -150 -120 -90 -60 -30 0 30 60 90 120 150 ф.град

Рис. 5. Распределения температур по поверхности вала при различных скоростях его вращения: 1 — 0,1п рад/с; 2 — 0,3п рад/с; 3 — 0,5п рад/с; 4 — 0,7п рад/с; 5 — п рад/с; 6 — 1,4п рад/с; 7 — 2п рад/с

На незначительном расстоянии от подшипника скольжения (2-4 мм) по осевой переменной распределение температуры в вале по окружности становится однородным. В связи с этим для исследования однородности температуры по окружной переменной в системе подшипников скольжения рассматривался один подшипник с исходными данными, принятыми выше. В проведенных расчетах при различных скоростях вращения функция удельной интенсивности тепловыделения оставались неизменными в силу введения условия pRlQ = pV = const. Максимальное и минимальное значения температуры в зоне трения вала достигаются в точке выхода поверхности вала из контакта и в точке входа в контакт. Расчетные зависимости температур поверхности вала от угловой координаты приведены на рис. 5. Несмотря на повышение коэффициента теплообмена с увеличением скорости вращения после выхода из контакта поверхность вала охлаждается меньше и минимальная температура контакта повышается. Кроме того, за счет сокращения времени контакта максимальная температура на поверхности вала снижается. Таким образом, увеличение скорости приводит к однородности распределения температуры по окружности.

Т/с 1 1.6 1,2 0,8 0,4

3 15 27 39 51 О, об/мин

Рис. 6. Зависимость разности максимальных и минимальных значений температур на поверхности вала от скорости вращения вала в момент времени £ =1 мин

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

4. Выводы

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

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

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

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

Обозначения: Qk — удельная интенсивность тепловыделения в зоне контакта к-го подшипника, и — температура вала, Тк — температура к-го подшипника, То — температура окружающей среды, к, N — номер и количество подшипников, Н1 — радиус вала, Н2к,Нзк — внутренний и внешний радиусы втулки к-го подшипника, Н^к — внешний радиус обоймы к-го подшипника, Зк — длина к-го подшипника, Ь — длина вала, т,(р,г — цилиндрические координаты, 2к-1,2к — осевые координаты начала и конца расположения к-го подшипника, 2<^о — угол контакта, О — угловая скорость, Ь — время, Ьт — время испытания, С1 — объемная теплоемкость материала вала, С2к,С33к — объемная теплоемкость материалов втулки и обоймы к-го подшипника соответственно, а1,а2,а33 — коэффициенты теплоотдачи с поверхности вала, втулки и обоймы соответственно, Л1 — коэффициент теплопроводности материала вала, Л2к, Лзк — коэффициенты теплопроводности материалов втулки и обоймы к-го подшипника соответственно, Лвоз — коэффициент теплопроводности окружающей среды, Не — критерий Рейнольдса, Сг — критерий Грасгофа, р — давление, V — линейная скорость.

ЛИТЕРАТУРА

1. Старостин Н. П., Тихонов А. Г., Моров В. А., Кондаков А. С. Расчет триботехнических параметров в опорах скольжения. Якутск: Изд-во ЯНЦ СО РАН, 1999.

2. Черский И. Н., Богатин О. Б., Старостин Н. П. Восстановление моментов трения в системе несмазываемых подшипников по замерам температуры // Трение и износ. 1986. Т. VII, № 5. С. 878-887.

3. Богатин О. Б., Старостин Н. П., Черский И. Н. и др. Экспериментальная оценка эффективности восстановления моментов трения в системе несмазываемых подшипников по замерам температур // Трение и износ. 1991. Т. 12, № 3. С. 442-445.

4. Старостин Н. П. Математическое моделирование теплового режима и температурная диагностика трения в системе цилиндрических подшипников скольжения // Мат. заметки ЯГУ. 1997. Т. 4, № 2. С. 161-170.

5. Васильева М. А., Кондаков А. С., Старостин Н. П. Исследование применимости упрощенных моделей тепловых процессов в радиальных подшипниках скольжения на основе численных экспериментов // Мат. заметки ЯГУ. 2008. Т. 15, № 2. С. 84-91.

6. Старостин Н. П., Кондаков А. С., Васильева М. А. Тепловая диагностика трения в радиальных подшипниках скольжения с учетом скорости и характера движения вала // Трение и износ. 2012. Т. 33, № 5. С. 454-464.

7. Кондаков А. С., Старостин Н. П., Васильева М. А. Трехмерная граничная обратная задача тепловой диагностики трения в подшипниках скольжения // Мат. заметки ЯГУ. 2012. Т. 19, № 2. С. 187-195.

8. Черский И. Н., Богатин О. Б., Борисов А. З. Анализ температурного поля полимерного подшипника скольжения в нестационарный период трения // Трение и износ. 1981. Т. 2, № 2. С. 231-238.

9. Самарский А. А. Теория разностных схем. М.: Наука, 1983.

10. Butt M. M., Taj M. S. A. Numerical methods for heat equation with variable coefficients // Int. J. Computer Math. 2009. V. 86, N 9. P. 1612-1623.

11. Shih T.-M., Sung Ch.-H., Yang B. A numerical method for solving nonlinear heat transfer equations // Numerical Heat Transfer. Part B: Fundamentals. 2008. V. 54, N 4. P. 338-353.

12. Sankar M., Park J., Kim D., Do Y. Numerical study of natural convection in a vertical porous annulus with an internal heat source: Effect of discrete heating // Numer. Heat Transfer. Part A: Applications. 2013. V. 63, N 9. P. 687-712.

13. Fu W.-Sh., Tong B.-H. Numerical investigation of heat transfer of a heated channel with an oscillating cylinder // Numer. Heat Transfer. Part A: Applications. 2003. V. 43, N 6. P. 639-658.

14. Liu Ch.-Sh. An iterative method to recover the heat conductivity function of a nonlinear heat conduction equation // Numer. Heat Transfer. Part B: Fundamentals. 2014. V. 65, N 1. P. 80-101.

15. Nabongo D., Boni T. K. Numerical quenching for a semilinear parabolic equation // Math. Modeling Anal. 2008. V. 13, N 4. P. 521-538.

16. Dropkin D., Karmi A. Natural-convection heat transfer from a horizontal cylinder rotating in air // Trans. ASME. 1957. V. 79, N 4. P. 741-749.

Статья поступила 4 ноября 2014 г.

Тихонов Роман Семенович, Старостин Николай Павлович Институт проблем нефти и газа СО РАН

ул. Октябрьская, 1, Якутск 677891, Республика Саха (Якутия) roman_tikhon@mail. ru, nikstar56@mail .ru