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

МОДЕЛИРОВАНИЕ ПРОЦЕССОВ

УДК 621.8

Т. А. Фомина

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

Ключевые слова: волновой шаговый двигатель, пневмогидрораспреде-литель, импульсы давлений, волна деформирования, поршневая группа, квазиодномерное движение, схема Лакса-Вендроффа.

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

Приведено описание математической модели расчета движения газа или жидкости (рабочего вещества) вдоль канала в пневмоги-дрораспределителе волнового шагового двигателя.

E-mail: fota@bk.ru

в зависимости от направления движения бегущей волны деформирования гибкого колеса. Массив с отверстиями, на котором закреплен пружинный пакет (см. рис. 1), находится в покое.

Рис. 1. Пружинный пакет (разрез)

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

Левая часть Правая часть

Рис. 2. Схема изменения взаимного положения каналов во времени

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

В настоящей работе приведена математическая модель расчета режимов течения рабочего вещества в системе, состоящей из двух каналов (I и II), непрерывно сдвигающихся относительно друг друга (рис. 2). При этом канал I может быть каналом высокого давления и тогда рабочее вещество движется из канала I в канал II (нагнетание). Если же канал I будет каналом низкого давления, то движение отработавшего рабочего вещества будет происходить из канала II в канал I (слив). При этом площадь поперечного сечения канала низкого давления должна быть подобрана так, чтобы обеспечить полный слив рабочего вещества из канала II.

Вычисление динамических характеристик рабочего вещества выполняется с использованием квазиодномерных уравнений газовой динамики в модифицированной форме. Они предназначены для расчета квазиодномерных движений, когда, строго говоря, движение не является одномерным вдоль оси канала, однако кинетической энергией и импульсом движения перпендикулярно этой оси можно в первом приближении пренебречь. Уравнения учитывают изменение площади сечения по длине канала [2, 3] и имеют вид

ди + д (Ё(и)Б (х)) =0; дЬ Б (х) дх

д (pv) 1 д ((p + pv2) S (x))

dt

+

S (x)

dx

p

1 dS (x) S (x) dx

= 0,

где и = ^ р ^; ^(и) = ^ (е +-р) ; Б(х) - площадь поперечного сечения канала, х — продольная координата. Исследуемыми параметрами являются: р, у,р — плотность, скорость движения, давление рабочего вещества; е — полная энергия рабочего вещества на единицу объема.

Давление определяется из уравнения состояния р = / (Е,р), в котором Е — внутренняя энергия газа на единицу массы. В качестве уравнения состояния используется уравнение Менделеева-Клапейрона р = (7 — 1) (е — 1 ру2), так как рассматривается адиабатический процесс. ж

С помощью замены г = J Б (х) д,х приведенные уравнения запи-

0

сываются в дивергентной форме:

f + дЕШМ = 0, (1)

дЬ дг

д (ру) + д ((р + ру2) Б (г)) — ^БМ = 0 дЬ дг ¿г

Уравнения (1) дополняются граничными и начальными условиями, а также условиями сопряжения каналов.

На левой границе канала I задается условие наличия источника

постоянного давления: р(г;,Ь) = р0, р(г;,Ь) = р0, — (г;,Ь) = 0, где р0, р0 — значения давления и плотности в канале при Ь = 0. На правой границе канала II граничное условие задает условие жесткой стенки:

др ,(> = °-11

Начальные условия определяют профиль распределения характеристик рабочего вещества при Ь = 0: р(г, 0) = р0, V (г, 0) = 0. Для канала I (г; ^ г ^ 0) р0 = 10 атм, для канала II (0 ^ г ^ гр) р0 = 1 атм.

В точке г = 0 в каждый момент времени Ь устанавливаются условия сопряжения каналов I и II. Они зависят от площади открытого отверстия Б0. Если каналы не взаимодействуют, то на правой границе канала I и левой границе канала II принимается условие жесткой

(др (0,Ь) = 0, др (0,Ь) = 0, у(0,ь) = 0) .В противном случае, дг дг

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

Исследование системы дифференциальных уравнений (1) для всей конструкции целиком затруднено, что связано с наличием скачка пло-

(zp,t) = 0, /(zp,t) = 0, v(zp,t) = 0.

стенки

щади сечения между частями, зависящего не только от координаты г, но и от времени (см. рис. 2). Поэтому исходная конструкция разбивается на две независимые части — правую и левую, имеющие один общий узел (см. рис. 2). Левая часть описывает характеристики рабочего вещества в канале I, правая — в канале II.

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

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

Численная реализация алгоритма проводится в нормированных пе-

ременных:

z =

s (z ) =

S (z)

С ' 4 > с '

^шах ^шах

где Сшах — максимальная площадь сечения правого канала.

Под относительной площадью отверстия понимается нормирован-

ная площадь сечения между левой и правой частями С о (г) =

Со (г)

Сшах

(С0 (г) — площадь отверстия, единственного во всей конструкции сечения). На рис. 3 приведен закон изменения относительной площади отверстия по времени С0(г). Верхней кривой соответствует канал II, выполненный в форме цилиндра. В этом случае максимальное значение площади отверстия С0 (г) равно Сшах и максимальное значение относительной площади отверстия равно единице. С увеличением радиуса цилиндрической части канала II максимальное значение относительной площади отверстия уменьшается (нижние графики, см. рис. 3).

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

z

При численной реализации (1) левая и правая части равномерно разбиваются по оси г на отрезки Дг. Для каналов I и II они могут быть разными, но не сильно отличающимися друг от друга. Число точек по оси г для левой части обозначим N, для правой — N (см. рис. 2).

На этапе "предиктор" вычисляются параметры рабочего вещества между узловыми точками в промежуточный момент времени по соотношениям (2) (£ € [0,1]):

г+1/2 = 2 (ип+1 + и,) т дг (^с(гт) - ^ с(^)) , (2) Мт/2 = 2 (Мт + (НГ\ Т

¿At

Т((p + pv2)n+i s (zi+i) - (p + pv2)n S (*)) ±

±2Az (-РП+1 + РП) (S (z»+i) - S (*)) ,

(pv)n+ä / , 2 n+ä _ A+1/2 -n+5 _ / i\ „n+5 1 -n+5 l„.n+s *

: pn+. = (7_ i) _i pH-« („„+ у

, Pi+1/2 = (Y X) lei+1/2 2Pi+1/2 ^i+1/2j

^+1/2 = рп+й , Рг+1/2 = (7 - 1) ег+1/2 - 2 Р - - I V. Р г+1/2 V 2

На этапе "корректор" уточняются значения в узловых точках по формулам (3) в момент времени г + Дг с использованием дополнительной процедуры усреднения в смежном узле и условий на границах:

иП+1 = иП Т Д (^Т++1/2С (гг+1/2\ - (гг-1/2\) , (3)

МГ+1 = МП Т Д| ((Р + Р™2)^(гг+1/2)

(р + pv2 )П+/2С (гг-1/2^ ±

n+ä е/— 1/2S

±A (pn+1/2 + рП—1/2) (S (zi+1) - S (zi)),

"Az V i+1/2 i—1/2 vn+1 = , РП+1 = (y - 1) (V - 1РП+1 (vn+1)2

В соотношениях (2) и (3) верхний знак используется при вычислении правой части конструкции, нижний — для левой. Различие появляется из-за выбора направления нумерации узлов сеток (см. рис. 2).

В крайнем правом узле (г = канала II выполняются условия жесткой стенки: 1 = 0, р£1 = р^; 1 = е^; 1 = Р^-г

В крайнем левом узле (г = N1) канала I выполняются условия наличия источника постоянного давления: ^^ = , р^1 = Р0;

п + 1 1 ( п + 1 N 2 п + 1

ел+г = ео + 2Ро ) ; Рд+ = Ро, где ро, ро, ео — значения давления,

плотности и удельной полной энергии в канале I при £ = 0.

Характеристики рабочего вещества в общем для каналов I и II

узле (г = 1) зависят от площади входного отверстия (£о). Если

£о (£п+1) = 0, т.е. каналы I и II не взаимодействуют между собой,

то в узле г = 1 для обоих каналов задается условие жесткой стенки:

2Д+ _

<+1 = 0, РП+1 = (7 - 1) еГ1, РП+1 = РП Т 2д=т (ру)Г++1/2 £ (*+1/2);

2Д£ _

еП+1 = еП Т ((е + р) у)™+1/2 £ (^г+1/2). В указанных выражениях верхний знак используется при вычислении правой части конструкции, нижний — при вычислении левой.

Если £о (¿п+1) = 0, тогда применяется интерполяция между соответствующими параметрами в первых ячейках слева и справа. Она проводится на физической оси х. Используя переменные для канала I х = х (г2), для канала II хр = х (г2) и обозначая индексом I и II значения соответствующих характеристик в каналах I и II, получаем:

/ N1+1 (ру)2,П + (ру)2Д п+1 р2,н + р2Д

М1 =--— > р1+ = +——,

х1 + хр х1 + хр

еп+1 + еп+1 / \1+1

п+1 _ е2,н + е2Д „1+1 _ (р^)1

е1 — 1 ил —

xl + Pn

1

Р1+1 = (7 - 1)(е?+1 - 1 р?+1 «+1)2) . Для расчета первой ячейки вводится усредненное давление:

1+1 Р1+1£ о + р2 (£под — £ о)

Рп+ =-=—---— для канала I;

S

под

n+1 P?+1S 0 + p2 + 1(Smin — S 0) тт

Р2+1 =-=—--1 — для канала II,

__S min

где Smin — минимальная приведенная к Smax площадь сечения канала II;, Sncw — приведенная к Smax площадь сечения канала I.

Для устранения осцилляций решений вводится линейная искусственная вязкость и проводится усреднение давления по трем точкам. Таким образом, решение получается достаточно гладким.

Кроме того, известно, что при увеличении коэффициента ö при At в соотношении (2) на этапе "предиктор" осцилляции также уменьшаются, но точность решения при этом снижается. Особенность предлагаемого алгоритма, отличающая его от существующих модификаций схемы Лакса-Вендроффа, состоит в том, что в разработанном программном комплексе пользователю предоставляется возможность интерактивного регулирования коэффициента ö для достижения необходимой гладкости решения.

Расчет значения шага по времени At для перехода на следующий этап проводится из условий устойчивости Куранта (4) совместно для левой и правой частей:

At = min (Atp, Ati), (4)

где

Atl = min

ie[1;N;-1]

|AXi|

2 + 1 Ji+1/2

I rtn+1 + Ci+1/2,

Atp = min

ie[1;Np-1]

Axi

2+1 Ji+1/2

i /rnn+1 + Ci+1/2,

C n + 1 Ci+1/2

\

7РП++1/2 Pn+1 Pi+1/2

скорость звука в ячейке.

С помощью построенной модели исследованы зависимости давления от скорости нагнетания и скорости слива рабочего вещества и геометрических характеристик правой части: большого радиуса конусной части, длин конусной и цилиндрической частей. Анализируемым показателем является изменение давления на правой стенке канала II (в полости поршневой группы). В начальный момент времени в левой области задается давление, равное ро = 10 атм, в правой — ро = 1 атм. Плотность р и энергия е в начальный момент вычисляются

соответственно по формулам: ро =

Ро

ео =

Ро

где y — по-

ЛСот' о 7 — Г казатель адиабаты (7 = 1,4); Л — универсальная газовая постоянная (Л = 8,3 Дж/(моль-К)); Со — абсолютная температура рабочего вещества (Со = 293 К); т — молярная масса (для воздуха т = 28,8 г/моль). Для канала I ро = 14,45 кг/м3, ео = 2 533 125 Дж/м3; для канала II ро = 1,445 кг/м3, ео = 253312,5 Дж/м3.

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

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

О 1 2 3 4 5 6 7 8 9 t, мс

Рис.4. Зависимость давления в поршневой группе от периода вращения внутренней части распределителя

Ртах? а™

2 4 6 8 10 12 14 16 18 *,мс

Рис. 5. Зависимость максимального давления в полости поршневой группы от периода вращения внутренней части распределителя

При более подробном анализе зависимости максимального давления в полости поршневой группы (ртах) от периода вращения внутренней части распределителя (£) (рис. 5) можно сделать вывод, что потери давления не превышают 10% при £ ^ 6 мс, что соответствует частоте вращения внутренней части распределителя, равной 10000 об/мин.

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

На рис. 6 показана гистограмма формируемых импульсов давления для случая, когда сечение отверстия для слива является круговым. Пе-

О 1 2 3 и мс

Рис. 6. Гистограмма расчетных импульсов давления в случае неполного слива давления

Отверстия

Отверстия для нагнетания

Рис. 7. Изменение внешнего вида развертки внутренней части распределителя в зависимости от формы отверстия для слива канала I

риод вращения внутренней части распределителя равен 4 мс, давление в полости низкого давления распределителя — 1 атм, рабочее давление источника — 10атм. Из рис.6 следует, что при заданных параметрах давление в полости поршневой группы не достигает максимального значения, что приводит к снижению КПД примерно на 20 %. Кроме того, минимальное давление рабочего вещества при сливе приблизительно равно 2 атм (в 2 раза выше необходимого уровня), что также приводит к снижению КПД. Максимальное давление можно увеличить, снижая скорость вращения внутренней части распределителя (см. рис. 5). Также можно снизить нижнюю границу давления.

Второй способ уменьшения минимального давления заключается в изменении сечения канала II для слива с круглого на прямоугольное со скругленными сторонами (рис. 7), при этом, если ширина прямоугольной части отверстия больше 70 % радиуса отверстия для нагнетания, то давление спускается практически полностью (рис. 8).

Pmin? а™

Рис.8. Зависимость минимального давления в полости от ширины прямоугольной части отверстия для слива

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

СПИСОК ЛИТЕРАТУРЫ

1. Клеников С. С., Майков А. И. Разработка математической модели расчета пружинного пакета волнового редуктора методом конечных элементов (МКЭ). - М.: Изв. вузов. Машиностроение. - 2008. - № 11. - С. 25-30.

2. Рихтмайер Р.,Мортон К. Разностные методы решения краевых задач. -М.: Мир, 1972.-421 с.

3. Р о у ч П. Вычислительная гидродинамика. - М.: Мир, 1980. - 618 с.

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