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

ВОЗМОЖНОСТИ МОДЕЛИРОВАНИЯ ТЕЧЕНИЙ ЖИДКОСТИ И ГАЗА С ПОМОЩЬЮ СОВРЕМЕННЫХ ПРОГРАММНЫХ ПРОДУКТОВ

А.Н. Кочевский, канд. техн. наук

Сумский государственный университет

ВВЕДЕНИЕ

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

Эти факторы обуславливают растущий интерес к средствам математического моделирования течений жидкостей и газов, позволяющих прогнозировать характеристики течений и параметры работы устройств на стадии их проектирования, до изготовления в металле. Раздел науки, решающий проблему моделирования течений с тепломассопереносом в различных технических и природных объектах, называется вычислительной гидродинамикой (ВГД), в англоязычной литературе - CFD (Computational Fluid Dynamics).

По мере роста вычислительной мощности компьютеров, которые становились доступными по цене все большому числу пользователей, с 70-х годов XX века началось бурное развитие коммерческих программ вычислительной гидродинамики. До начала 90-х годов ХХ века эти программы устанавливались лишь на мощные компьютеры класса “рабочие станции”. В 90-х годах ХХ века дешевые персональные компьютеры догнали по мощности рабочие станции, и появились программные продукты в области ВГД, предназначенные для персональных компьютеров.

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

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

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

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

МОДЕЛИРОВАНИЕ НЕСЖИМАЕМЫХ ТЕЧЕНИЙ ЖИДКОСТИ И ГАЗА

1. Ламинарное течение

Расчет течения жидкости или газа в современных программных продуктах выполняется путем численного решения системы уравнений, описывающих наиболее общий случай движения жидкой среды. Таковыми являются уравнения Навье - Стокса (1) и неразрывности (2):

,(1)

.(2)

Здесь использована сокращенная запись уравнений, i,j = 1,…,3, предполагается суммирование по одинаковым индексам, x1, x2, x3 - оси координат, t - время. Полная запись этих уравнений в криволинейной системе координат приведена, напр., в [3]. Член fi выражает действие массовых сил.

В этой системе из 4 уравнений независимыми искомыми параметрами являются 3 компоненты скорости u1, u2, u3 и давление p. Плотность с жидкости, а также газа при скоростях примерно до 0,3 числа Маха можно полагать величиной постоянной.

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

2. Турбулентное течение

Течения в технике, как правило, являются турбулентными. Непосредственное моделирование турбулентных течений путем численного решения уравнений Навье - Стокса, записанных для мгновенных скоростей, все еще является крайне затруднительным, а, кроме того, интерес представляют, как правило, не мгновенные, а осредненные по времени значения скоростей. Таким образом, для анализа турбулентных течений вместо уравнений (1) используют уравнения Рейнольдса (3):

,(3)

где , , - осредненные по времени значения скоростей;

, , - пульсационные составляющие скоростей.

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

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

1. Течение с перемешиванием нескольких несжимаемых сред

При рассмотрении течения смеси двух или нескольких жидкостей (или газов - при малых скоростях течения) с различными плотностями плотность смеси сm уже не является постоянной величиной и зависит от концентрации. Для вычисления концентрации некоторого компонента смеси в каждой точке пространства система уравнений ММ дополняется еще одним дифференциальным уравнением в частных производных (ДУЧП) - уравнением переноса концентрации C:

,(4)

где Rr - скорость процесса растворения или химической реакции, если таковые происходят.

Величина Jj выражает интенсивность молекулярной диффузии вдоль координатного направления xj и связана с концентрацией C согласно второму закону Фика:

,

где D - коэффициент диффузии.

При смешивании n жидких компонентов, плотности которых составляют с1, с2, …, сn, для моделирования процесса перемешивания в систему уравнений ММ необходимо включить n уравнений вида (4) - по одному для каждого компонента. Допустим, к примеру, что в некоторую емкость, заполненную жидкостью 1, по нескольким трубам подводят другие жидкости (2, 3, …, n), откачивая при этом жидкость 1. В начальный момент времени концентрация жидкости 1 в емкости полагается равной единице, концентрация остальных компонентов - нулю. В трубах, подводящих жидкости к емкости, наоборот, единице равна концентрация соответствующей жидкости. В результате выполнения расчета в каждый момент времени в каждой точке пространства будет найдено значение концентрации C1, C2, …, Cn каждого из компонентов.

Плотность смеси вычисляется следующим образом:

,(5)

где C1 + C2 + … +Cn = 1.

Вязкость смеси мm вычисляется следующим образом:

мm = C1 м1 + C2 м2 + … + Cn мn,(6)

где м1, м2, …, мn - вязкость отдельных компонентов. При расчете течений с теплопереносом аналогичным образом вычисляются удельная теплопроводность и удельные теплоемкости смеси.

Уравнение переноса концентрации (4), записанное для осредненных по времени величин для моделирования турбулентных течений, имеет вид

.(7)

2. Течение с перемешиванием нескольких сжимаемых сред

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

МОДЕЛИРОВАНИЕ ТЕЧЕНИЙ С ТЕПЛОПЕРЕНОСОМ И СЖИМАЕМОСТЬЮ

1. Изменение плотности и температуры неподвижной среды при ее принудительном сжатии

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

p = с R T,(8)

где R = 287 Дж / (кг К) - универсальная газовая постоянная.

Уравнение состояния (8) справедливо для совершенного газа, т.е. газа, в котором взаимодействие между молекулами совершается путем упругих столкновений, а линейный размер молекулы мал по сравнению со средним межмолекулярным расстоянием. При очень низких температурах и/или высоких давлениях использование этого уравнения может приводить к существенным расхождениям с экспериментальными данными для реального газа. Более точные результаты в этом случае дает, например, уравнение состояния, предложенное Ван-дер-Ваальсом:

,(9)

Где

,,

pc и Tc - давление и температура, соответствующие фазовому переходу.

2. Моделирование течений с теплопереносом в случае несжимаемой или слабо сжимаемой среды

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

В течениях несжимаемой и слабо сжимаемой среды (при скоростях потока менее 0,3 числа Маха) плотность среды зависит лишь от температуры, а влияние перепадов давления на изменение плотности и температуры незначительно. С вычислительной точки зрения необходимость учета изменения температуры требует дополнить систему (1) - (2) еще одним ДУЧП - уравнением энергии:

,(10)

где член fi выражает действие массовых сил.

В это уравнение не включены внутренние источники тепла, связанные, например, с химическими реакциями. Полная энтальпия H связана с полной энергией E, внутренней энергией e и статической или удельной энтальпией h следующими соотношениями:

.

Для совершенного газа полагают, что статическая энтальпия h пропорциональна температуре T:

h = cp T,

где cp - удельная теплоемкость при постоянном давлении.

Величина Qj выражает поток энергии, переносимой путем теплопроводности вдоль координатного направления xj, и связана с температурой T согласно закону Фурье:

,(11)

где л - коэффициент теплопроводности.

Для течений во вращающейся системе координат вместо полной энтальпии H в уравнении (10) следует использовать ротальпию

,

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

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

Уравнение энергии (10) позволяет, например, моделировать процесс нагрева холодной жидкости, текущей в полости с горячими стенками. Постановка задачи требует указать температуру жидкости в начальный момент времени, в т.ч. температуру во входном и выходном сечениях. В качестве граничных условий задается также температура стенок полости (если она поддерживается фиксированной в течение рассматриваемого процесса) либо поток энергии Qj через стенку (если фиксированной является интенсивность подвода тепла).

Уравнение энергии (10), записанное для осредненных по времени величин для моделирования турбулентных течений, имеет вид

.(12)

Отметим, что запись уравнений (2), (3) и (12) в такой форме пригодна для моделирования течений несжимаемой среды. Для течений сжимаемой среды плотность тоже подвержена флуктуациям, и при выполнении осреднения по Рейнольдсу члены, содержащие плотность, представляют собой новые неизвестные величины. Чтобы избежать этого, выполняют осреднение по Фавру, то есть мгновенные значения переменных осредняют не только по времени, но и по массе. Эта процедура и получаемые уравнения описаны, напр., в [4].

3. Моделирование течений, вызванных естественной конвекцией

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

4. Моделирование течений в случае сжимаемой среды

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

.(13)

Для течений с низкими скоростями это уравнение сводится к виду (10). Как показано в работе [4], анализ порядков слагаемых позволяет в этом случае отбросить члены с вязкостью и давлением. Но при моделировании течений с высокими скоростями (свыше 0,3 числа Маха) эти члены играют существенную роль.

При рассмотрении течения сжимаемой среды (т.е. среды, в которой плотность существенно зависит от перепадов давления), вводя в систему уравнение энергии в виде (13), мы вводим два новых независимых искомых параметра - температуру и плотность. Для замыкания системы используется уравнение состояния (8), устанавливающее связь между этими параметрами.

При рассмотрении течения сжимаемой среды давление входит в систему уравнений не только в виде первых производных, но и явным образом - в уравнение состояния (8). Как следствие, постановка задачи требует указать распределение давления (либо плотности) в начальный момент времени.

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

Можно показать, что для течений сжимаемой среды численное поведение системы уравнений ММ с включенным в нее уравнением (13) существенно зависит от числа Маха M, определяемого как отношение средней скорости потока в некотором сечении к скорости звука a в данной среде. Скорость звука определяется следующим выражением:

.

Скорость звука - это скорость, с которой в сжимаемой среде распространяются акустические волны (волны сжатия малой амплитуды). Для идеального газа это выражение сводится к виду

,

где г = cp / cv

- отношение удельных теплоемкостей при постоянном давлении и постоянном объеме.

Если скорости в потоке превышают скорость звука (M > 1), то система уравнений ММ проявляет гиперболическое поведение, что не требует задания граничных условий на выходе из расчетной области. При этом в получаемом решении наблюдаются так называемые ударные волны - поверхности, на которых происходит скачкообразное изменение параметров расчета. Это явление полностью соответствует экспериментальным наблюдениям. Условия на ударной волне, при которых в направлении нормали выполняются законы сохранения массы, импульса и энергии, называются условиями Ренкина - Гюгонио и сформулированы, напр., в работах [5, 6].

5. Задание зависимости свойств среды от температуры

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

,

где T - температура в некоторой точке пространства; с(T) - соответствующее значение плотности; Tn, Tn+1 и сn, сn+1 - ближайшие значения температуры, для которых пользователь указал плотность, и соответствующие указанные пользователем значения плотности.

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

6. Моделирование теплопереноса в твердой среде

Сопряженный теплообмен. Приведенное выше уравнение энергии (10) может использоваться для моделирования теплопереноса не только в жидкой, но и в твердой среде. Поскольку течение в этом случае отсутствует (т.е. скорость u = 0), а полная энтальпия сводится к виду H = c T, где c - удельная теплоемкость материала твердого тела, уравнение энергии (10) принимает вид

,(14)

где ST - источник энергии внутри твердого тела (если таковой имеется).

Уравнение (14) - единственное уравнение, требуемое для моделирования теплопереноса в твердой среде. Единственным неизвестным в этом уравнении является температура (т.е. распределение температуры внутри твердого тела). В качестве граничных условий, требуемых для решения этого уравнения, необходимо распределение температуры на поверхности расчетного объема, соответствующего твердому телу. Вместо температуры на части поверхности можно задать тепловой поток Qj, определяемый формулой (11).

В инженерной практике нередко встречается необходимость решения задач на сопряженный теплообмен, подразумевающих совместный расчет течения жидкости или газа при наличии теплопереноса и расчет теплопереноса в охватывающей их толстостенной твердотельной оболочке. При этом течение жидкости или газа моделируется, напр., уравнениями (1) - (2) и (10), а теплоперенос в твердом теле - уравнением (14). На поверхности раздела между жидкой и твердой фазами ставится условие равенства температуры и тепловых потоков со стороны жидкой и твердой фаз. Один из примеров задачи на сопряженный теплообмен - расчет поля температур на лопастях газовой турбины, принудительно охлаждаемых изнутри холодной жидкостью. Такого рода задача решалась, например, в работе [7].

МОДЕЛИРОВАНИЕ МНОГОФАЗНЫХ ТЕЧЕНИЙ

1. Моделирование многофазных течений

В природе и технике встречаются многофазные течения весьма различной природы. Их обзор приведен, напр., в [8]. Будем условно различать следующие типы многофазных течений. Первый случай - рассматриваемый объем полностью заполнен веществом одной фазы (напр., жидкости), а вещество другой фазы встречается в этом объеме в виде дискретных частичек (твердой фазы) или пузырьков (газообразной фазы), причем объемная доля вещества другой фазы невелика (до 10% общего объема). Второй случай - рассматриваемый объем частично заполнен жидкостью, а частично - газом, которые не смешиваются между собой и отделены друг от друга свободной поверхностью. Третий случай - самый сложный - вещества различных фаз могут смешиваться между собой (растворяться / выделяться из раствора), причем объемная доля вещества другой фазы велика (свыше 10% общего объема). Для моделирования этих типов многофазных течений используются различные подходы, которые мы сейчас и рассмотрим.

2. Модель дисперсных частиц

Этот подход используется для моделирования двухфазных течений, в которых вещество одной из фаз представлено в виде дисперсных частиц, причем объемная доля, занимаемая этими частицами, невелика (до 10% общего объема). Примерами таких течений являются разбрызганные в потоке воздуха капли воды, воздушные пузырьки в потоке жидкости, а также твердые частицы в потоке воздуха или воды. Вещество, образующее основную фазу, полагается сплошной средой, и его течение моделируется уравнениями Навье - Стокса (1) (или Рейнольдса (3)) и неразрывности (2). Вещество, присутствующее в потоке в виде дискретных частиц, не образует сплошную среду, отдельные частицы взаимодействуют с потоком основной фазы и друг с другом дискретно. Для моделирования движения частиц рассеянной фазы используется подход Лагранжа, т.е. отслеживается движение отдельно взятых частиц рассеянной фазы под действием сил со стороны потока основной фазы.

Течения, для моделирования которых требуется эта модель, имеют место, напр., в распылительных сушилках, циклонных сепараторах, тканевых фильтрах, при эрозии клапанов,при впрыске топлива в двигатель. Задачи такого рода нередко встречаются в химической промышленности, см., напр. [9].

Для простоты полагают, что частицы рассеянной фазы имеют шарообразную форму. Силы, действующие на эту частицу, обусловлены разностью скорости частицы и скорости потока основной фазы, а также вытеснением этой частицей среды основной фазы. Уравнение движения такой частицы было выведено в работе [10] и имеет вид:

.(15)

Здесь mp - масса частицы, d - диаметр частицы, v - скорость, м - динамическая вязкость вещества основной фазы, Ccor - его коэффициент вязкого сопротивления; - угловая скорость вращения, - радиус вектор (при рассмотрении движения в относительной системе отсчета). Индекс p (particle) относится к частице, индекс f (fluid) - к веществу основной фазы.

Левая часть уравнения (15) представляет собой сумму всех сил, действующих на частицу, выраженную через массу и ускорение этой частицы. Первый член в правой части выражает торможение частицы в результате вязкого трения о поток основной фазы согласно закону Стокса. Второй член - сила, приложенная к частице, вследствие перепада давления в основной фазе, окружающей частицу, вызванного ускорением потока основной фазы. Третий член - сила, требуемая для ускорения веса основной фазы в объеме, вытесненном частицей. Эти два члена нужно учитывать, когда плотность основной фазы превосходит плотность частиц, например, при рассмотрении пузырьков воздуха в потоке жидкости. Четвертый член (Fe) - внешняя сила, непосредственно действующая на частицу, например, сила тяжести или сила электрического поля. Последние два члена - центробежная сила и сила Кориолиса, которые имеют место лишь при рассмотрении движения в относительной системе отсчета. Кроме того, иногда в правой части (15) бывает необходимо учесть некоторые дополнительные силы (напр., при наличии в потоке значительного перепада температур).

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

Коэффициент вязкого сопротивления Ccor при умеренных числах Рейнольдса 0,01 < Rep < 260 можно вычислить, напр., по формуле

где Rep = сf | vf - vp | d / м, б = log Rep.

Отметим также возможность моделировать тепло- и массообмен между дисперсными частицами и основным потоком, имеющуюся в современных программных продуктах, например, испарение капли жидкости при достаточно низком давлении газа в окружающем потоке или достаточно высокой температуре. Алгоритмы, реализованные во Fluent, CFX, STAR-CD, позволяют моделировать и воздействие на поток вещества основной фазы со стороны движущихся в нем дискретных частиц. В первом приближении плотность и вязкость вещества основной фазы, а также некоторые другие величины умножаются на (1 - бp), где бp - удельный объем, занятый дискретными частицами. Далее на каждом шаге по времени вычисляются изменения массы, импульса и энергии дискретных частиц, и эти изменения добавляются соответственно в уравнения сохранения массы (2), импульса (1) и энергии (10) для потока основной фазы в виде источниковых членов. Таким образом, расчет течения основной фазы и расчет движения дисперсных частиц выполняется совместно.

Если поток вещества основной фазы является турбулентным, то траектория движения дисперсных частиц не является детерминированной, поскольку зависит от интенсивности и направления турбулентных пульсаций. Способ моделирования воздействия турбулентных пульсаций основного потока на движение дисперсных частиц был предложен, в частности, в работе [11].

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

3. Модель течений со свободной поверхностью

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

Согласно данному подходу, математическая модель для аппроксимации свободной поверхности дополняется уравнением переноса функции заполнения F, выражающей “концентрацию жидкости в газе” (при рассмотрении течения жидкости с газом). Отсюда следует название модели течения - модель VOF (Volume Of Fluid, объем, занятый жидкостью).

.(16)

В области, занятой жидкостью, F = 1, в области, занятой газом, F = 0. Лишь в ячейках, через которые проходит свободная поверхность, 0 < F < 1. В качестве начального условия задается исходное положение свободной поверхности.

Алгоритм численного расчета при использовании такой модели описан в [12, 13].

Примеры задач, для которых применим данный подход: заполнение топливного танкера, плескание или кипение в емкости со свободной поверхностью, моделирование обтекания морских судов. В работе [14], в частности, программный продукт FlowVision с моделью VOF использовалась для моделирования течения воды у колеса автомобиля, въезжающего в лужу в режиме аквапланирования.

4. Многофазная модель смешения (Multiphase Mixture Model)

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

Для моделирования течения двух или нескольких фаз в этой модели используется одно уравнение неразрывности, один набор уравнений движения и одно уравнение энергии, записанные относительно осредненных по массе значений скорости и плотности смеси. Так, уравнение неразрывности в этой модели имеет вид:

,(17)

где сm - плотность смеси (формула (5)), um - осредненная по массе скорость, umj - проекция скорости на ось xj, по индексу j предполагается суммирование, а член m, по умолчанию равный нулю, представляет массоперенос вследствие кавитации и/или других физических эффектов.

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

Уравнение движения в проекции на ось xi в этой модели имеет вид:

(18)

где мm - осредненная по массе вязкость (формула (6)), uk - скорость вещества k-й вторичной фазы, umj - проекция этой скорости на ось xi, (uk - um) - скорость проскальзывания вещества k-й вторичной фазы относительно осредненной по массе скорости um. Уравнение (18) отличается от (1) наличием последнего члена, моделирующего взаимное проскальзывание фаз.

Аналогичным образом записывается уравнение энергии.

Уравнение неразрывности для отдельно взятой k-й вторичной фазы можно представить в виде

,(19)

где сk - плотность вещества k-й фазы. Из этого уравнения можно определить объемную долю Fk, занимаемую веществом k-й фазы в некоторой ячейке пространства.

5. Многофазная модель Эйлера

Данная модель является наиболее общей и наиболее сложной среди моделей многофазного течения. Вещество в каждой из фаз полагается сплошной средой, и движение вещества каждой из фаз моделируется собственной системой уравнений Навье - Стокса (Рейнольдса), неразрывности и энергии.

Согласно этой модели, уравнения движения, записанные для каждой фазы, решаются совместно. Алгоритм расчета таких течений был предложен, в частности, в работе [15], и реализован в CFX, Fluent и STAR-CD. Данная модель является наиболее требовательной к вычислительным ресурсам компьютера - и к размеру оперативной памяти, и к быстродействию процессора.

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

6. Дополнительные рекомендации по выбору модели многофазного течения

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

- Отношение в массы вещества дисперсной фазы (d) к массе вещества несущей фазы (c):

,(20)

где Fd и Fc - объемные доли, г - отношение плотности дисперсной и несущей фазы,

г = сd / сc;

это отношение может составлять свыше 1000 для твердых частиц в потоке газа, около 1 для твердых частиц в потоке жидкости, и менее 0.001 для частиц газа в потоке жидкости.

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

- Число Стокса St:

,(21)

где td - время, характеризующее движение частиц,

td = (сd dd2) / (18 мc),

dd - диаметр частицы, мc - вязкость вещества несущей фазы,

tc = Lc / Uc

- время, характеризующее течение несущей фазы, Lc - характерная длина, Uc - характерная скорость.

При St << 1,0, частицы дисперсной фазы почти не отклоняются от линий тока несущей фазы, и можно использовать любую модель течения (как правило, модель смешения - как наименее ресурсоемкую). При St > 1,0, траектория частиц дисперсной фазы совершенно не совпадает с линиями тока несущей фазы, и модель смешения в этом случае непригодна: нужно использовать либо модель дисперсных частиц, либо многофазную модель Эйлера.

7. Моделирование течений с кавитацией

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

В современных программных продуктах встречаются различные подходы к моделированию кавитации. Один из наиболее перспективных подходов предложен, напр., в [18]. Для моделирования кавитации используется, как правило, модель течения со свободной поверхностью или многофазная модель смешения. Многофазную модель Эйлера использовать нецелесообразно вследствие низких значений параметров в и St (формулы (20) и (21)) в таких течениях.

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

ПРОЧИЕ ВОЗМОЖНОСТИ МОДЕЛИРОВАНИЯ ТЕЧЕНИЙ

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

1. Течения неньютоновской жидкости

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

,

м - динамический коэффициент вязкости, не зависящий от D.

В неньютоновских жидкостях, к которым относятся, например, парафин, воск, мед, смола, коэффициент вязкости м проявляет зависимость от D.

Зависимость м (D) должна быть задана в качестве исходных данных для расчета течения. К настоящему времени был предложен целый ряд приближенных законов, описывающих эту зависимость для различных неньютоновских жидкостей. Более подробные сведения об этом имеются, например, в [19].

2. Течения сквозь пористую среду

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

Для анализа таких течений, как правило, множество мелких твердых элементов, препятствующих течению, условно считают сплошным равномерно распределенным сопротивлением. Это сопротивление моделируют, вводя в правую часть уравнений движения (1) дополнительный источниковый член Si следующего вида:

,(22)

где Dij и Cij - матрицы размером 3 х 3, преопределенные на стадии постановки задачи.

В случае однородной (изотропной) пористой среды диагональные элементы матрицы Dij равны 1 / б, диагональные элементы матрицы Cij равны C2, а все остальные элементы этих матриц равны нулю. Здесь б - проницаемость, а C2 - коэффициент внутреннего сопротивления. Тогда Si принимает вид

(23)

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

.

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

.

Теплопроводность пористой среды лeff можно рассчитать как среднюю по объему теплопроводность вещества жидкой лf и твердой лs фаз в пористой среде:

лeff = г лf + (1 - г) лs,

где г - пористость среды, т.е. доля объема, занятая жидкой средой.

3. Течения с радиационным теплопереносом

Как известно, существует 3 способа передачи тепла на расстояние - конвективный, диффузионный и радиационный перенос. Уравнение энергии в виде (10) или (13) позволяет моделировать первые два вида теплопереноса. Радиационным теплопереносом во многих задачах, встречающихся в технике, можно пренебречь. Радиационный теплоперенос может играть доминирующую роль при наличии очень больших перепадов температур. Тепловой поток Qrad, возникающий вследствие радиационного теплопереноса, напр., от горячей стенки к холодной, выражается формулой

Qrad = у (Tmax4 - Tmin4),

где Tmax и Tmin - температура горячей и холодной стенки;

у = 5,67 * 108 Вт / (м2 К4) - константа Стефана-Больцмана.

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

Жидкость или газ, протекающий в области с радиационным теплопереносом, нагревается от горячих стенок путем конвекции и диффузии. Кроме того, если оптическая прозрачность жидкой среды отличается от нуля, энергия радиационного теплопереноса передается жидкой среде и непосредственно. Для моделирования этого процесса в уравнение энергии (10) или (13) добавляется член, вид которого определяется принятой моделью радиационного теплопереноса. Обзор современных моделей радиационного теплопереноса представлен, напр., в книге [20].

Наглядный пример задачи на радиационный теплоперенос входит в комплект демо-версии FlowVision.

4. Течения с химическими реакциями и горением

В математическую модель, помимо прочего, входят и уравнения переноса концентрации (4), позволяющие определить концентрацию отдельных компонентов смеси и уравнения энергии (10) или (13). Кроме того, добавляются уравнения, описывающие собственно процесс химической реакции или горения. Для различных типовых реакций (в особенности для горения смеси метана с воздухом) соответствующие модели разработаны и в той или иной мере реализованы в современных программных продуктах. Кроме того, интерфейс этих программных продуктов содержит возможности пользовательского программирования, позволяя достаточно опытному пользователю создать подходящую модель для некоторой специфической реакции.

Отметим, что в течениях с горением, как правило, наблюдаются очень сильные турбулентные пульсации, существенно влияющие на картину течения. Для учета этого явления, напр., во FlowVision, k - е модель турбулентности дополняется специальным уравнением для моделирования пульсаций.

5. Моделирование акустики

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

Отметим, что для задач моделирования акустического шума при расчете течения в качестве модели турбулентности важно использовать модель LES или DNS [2], которая обеспечивает наилучшее разрешение происходяего в потоке пульсаций давления, поскольку эти пульсации и являются источником акустического шума.

ВЫВОД

Проведен обзор моделей течения жидкостей и газов, реализованных в ведущих программных продуктах ВГД, предназначенных для моделирования течений многокомпонентных и многофазных сред, сжимаемых течений, течений с теплопереносом, кавитацией и прочими явлениями. Обзор литературы свидетельствует, что перечисленные программные продукты (CFX, Fluent, STAR-CD и др.) позволяют адекватно моделировать сложные физические эффекты различной природы, в том числе для задач, в которых проведение физического моделирования крайне затруднительно.

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

Aksenov A.A., Dyadkin A.A., Gudzovsky A.V. Numerical Simulation of Car Tire Aquaplaning // Computational Fluid Dynamics '96, J.-A. Desideri, C.Hirsch, P.Le Tallec, M.Pandolfi, J.Periaux edts. - John Wiley&Sons, 1996. - P. 815-820.

Cebeci T., Smith A. M.O. Analysis of Turbulent Boundary Layers. - Academic Press, New York, 1984.

Crowe C., Sommerfield M., Yutaka Tsuji. Multiphase Flows with Droplets and Particles. - CRC Press, 1998.

Dukowicz J.K. A Particle-Fluid Numerical Model for Liquid Sprays // Journal of Computational Physics. - Vol. 35, 1980. - P. 229-253.

Harlow F.H., Welch J. E. Numerical Calculation of Time-Dependent Viscous Incompressible Flows of Fluid With Free Surface // Phys. Fluids. - No. 8, 1965. - P. 2182-2187.

Hinze, J.O. Turbulence. - McGraw-Hill, New York, 1975.

Hirt C.W., Nicholls B. D. Volume of Fluid (VOF) method for dynamical free boundaries // J. Comput. Phys. - No. 39, 1981. - P. 201-225.

Siegel R., Howell J.R. Thermal Radiation Heat Transfer. - Hemisphere Publishing Corporation, Washington D.C., 1992.

Singhal A. K., Li H. Y., Athavale M. M., Jiang Y. Mathematical Basis and Validation of the Full Cavitation Model. // ASME FEDSM'01 - New Orleans, Louisiana, 2001.

Tanner R.I. Engineering Rheology. - Clarendon Press, Oxford, rev. edition, 1988.

Vasquez S.A., Ivanov V.A.A Phase Coupled Method for Solving Multiphase Problems on Unstructured Meshes // In Proceedings of ASME FEDSM'00: ASME 2000 Fluids Engineering Division Summer Meeting. - Boston, June 2000.

Абрамович Г. Н. Прикладная газовая динамика. - 3-е изд. // Гл. ред. физ-мат. лит-ры изд-ва “Наука”: М., 1969. - 824 с.

Дейч М.Е., Филиппов Г. А. Газодинамика двухфазных сред. - М.: Энергия, 1968. - 424 с.

Кнэпп Р., Дейли Дж., Хэммит Ф. Кавитация. - М.: Мир, 1974. - 688 с.

Кочевский А.Н. Расчет внутренних течений жидкости в каналах с помощью программного продукта FlowVision // Вестник СумГУ. - Сумы, 2004. - № 2 (61). - С. 25-36.

Кочевский А.Н., Неня В. Г. Современный подход к моделированию и расчету течений жидкости в лопастных гидромашинах // Вестник СумГУ. - Сумы, 2003. - № 13 (59). - С. 195-210.

Кочевський О.М. Оптимізація геометричних параметрів відвідних пристроїв насосів високої швидкохідності з лопатевою системою типу НР / Дис... канд. техн. наук. - Суми: СумДУ, 2001. - 195 с.

Парьохін О.В., Склабінський В. І., Кочевський О. М. Розрахунок гідродинаміки газового потоку у внутрішньому конусі вихрового гранулятора // Вісник Сумського національного аграрного університету. - Суми, 2002. - Вип. 8 - С. 64-69.

Приходько О.А., Сьомін Д.О. Технічна аеромеханіка: Навчальний посібник. - Луганськ: Вид-во Східноукраїнського нац. ун-ту ім. В. Даля, 2002. - 170 с.

Слитенко А.Ф., Рассохин Е.В. Повышение экономичности газотурбинных установок за счет уменьшения расхода охлаждающего воздуха // Совершенствование турбоустановок методами математического и физического моделирования: Сб. научн. трудов. - Харьков: Ин-т проблем машиностроения НАН Украины, 2003. - Т. 1. - С. 143-148.