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

УДК 637.5.037+66-973

Моделирование промерзания биомассы при криозаморозке

Канд. техн. наук А.В. Зайцев, zai_@inbox.ru, С.Е. Кублицкий, stanleykub@mail.ru, д-р техн. наук В.В. Пеленко, pelenko1@rambler.ru

Университет ИТМО 191002, Россия, Санкт-Петербург, ул. Ломоносова, д.9

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

Ключевые слова: криозаморозка; тепловой баланс; уравнение теплопроводности; граничные условия; теплофизические характеристики; численное решение.

DOI: 10.17586/2310-1164-2017-10-1-9-19

Simulation of biomass freezing at cryofreezing

Ph.D. Andrei V. Zaitsev, zai_@inbox.ru, Stanislav E. Kublitskiy, stanleykub@mail.ru, D. Sc. Valery V. Pelenko, pelenko1@rambler.ru ITMO University 9, Lomonosov str., St. Petersburg, 191002, Russia

The task of modeling thermal processes during cryogenic freezing of the objects of biological origin (in particular, the foodstuff) is formulated and solved by numerical methods. Mathematical model of cryofreezing for biological objects in the form of heat balance equation (heat conductivity equation) is presented and its numerical realization is made by the method of final differences. To ensure the information value of the solution propsed, its analysis, as well as the possibility of greater variation and changes in the range of the task parameters in the research process of cryogenic freezing, our own calculation program of the object processes is written in Fortran. The following parameters are chosen as the boundary conditions: constant temperature on the heat removal surface which is the boundary condition of the type I, and constant density heat flow on this surface which is the boundary condition of the type II. The dependence of thermophysical properties of biological tissues from the target temperature contributes significant nonlinearity to the thermal conductivity equation and makes it difficult to solve the problem by analytical methods. Three stages of cryogenic freezing for the objects of different geometric shapes are considered: cooling the material to the cryoscopic temperature, freezing the liquid phase, cooling the ice to the desired temperature. The application of the technique proposed allows forecasting the cooling regime for the purpose of the necessary temperature achievement in the specified points of the object in the given time. Thus, at the relatively lower energy consumption, such technologies allow to obtain products and high quality objects while preserving most of their native properties, biologically active elements, and to implement energy-saving processes.

Keywords: cryogenic freezing; heat balance; thermal conductivity equation; boundary conditions; thermophysical properties; numerical solution.

Введение

Криовоздействие на биоткань заключается в обработке биологических объектов при использовании криогенных температур (ниже 120 К). В настоящее время такая технология представляет собой перспективное направление развития медицины, биологии, пищевой промышленности и других отраслей экономики [1]. Актуальность вопроса расширяется и глубже раскрывается особенно на фоне растущего спроса на технологии повышения операционной эффективности технологических процессов в медицине и пищевой промышленности [2-5].

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

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

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

Объекты исследования

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

Охлаждаемый объект представляет собой некое полуограниченное (рисунок 1, а) или ограниченное (рисунок 1, б) тело, от некоторой поверхности которого средствами криогенных технологий отводится теплота при постоянной температуре ТО (граничное условие I рода) или постоянной плотности теплового потока qo (граничное условие II рода). Возможно и более сложное задание этих граничных условий (ГУ), когда они являются функцией времени: ТО = То(т) или qo = qo(т).

Рисунок 1 - Расчетная схема объекта:

- поверхность теплоотвода;

- поверхность теплообмена с окружающей средой;

- граница раздела фаз

В процессе охлаждения часть тела претерпевает фазовый переход из условного жидкого состояния в твердое, т.е. происходит процесс льдообразования.

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

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

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

Граничные и начальные условия. Рассмотрим математическую формулировку задачи для схемы, приведенной на рисунке 1, а.

Необходимо рассчитать температурное поле T (r, z, т) T2 (r, z, т). В начальный момент времени

т = о температура во всех точках объекта одинакова и равна T (r, z,0) = T0 . В определенный момент расчетная зона разделяется на две части с различным фазовым состоянием и соответствующим распределением температур: T1 (r, z, т) - для твердой фазы (лед) и T2 (r, z, т) - для жидкой фазы.

Для определенности рассмотрим следующие граничные условия (рисунок 2):

- граничное условие I рода на поверхности теплоотвода

T (r, z, т) = Тх = const, при z = 0; 0 < r < r0 ; (l)

- условие теплоизоляции (II рода) на внешней поверхности

dT (r, z, т) z = 0; r > r0 ,пЛ

= 0, при ' 0; (2)

dz

- условие неразрывности теплового потока на границе раздела фаз Fi

. dT!(r,z,т) = . dT2(r,z,т) при r,z e Fx (3)

Aj---= -;-> при 1 j (3)

an an

где Aj, A2 - коэффициенты теплопроводности соответствующих областей, ВтДм-К);

n - нормаль к поверхности Fi;

- условие симметрии относительно оси о-z, т.е. отсутствие теплового потока через эту границу

= 0 , при r = 0 ; (4)

dr

- условие постоянства температуры на поверхности F2 вдалеке от зоны охлаждения

T (r, z, т) = T0 = const, при г, z e F2 r, z e F . (5)

А

Жидкость

Т2 (г, г, т) чч ^

Лед X -

Ч \

Т1 (г,г,т) Fl \

оТТЕ

Рисунок 2 - Граничные условия: Го - радиус поверхности теплоотвода; Fl - поверхность фазового перехода; F2 - удаленная поверхность прекращения

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

Уравнение теплопроводности. Запишем уравнение теплопроводности для сформулированной задачи [9-12]

dT ( г, г, т) срр =х

d 2Т (г , т) d 2Т ( г, г, т)

dr2

+ -

dг2

+ Чу

(6)

где ср - удельная теплоемкость, ДжДкг-К); р - плотность, кг/м3;

Чу - объемная плотность тепловыделений при фазовом переходе (теплота кристаллизации), Вт/м3.

Давление всей системы принимаем постоянным, равным атмосферному. Теплофизические свойства объекта охлаждения, рассчитываемые по рекомендациям в существующих литературных источниках [13], приведены в таблице 1.

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

Таблица 1 - Теплофизические свойства биологических тканей в зависимости от температуры

z

Свойство Значение

теплопроводность X, Вт/(м-К) 0,5 при Т > 273 К; 15,98-0,05б7Т при 273 К > Т > 251 К; 1005Т-115 при Т < 251 К

объемная теплоемкость С, кДж/(м3-К) 3600 при Т > 273 К; 15,44 при 273 К > Т > 251 К; 3,98Т при Т < 251 К

скрытая теплота г, МДж/м3 300

теплоприток от прохождения (перфузии) крови WbCb, кВт/(мз-К) 40

Методы исследования

Как отмечалось ранее, для максимальной информативности модели и эффективности анализа, а также возможности варьирования большинства параметров задачи в процессе исследований, была написана собственная гибкая программа на языке программирования Фортран.

Применение метода конечных разностей. В описанной выше постановке наиболее эффективным способом решения задачи представляется численное решение приведенных уравнений методом конечных разностей [14].

Построим сеточную модель расчетной области. Обозначим индексом i номер текущего узла по координате г, индексом ] - по координате z, индексом к - номер текущего временного слоя (шага по времени т).

Перейдем от поля температур Т (г, т) к сеточной функции Т,

Для упрощения вида разностного уравнения введем обозначения:

ГТ1 ГТ~Ч ГТ1 ГТ~Ч ГТ1 ГТ~1 ГТ1 ГТ1 ГТ1 ГТ1 ГТ1 гт~

Т -1,j,k = Т -1; Т +1,j,k = Т+1 ; 5j-1,к = 1j-1; 5j+1,к = 1j+1; 5j,k = 1 ; 5j,k+1 = 1 .

Обозначим шаги по координатам и по времени соответственно Ат, Аг и Ат .

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

д2Т _Т+1 - 2Т + д2Т_1+! - 21 +1 -

дт г

Ат2

дг2

Аг2

Соседние узлы могут принадлежать подобластям с различающимися теплофизическими

д2Т

свойствами. Тогда, например, для X- следует записать

дт 2

X

д2Т_ х,+1 {1+1 - Т) + х,-11_1 - Т)

дт2

Ат2

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

Т' - Т_ а+1 (Тм - Т) + 0,-1 (1,-1 - Т) , — — -1) + — (— -1) qy (т)

Ат

Ат2

Аг 2

(т)р(т)

(7)

где а = (срр) - коэффициент температуропроводности, м2/с.

Объемную плотность источника теплоты в процессе фазового перехода определяем с учетом скрытой теплоты при замерзании вещества qv = т / дт, Вт/м3. Преобразуем уравнение (7)

Т ' = Т +

а

7+1

(Т+1 - Т)+а-1 (Т-1 - Т)

Ат2

Ат +

а--+1 (1-+1-1)+а--1 (Т--1-1) Аг 2

Ат + -

СрР

(8)

Это уравнение годится для любого фазового состояния; при переходе в иное фазовое состояние следует соответственно изменить значения коэффициента температуропроводности, теплоемкости и плотности и учесть последнее слагаемое.

Введем нумерацию для обозначения различных геометрических элементов рассматриваемой области в разрабатываемой программе Фортран (рисунок 3).

Перепишем граничные условия (1)-(4) в конечных разностях.

Т-=0 = Тх = const, при 0 < т < т0;

Т-=1 - Т-=0

Аг

= 0, при т > т0;

X, Тм^Т=Х21-5=1; X,Т-+■-1

Ат

Ат

Аг

Т=1 - Т=0 Ат

= X

Т - Т

--1

Аг

, при Т = 0;

= 0.

(9) (10)

(11)

(12)

с

т

Рисунок 3 - Нумерация обозначений геометрии расчетной области в программе: 1 - внутренняя область (тело); 2 - граница теплоотвода; 3 - теплоизолированная граница; 4 - теплоизолированная граница (ось симметрии); 5 - граница раздела фаз; 6 - условная внешняя граница

Здесь на границе фазового перехода вместо производных по нормали используются их проекции на оси г и z. Граница фазового перехода определяется динамически из условия Т = 0°С. Записи i = 0, ] = 0 в сеточной модели означают узлы на границах г = 0 и z = 0, соответственно; i = 1, = 1 - соседние узлы внутри тела.

Особенности расчета фазового перехода. Рассмотрим процесс перехода элемента объекта из «жидкого» состояния в «твердое», т.е. процесс замерзания. Будем условно называть состояние элемента «жидкость» и «лед». Выделяем 3 этапа: охлаждение жидкости, замерзание, охлаждение льда.

Считаем, что процесс замерзания начинается при достижении элементом температуры 0°С. В результате применения уравнения (8) можно определить продолжительность каждого этапа. 1. Охлаждение жидкости от Т до Т'=0°С.

0 = Т +

'о+1 (Т+1 - Т) + Ом (Тч - Т) + (Т+1 - Т) + а]_1 (Т-_! - Т)

Аг2

Аг2

Ат1,

откуда

Ат1 = --

Т

а+1 (Т+1 - т )+а,-1 (Т-1 - Т ) + О-+1 (Т+1- Т)+а- (Т-1- Т)

Аг

Аг2

(13)

2. Замерзание жидкости при Т'=Т=0°С.

0 =

а,+1Т+1 + аг-1Тг-1

Аг 2

+

а.! +Т+1 + аУ-Т-1

Аг 2

Ат2 +-

срР

откуда

Ат0 = --

г \ е

( срР)

)+1т,+1+а-т-1, а+1Т]+1+а]-1Т]-1

а.

Аг2

+

Аг2

3. Охлаждение льда от Т=0°С до Т'

Т ' =

а,+Т+1 + а,-1Т,-1 + а1+1Т]+1 + а]-1Т]-1

Аг

2

Аг

2

Ато

откуда

г

Т'

Ат3 =-^-^^ . (15)

д+Д+1 + а-1Тг-1 + а)+1Т)+1 + а] -1Т]-1 Аг2 А^2

Таким образом, применение уравнений (13)-(15) позволяет разбить временной шаг Дт на определенные части и избежать необходимости применения граничных условий (11), учитывая однако скачкообразное изменение теплофизических свойств на этой границе.

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

Обсуждение результатов

Особенности расчета в области «жидкого» состояния. За границей раздела фаз 5 (рисунок 3) начинается область жидкого состояния, которая простирается до границы 4 и далее, где температура остается неизменной и равной Т0. В математическом смысле эта граница удалена в бесконечность. В физическом смысле можно считать, что граница достигнута, когда разница температур между двумя соседними вдоль координаты элементами (узлами разностной сетки) становится меньше заданного минимального значения, например ДТтп = 0,001°С. Тогда граничное условие (5) в разностной форме:

Т . - Т < АТ . ; Т. , - Т <АТ. . (16)

1+1 тт' ] +1 тт

Система разностных уравнений (8)-(10); (12)-(1б) является замкнутой и позволяет численным методом получить поле температур для сформулированной задачи. В соответствии с этими уравнениями разработан алгоритм и компьютерная программа, которые не являются предметом рассмотрения данной статьи.

Рисунок 4 иллюстрирует передвижение границы раздела твердой и жидкой фаз для следующих условий: шаг сетки по обеим координатам 1 мм; шаг по времени 0,1 с; радиус поверхности теплоотвода (криовоздействия) 4 мм; температура криовоздействия соответствует азотному уровню (-196°С); начальная температура биомассы 20°С.

г, мм 20

15

10

5

0

Рисунок 4 - Расчетное положение границы раздела фаз

Рисунок 5 отражает одну из важнейших характеристик процесса - скорость промерзания до заданной температуры на заданную глубину. Данное семейство кривых хорошо согласуется с известными экспериментальными данными, приведенными например в [15-16], что подтверждает адекватность представленной модели.

Рисунок 5 - Динамика изменения температуры биообъекта на различном расстоянии от места

криовоздействия

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

Вместо граничного условия (5) на внешней границе 6 используется условие теплоотдачи при теплообмене на внешней поверхности (ГУ III рода).

- Расчетная схема объекта: поверхность теплоотвода;

поверхность теплообмена с окружающей средой;

граница раздела фаз;

теплоизоляция

Рисунок 6 -

Рассмотрим эллипсовидную форму сечения. Поверхность теплоотдачи описывается уравнением

г2 и - с)2 , Л

^ + = 1, (17)

а Ь

где а и Ь - большая и малая полуоси, соответственно; с - сдвиг большой полуоси эллипса вверх относительно оси z = о.

Тогда граничное условие на внешней границе

¿Т(г % т)

X-—— = а(Т - Гос), при выполнении (17).

¿п

В конечных разностях граничное условие на внешней границе:

Г - Т Х ^ = а(Т - Тос).

Окончательно получаем, учитывая, что г = (г -1)Дг и z = (7' -1)^:

„ ХТ , + аАгТос 4 — I 2 г7 ,4.12 _

Г= _м-ос, при = —т——-Г(/- 1)Аг 1 + с. (18)

аАг + Х а (j -1)^ ^ ^

Кроме того, из системы исключается уравнение (16).

Альтернативным вариантом моделирования объекта с неправильной формой тела является задание массива координат (г, z) и численная аппроксимация по этим точкам кривой, описывающей теплоотдающую поверхность тела.

Рассмотрим следующий пример: Дг = 4 мм; Дz = 4 мм; Дт = 0,2 с; а = 10 см; Ь = 8 см; с = 1 см; Го = 4 см; Т 0 = 20°С; Тх = -196°С; Тос = 20°С; а = 10 Вт/(м2-К).

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

z, см

Рисунок 7 - Динамика передвижения изотерм процесса промерзания биомассы:

„ууууууу. - граница между твердой и жидкой фазами (X = о°С); ---- граница распространения температурного воздействия (Т = 20°С)

Выводы

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

биообъектов, в частности продуктов питания [17]: скорость оценки результата заморозки на порядок выше скорости самого процесса заморозки, поэтому включение системы оценки результата криовоздействия в реальном времени позволяет автоматизировать управление процессом заморозки, минимизируя энергозатраты, стоимость установки, незапланированное разрушение клеточной структуры объекта.

Работа выполнена при финансовой поддержке РФФИ (проект №16-38-00728 мол_а)

Литература

1. Зайцев А.В. Криогеника в начале XXI века // Научный журнал НИУ ИТМО. Серия: Холодильная техника и кондиционирование. 2014. № 1.

2. Boonsumrej S., Chaiwanichsiri S., Tantratian S., Suzuki T., Takai R. Effects of freezing and thawing on the quality changes of tiger shrimp (Penaeus monodon) frozen by air-blast and cryogenic freezing. Journal of Food Engineering, 2007, no. 80(1), pp. 292-299.

3. Ali S., Zhang W., Rajput N., Khan M.A., Li C.-B., Zhou G.-H. Effect of multiple freeze-thaw cycles on the quality of chicken breast meat. Food Chemistry, 2015, no.173, pp. 808-814.

4. Alizadeh E., Chapleau N., De Lamballerie M., LeBail A. Effects of freezing and thawing processes on the quality of Atlantic salmon (Salmo salar) fillets. Journal of Food Science, 2007, no. 72(5), pp. E279-E284.

5. Chun H.H., Choi E.J., Han A.R., Chung Y.B., Kim J.S., Park S.H. Changes in quality of hanwoo bottom round under different freezing and thawing conditions. Journal of the Korean Society of Food Science and Nutrition. 2016, no. 45(2), pp. 230-238.

6. Шлапак С.В., Шокина Ю.В. Влияние режимов замораживания на показатели пищевой и биологической ценности рыбных пельменей, изготовленных с использованием ската звездчатого // Материалы международной научно-практической конференции «Современные эколого-биологические и химические исследования, техника и технология производств». Мурманск: Издательство МГТУ, 2015. С. 202-209.

7. Быков Р.В. Криозамораживание продуктов // Мясные технологии. 2010. № 3. С. 30.

8. Соболев С.Л. Уравнения математической физики. М.: Наука, 1966. 444 с.

9. Архаров А.М., Архаров И.А., Антонов А.Н. и др. Машины низкотемпературной техники - Криогенные машины и инструменты. М.: Изд-во МГТУ им. Н.Э. Баумана, 2011. 582 с.

10. Rabin Y., Korin E. An efficient numerical solution for the multidimensional solidification (or melting) problem using a microcomputer. Int. J. Heat Mass Transfer. 1993, V. 36, no. 3, pp. 673-683.

11. Пушкарева А.Е. Методы математического моделирования в оптике биоткани. СПб.: Изд-во Университета ИТМО, 2008. 103 с.

12. Rabin Y., Shitzer A. Numerical Solution of the Multidimensional Freezing Problem During Cryosurgery. Journal of Biomechanical Enguneering, 1998, V. 120, pp. 32-37.

13. Tanaka D., Shimada K., Rabin Y. Two-Phase Computerized Planning of Cryosurgery Using Bubble-Packing and Force-Field Analogy. Journal of Biomechanical Engineering. 2006, V. 128, pp. 49-58.

14. Андерсон Д., Таннехилл Дж., Плетчер Р. Вычислительная гидромеханика и теплообмен. В 2-х т. М.: Мир, 1990. 728 с.

15. Цыганов Д.И. Криомедицина: процессы и аппараты. М.: Сайнс-Пресс, 2011. 304 с.

16. Hewitt P.M., Zhao J., Akhter J., Morris D. L. A Comparative laboratory study of liquid nitrogen and argon gas cryosurgery systems. Cryobiology. 1997, V. 35, pp. 303-308.

17. Быков Р. LINDE GAS RUS: криозаморозка продуктов // Империя холода. 2009, январь, С. 43-44.

References

1. Zaitsev A.V. Kriogenika v nachale XXI veka [Cryogenics at the beginning of the 21st century]. Scientific journal NRUITMO. Series: Refrigeration and Air Conditioning. 2014, no. 1.

2. Boonsumrej S., Chaiwanichsiri S., Tantratian S., Suzuki T., Takai R. Effects of freezing and thawing on the quality changes of tiger shrimp (Penaeus monodon) frozen by air-blast and cryogenic freezing. Journal of Food Engineering, 2007, no. 80(1), pp. 292-299.

3. Ali S., Zhang W., Rajput N., Khan M.A., Li C.-B., Zhou G.-H. Effect of multiple freeze-thaw cycles on the quality of chicken breast meat. Food Chemistry, 2015, no.173, pp. 808-814.

4. Alizadeh E., Chapleau N., De Lamballerie M., LeBail A. Effects of freezing and thawing processes on the quality of Atlantic salmon (Salmo salar) fillets. Journal of Food Science, 2007, no. 72 (5), pp. E279-E284.

5. Chun H.H., Choi E.J., Han A.R., Chung Y.B., Kim J.S., Park S.H. Changes in quality of hanwoo bottom round under different freezing and thawing conditions. Journal of the Korean Society of Food Science and Nutrition. 2016, no. 45(2), pp. 230-238.

6. Shlapak S.V., Shokina Yu.V. Vliyanie rezhimov zamorazhivaniya na pokazateli pishchevoi i biologicheskoi

tsennosti rybnykh pel'menei, izgotovlennykh s ispol'zovaniem skata zvezdchatogo [Influence of the modes of freezing on indicators of nutrition and biological value of the fish pelmeni made with use of a slope star-shaped]. Proceedings of the international scientific-practical conference "Modern ecological and biological and chemical research, engineering and technology industries". Murmansk, MGTU Publ., 2015, pp. 202-209.

7. Bykov R.V. Kriozamorazhivanie produktov [Cryofreezing of products]. Meat technology, 2010, no. 3, P. 30.

8. Sobolev S.L. Uravneniya matematicheskoi fiziki [Equations of mathematical physics]. Moscow, Nauka Publ., 1966, 444 p.

9. Arkharov A.M., Arkharov I.A., Antonov A.N. i dr. Mashiny nizkotemperaturnoi tekhniki - Kriogennye mashiny i instrument [Machinery of the low-temperature equipment - Cryogenic machinery and tools]. Moscow, MGTU im. N.E. Baumana Publ., 2011, 582 p.

10. Rabin Y., Korin E. An efficient numerical solution for the multidimensional solidification (or melting) problem using a microcomputer. Int. J. Heat Mass Transfer. 1993, V. 36, no. 3, pp. 673-683.

11. Pushkareva A.E. Metody matematicheskogo modelirovaniya v optike biotkani [Methods of mathematical modeling in biofabric optics]. St. Petersburg, ITMO University Publ., 2008. 103 p.

12. Rabin Y., Shitzer A. Numerical Solution of the Multidimensional Freezing Problem During Cryosurgery. Journal of Biomechanical Enguneering, 1998, V. 120, pp. 32-37.

13. Tanaka D., Shimada K., Rabin Y. Two-Phase Computerized Planning of Cryosurgery Using Bubble-Packing and Force-Field Analogy. Journal of Biomechanical Engineering. 2006, V. 128, pp. 49-58.

14. Anderson D., Tannekhill Dzh., Pletcher R. Vychislitel'naya gidromekhanika i teploobmen [Computing hydromechanics and heat exchange]. Moscow, Mir Publ., 1990. 728 p.

15. Tsyganov D.I. Kriomeditsina: protsessy i apparaty [Cryomedicine: processes and devices]. Moscow, Sains-Press Publ., 2011, 304 p.

16. Hewitt P.M., Zhao J., Akhter J., Morris D. L. A Comparative laboratory study of liquid nitrogen and argon gas cryosurgery systems. Cryobiology. 1997, V. 35, pp. 303-308.

17. R. Bykov. LINDE GAS RUS: cryofreezing of products. Empire of cold. 2009, January, pp. 43-44.

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