CC BY
18
УДК 66.047.25:664
Аналитическая оценка длительности предварительного замораживания в технологии вакуумной сублимационной сушки термолабильных материалов
Д-р техн. наук Г. В. СЕМЕНОВ, д-р техн. наук К. П. ВЕНГЕР, ХУССЕЙН МОХАМЕД МАГРАБИЕ СЛАМА
Московский государственный университет прикладной биотехнологии 109316, г. Москва, ул. Талалихина, 33
An analytical model for calculating the duration of the low-temperature freezing air from turboholodilnoy machine technology to vacuum freeze-drying. Used to model the heat balance integral method with a quadratic polynomial of the temperature distribution in the solution of the heat equation. The relative error of calculations in comparison with the experimental data was about 4 % for symmetric and 2 % for the asymmetric heat.
Key words: mathematical model, freezing duration, Stephen's problem, polynome, vacuum freeze drying, low temperature air, turborefrigeration machine.
Ключевые слова: математическая модель, продолжительность замораживания, задача Стефана, полином, вакуум-сублимационная сушка, низкотемпературный воздух, турбохолодильная машина.
Сублимационная сушка получила широкое применение в фармацевтических производствах и отраслях прикладной биотехнологии для достижения сохранности комплекса нативных свойств обезвоженных материалов. Уровень сохранности качества определяется как собственно параметрами процесса сушки, так и условиями и режимами предварительного замораживания. При этом именно роль замораживания в формировании показателей качества является доминирующей. Высокоинтенсивный теплообмен, происходящий на границе замораживаемый материал — теплоотводящая среда, а следовательно, достаточно высокие скорости перестройки температурных полей, формирования и продвижения фронта фазового перехода обеспечивают мелкокристаллическую структуру льда и равномерное распределение в объеме всех компонентов замораживаемого объекта. По этим соображениям принципиально важной является разработка максимально точных методов аналитической оценки взаимосвязи режимов теплообмена, теплофизических свойств объекта и времени замораживания до заданной конечной температуры.
Сегодня наряду с распространенными вариантами замораживания, базирующимися на применении аммиачных или фреоновых холодильных машин, перспективным является использование воздушной турбохолодиль-пой машины (ВТХМ), детандер которой одновременно обеспечивает низкую температуру от минус 60 до ми-
нус 120 °С и скорость потока воздуха от 5 до 25 м/с. К преимуществам данной системы хладоснабжен ия относится использование естественного и, следовательно, экологически безопасного и дешевого хладагента — атмосферного воздуха [1].
На кафедре «Холодильная техника» МГУПБ проводятся исследования, связанные с использованием ВТХМ в сублимационной сушке.
Одним из этапов данных исследований является разработка математической модели расчета продолжительности замораживания пищевых продуктов с использованием низкотемпературной воздушной системы хладо-снабжения при несимметричных условиях теплообмена, характерных для сублимационной установки.
Следует отметить, что для таких условий теплообмена известна математическая модель расчета времени замораживания пищевых продуктов в трехзонном азотном аппарате [2].
В данной модели использованы приближенные методы решения основной задачи теплопроводности и на базе формулы Р. Планка получены уравнения для определения продолжительности процесса замораживания на момент достижения заданной конечной температуры в точке встречи фронтов кристаллизации с учетом их несимметрии. Адекватность данной модели расчета составила 17 %.
В предлагаемой модели решение основывается на интегральном методе теплового баланса Т. Р. Гудмана [3]. При этом используется допущение, что распределение температуры в слоях продукта аппроксимируется квадратичным полиномом, который удовлетворяет граничным условиям и уравнению теплопроводности. Применение таких полиномов для решения задачи теплопроводности рассматривалось в работах А. М. Бражникова, Т. Р. Гудмана и ряда других авторов [3,4].
В расчетах использовалась компьютерная программа Maple Solver 7.
Процесс замораживания при несимметричных условиях разделяется на три характерные стадии, которые рассматривались последовательно для объекта, имеющего форму пластины:
1) стадия охлаждения: от начальной температуры продукта íHa4 до достижения криоскопической температуры tKр на любой поверхности продукта;
2) стадия замораживания: до достижения криоскопической температуры tKp в термическом центре продукта;
3) стадия домораживания: до достижения заданной конечной температуры ¿-,ад в термическом центре продукта.
Стадия охлаждения
С использованием на противоположных поверхностях продукта различных значений коэффициентов теплоотдачи а-| ф а2 и равных температур охлаждающей среды
cpl
tCp2 = tcр составлена расчетная схема стадии охла-
ждения (рис. 1).
х= О
Рис. 1. Расчетная схема стадии охлаждения для пластины при несимметричных условиях теплообмена:
1,2 — возмущенные зоны; 3 — невозмущенная зона;
Ь,, Ь\, Ьг — толщина, соответственно, всей пластины и возмущенных слоев, м;
Ь\(т), ¿2(т) —координаты температурных фронтов
Аналитическое описание стадии охлаждения связано с решением уравнения теплопроводности для левого и правого возмущенных слоев (см. рис. 1), где при несимметричных условиях теплообмена образуется различная толщина Ь\ (г) и Ь2 (т),
д2і (х, г) dt
дх2 дт
(1)
при
О < т ^ ті, 0 < х < Ь\(т) и Ь — Ь2(т) ^ х ^ Ь,
где т — продолжительность, с; х — координата, м;
а — коэффициент температуропроводности продукта, м2/с;
£(ж, т) — изменение температуры по времени и координате.
Для решения приняты следующие начальные и граничные условия:
— начальное условие (т = 0) при Ьі(0) = 1-2(0) = 0
ti (z,0) = t2 (ж, 0) = t„ач!
(2)
— на температурных фронтах при 0 < т < тх,х = Ьі(т) и х = L — Ь2(т)
h (х, т) = t
ді{х,т)
нач И — U,
ОХ
(3)
— граничные условия на поверхностях пластины (х — 0 и х = L)
, dt\ (х,т) . . ,
А---^------|х=о = <*1 [¿1 (х, т) |х—о - ícp], (4)
дх
dt2 (х,т) дх
Ix=L = «2 \р2 {х, т) \X=L - ícp ] , (5)
где Л — коэффициент теплопроводности незамороженного продукта, Вт/(м-К).
Решение поставленной задачи проводилось с учетом безразмерных переменных
х L(r) та aL
е=Г ч=—• Fo = e = T2' Bl = T’
т(Є.0) =
і(х,т) - tK
^кр tCp
и распределения температуры в возмущенных слоях Т1^,в)=с1(в)+с2{в)^ + с3(в)Є
при 0^£^гц(е),
(6)
Г2(£,0) = а, (9) + а2(9)^ + а3(в)Є
Стадия замораживания
(7)
где £, r¡ — координаты температурных фронтов;
Fo, Bi — числа Фурье и Био;
Г(£, 9) — безразмерная температура; сі(0), с2(9), сз(9), (9), а2(0), а3(0) — коэффициен-
ты, зависящие от ВЦ, Bi2, tHa4, т/і(9) и r/2(0).
Значения коэффициентов сі, с2, с3, аь а2, а3 определялись матричным способом с помощью компьютерной программы Maple Solver 7.
Путем интегрирования уравнения теплопроводности по величине £ в пределах от нуля до 771 и с учетом преобразований получено окончательное выражение, позволяющее численно определить значение r¡i(6), а также продолжительность стадии охлаждения 9,
А\ (ln (BijTji + 2) — ln 2) +
+А2 ^(B¡x?7i + 2)2 — 22^ = 9Аз, (8)
где А\, А2 и Аз — коэффициенты, зависящие от Вії и tH¡í4,
л а^нач (4 2£нач 1
4ві7‘із~віГ
А 2 —
1 (2г„ач - і Вії
Аз — —2В11 (¿нач + 1) ■
Идентичное уравнение получено и для второго возмущенного слоя 1 — т)2\
УІ (ln (Ві2?72 + 2) - ln 2) + +У2 ((Ві2772 + 2)2 - 22) = 9Y3.
(9)
Стадия замораживания начинается с момента образования фронта кристаллизации на левой поверхности пластины за время т\. При этом в продукте образуется несколько зон: замерзшая зона, в которой температура ниже криоскопической, незамерзшие зоны и невозмущенная зона, в которой температура равна начальной. На этапе замораживания существует одновременное перемещение двух температурных фронтов Ь\{т) и Ь2(т) и двух фронтов кристаллизации ¿:(г) и <52(т).
Предполагается, что стадия замораживания разделится на два этапа:
— первый этап связан с достижением правой поверхностью пластины криоскопической температуры продукта Т2(1,92) = Ткр;
— второй этап характеризуется исчезновением невозмущенной зоны в момент встречи двух температурных фронтов Ь\(т3) = Ь2(т3) и последующей встречей фронтов кристаллизации £1(^3) = е2{93) = ео-
Первый этап стадии замораживания
Расчетная схема первого этапа стадии замораживания показана на рис. 2. Значения теплофизических характеристик замороженного продукта отличаются от значений этих характеристик в незамерзшей зоне, поэтому в дальнейшем для замерзшей зоны используется индекс «зам».
Постановка задачи первого этапа представляет собой систему уравнений Фурье для зон 1 и 4, а физическая модель правой стороны пластины (зоны 2и 3) не меняется.
В ходе решения использованы следующие безразмерные комплексы:
Совместное решение уравнений (8) и (9) позволило определить значения коэффициента асимметрии теплоотдачи
Bij =
ot\L
B¡2 —
a2L
к =
Bi,
B¡9
a2
(10)
— при симметричных условиях теплообмена (к = 1):
*1=Ль У2 = А2, Уз = А3;
— при несимметричных условиях теплообмена (А; > 1):
. =21 = =
Аг А3 В12-
Стадию охлаждения можно считать законченной, когда криоскопическая температура достигается на левой поверхности пластины (£ = 0), причем ранее, чем на ее правой поверхности (£ = 1), где а\ > а2.
Температура в промежутке между двумя температурными фронтами постоянна и равна начальной.
і(х,т) tKр Т {(; 9) — ^tp
^нач ^кр
е+ нач ¿кр) с. Cp3aM(íKp ícpl) л ^зам
olí — -------------------, М4 = ------------------------, Л = --------,
г г а„ез
где 511, 51^ — число Стефана;
г — удельная теплота льдообразования в продукте, Дж/кг;
ср и ср зам — удельная теплоемкость незамороженного и замороженного продукта при постоянном давлении, Дж/(кгК).
дТ4(Ї.Є)
Вії
Рис. 2. Расчетная схема первого этапа стадии замораживания: 1,2 — незамерзшие зоны; 3 — невозмущенная зона;
4 — замерзшая зона;
¿і (г) — координата фронта кристаллизации
Безразмерное уравнение теплопроводности для зон 7 и 4 имеет вид:
(П)
Т4(Єі, 0) — Ті(єі, 0) — Ткр.
Віієі(0)
— при О<£^£і(0), (15)
э27Ж,0) 9Т, Э2Т4(£,0) дТ4
дЕ? дв ' д£2
при 01 < 9 ^ 02, єі(0) < £ < »7і(0) и 0 < С < еі(0).
Начальное условие в безразмерном виде:
Ті (»7ь 01 ) = Тнач, 6(00=0 при 0 = 0і. (12)
Граничное условие в безразмерном виде при £ = 0 представляет собой уравнение, аналогичное (4); на температурном фронте при £ = »7і (0) — уравнение (3).
На фронте кристаллизации £і(0)
(13)
Ті (?, 0) - Тнач (Тнач Ткр) ( ^ I )
при Єі(0) < ^ 771 (0). (16)
В ходе решения данного этапа задачи, с учетом дифференциальных соотношений Л. С. Лейбензона [5], получены уравнения, позволяющие определить значения 771(0) и £1(0) в любой момент времени,
(Т -Т ) + Т - Т —= 77і(0) — £!(0)) ' 4 кр' + "ач ¿9 1к^9
-(2Тнач+Ткр) - -^.у (17)
Ві,
- 281,
7нач1 ТКрі ________ СІ Єї
1 + ВІ!£і(0) 17,і(0)-£і(0) ¿0
при 0Х < 0 ^ 02, £ = Єі- (18)
Окончательное решение системы уравнений (17) и (18) позволило получить зависимость координат температурного фронта »7і(0) и фронта кристаллизации єі(0) от величины Бо (рис. 3).
Задача сводится к решению уравнения теплопроводности с подвижными границами. Условие теплового баланса на фронте кристаллизации (условие И. Стефана):
(14)
Решение поставленной задачи выполнялось с учетом линейного распределения температуры в замороженном слое и в виде полинома второй степени в незамороженном слое:
Рис. 3. Зависимость координат температурного фронта »71(6) и фронта кристаллизации е\(в) от значений Ро
Первый этап стадии замораживания заканчивается, когда температура на правой поверхности пластины становится равной криоскопической: Т2(1,02) = Ткр.
На данном этапе рассматривается одновременное перемещение нескольких температурных фронтов и фронтов кристаллизации. При этом пластина разделяется на пять зон: незамерзшие зоны 1 и 2 толщиной (щ{9) - ¿і(0)) и (<Ь(0) — 1і{9)), замерзшие зоны 4 и 5 толщиной <5і(0) и (1 — <Ь(0)), невозмущенная зона.?толщиной (т^(0) - г?і(0))-Расчетная схема для второго этапа замораживания показана на рис. 4. Для определения продолжительности второго этапа стадии замораживания решаются дифференциальные уравнения теплопроводности для зон 2 и 5.
Ш,в) = Гнач - (Твач - Гкр)
(21)
при 92<9 ^ 03, %(0) < £ ^ £г(0).
В ходе решения получены окончательные выражения для температурного фронта (22) и фронта кристаллизации (23), решение которых также позволило построить зависимость их координат (г]2(9) и є2(0)) от значений Фурье (Бо) — рис. 5:
( 2(Тилч2 - Ткр2)\ ^ СІ7Ь , гг ¿£2 _
V т(в)-є2(в) ) нач2с19 + крсі0 -
1,2 — незамерзшие зоны; 3 — невозмущенная зона;
4, 5 — замерзшие зоны
Постановка задачи на втором этапе также включает в себя начальные и граничные условия, условие Стефана и ряд введенных безразмерных переменных:
в = = f Т>~,4
[ ‘'нач ^кр
Sti = St2, St4 = Sts. (19)
Ш,в) =
Азаі
t(x,r) -
t ' Lcp
Методика решения данного этапа аналогична вышерассмотренным, однако в ходе вычислений и преобразований принимаются:
— температурный градиент для замороженного слоя (зона 5)
дШ,9) В12
<9£ Ві2є2(0) + 1
при 02 < 9 < 03,є2(0) ^ ^ 1;
(20)
Рис. 5. Зависимость координат температурного фронта ^(0) и фронта кристаллизации єг(б) от значений Го
В момент окончания второго этапа замораживания температура в центре продукта Т{щ, 03) становится ниже, чем начальная температура продукта, и продолжает понижаться до тех пор, пока в термическом центре будет достигнуто значение криоскопической температуры. В данный момент начинается фазовый переход (льдообразование), а при встрече фронтов кристаллизации £і(#з) = £г(03) = £о фазовый переход и вся стадия замораживания заканчиваются.
Стадия домораживания
Математическая модель указанной стадии идентична модели стадии охлаждения, однако отличается от нее тем, что температура в термическом центре пластины изменяется. Расчетная схема стадии домораживания показана на рис. 6.
Рис. 6. Расчетная схема стадии домораживания
Постановка задачи на этом этапе представляет собой систему уравнений Фурье для зон 6 и 7, соответствующие им начальные и граничные условия, а также используемые для решения безразмерные комплексы
Ио =
Г(£,0) = *зад. (24)
^•кр ^ср
/Віі№о)+3)\ ВііЄо + 2 К
Вії (Т’ц(Ро) + 1) і 2
го(Віі£о + 2)
£2; (25)
- Ві2є§ + 2Тц(Ро)є0 - Ві2Тц(Ро) - 2Тц(Ро)+
+2Ві2Тц(Ро)є(]^ I ^Ві2£о — 2єо + Віг + 2 — ВігЄо^-
+
2ВІ2Єо(Тц(Ро) + 1)
Ві2£о — 2є0 + Ві2 + 2 — 2Ві2£о
Ві2(Гц(Ро) + 1)
Ві2Єд — 2є0 + Ві2 + 2 — 2Ві2£0
е-
В итоге получено следующее уравнение:
1п(Тц + 1)=1п(Ткр + 1) + С1(0-04),
(27)
где Тц — конечная температура в термическом центре пластины;
С] — безразмерный комплекс,
<?!=-
6Вц
бєо + 2ВЦе0
(28)
Стадия домораживания заканчивается, когда температура в термическом центре продукта достигает заданного значения Т(е0, в5) = Тзад.
Для определения продолжительности этапа домораживания до конечной среднеобъемной температуры Т„ использовались выражения (25), (26), (27) с учетом, что
Т5 =
Тц + Г(£,Ро)
(29)
где
Т(€, Ро) =
Т7(0,Ро) +Т7(1,Ро)
Алгоритм решения поставленной задачи идентичен алгоритму решения стадии охлаждения, однако при этом распределение температуры задавалось полиномами вида
Г.«.*
Проверка адекватности предложенной математической модели проводилась в сравнении с экспериментальными данными по замораживанию широкого ассортимента штучных пищевых продуктов низкотемпературным воздухом от ВТХМ от минус 60 до минус 120 °С. При этом использовались данные, приведенные в работе А. В. Бобкова [1], а также данные, полученные авторами данной статьи при замораживании в этой же установке с последующей сублимационной сушкой пастообразной клейковины (таблица). Основныеусредненныетеплофи-зические характеристики продукции рассматриваемого ассортимента принимались согласно данным справочной литературы [6].
В расчетах варьировались следующие величины:
— температура охлаждающей среды в интервале от минус 60 до минус 120 °С;
— начальная температура продукта на уровне 18-200 С;
— конечная среднеобъемная температура на уровне минус 18 °С;
— толщина продукта в интервале от 6 до 32 мм;
— коэффициент несимметрии к в интервале от 1 до 1,4.
Наименование продукта Толщина продукта L, м Начальная температура продукта ¿нач, °С Температура воздуха £ср, °С Коэффициент несимметрии к Продолжительность Т, ч Погрешность Е, %
Эксперимент Расчет
Бифштексы 0,024 20 -120 1 0,32 0,299 6,56
Котлеты 0,016 20 -120 1 0,197 0,192 2,54
Биточки 0,032 20 -120 1 0,473 0,448 5,28
Шницели 0,008 20 -100 1 0,13 0,135 3,85
Клейковина 0,01 18 -60 1 0,158 0,153 3,26
0,006 18 -60 1 0,138 0,133 3,62
0,012 18 -60 1 0,180 0,173 3,88
0,01 18 -60 1,1 0,159 0,156 1,88
1,2 0,165 0,163 1,21
1,3 0,173 0,170 1,73
1,4 0,174 0,171 1,72
Полученные результаты доказывают, что предлагаемая аналитическая модель обеспечивает достаточно высокую точность вычислений продолжительности замораживания для рассматриваемого ассортимента пищевых продуктов. Так, при симметричных условиях теплообмена к = 1 величина относительной погрешности — от 2 до 7 %, и ее среднее значение составило Е = 4,14 %, а для несимметричных условий к = 1,1—1,4 значения погрешности — на уровне Е = 1,2—1,9 %.
Таким образом, разработанная математическая модель с компьютерной программой расчетов позволяет определить продолжительность как всего процесса замораживания, так и отдельных его этапов, характер распределения температуры в любой момент времени для объектов в форме пластины при симметричных к = 1 и несимметричных 1 < к < 2 условиях теплообмена.
Полученные результаты аналитических расчетов использовались в дальнейших исследованиях, связанных с оценкой условий и режимов предварительного замораживания пищевых продуктов низкотемпературным воздухом от воздушной турбохолодильной машины в технологии их вакуумной сублимационной сушки.
Список литературы
1. Бобков А. В. Разработка проточной системы хла-доснабжения туннельного скороморозильного аппарата с использованием низкотемпературного воздуха от тур-борефрижераторной установки: Автореф. дис. ... канд. техн. наук. — М.: Франтэра, 2004.
2. Аксельрод И. Л., Венгер К. П., Антонов А. А. Расчет продолжительности быстрого замораживания пищевых продуктов в туннельном аппарате с проточной азотной системой хладоснабжения // Вестник МАХ. 2004. № 3.
3. Goodman T. R. Application of integral methods to transient nonliner heat transfer, in: advances in Heat Transfer, Academic Press, San Diego, Cal., USA, 1964. Vol. 1.
4. Бражников А. М. Теория термической обработки мясопродуктов. — М.: Агропромиздат, 1987.
5. Лейбензон Л. Собрание трудов АН СССР. — М.: Изд-во АН СССР, 1955. Т. 4.
6. Антонов А. А., Венгер К. П. Азотные системы хладоснабжения для производства быстрозамороженных пищевых продуктов. — Рязань: Узорочье, 2002.
