УДК 66.01/07(075.8):532.533
МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ МЕХАНИКИ МНОГОФАЗНЫХ СРЕД С ДИСПЕРСНЫМИ ЧАСТИЦАМИ С УЧЕТОМ МЕЖФАЗНОГО ВЗАИМОДЕЙСТВИЯ
А.К. Некрасов,
канд. физ.-мат. наук, доцент, Московский государственный
машиностроительный университет (МАМИ)
E-mail: [email protected]
Е.И. Некрасова,
канд. техн. наук, доцент,
Московский финансово-юридический университет МФЮА E-mail: [email protected]
Аннотация. В статье приведены математическая модель и результаты численного моделирования движения гетерогенной вязкой несжимаемой среды с твердыми сферическими дисперсными частицами в условиях нестационарной свободной тепловой конвекции, вызываемой тепловым воздействием на боковой стенке двумерной полости и вводом в нее твердых сферических частиц. Целями являются: учет влияния межфазного взаимодействия частиц и несущей среды на скорости движения несущей среды и дисперсных частиц, на траектории движения и на распределение их в объеме полости.
Ключевые слова: гидродинамика, математическое моделирование, дисперсные потоки, метод конечных разностей, нестационарный теплообмен.
Abstract. Offered results of modeling of movement of the heterogeneous viscous incompressible environment with firm spherical disperse particles in conditions of non-stationary free thermal convection, to caused constants on time thermal influence on one of lateral walls of a bidimentional cavity and input in settlement area of firm spherical particles are resulted. The purposes are: revealing of influence of interphase interaction of disperse particles and the carrying environment for absolute and relative speeds of movement of the carrying environment and disperse particles, construction of trajectories of their movement and distribution of particles in volume of the heterogeneous environment.
Keywords: hydrodynamics, numerical modeling, disperse flows, method of finite differences, non-stationary heat transfer.
Движение неоднородных многофазных сред с различного рода частицами при свободной тепловой конвекции широко распространено в природе, технике и различных технологиях. Во многих случаях движение вызывается тепловыми и механическими воздействиями на дисперсные среды, которые усложняют гидродинамику и могут оказывать существенное влияние на интенсивность протекающих в этих гетерогенных средах тепло- и массообменных процессов и химических превращений [1]. При исследовании таких многофазных сред одной из важнейших является проблема нестационарного динамического взаимодействия между дисперсной фазой и несущей вязкой несжимаемой средой при их относительном движении в неоднородных по пространству и переменных во времени полях скорости, температуры и концентрации.
В подавляющем большинстве работ, посвященных математическому моделированию механики дисперсных потоков, используются средние по объему однородные или известные из опытов неоднородные объемные или массовые концентрации дисперсных частиц [2]. Реальные же распределения дисперсных частиц в объеме несущей среды, а, следовательно, и их концентрация, непосредственно определяются их динамическим взаимодействием при относительном движении. Поэтому они, как правило, не являются однородными и должны определяться из совместного расчета относительного движения несущей среды и дисперсных частиц.
Для решения подобных задач плодотворным является метод [3-4], согласно которому векторное уравнение движения дисперсной частицы в лагранжевой системе координат решается совместно с уравнениями движения сплошной несущей среды в эйлеровой системе координат. Связь частицы со сплошным потоком учитывается через локальные компоненты напряжения, что составляет оригинальность данного подхода к изучению динамики дисперсной фазы.
В настоящей работе данный подход получил дальнейшее развитие применительно к математическому моделированию свободной тепловой конвекции гетерогенной среды с дисперсными частицами при учете межфазного взаимодействия на динамику несущей среды при постоянных по времени тепловых воздействиях на стенках.
Рассматривается движение гетерогенной вязкой несжимаемой среды с твердыми сферическими монодисперсными частицами, обусловленное нестационарной свободной тепловой конвекцией, вызываемой
постоянным по времени тепловым воздействием на одной из боковых стенок двумерной полости и движением твердых сферических частиц, вводимых в полость. Целями являются: определение абсолютной и относительной скоростей движения несущей среды и дисперсных частиц, вычисление и построение траекторий их движения и определение распределения частиц в объеме гетерогенной среды в различных условиях. Силы присоединенных масс, Бассэ, термофореза, а также взаимодействие частиц между собой и стенками полости не учитывались. Частицы, коснувшиеся стенок полости, исключаются из расчета и удаляются из полости. Теплообмен между дисперсными частицами и несущей средой не учитывался. Расчетная область (рис. 1) представляет собой длинную горизонтальную полость с квадратным сечением размером Ь х Ь, которая заполнена несущей вязкой несжимаемой средой плотностью р. Сверху, с заданным темпом по времени и нулевой начальной скоростью, в полость вводятся твердые сферические дисперсные частицы диаметром dp и плотностьюр . Считаем, что ввод частиц в виде плоской по длине полости струи обеспечивает минимальные градиенты параметров движения двухфазной среды по длине полости. Это позволит получить двумерное движение дисперсной среды в полости и рассматривать задачу в двумерной постановке для квадратного поперечного сечения.
Рис. 1. Расчетная область с тепловыми граничными условиями
В начальный момент времени несущая среда неподвижна и имеет однородную температуру Т. В течение всего процесса левая стенка по-
лости поддерживается при постоянной по времени температуре Т, > Т, правая стенка - при постоянной начальной температуре Т0 . Нижняя и верхняя стенки полости теплоизолированные. Для вязкой несжимаемой несущей жидкости на стенках выполняются условия прилипания.
Неоднородные тепловые граничные условия, приведенные на рис. 1, записаны в относительных температурах: $ = Т - Т0, 3№ = Т - Т0.
Для нахождения поля скорости движения несущей вязкой несжимаемой жидкости, абсолютной и относительной скоростей и построения траекторий движения дисперсных частиц в гетерогенной среде применен метод, который основан на совместном решении уравнений для несущей среды, записанных в эйлеровой системе координат, и уравнений для дисперсных частиц, записанных в лагранжевой системе координат. Этот метод подробно изложен в работе [3] и применен для моделирования относительного движения фаз в дисперсных потоках, как при вынужденной [4], так и при естественной конвекции [5].
Описание ламинарного движения вязкой несжимаемой несущей среды при свободной тепловой конвекции выполнено в приближении Обербека-Буссинеска [6] на основе системы нестационарных уравнений Навье-Стокса, неразрывности и энергии. Влияние частиц на течение несущей жидкости учитывается источниковыми членами в уравнениях Навье-Стокса, определяемыми точечными силами межфазного взаимодействия [7], которые локально распределены в расчетной области и равны силам гидродинамического сопротивления, оказываемым несущей жидкостью на частицы при их относительном движении, отнесенным к единице объема.
В безразмерных переменных относительно соответствующих масштабов (для пространства - Ь; для времени - Ь2/; для скорости - у/Ц для давления - ру2/Ь2; для силы межфазного взаимодействия - ру2/Ь3; для температуры - 3 w) векторная запись этой системы уравнений имеет вид:
дУд + (VУ)¥= -Vр + АУ + пОтв - В, (1)
сИ\У = 0, (2)
двд + (^ в = Рг1 Ав, (3)
где: VI, р - соответственно безразмерные вектор скорости несущей среды, время, давление; От = Оа. 3„ вт- критерий Грасгофа; Оа = gL3/v2 -
критерий Галилея; V - коэффициент кинематической вязкости, м2/с; g - ускорение свободного падения, м/с2; g- коэффициент термического расширения, 1/К; Вр - безразмерная сила межфазного взаимодействия, отнесенная к единице объема; Рг = v/a - критерий Прандтля; а - коэффициент температуропроводности, м/с2; 0 = - безразмерная температура; п - единичный вектор направления силы тяжести.
Краевыми условиями для системы (1)-(3) в двумерной квадратной расчетной области принимались:
t = 0: V = 0, 9 = 0, (4)
у = 0; 0 < х < 1, 1: д9/дх = 0, V = 0, (5)
X = 0, 0 < у < 1: 9 = 1, V = 0, (6)
X =1, 0 < у < 1: 9 = 0, V = 0, (7)
здесь X, у - безразмерные координаты по пространству.
Скорости Vp¡ для Ж дисперсных сферических частиц, движущихся в неоднородном поле скорости несущей среды, находятся из решения уравнений движения частиц в двухфазной среде, записанных в лагран-жевой системе координат для условий действия различных сил
т . dV /dx = Г + ГО + ¡ = 1,...,Ы, (8)
p¡ p¡ р О А ' ' р' ^ 7
где: т. - масса ¡-дисперсной частицы, кг; Г,, ГО, Г - векторы, соответственно, действующих на частицу сил: гидродинамического сопротивления, тяжести и Архимеда, Н.
Сила гидродинамического сопротивления для дисперсных частиц вычисляется по формуле
Гр = - 0,5Ср пгр2рм2е, (9)
где: Ср - коэффициент сопротивления; гр = dp/2 - радиус дисперсной частицы, м; м = \Vотн \ - модуль вектора относительной скорости частицы, м/с; е - единичный вектор направления относительной скорости частицы.
Входящий в соотношение (9) коэффициент сопротивления Ср для твердой сферической частицы, движущейся в несущей жидко-
сти с локальной относительной скоростью при числе Рейнольдса Re = wd/v < 1 определялся по соотношению Cp=24/Re [1].
Вектор скорости дисперсной частицы V представляется в виде суммы векторов скорости сплошной среды V и относительной скорости частицы V
отн
V = V + V .
р отн
Вектор относительной скорости V может быть выражен через мо-
отн
дуль вектора относительной скорости w и его угол направления а
V = (w cos а, w sin а) = we.
отн
Угол а равен углу поворота от орта i оси координат х до вектора е.
В работах [3; 4] показано, что для определения скорости частицы векторное уравнение (8) приводится к двум скалярным нелинейным уравнениям с начальными условиями относительно неизвестных w и а.
Для определения текущих координат и движущейся в пространстве расчетной области дисперсной частицы решались уравнения, также записанные в безразмерном виде
dx/dt = U + w cos а, dy/dt = V + w sin а
с начальными условиями при ? = 0: х = х0 и у = у
Решение сформулированной задачи получено численно методом конечных разностей. Система уравнений для несущей среды (1)-(7), записанная в переменных вихрь - функция тока - температура, решалась методом переменных направлений [6; 8]. Системы из двух скалярных нелинейных уравнений, для модуля вектора относительной скорости частицы м и угла направления вектора а, решались методом Ньютона.
Математическое моделирование проведено для газовзвеси при значении числа Прандтля несущей среды Рг = 0,7. Твердая фаза с относительной плотностью частиц р/р= 0,0005, что соответствует частицам кварцевого песка, подается для обеспечения двумерного течения с однородным по длине и ширине верхней поверхности полости распределением частиц. Подача частиц в полость производится синхронно по всей поверхности подачи с шагом по времени равным прД?. Здесь п - число шагов интегрирования по времени при численном решении
задачи между подачей в полость очередной порции частиц; At - шаг интегрирования по времени. Безразмерный шаг интегрирования по времени в расчетах, результаты которых приводятся ниже, принимался постоянным и равным 0,00002.
На рис.2 приведены расчетные поля температуры (а), векторов скорости несущей среды (б) и распределения дисперсных частиц с векторами скорости (в) при Ог = 105, Рг = 0,7, dр = 10 мкм, пр = 10. Из рисунка видно, что учет в расчете обратного влияния дисперсных частиц на несущую среду привел к существенной перестройке течения. При сравнении приведенного на рисунке распределения температуры с распределениями температуры при свободной термогравитационной конвекции в квадратной полости вязкой несжимаемой однофазной жидкости, приведенными в работах [6-8], можно видеть, что падающие дисперсные частицы сужают область восходящего течения у нагретой левой стенки полости. В центральной части полости на месте застойной зоны с практически нулевой скоростью, которая характерна для случая однофазной среды, возникают дополнительные нисходящее и восходящее течения. Они обусловлены воздействием падающих частиц наиболее удаленных от нагретой левой стенки. Эти течения хорошо видны на поле векторов скорости несущей среды (рис. 2б). На рис. 2в видно, что часть частиц, находящиеся в начале движения вблизи нагретой стенки, захватываются восходящим течением несущей среды и витают некоторое время у нагретой стенки.
Как показало моделирование, в характерном сечении в движении участвуют дисперсные частицы в количестве N = 430 ± 20. Отклонение количества частиц от номинала объясняется различием между количеством вводимых в полость частиц и количеством осаждающихся частиц, т. к. при учете межфазного воздействия частиц на несущий поток седиментация частиц носит периодический характер относительно некоторого среднего значения. Поэтому течение двухфазной среды не достигает стационарного состояния, а происходят плавные колебания параметров течения. Средняя объемная концентрация частиц, определяемая для всей полости в рассматриваемом случае, ар = 0,003. Но на рисунке видно, что в полости имеется обширная зона свободная от частиц. А в зоне полости, в которой сосредоточены дисперсные частицы, их распределение характеризуется большой неоднородностью, следовательно, локальная объемная концентрация частиц будет значительно больше средней по объему полости.
•-•чщГ1
Л к
II.
11 Шг/ .^Ж^ЛШк^
""'ЩЩД'
■""щит1
"'ЩШЧ' ""¿тт'
«47/////П''
а б в
Рис. 2. Расчетные поля температуры (а), векторов скорости несущей среды (б) и распределения дисперсных частиц (в) при Gr = 105, Pr = 0,71, dp = 10 мкм
Количество дисперсных частиц, входящих в образующуюся двухфазную систему, зависит от темпа подачи частиц в полость (величина шага по времени между вводом частиц), диаметра и плотности частиц, физических свойств несущей среды (критерия Рт), геометрических параметров полости и тепловых условий на стенках полости (критерия От).
Частицы, относительная скорость движения которых соизмерима со скоростью движения несущей среды, начинают витать в полости. Витание наиболее интенсивно происходит вблизи нагретой стенки, у которой скорость движения несущей среды имеет наибольшие значения.
Увеличение числа Грасгофа до 106 (за счет увеличения 3 ^), при сохранении значений прочих параметров, приводит к повышению скорости движения несущей среды при свободной конвекции. Снижается влияние дисперсных частиц. На рис.3 приведены мгновенные поля распределения температуры (а) и векторов скорости (б) несущей среды,
а б в
Рис. 3. Расчетные поля температуры (а), векторов скорости несущей среды (б) и распределение дисперсных частиц (в) при Gr = 106,
Pr = 0,7, d = 10 мкм
' ' p
которые мало отличаются от их распределений для случая отсутствия влияния частиц. Число дисперсных частиц в двухфазной системе в этом случае увеличивается до N = 620 ± 30 и средняя по полости объемная концентрация частиц также возрастает до ар = 0,004. Видно, что число витающих частиц в полости увеличилось. Это связано с тем, что скорость движения несущей среды при этом числе Грасгофа становится соизмеримой с относительной скоростью движения дисперсных частиц.
Так, при сохранении темпа подачи дисперсных частиц и прочих параметров, увеличение диаметра частиц до ^ = 50 мкм число частиц в двухфазной системе из-за быстрой седиментации резко уменьшается до N = 22 ± 1. Но средняя объемная концентрация, из-за увеличения диаметра частиц, уменьшается незначительно (ар = 0,0037).
На рис. 4 представлены некоторые результаты для движения дисперсных частиц диаметром dp = 50 мкм при Ог = 106, Рг = 0,7, пр = 5. С увеличением массы частиц возрастает скорость их осаждения, которая значительно превышает скорость свободноконвективного движения несущей среды. Видно, что характер движения дисперсной среды меняется (рис. 4а).
а б в
Рис. 4. Расчетные поля температуры (а), векторов скорости несущей среды (б) и распределение дисперсных частиц (в) при Gr = 106, Pr = 0,7, dp = 50 мкм и п = 5
Межфазное взаимодействие при этих значениях параметров у горячей стенки полости приводит к движению дисперсной среды в обратном направлении по сравнению со свободноконвективным движением. В центре полости появляется восходящий поток среды. У холодной стенки нисходящий поток среды сохраняется, но скорости движения в нем по сравнению с первыми двумя потоками значительно меньше.
Число частиц в полости в квазиустановившемся режиме в этом случае Np = 38 ± 2, а средняя объемная концентрация частиц ар = 0,0064.
Проведенное численное исследование выявило существенные неоднородности распределения частиц в объеме гетерогенной среды, обусловленные гидродинамикой, формируемой силовым межфазным взаимодействием дисперсных частиц и несущей сплошной средой.
ЛИТЕРАТУРА
1. Нигматулин Р.И. Основы механики гетерогенных сред. - М., 1978.
2. Волков К.Н., Емельянов В.Н. Течения газа с частицами. - М., 2008.
3. L.P. Kholpanov, B.R. Ismailov, P.Vlasak. Modelling of Multiphase Flow Containing Bubles, Drops and solid Particles // Eng. Mech., 2005, Vol. 12, No. 6.
4. Холпанов Л.П., Ибятов Р.И.. Математическое моделирование динамики дисперсной фазы // ТОХТ. 2005. Т. 39, № 2.
5. Холпанов Л.П., Некрасова Е.И., Некрасов А.К. Математическое моделирование динамики дисперсной фазы при свободной гравитационной конвекции вязкой несжимаемой жидкости в квадратной полости // ИФЖ. 2008. Т. 81, № 1.
6. Полежаев В.И., Бунэ А.В., Верезуб Н.И. и др. Математическое моделирование конвективного тепломассообмена на основе уравнений Навье - Стокса. - М., 1987.
7. P.G. Saffman. He Setting Speed of Free and Fixed Suspensions // Studies in Applied Mathematics. 1973. Vol. LII, №2.
8. G. De. Vahl Davis. Natural convection of air in a square cavity: a bench mark numerical solution // Int. J. for Numerical Methods in Fluids. 1983. Vol. 3, №3.