Расчет течения степенной жидкости в одношнековом экструдере Текст научной статьи по специальности «Физика»

2017 Математика и механика № 49

УДК 532.516.5

Б01 10.17223/19988621/49/8

М.П. Бессонова, М.А. Пономарева, В.А. Якутенок

РАСЧЕТ ТЕЧЕНИЯ СТЕПЕННОЙ ЖИДКОСТИ В ОДНОШНЕКОВОМ ЭКСТРУДЕРЕ1

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

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

При формовании изделий из высоковязких полимерных композиций одним из важнейших элементов технологической оснастки является экструдер, с помощью которого осуществляется транспорт перерабатываемой массы в пресс-форму. Кроме этого, экструдер выполняет роль шнекового смесителя. Экструзионной переработке полимеров посвящена обширная литература. Наиболее полными можно считать работы [1-4]. Помимо этого постоянно появляются исследования более частного характера [5-8]. Внимательный анализ состояния исследований показывает, что в настоящее время нет устоявшейся методики расчета процесса экструзии. Используются либо подходы со значительными, иногда неоправданными упрощениями, либо исследования носят общетеоретический характер и используют трудно воспроизводимые методы. В настоящей работе используется подход, основанный на применении многократно проверенных физических предположений о характере сложного движения экструдируемого материала и в то же время учитывающий максимальное количество факторов, определяющих процесс. Такой подход представляется наиболее перспективным для создания модели процесса экструзии, адекватно описывающей поля давлений, напряжений, вязкости и температур и позволяющей создать эффективные прикладные программы, работоспособные в максимально широком диапазоне определяющих параметров.

Математическая постановка задачи

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

1 Работа выполнена при финансовой поддержке Минобрнауки РФ в рамках государственного задания №9.9625.2017/БЧ.

экструдируемой массы. Будем рассматривать обращенное движение: шнек считается неподвижным, а движется корпус экструдера. Данное предположение часто используется при моделировании течения в экструдере [4].

Систему координат расположим как на рис. 1. Здесь: $ - ширина канала шнека; Н - его глубина; х - координата, направленная вдоль канала, у - поперек канала; 2 -по высоте канала. Все остальные обозначения (для составляющих скорости, давления, напряжений и т.д.) являются общепринятыми.

Предполагается, что боковые стенки не оказывают влияния на течение ($ / Н > 3 [9]) и нет утечек через гребни шнека. По направлениям х и у действуют градиенты ¡р / (х

и ¡р / (у. Тогда: у2 = 0, уу = уу (2), ух = Ух (2),

Рис. 1. Область решения Fig. 1. Solution region

Р=Р (х, у ).

Уравнения движения при этих предположениях приобретают форму

d т x

dp dx '

d т

yz

¡2 дх ¡2 Уравнение неразрывности вырождается, так как

dp ¿У '

(i)

dv dvy dv

= 0, —y— = = 0. dx dy dz

Второй инвариант тензора скоростей деформаций запишется в виде

Л Ч '2

2 2

dv

y

V

dz

у

âvx dz

(2)

Выражение для вязкости запишется следующим образом:

n = в

n-1

12 Л 2

= в

^Л 2 V dz У

n-1

dvx

dz

С учетом вида компонент тензора напряжений

т=п—-, тЛ!7 = П"

х2 ¡2 у2

и уравнений (1) получаем

d_ <

dz

в

d J

dz

в

V dz У

V dz у

V dz

dvy

dz

dv, л2^2

+ V dz

n-1

2 dvx

dz

n-1

fdVy\2 / л.. X2

2 dvy_

dz

dp dx '

dp

dp dp

Введем обозначения — = A1, — = A2. В результате интегрирования первого

dx dy

уравнения системы (1) имеем

I dTxz =| AJdZ,

Т xz (Z) = {Т xZ (Z1) - A1 Z1}+ A1z. В соответствии с [1], постоянную т xz (z1) - A1 z1 выберем равной

т xz (z1) - A1 z1 = -AlC1H, где C1 = const,

т.е. т xz (Z1) = A1( Z1 - CH).

Преобразовав соответствующим образом второе уравнение (1), получим сис-

тему

Тxz (z) = 4(Z - CH),

т yz (z) = A,( z - C2 H).

Введем следующие безразмерные переменные:

e z vx vy A2 dpldy

£ = — (0 < £ < 1), v, =—, v2 = —, a = — = , , Hv 71 VB VB A1 dp dx

тогда

txz = B

dvx

(dv.

. +1 — dz ) ^ dz

n-1

2 ^ = A1 (z - CH), dz

(3)

Bvn [( di Y +( dv2 Y 2 dv1

Hn- 1H _f d £ ) +1 d£ J _ ~d£

n-1

Bvn

AH

dv1 ^ ( dv2 ) + 1

2 di = £-C1.

d f 1

Здесь VB - скорость корпуса в обращенном движении.

Следуя [1], введем обозначение а =

Vb ( B

H I HA

а = -

Bvn

Hn+l a

тогда

dv1 V ( dv2 ) + f "df

2 di=£-C„

d £ 1

Для второго уравнения (3) соответственно будем иметь

n A1

A2

dvj ^ ( dv2

df J +1df

2 ^ = £-C2 d £

n-1

n-1

n-1

п—1

а а

' ёу, V (ёу2 л2 11 + 2

ёI) I ёI

2

ёI

= 1— С2.

(5)

ёу, ёу2

Из (4) и (5) можно получить соотношения для —1 и —2 . Для этого поделим (4)

ё I ё I

на (5):

ё I ) I — С, ёу, 1 I — С, (у2

(у2 \ 1 — С2 ё| а I — С2

ё I.)

Подставим полученное выражение в (5)

а а

А. (I—с, )2 I а2 (I —С2 )2

■ +1

2 (ё^2 ё I

= 1— С2

ёу2 ё I

(I — С2

(I —С, )2 +а2 (I —С2)2 а2 (I —С2 )2

—[(I — С, )2 + а2 (I — С2 )2 ] [(I —С2 )1—п ап а1'7 [ и У 2' ]

= (; — С2)а-1 [(I — С, )2 + а2 (I — С2)2]

1—п

1—п

2п

Пусть у = [(; — С, )2 + а2 (I — С2 )2 ]"

1—п

тогда

Аналогично, имеем

ё2=- С2).

ё I а

ёу2 I — С2 ёу, ■ = а-

ё I I — С, ё I Подставляя данное выражение в (4), получим

а

а 2(! — С2)2 [ (I—С, )2

+1

п—1

Ч * Г =I-С-

(6)

ёу,

(I —С,)

а

а2 (I — С2 )2 +(! — С, )2

(I—С, )2

п—1

Таким образом,

dt - C,).

d I а

(7)

Граничные условия задачи заключаются в следующем: на нижней границе (£ = 0) vx = 0, vy = 0 ; на верхней границе (£ = 1) v1 = cos ф, v2 = sin ф. Используя

(6), (7) и граничные условия на нижней границе, получим

I

VI) = - Jy(|- q)d I,

,, ik, i.

ГУ J

(S)

п УБ ( В \ п Л2 Е

Здесь учтено, что а = -I и а =от Е , т.е. от г не зависят, так как

НI ЯЛ! J Л! ь

Р = Р( х, у).

Обозначим расход на единицу ширины канала 0, следовательно, в безраз-

0 1

мерном виде ^ = —— = I У1йЕ . Тогда, используя (8),

УБЯ 1

q =

11

rv "

s

JV(|- Q)d I

d |.

Последнюю формулу можно упростить с помощью следующего соотношения:

1 Г х 1 1 Г1 1 1

I I /(у)йу йх = | I /(у)йх йу = | /(у)(1 - у)йу,

тогда

q =-1J y(i-C )(1 -I)d I.

а

Аналогично для безразмерного расхода в поперечном (оси у) направлении д2 будем иметь

q2 = - JV(I- -2)(1 -|)d I.

rv *

Выпишем полученную систему уравнений

1J y(|- Cj)(1 -|)d | = q„ а 0

1

-JV(i- C2)(1 -|)d|= 0,

rv

(9)

-||у(| —С,)ё | = о,

а о

1

-Гу^—С2)ё | = о,

г/ ^

(11)

(12)

у = [(| — С,)2 + а2(| — С2)2] 2п .

Здесь

1 Н ( Н

-А,

а _А2 Н (Н

А1 Г =

Н

1

Н АЧ л

--М

В Ап V А1 )

а ¥В V В ) а А, ¥В V В 1) ¥В Зная неизвестные а, а, С,, С2, можно найти А,, А2.

п = Н_ = ^

V ВА,п—1)

Метод решения

Для решения системы нелинейных уравнений (9) - (12) используется метод Ньютона. Перепишем систему, введя новые, более удобные для последующего решения, обозначения неизвестных. Пусть у, = 1/ а , у2 = а / а , тогда а = у2 / у,,

у3 = С,, у4 = С2 и система примет вид:

у, }у(| —уз)(1—|)ё I — д, = о,

0

1

у 2 {уй— у4)(1 — Оё | = о,

0

1

у, {у(| — уз)ё ^ °

0

1

у2 {у(| — у4)ё^ 0,

у =

(I— уз)2 + 4(| — у4)2 у,2

1—п 2п

В общем виде систему уравнений можно представить в форме

/ (у}) = о, (/ (у,,у2,уз,у4) = о), г,) = 1 -4.

Приращения (дифференциалы) функций / записываются в виде

д/г

ду} ;

например:

д/ д/ д/ д/ / ёу1 + ^ (у2 ¡уз + ^ (у4.

су, ду2 суз ду4

В конечно-разностном виде имеем

f (k+i) _ f ) _

(k) _

f I ((

,,(k )ï

Так как должно быть /<-к+1) = 0, у(к+1)

f f Ik) df

(k) _ Al,(k+1)

_ y (k ) =ду

то

дУ,

Ау^+1) =-йк) (),у2*),у?),)).

В результате решения полученной системы будет найдено значение Ау^

(к+1)

следовательно, и y

(к+1) _ У( к )

= y,) _Ду

,( к+1)

и т. д., пока

\fi~к) | >е (8- малое положи-

тельное число). Вычисление производных осуществляется по формулам

/ (*+1) _ / (*)

С f I(к) df

y(k+1)

-y(k) '

V 71 /

Таким образом, алгоритм метода решения состоит из следующих действий: I. Задаются начальные значения: к = 0, у* (т.е. у\ , у2 , уЗ, у4 ).

II. Вычисляется матрица Якоби

% =—, г,1 = 1 - 4

следующим образом:

1. Вычисляется /« = /« (у(к),у2к),у3к),у4к)), г = 1,4.

1а. Если условие тах|/к^ <Е выполняется, то решение найдено, иначе происходит переход к следующему пункту.

2. Вычисляется

f (к ) _ ij

/Зк) (у(к), у2к), уЗ к) + к, у4' /г{4 ) ((), у2к), уЗ к), у4к) +

где к - малое положительное число (к = 0.001). 3. Вычисляется матрица Якоби

f к) (У1(к) + h,(2к),(3к),(4к)),

) (к),(2к) + h,(3к),yf)),

.(к) )

h),

i, j=1,4,

а(к ) = i,

\ду] j

\( к )

f (к ) _ f к iji

i, j = 1,4.

III. Решается система

и находится Ду(к ).

) Д((к) =

af ду(

= _ f

( к )

а

IV. Вычисляется у(к+1) = у(к) + Ду(к).

V. Переходим к пункту I.

В результате решения системы уравнений получаются значения следующих величин: у, = 1/а , у2 = а / а , у3 = С,, у4 = С2.

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

Вычисление температуры расплава вдоль канала экструдера

Предположения относительно гидродинамики прежние (сложный сдвиг):

У2 = ^ Уу = Уу (2X УХ = УХ (2X Р = Р(^ у).

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

H

Th

, VB cos ф

0 *

Рис. 2. Область решения для расчета температуры Fig. 2. Solution region for the temperature calculation

В предположении пренебрежимо малой теплопроводности в направлении х (д2Т / дх2 = о) и того, что Т = Т(х, 2) (т.е. по у происходит выравнивание температуры за счет циркуляции), уравнение переноса тепла с учетом вязкой диссипации можно записать в виде:

дт „ д2Т _ pevx — = Х — + ф, дх dz

(13)

^ J2 I2 ф = л—, — 2 2

dVx

dz

2 +,' Ь

dz

Эффективная вязкость тогда

I.

ф = Л^ = BI

2 I 2

n—1

л=в I i ^ "

2 1 Ii) = в < ^

n+1 2

= в

dVx dz

( dvy

y

+I -

^ dz

n+1 2

(14)

Функция Ф представляет собой количество тепла, выделяемого в потоке за счет вязкой диссипации. Условие на входе в канал шнека экструдера для уравнения (13):

х = о, Т (о, 2) = ТН, 1 < 2 < Н.

z

w

w

Граничные условия на стенках канала имеют вид

Т(х,о) = Т№, Т (х, Н) = Т№. Выберем в качестве безразмерных переменных величины

Ух х г 2 уу

У, =—,п= —,|=—,у2 = —. 1 Ув Н 2 ГВ

Уравнение (13) в безразмерных переменных с учетом (14) примет вид

п+1

рсУву дТ Хд2Т + ВУвп

-У,-= —---— + В-

Ь 1 дп Н2 д|2 Нп

ёу,

( ёу2

+ 1 -

V ё |

В результате умножения полученного уравнения на Н2 /Х будем иметь

п+1

рсУвН2 у дТ = + В УВ хь 1 дп д|2 хнп—1

2

п +1

ёу1

ёу2 ё I

(15)

Введем следующие обозначения. Пусть вг = рсУВН /ХЬ - число Гретца, Вг = ВУ^+1 /ХНп—1ДТ - число Бринкмана (ДТ = Тш — ТН ), тогда ВгТ = ВУВп+1 / ХНп—1 (°С), т.е. ВгТ = ВгДТ. С учетом введенных обозначений, уравнение (15) примет вид

п +1

^ дТ д2Т „ вгу, — =--+ ВгТ

1 дп д|2 Т

ёу1

( ёу2

+ 1-

V ё I

-|( п+1)/2

Г 2 2 ~|(

и пусть для сокращения записи Ф = ВгТ I (ёу, /ё|) +(ёу2/ё|) I . Тогда приходим к уравнению

дТ д2Т

(16)

^ дТ д2Т _ вгу, — = —- + Ф.

1 дп д|2

Для решения уравнения теплопроводности (16) применяется маршевый метод, использующий алгоритм прогонки.

Результаты

Для расчета течения степенной жидкости в одношнековом экструдере были выбраны следующие геометрические, реологические и теплофизические параметры: длина экструдера Ь = 25о мм, ширина канала шнека $ = 3о мм, глубина канала шнека Н = 1о мм, угол наклона нарезки ф = 17°, число заходов г = 1, коэффициент консистенции В = 1-Ю4Па• с , коэффициент теплопроводности Х = о.1 Вт/(м • К), коэффициент теплоемкости с = 8оо Дж/(кг • К), плотность материала р = 1ооо кг/м3. Начальная температура поступающего в экструдер материала ТН = 2о °С . Температура стенок Т№ = 1оо °С. Скорость поверхности корпуса экструдера в обращенном движении принималась равной УВ = о.о7 м/с. Производительность экструдера 3о кг/ч.

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

Рис. 3. Профили продольной v1 (а), поперечной v2 (b) и осевой v (с) составляющих скорости для различных значений показателя нелинейности n: 1 - 1.0, 2 - 0.6, 3 - 0.4

Fig. 3. Profiles of the (а) longitudinal v1, (b) transverse v2, and (с) axial vt velocity components at different power-law indices: n = (1) 1, (2) 0.6, and (3) 0.4

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

p, ат 60

Т °с

2 3

0.8 х, м

Рис. 4. Изменение давления (a) и средней температуры (b) по длине канала шнека экструдера для различных значений показателя нелинейности, n: 1 - 1.0, 2 - 0.6, 3 - 0.4 Fig. 4. Variation in (a) pressure and (b) average temperature along the screw extruder channel at different power-law indices: n = (1) 1, (2) 0.6, and (3) 0.4

Заключение

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

ЛИТЕРАТУРА

1. Янков В. И., Боярченко В. И., Первадчук В.П. Переработка волокнообразующих полимеров. М.; Ижевск: НИЦ «Регулярная и хаотическая динамика». 2005. Т. 2. 998 с.

2. Торнер Р.В., Акутин М.С. Оборудование заводов по переработке пластмасс. М.: Химия, 1986. 400 с.

3. Тадмор З., Гогос К. Теоретические основы переработки полимеров. М.: Химия, 1984. 632 с.

4. Раувендааль К. Экструзия полимеров / пер. с англ. под ред. А.Я. Малкина. СПб.: Профессия, 2008. 768 с.

5. Субботин Е.В., Щербинин А.Г., Труфанова Н.М. Численное исследование процессов течения полимеров в условиях фазового перехода в винтовых каналах экструдеров при производстве пластмассовой изоляции // Известия ТПУ. 2012. Т. 320. № 4. С. 171-177.

6. Sobhani H., Ghoreishy M. H. R., Razavi-Nouri M., Anderson P. D., Meijer H. H. E. Modelling of polymer fluid flow and residence time distribution in twin screw extruder using fictitious domain method // Plastics, Rubber and Composites. 2011. V. 40. Iss. 8. P. 387-396.

7. Althaus T. O., Windhab E. J. Extrusion of highly unsaturated wet powders: Stress fields in extruder barrels // Powder Technology. 2011. V. 211. Iss. 1. P. 10-18.

8. Baron R., Vauchel P., Kaas R., Arhaliass A., Legrand J. Dynamical modelling of a reactive extrusion process: Focus on residence time distribution in a fully intermeshing co-rotating twin-screw extruder and application to an alginate extraction process // Chemical Engineering Science. 2010. V. 65. Iss. 10. P. 3313-3321.

9. Пономарева М.А., Филина М.П., Якутенок В.А. Циркуляционное течение высоковязкой неньютоновской жидкости в канале одношнекового экструдера // Вестник Томского государственного университета. Математика и механика. 2016. № 2 (40). С. 97-107.

Статья поступила 01.09.2017 г.

Bessonova M.P., Ponomareva M.A., Yakutenok V.A. (2017) CALCULATION OF A POWER-LAW FLUID FLOW IN A SINGLE-SCREW EXTRUDER. Tomsk State University Journal of Mathematics and Mechanics. 49. pp. 81-93

DOI 10.17223/19988621/49/8

This paper is devoted to the creation of a physical and mathematical model of the fluid with a complex rheological behavior flowing in a single-screw extruder channel taking into account the presence of a circulation flow. The power-law model is used to describe the rheological behavior of the extruding polymer. The reversed motion is considered: the screw is assumed to be fixed and the extruder cylinder is movable. Providing a predetermined flow rate, four nonlinear equations describing the flow are obtained. The system of equations is solved by Newton's method. The calculation procedure for the polymer melt temperature along the extruder channel is presented. The heat conduction equation is written in a two-dimensional form taking into account the heat generation due to the viscous dissipation of energy and it is solved by marching method. The power-law fluid flow characteristics are calculated at different nonlinearity degree.

Keywords: single-screw extruder, non-Newtonian fluid, non-isothermal flow, mathematical model, circulation flow.

BESSONOVA Mariya Petrovna (Tomsk State University, Tomsk, Russian Federation) E-mail: bessonova.mp@mail.ru

PONOMAREVA Mariya Andreevna (Candidate of Physics and Mathematics, Assoc. Prof., Tomsk State University, Tomsk, Russian Federation) E-mail: pma@ftf.tsu.ru

YAKUTENOK Vladimir Al'bertovich (Doctor of Physics and Mathematics, Prof., Tomsk State University, Tomsk, Russian Federation) E-mail: yva@ftf.tsu.ru

REFERENCES

1. Yankov V.I., Boyarchenko V.I., Pervadchuk V.P. (2005) Pererabotka voloknoobrazuyushchikh polimerov [Fiber-forming polymer processing]. Moscow, Izhevsk: R&C Dynamics.

2. Torner R.V., Akutin M.S. (1986) Oborudovanie zavodov po pererabotke plastmass [Equipment for plastic processing plants]. Moscow: Khimiya.

3. Tadmor Z., Gogos C.G. (1979) Principles of Polymer Processing. New York: John Wiley & Sons.

4. Rauwendaal C. (1986) Polymer Extrusion. Munich: Hanser Publishers.

5. Subbotin E.V., Shcherbinin A.G., Trufanova N.M. (2012) Chislennoe issledovanie protsessov techeniya polimerov v usloviyakh fazovogo perekhoda v vintovykh kanalakh ekstruderov pri

proizvodstve plastmassovoy izolyatsii [Numerical study of a polymer flow in the phase transition conditions in extruder screw channels in the plastic insulation production]. Bulletin of the Tomsk Polytechnic University. 320(4). pp. 171-177.

6. Sobhani H., Ghoreishy M.H.R., Razavi-Nouri M., Anderson P.D., Meijer H.H.E. (2011) Modelling of polymer fluid flow and residence time distribution in twin screw extruder using fictitious domain method. Plastics, Rubber and Composites. 40(8). pp. 387-396. DOI: 10.1179/1743289810Y.0000000032.

7. Althaus T.O., Windhab E.J. (2011) Extrusion of highly unsaturated wet powders: Stress fields in extruder barrels. Powder Technology. 211(1). pp. 10-18. DOI: 10.1016/j.powtec.2011. 03.002.

8. Baron R., Vauchel P., Kaas R., Arhaliass A., Legrand J. (2010) Dynamical modelling of a reactive extrusion process: Focus on residence time distribution in a fully intermeshing co-rotating twin-screw extruder and application to an alginate extraction process. Chemical Engineering Science. 65(10). pp. 3313-3321. DOI: 10.1016/j.ces.2010.02.019.

9. Ponomareva M.A., Filina M.P., Yakutenok V.A. (2016) Circulatory high-viscosity non-Newtonian fluid flow in a single-screw extruder channel. Tomsk State University Journal of Mathematics and Mechanics. 2(40). pp. 97-107. DOI: 10.17223/19988621/40/10.