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

УДК 66.063.8

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

© 2015 Н А. Газизуллин

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

Поступила в редакцию 15.02.2015

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

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

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

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

Газизуллин Назым Абдуллович, кандидат технических наук, доцент кафедры высшей математики. E-mail: gnazym @gmail. com

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

Рис. 1. Турбинная мешалка

В данной работе проведено численное моделирование неизотермического ламинарного течения несжимаемой неньютоновской жидкости, описываемой реологической моделью вяз-коупругой Олдройд-Б жидкости, в аппарате с турбинной мешалкой Раштона. Исходными уравнениями, описывающими неизотермическое течение жидкости в аппарате, будут уравнения движения, энергии и неразрывности в виде [4]

д -— (ру)+ V • руу = F + V • (- pI + а), д (1)

дT

— + - •VT = aV ^,

дt

V-v = 0,

(2)

(3)

где р - плотность жидкости, р - давление, а -тензор-девиатор напряжений, I - единичный тензор.

Реологическая модель вязкоупругой Олд-ройд-Б жидкости записывается в виде соотношения [5]

С учетом (5) уравнение движения (1) примет вид

V Г V

а + ^а = 2л| Б + Х2Б I,

(4)

где п - сдвиговая вязкость, А4 - время релаксации, А,2 - время ретардации, а тензор скоростей деформаций определяется как

Б = I

2

+ )Т ].

а = т + 2л2 Б,

(5)

где вязкоупругая составляющая тензора напряжений т удовлетворяет уравнению

т + ^т = 2^ Б.

(6)

В формулах (5) и (6) постоянные п и п2, представляющие собой соответственно неньютоновский и ньютоновский вклады в общую вязкость, определяются соотношениями

Л = Л1 +Л2> X 2 =(1 -Э)ЛХ 1, Р = Л

Л

д

—(ру )+V • (руу - л^ ) =

дХ

= Р-Vp + V • (т - 2рлБ) .

(7)

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

дИ дИ

— = -и--ъ w,

дХ дг

(8)

Соотношение (4) содержит также верхнюю конвективную производную, определяемую в виде [5]

V да ^ / ^

а =--ъ V • Vа - (Vv) • а - а • Vv,

дХ

где верхний индекс Т обозначает операцию транспонирования.

Введем подвижную, связанную с вращающейся мешалкой, цилиндрическую систему координат г,у,г. Проекции вектора скорости обозначим соответственно через u,v,w. Проекции объемной силы в подвижной системе координат содержат ускорение Кориолиса и центробежное ускорение, а также ускорение силы тяжести и имеют вид

р = р(ю2 г + 2шv), РФ=-2Рюи, Р = -Р£.

Для решения системы уравнений (1) - (4) используем метод расщепления напряжений [6, 7], согласно которому тензор напряжений расщепляется на вязкоупругую составляющую и вязкую составляющую с помощью равенства

где И - высота свободной поверхности над мешалкой, свидетельствует о том, что скорость движения свободной поверхности в направлении нормали к ней совпадает с нормальной составляющей скорости движения жидкости [8]. Динамическое условие является условием равенства сил, действующих на свободную поверхность [9]. Будем предполагать, что силы поверхностного натяжения незначительны, так как радиус кривизны поверхности жидкости в аппарате достаточно велик. Тогда динамическое условие может быть записано в виде

п

(- р1 + а) = -р0 п • I,

(9)

где р0 - атмосферное давление, а вектор нормали к свободной поверхности жидкости

п = ^,0;1

дг

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

и = 0, V = -юг, w = 0,

а на поверхности вала и мешалки соответственно

и = 0, V = 0, w = 0,

где ю - угловая скорость вращения вала и мешалки.

На оси вращения потока под мешалкой примем

дw

и = 0, V = 0, — = 0.

дг

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

V

наружными нагревательными элементами, поддерживающими постоянную температуру 71. Тогда на боковой стенке аппарата будет справедливо граничное условие Т=Т\. На свободной поверхности жидкости и на дне аппарата примем Т=Т0, где Т0 - температура окружающей среды. На поверхности вала и мешалки примем адиабатическое условие ЭТ/Эп=0. Поскольку форма свободной поверхности жидкости неизвестна и должна быть найдена в результате расчетов, то перейдем от физической области течения к расчетной области с известными границами. Для этого физическую область поделим на две подобласти, нижнюю и верхнюю, горизонтальным сечением, проведенным через верхнюю поверхность мешалки. Введем безразмерные координаты и функции

* иг г =

* г *

Г =—, ф = ф ,

Ь ' ь

*

2 =•(

г Ь

г - Н - Ь

^ * V V = —.

и'

*

с =

к

*

Р =

Ьс ци'

при г < Н + Ь при г > Н + Ь ь(р - ро )

0 =

ци

Т - Т

т - т

Т1 7 о

где Н - высота расположения мешалки над дном аппарата; Ь - высота лопасти мешалки, а в качестве характерной длины Ь и характерной скорости потока и выберем соответственно диаметр мешалки ё и окружную скорость (в неподвижной системе координат) конца лопасти. Заметим,

что свободной поверхности жидкости соответст-

* 1

вует при этом значение 7 =1.

После преобразования координат уравнение неразрывности (3) сохранит вид

Ч-У = 0.

При этом проекции вектора скорости в расчетной области будут определяться следующим образом:

и = ум , V = у V , Ж = м> - иг

ду

дг *

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

д

1

д

Л(лф)+

дг Уккк д хк

1

д

"г'ъ

(

к1к2 к3

п 2 3 Ли, Ф

ккк д X.

ккк г дФ

д X

+ ^

к У

(10)

где приняты следующие обозначения:

где у=1 для нижней подобласти и у=к* для верхней подобласти; к =к/ё, а и .V - компоненты вектора скорости в безразмерной физической области, определяемой преобразованием

* г * * г

г =-, ф =Ф, г =-.

а а

к = 1, к = г , к = 1,

* * *

= г , х2 = ф , х3 = г , и = и, и2 = V, и3 = Ж, Г = 1; Л = л:Б3е,

где к - индекс суммирования; Ке=рпа2/п - центробежное число Рейнольдса; п - число оборотов мешалки в единицу времени, а 8Ф - член типа источника, который определяется соответственно искомой функции Ф.

Уравнение (6) в расчетной области в проекциях на оси координат также может быть представлено в виде обобщенного уравнения (10). При этом следует считать Г=0; Л=Же; We=A,1ю - число Вайссенберга. Уравнение (2) в безразмерной форме также приводится к виду (10), где Г=1; Л=лЯеРг; Рг=^/(ра) - число Пран-дтля; а - коэффициент температуропроводности.

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

Проинтегрируем уравнение (10) по контрольному объему и временному интервалу Дг. В результате, с учетом уравнения неразрывности, получим соотношение, связывающее значения искомой функции Ф в узловой точке Р с ее значениями в центрах Е, Ж, N 8, Т, В соседних ячеек, в виде

к

к

к

ap Ф р = aE Ф e + aw Фш + aN Ф N + + as ф + ar ф + aB Фв + Sp AV,

(11)

где - узловое значение источникового члена; А У - объем ячейки.

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

zz - Р + Po - 2°r

8h 8r

- +

+ (^rr - P + Po )| — I = 0

8h

(12)

где ат огг - компоненты тензора напряжений.

Будем полагать, что краевые углы смачивания вала и боковой стенки аппарата прямые. Это позволяет считать, что в любой момент времени в непосредственной близости от вала и стенки будут справедливы равенства

h (t *, г/ )= h (t *, r* + Ar *), h* (t *, R )= h* (t \ R*-Ar *),

(13)

(14)

где г* - радиус вала; Я - радиус аппарата. Перепишем условие (12) в безразмерном виде

8h

8r *

* * / * * \ CTzz - Р +(CTrr - Р )

ah*

8r *

2al

(15)

V - V(*) =^(r*2 - г*2 )5h(k),

(16)

где У0 - объем жидкости над мешалкой с невозмущенной свободной поверхностью; у (к) - объем жидкости над мешалкой на к-ой итерации, который вычислялся на каждой итерации путем

1**(к)

численного интегрирования по значениям п^ . Следует отметить, что в соответствии с формулой (16) поправка 5П (к) не влияет на форму свободной поверхности, а лишь корректирует ее по высоте. Окончательно с учетом поправки 5П (к) скорректированные значения высоты свободной поверхности П* ) могут быть найдены как

h*k) = h*(k) +5h(k).

(17)

Расчет поля течения проводился на основе алгоритма SIMPLE [10], в котором используется дискретизация уравнений по методу контрольных объемов на сетках с расположением узлов в шахматном порядке.В алгоритме SIMPLE должно задаваться начальное предполагаемое поле

давления. Это поле давления в соответствии с

**

рекомендациями [10] полагалось равным p =0 для того, чтобы погрешности округления при расчете градиентов давления не достигали слишком больших значений.

Вначале по уравнениям (10) вычислялись компоненты тензора напряжений и предварительные значения компонент u ,v ,w вектора скорости, не удовлетворяющие уравнению неразрывности. С учетом приближенного решения для скорости находилась поправка к давлению dp из уравнения

Решение уравнения (15) проводилось на каждой итерации методом Рунге-Кутта [14] при

* 7*/* * * * \

начальном условии п0 = п (X - Д , г* + Дг ) с учетом условий (13) и (14) при фиксированном значении производной дП /дг в правой части (15), взятом с предыдущей итерации. При вычисленных на к-ой итерации значениях высоты

свободной поверхности П* *к) может не выполняться условие постоянства объема жидкости в аппарате (при этом достаточно учитывать объем жидкости над мешалкой). Таким образом, возникает необходимость введения некоторой поправки 5П (к) к значениям П* *к) на каждом итерационном шаге. Эта поправка находилась из соотношения

apbpp = aE ЪРе + aw bpw + aN 5Pn +

+ as bps + ar bpT + aB bpB + bp,

(18)

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

Затем с учетом поправок др рассчитывались скорректированные значения компонент скорости на гранях е, w, п, *, X, Ь ячеек, удовлетворяющие уравнению неразрывности, по поправочным формулам, температуры и давления

р* = р** +абр

где а - параметр релаксации. В расчетах принималось а=0,8.

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

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

2

2

Дискретные уравнения (11) и (18) решались методом прогонки [15]. В качестве критерия сходимости рассматривалась сумма модулей невязок по всем контрольным объемам при решении уравнений (11). Расчеты проводились по этому критерию с точностью до 10-6 на равномерной сетке. В расчетах принималось

Н0 = Б; а/Б = 0,3; а /Б = 0,05; Ь / а = 0,2; Н / Н0 = 0,4,

где Б - диаметр аппарата; Н0 - высота невозмущенной поверхности жидкости над дном аппарата; а - диаметр вала.

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

Рис. 2. Радиально-осевая циркуляция в аппарате при Яе = 200, We = 5

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

жидкость выдавливается из аппарата и поднимается вдоль вала, создавая эффект поднятия. На рис. 3, 4 представлены результаты расчетов изотерм В=еотг в меридиональном сечении аппарата.

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

Рис. 3. Картина изотерм в аппарате при Яе = 50, Рг = 0,05

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

Рис. 4. Картина изотерм в аппарате при Яе = 200, Рг = 0,05

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

1. Стренк, Ф. Перемешивание и аппараты с мешалками. - Л.: Химия, 1975. 384 с.

2. Брагинский, Л.Н. Перемешивание в жидких средах / Л.Н Брагинский, В.И. Бегачев, В.М. Барабаш. - Л.: Химия, 1984. 336 с.

3. Манусов Е.Б. Расчет реакторов объемного типа / Е.Б Манусов, ЕА. Буянов. - М.: Машиностроение, 1978. 112 с.

4. Астарита, Дж. Основы гидромеханики неньютоновских жидкостей / Дж. Астарита, Дж. Маруч-чи. - М.: Мир, 1978. 312 с.

5. Larson, R.G. The structure and rheology of complex fluids. - New York: Oxford University Press, 1999. 663 pp.

6. Xue, S.C. Three-dimensional numerical simulations of viscoelastic flows through planar contractions / S.C. Xue, N. Phan-Thien, R.I. Tanner // Journal of Non-Newtonian Fluid Mechanics. 1998. vol. 74. P. 195-245.

7. Phillips, T.N. Viscoelastic flow through a planar contraction using a semi-Lagrangian finite volume method / T.N. Phillips, A.J. Williams // Journal of Non-Newtonian Fluid Mechanics. 1999. vol. 87. P. 215-246.

8. Лаврентьев М.А. Проблемы гидродинамики и их математические модели / МА. Лаврентьев, Б.В. Шабат. - М.: Наука, 1977. 408 с.

9. Ландау Л.Д. Гидродинамика / Л.Д. Ландау, Е.М. Лифшиц. - М.: Физматлит, 2003. 732 с.

10. Патанкар, С. Численные методы решения задач теплообмена и динамики жидкости. - М.: Энерго-атомиздат, 1984. 152 с.

11. Harlow, F.N. Numerical calculation of time-dependent viscous incompressible flow of fluid with free surface / F.N. Harlow, J.E. Welch // Physics of Fluids. 1965. Vol. 8. P. 2182-2189.

12. Газизуллин Н.А. Расчет скорости циркуляции жидкости со свободной поверхностью в аппарате с мешалкой // Известия Самарского научного центра Российской академии наук. 2012. Т.14, № 1(2). С. 356-359.

13. Газизуллин Н.А. Расчет глубины воронки в аппарате с мешалкой // Химическая промышленность сегодня. 2013. №11. С. 51-56.

14. Холл, Д. Современные численные методы решения обыкновенных дифференциальных уравнений / Д. Холл, Д. Уатт. - М.: Мир, 1979. 312 с.

15. Марчук, Г.И. Методы вычислительной математики. - М.: Наука, 1989. 608 с.

HEAT EXCHANGE OF VISCOELASTIC FLOW IN THE DEVICE WITH TURBINE MIXER

© 2015 N.A. Gazizullin Kazan National Research Technological University

By means of iterative procedure on the basis of SIMPLE algorithm the numerical modeling of convective heat exchange of viscoelastic flow in the device with turbine mixer is carried out. Results of calculations are presented in the form of streamline contours of secondary fluid circulation and isotherms.

Keywords: device with mixer, mixing, viscoelastic fluid, heat exchange, circulation, streamlines

Nazym Gazizullin, Candidate of Technical Sciences, Associate Professor at the Higher Mathematics Department. E-mail: gnazym@gmail.com