ЖУРНАЛ ВЫЧИСЛИТЕЛЬНОЙ МАТЕМАТИКИ И МАТЕМАТИЧЕСКОЙ ФИЗИКИ, 2015, том 55, № 10, с. 1756-1761
УДК 519.634
ИНТЕРПОЛЯЦИЯ РИ-ЧОУ В СИЛЬНЫХ ЦЕНТРОБЕЖНЫХ ПОЛЯХ1)
© 2015 г. С. В. Боговалов, И. В. Тронин
(115409 Москва, Каширское шоссе, 31, МИФИ) e-mail: ss433@mail.ru; ivtronin@mephi.ru Поступила в редакцию 03.09.2014 г.
Переработанный вариант 03.03.2015 г.
Получены интерполяционные формулы Ри-Чоу из уравнений Навье—Стокса и уравнения непрерывности. Произведено обобщение этих формул на случай динамики газа в сильных центробежных полях газовых центрифуг, достигающих величины 106 g. Библ. 8.
Ключевые слова: уравнение Навье—Стокса, интерполяционная формула Ри-Чоу, сильные центробежные поля газовых центрифуг, разностный метод численного решения.
DOI: 10.7868/S0044466915100087
1. ВВЕДЕНИЕ
Газовые центрифуги (далее ГЦ) используются в промышленных масштабах для разделения изотопов урана. Центробежное ускорение в них достигает величины порядка 106 g. Одним из наиболее перспективных средств для получения информации о динамике газа в ГЦ является численное моделирование газовых течений в этих устройствах. Поэтому во всем мире сейчас разрабатываются методы расчета газовой динамики в условиях сильных центробежных полей. Они возникают при расчетах газового течения во вращающейся системе координат. Необходимость в проведении расчетов в этой системе координат связана с задачами, которые возникают при конструировании газовых центрифуг.
Одним из секретов центробежного разделения изотопов является использование медленного осевого циркуляционного течения в рабочей камере ГЦ для кардинального увеличения эффективности разделения. Основной задачей всех газодинамических расчетов при проектировании ГЦ является определение параметров, при которых это циркуляционное течение вдоль оси вращения обеспечит максимальную разделительную способность. Скорость газа в этом течении достигает величин порядка нескольких мм/с вблизи стенки ротора. Скорость ротора составляет около 600—700 м/с. Выделять циркуляционное течение со скоростями в несколько мм/с на фоне вращения газа в 600—700 м/с проще всего во вращающейся системе отсчета, где твердотельная скорость газа равна нулю. По этой причине начиная с пионерских работ [1], [2] и заканчивая совсем недавними публикациями [3], [4] все численные расчеты газовых течений в ГЦ выполнялись во вращающейся системе отсчета.
До недавнего времени расчеты циркуляционного течения выполнялись в осесимметричном приближении. Современное состояние дел диктует необходимость в расчетах в полностью 3D геометрии. Такие расчеты необходимо проводить на сетках с ячейками произвольной формы, чтобы можно было точно воспроизводить геометрию устройства. Для дискретизации уравнений гидродинамики на произвольных сетках обычно используется метод контрольных объемов. Он требует расчета потоков массы, импульса и энергии через поверхности контрольных объемов. Ключевым элементом любой численной схемы является метод интерполяции газодинамических переменных от своих узловых значений к значениям на поверхности контрольных объемов. В каждом методе для интерполяции используются градиенты давления, которые достигают огромных величин. Достаточно сказать, что в модельной ГЦ Игуасу давление меняется на 5, 6 порядков при изменении радиуса, всего на 1 см. Прямое использование известных методов интерполяции приводит к возникновению нефизичных массовых потоков, которые делают невозможным выполнение расчетов с необходимой точностью. Поэтому одной из основных задач при
1) Работа выполнена при финансовой поддержке МОНРФ (№ 3.726.2014/К).
1756
разработке специализированных численных схем для моделирования течения газа в ГЦ является разработка интерполяционной схемы для расчета массовых потоков на границах контрольных объемов.
Одним из популярных на сегодня методов расчета массовых потоков на поверхности контрольных объемов, использующим узловые значения переменных, был метод, предложенный в начало 80-ых годов прошлого века в [5]. В настоящее время этот метод интерполяции широко используется в большинстве численных схем для расчета гидро- и газодинамических течений. Вместе с тем, этот метод интерполяции имеет ряд недостатков, с которыми до сих пор пытаются бороться. Немаловажным методическим недостатком метода является то, что он обычно формулируется на уровне системы уравнений, получаемых после дискретизации уравнений Навье— Стокса, уравнения непрерывности и уравнения для энергии. При решении стандартных задач гидро- и газодинамики это обычно не создает проблем, поскольку уже существуют апробированные схемы, про которые известно, что они дают приемлемые результаты. Но когда возникают нестандартные задачи, а к ним можно отнести динамику газа, в газовых центрифугах, методические проблемы трансформируются в принципиальные технические вопросы, которые трудно решить, не прояснив основы самого метода.
В данной статье по сути впервые интерполяционные соотношения, предложенные Ри-Чоу, выведены из первых принципов и на их основе получены новые интерполяционные соотношения, которые позволяют проводить моделирование динамики газа в условиях сверхсильных центробежных полей, достигающих 106 §.
В данной работе мы приводим только вывод основных интерполяционных соотношений без их верификации. Связано это с тем, что, в силу ряда обстоятельств, верификация предложенной интерполяции была опубликована ранее в двух работах. В [6] приведены результаты верификации для задач стационарного течения газа в газовой центрифуге, а в [7] приводятся результаты верификации для нестационарных задач.
2. ИНТЕРПОЛЯЦИОННАЯ СХЕМА РИ-ЧОУ В СИЛЬНЫХ ЦЕНТРОБЕЖНЫХ ПОЛЯХ
Движение газа в сильном центробежном поле определяется уравнением Навье—Стокса, записанном для вращающейся системы отсчета в виде
+ ^ + 2 ре ^ = -дХ + р + рП2 П, (1)
от дхк 0X1 дхк
уравнением непрерывности
0Р + дРУ = о
д? дхк
и уравнением, определяющим изменение плотности энергии газа, которое мы здесь приводить не будем, поскольку оно нам в дальнейшем не понадобится (см. [8]). Здесь р — плотность газа, у- — компоненты скорости, р — давление газа, П — угловая скорость вращения ротора ГЦ, г — компоненты радиуса-вектора, спроектированные на плоскость, перпендикулярную оси вращения, т к — тензор вязких напряжений.
Для удобства введем некоторые обозначения. Член вида 2ре ¡цкПцУк представим как уккр vk, где ук = е цкПц. Применив операцию дивергенции к уравнению (1), получим уравнение вида
афу + фуу + ^ =-°£ + ^Тк + Арп2 (2)
-л, -л л л 11кг к ") ^ ^ ^ г I V/
д? дхк дхкдхк дхк дх2 дх1дхк дхк
Дискретизируем производную по времени в первом члене этого уравнения в соответствии с формулой
д дРУ к = 1 |"дРу- - др^^Л (3)
д? дхк А ?V дхк дгк ;
n n- 1 / 1ТТ
где р vt и р vi — плотности потока, массы на временном слое n и n — 1. Используя уравнение непрерывности и уравнение (3), уравнение (2) приводится к виду
d-+ '-it-11 = a/-2- - Pt - дхр^+дх^(рv)+ддрР).
dt дхк Vдх2 dxjdxk dXj dXj дхдхк J
Это уравнение имеет вид уравнения непрерывности. Плотность потока массы в нем на временном слое n равна
f \П , \n -1 ./ др" 2 д(%1к)П , ,П д(р Vj vk)n\
(pVk) = (pVk) + Atl - д-- + р О /* + -д^ - Ykj(p Vj) - VF ' k- J. (4)
v дxk дХ( дх1 J
Это выражение можно интерпретировать как плотность потока массы на поверхности контрольного объема на временном слое n. На самом деле это выражение задает уравнение на плотность потока массы (р v^)n на временном слое n, в которое входит плотность потока массы на предыдущем временном слое и значения переменных на узловых точках. В дальнейшем будем рассматривать это выражение как плотность потока, массы в точках интегрирования (IP точках) на поверхности контрольных объемов. Аналогичное по форме уравнение можно получить для узловых значений плотности потока массы непосредственно из уравнения Навье—Стокса.
3. УПРОЩЕНИЯ РИ-ЧОУ
В оригинальном методе Ри-Чоу используется упрощенное представление тензора вязких напряжений в виде
дГ«к = ц ^. дх дх2
Пространственная дискретизация вязкой силы в этом случае будет иметь вид
^ = ^XWn bVk(nb) - Vkc? W
nb
nb
где — веса, с которыми входят значения скорости в соседних пЬ точках к центральной точке с. Точкой с может быть или узел или ¡р точка на поверхности контрольного объема.
Основная идея метода интерполяции Ри-Чоу состоит в том, чтобы использовать средние по узловым точкам значения потока массы для оценки конвективного члена и члена, содержащего
(„ь) на ¡р точке. Пусть уравнение (4) записано для узлов. Значение плотности потока массы, определенное по узловым значениям, равно
<(р^)„)р = <(р^)„-% - д(<дХ) + (рп2гк)р +
k jp
(
+
<W Xwnbvk(nb)
v
J
)jp - <рсТkc>jp < р Xwnb)jl
р n
-1 <Yk,fiv),p - <^
nb nb
Здесь также использовано приближенное равенство
(5)
((PcVtcf-X wnb>jp = ^(рcVkc)n> ip (-1- X wnb>jp .
nb
nb
Усреднение <...) в этих уравнениях производится по некоторой процедуре, конкретный вид которой сейчас не существенен.
c
c
После вычитания из уравнения (4) уравнения (5) получаем
(Р Ук)I - <(рУк)% = (р Ук1 - <(р Ук)--1).р + А?(-((|Р);. - < др).Р) + (РП2Гк)р - <рП2Гк).р -
- И
\
(Рv.YlpZ — - <(Pvk)">ip<-LZWnb))fp| -Yk((Pv)ПР- <(PV)%)
nb
nb
В соответствии с методикой Ри-Чоу мы полагали при вычитании, что
<д Р V Vk> =
< ас, >ip
дХ:
и
(6)
Ц Z wnbv:(nb)
>,p = HZ WnbV( nb)
nb nb
Для того чтобы упростить выражение (6), используем также формулу
Z - = ^ Z Wnb)p •
Pip
nb nb
С учетом этого, выражение (6) можно переписать в виде
(Р Vk)Пр - <(Р Vk)n>,p = (Р Vk)nt- 1 - <(Р vk)n-:>,p + л(-(Щ),р - <дХ->:p + (Р^2rk\
<pQ2rk>,p - иZ — ((Р Vk)Пр - <(Р Vk)%) - Yk,((Р^)Пр - <(Р v)n>ip) Р p
(7)
nb
4. МОДИФИЦИРОВАННАЯ ИНТЕРПОЛЯЦИЯ РИ-ЧОУ
Следующий шаг состоит в том, чтобы выразить плотность потока, массы на 1р-й точке в явном виде. Пусть
= (РУк)1- <(РУк)n)iP,
Ft = (- дХ- + Р^ Л •
( dx, ) ,p
Тогда уравнение (7) можно переписать в виде
(
t = t 1 + At (F- <F>n) - иZ— И -YikYk
nb
Перепишем это уравнение в другом виде:
(( ,
|| 1 + иУ ^
At HZ п.
Определим новую матрицу Aik как
Р p)
8 ik + Yik
Гк = ^Ат" + (F - < F:> n)
Для дальнейшего прочтения статьи необходимо приобрести полный текст. Статьи высылаются в формате PDF на указанную при оплате почту. Время доставки составляет менее 10 минут. Стоимость одной статьи — 150 рублей.