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

© «Вестник ИГЭУ» Вып. 4 2007 г.

УДК 004.942::[544.27:544.015.5]

О НОВОМ РЕШЕНИИ ПРОБЛЕМЫ РАСЧЕТА РАСПРЕДЕЛЕНИЙ ПРИ МОДЕЛИРОВАНИИ НЕКОТОРЫХ МНОГОФАЗНЫХ СИСТЕМ

ПЕКУНОВ В.В., канд. техн. наук

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

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

NEW SOLUTION OF ASSIGNMENT CALCULATION PROBLEM UNDER SOME MULTIPHASE SYSTEMS SIMULATION

V.V. PEKUNOV, Ph.D.

This paper is devoted to the problem of sufficient drop assignment search under numerical simulation of multiphase systems, which include drop phase. The author suggests a new approach towards piecewise-linear assignment calculation by means of interpolation regarding density and drop component concentration. This method was approved in numerical experiments. There was demonstrated adequate accuracy and low computing labouriousness of the approach suggested.

Key words: sufficient drop assignment, interpolation regarding density, average accumulated relative accuracy, analysis regarding trends.

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

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

2. Эйлерово-эйлеровы модели (модели сплошных сред). Выделим четыре основных подхода:

а) применение фиксированных распределений капель (дает недостаточную точность [1]);

б) подход А.Е.Алояна [2] с решением кинетического уравнения (имеет высокую трудоемкость);

в) атмосферные модели (весьма приближенные, имеющие много допущений о характерных распределениях капель (см., например, [3]));

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

Новый, достаточно точный и гибкий, обладающий невысокой трудоемкостью и не имеющий

существенных ограничений подход к решению сформулированной выше проблемы развит в рамках предложенного варианта эйлерово-эйлеровой модели капельной фазы [5].

Построение адекватных распределений капель. Предлагаем следующий подход:

1. Секционируем капельную среду по диапазонам диаметров D, применив кусочно-линейные распределения вида

п1 (й) = ар + Ь1

в каждой ьй секции. Пусть всего имеется Z секций.

2. Для корректности расчета будем рассматривать межъячеечный перенос лишь для физических характеристик ьго компонента (плотности рк и

концентрации Ы'к ), соответствующего ьй секции.

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

а) параметры линейного распределения (при полном заполнении секции, то есть для капель с диаметрами от di до di+1);

б) параметры равномерного распределения (при этом ai = 0) и координату начала Xi или конца yi незаполненного участка (при частичном заполнении секции). Подчеркнем, что в некоторых случаях (на краях распределений) знание положения начала/конца заполненного участка важнее, чем построение точного приближения верхней границы

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

Изложенный подход допускает расчет распределения с разрывами: в таких случаях некоторые секции просто будут пустыми, то есть П (й) = 0.

Определим функцию mode(i), определяющую тип процесса (конденсация (К), испарение (И) или стабилизация (С)) для каждой i-й секции:

K, если (8,- (di) > 0} л (8,- (d,+1) > 0}, И, если (8, (d,) < 0} л (8, (d,+1) < 0} л {( > 1}, mod e (() = < С, если (8, (d,)8, (d,+1) < 0} v v(секция пуста}v v[(8, (d, )< 0}л(8, (d,+1)< 0} л (/ = 1}],

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

8, (d) = Сп - С,ов (d)

определяет знак градиента концентрации, причем СП - концентрация пара в среде; С'пов (d) - концентрация пара на поверхности капли диаметром d в i-й секции.

Пустой считаем секцию, для которой выполняется условие

( <ЕЙ $ < ^ М i > •

К 6

где в - некоторая заданная малая величина; pk -

плотность вещества капель в i-й секции. Вторая и третья части условия характеризуют случай физически некорректного (для данной секции) сочетания значений плотности и концентрации, причем третья часть является менее строгой. Такое сочетание может возникнуть на нескольких первых итерациях расчета поступления капель в пустую ячейку.

Частично заполненной будем считать i-ю секцию, для которой выполняется условие

((/ = 1) л (mode (/) = K)} v ((, = Z) л (mode (/) = И)} v v((/ < Z) л (mode(() ф mode(, +1)} л (mode(() = И)} v v((/ > 1) л (mode(() ф mode((-1)} л (mode (() = K)} v v(forw(,)} v (back (,)},

все прочие непустые секции считаем полностью заполненными.

Логические функции forw(i) и back(i) выделяют два особых случая (промежуточные секции, содержащие «пограничные» капли):

(mod e(( -1) = K)л forw(,) = ( > 1} л ( < Z} л<! л (mod e (,) = K)л

л [секция (, +1) пуст]

[секция (, -1) пуста] л(п^ e (,) = И) л (mode (, +1) = И)

Ьаск() = ( > 1} л ( < Z} л

Прочие основные соотношения приведены в работе [5].

Рассмотрим вычисление параметров распределения в полностью заполненной секции [ж:; у], соответственно:

х, = di; У; = di+1-

Запишем интегральные выражения для плотности и концентрации секции: УI

| п! (О)dD = <,

J,

npi j n (D). D3dD = pk.

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

(vf - xf)

■Li-^ + b (y, - x, ) = Nk,

(v5-xf) + b (4) = 6pi

nP k

решением которой будут искомые значения а: и Ь:.

Выведем формулы для нахождения координаты начала ж: или конца у: незаполненного участка (при частичном заполнении секции). Для нахождения ж: выпишем те же интегральные выражения для плотности и концентрации относительно новой переменной ж^

-I +1

| п! (й^й = N,

d¡+1

6 рк | П (О)• D3dD = рк, */

из которых после преобразований получим нелинейное уравнение

(df+1 + xf)(

,+1

■x. ) = -f4pk_

" пРЛ'

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

УI = di+1.

Аналогично выводится нелинейное уравнение для нахождения ж::

(+у2)+у >=Ц| •

которое также имеет единственное решение в вещественных числах. Это уравнение решается, если в 1-й секции преобладает испарение. В данном случае

х1 = di.

Параметры распределения для частично заполненной секции следующие: е1 = 0;

N1

Ь =——.

УI - X

Апробация. Было проведено численное моделирование кинетики капельной фазы (процессов конденсации, испарения и поглощения газообразных загрязнителей) в пределах одной ячейки расчетной сетки. Использовались основные уравнения [5], записанные с применением изложенного выше подхода к поиску параметров распределения. Результаты сравнивались с данными, полученными

f

путем прямого моделирования соответствующего множества капель в соответствии с лагранжевым подходом.

Было проведено девять серий экспериментов, различающихся начальным распределением капель по диаметрам (нормальное r(d), перевернутое нормальное n(d) = max {Г(d)} - Г(d), линейное

ri(d)) и наличием загрязнителя SO2 (варианты: не рассматривается, поглощается каплями, уходит из капель). В каждой серии варьировались:

а) значение относительной влажности воздуха в полном [0,8; 1,2] и малом [0,99; 1,01] диапазонах, что позволило рассмотреть как режимы конденсации, так и испарения с поправкой Кехлера;

б) количество секций - от одной до десяти. Все распределения капель были заданы в диапазоне диаметров [10-6; 10-5], применялось равномерное разбиение на секции. Соответственно, ширина секции варьировалась в диапазоне [8,18210-7; 4,5-10"6] метра.

Эксперименты проводились при температуре 283 К, на интервале времени [0; 0,01], с шагом интегрирования по времени т = 0,0001. При расчетах с рассмотрением загрязнителя фоновая концентрация SO2 составляла 0,88-10"2 моль/м3 (200 млн-1). В опытах с поглощением начальная концентрация SO2 в каплях была равна 0,176 моль/м3, с уходом -1,76 моль/м3. При прямом моделировании рассчитывалось множество из 1000 капель.

В каждой серии вычислялась средняя накопленная относительная погрешность (в полном и малом диапазонах влажности) результатов, полученных при использовании нашего подхода, от результатов прямого моделирования капель (рис. 1). Каждая «строка» графиков соответствует одному из

распределений, каждый «столбец» - режиму моделирования загрязнителя.

Средние погрешности для Г-распределения ( рис. 1, б, в) лежат в диапазоне от 1 до 2% по плотности и от 1 до 10% по Э02. Средние погрешности для П-распределения (

рис. 1, д, е) лежат в диапазоне от 15 до 18% по плотности и от 17 до 35% по Э02. Более высокие погрешности для П-распределения объясняются тем, что основная масса капель сосредоточена в области больших диаметров, что повышает вероятность усиления погрешности при расчете межсекционных переходов и выходе капель за пределы рассматриваемого диапазона диаметров. Тем не менее полученные величины средних погрешностей для П-распределения вполне приемлемы.

Приведенные результаты отличаются наличием колебаний при различной ширине секций, вызванных погрешностями пересчета величин (в зависимости от удачности выбранного разбиения на секции) при межсекционных переходах. Было принято решение вести дальнейший анализ по трендам, в качестве которых априорно были выбраны линейные функции. Тренды строились с помощью метода наименьших квадратов (МНК), их адекватность проверялась по критерию Фишера (рис. 2).

Обратим особое внимание на точки с наибольшим отклонением от тренда (

рис. 1, б, в, з). Значимой величины (около 20%) такие отклонения достигают лишь (

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

Рис. 1. Средняя накопленная относительная погрешность в каждой серии экспериментов:- и---- погрешности по

плотности в полном и малом диапазонах влажности;........ и------ по Э02 в тех же диапазонах; а, б, в - нормальное распределение; г, д, е - перевернутое нормальное; ж, з, и - линейное; а, г, ж - отсутствие Э02; б, д, з - поглощение Э02 каплями; в, е, и - уходу Э02 из капель

0.01

0.005

0.004

0.003

0.002

0.001

0.006

0.004-

0.002" Г.

0.005 -

0 1 10 6 2 10 6 3-10 6 4-10 6 5 10 6

0 1 10 6 2-10 6 3 10 6 4 -10 6 5 ю"

0 МО 6 2 10 6 3-10 6 4 10 6 5 10 6

ж) з) и)

Рис. 2. Тренды средней накопленной относительной погрешности в каждой серии экспериментов: обозначения как на рис. 1

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

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

Заключение

Предложенный новый подход к расчету распределений капельных фаз отличается:

а) возможностью работы с произвольными распределениями капель, в том числе с разрывами;

б) низкими вычислительными трудозатратами (в 540 раз меньше по сравнению с лагранжевой моделью на 1000 капель);

в) приемлемой величиной относительной погрешности (до 2-18 % по плотности и 17-35 % по загрязнителям), которая остается ограниченной.

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

Работа была выполнена при финансовой поддержке Минобразования и науки (грант РНП.221.1.7280).

Список литературы

1. Seinfeld J.H., Pandis S.N. Atmospheric Chemistry and Physics, Wiley. - New York, 1998.

f. Aloyan A.E., Arutyunyan V.O., Louzan P.I. Numerical modeling of the gas-aerosol interaction in the atmosh-pere // Измерения, моделирование и информационные системы как средства снижения загрязнений на городском и региональном уровне: Тр. Междунар. науч. конф. «ENVIROMIS f00f». - Томск, f00f. - Т.1. - С. 158-164.

3. Zhang M., Lin W., Bretherton C.S., Hack J.J., Rasch P.J. A Modified Formulation of Fractional Stratiform Condensation Rate in the NCAR Community Atmospheric Model (CAMf) // J. Geophys. Res. - f003. - Vol. 108. -No. D1. - Р. ACL 10-1.

4. Kelleners P.H. Simulation of Inviscid Compressible Multiphase Flow with Condensation // Annual Research Briefs f003. - Center for Turbulence Research, NASA Ames Research Center, f003. - Р. 49-67.

f. Пекунов В.В., Ясинский Ф.Н. Математическая модель микроклимата в производственных помещениях с повышенной влажностью // Изв. вузов. Технология текстильной промышленности. - f006. - № f. - С. 1f8-133.

Пекунов Владимир Викторович,

ГОУВПО «Ивановский государственный энергетический университет имени В.И. Ленина», кандидат технических наук, доцент кафедры высокопроизводительных вычислительных систем, телефон (4932) 26-98-29.