CC BY
16
структуры и моделирование 2016. №3(39). С. 93-109
УДК 519.6
МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ СИЛЬНОГО СЖАТИЯ ОДНОМЕРНЫХ СЛОЁВ ГАЗА В КОНФИГУРАЦИИ Р. МИЗЕСА
Н.С. Новаковский
аспирант, e-mail: [email protected]
Уральский государственный университет путей сообщения (УрГУПС)
Аннотация. В статье изложены результаты численного исследования задачи сильного сжатия одномерных слоёв газа в конфигурации Р. Мизеса. Приведена математическая постановка задачи в виде двух характеристических задач Коши (ХЗК1 и ХЗК2) и одной задачи Коши для восстановления траектории движения сжимающего поршня в обратом направлении изменения времени. Приводятся алгоритмы численного решения этих двух ХЗК (полностью) при помощи стандартного метода характеристик с предварительным раскрытием особенности. Описан алгоритм восстановления траектории сжимающего поршня. Описан конечно-разностный метод «РОМБ» для численного решения системы уравнений газовой динамики (СУГД) в прямом направлении изменения времени. Представлены результаты численных исследований, в том числе сравнение с точным решением.
Ключевые слова: сильное сжатие газа, характеристическая задача Коши, центрированная волна (ЦВ) сжатия, метод характеристик, конечно-разностный метод «РОМБ».
Введение
Математическое описание процесса безударного изоэнтропического сжатия идеального газа до любого наперёд заданного значения плотности, в том числе до бесконечной плотности представляет интерес в связи с проблемой лазерного термоядерного синтеза [1,2] (более подробную библиографию см. в [3]).
В [3,7] разработана математическая теория безударного сильного сжатия идеального газа. В частности, для случая сжатия цилиндрически V =1 и сферически V = 2 симметричных слоев политропного газа с показателем 7 > 1 доказано, что непрерывная состыковка двух течений даёт решение задачи о безударном сильном сжатии до любой наперёд заданной плотности ненулевой массы газа.
В работе [10] предложен алгоритм расчёта безударного сильного сжатия в обратном направлении изменения времени одномерных слоёв первоначально однородного и покоящегося газа с р0 = 1 до любой наперёд заданной конечной постоянной плотности р* > 1 с последующим восстановлением закона движения поршня, сжимающего слой газа. В работе [11] представлены результаты расчётов вышеупомянутым алгоритмом ряда одномерных задач с различной симметрией. Максимально достигнутое значение плотности в расчётах ЦВ Ри-мана из работы [12] в 3 ■ 104 раз превышает исходную плотность. При этом среднее значение плотности в ЦВ равно 8 ■ 103.
Результаты расчётов сильного сжатия двумерных областей приводятся в работе [13]. В работах [14,15] подробно рассматриваются расчёты автомодельных двумерных и трёхмерных задач сильного сжатия. В работах [12,14,15] расчёты проведены в прямом направлении изменения времени, а движение сжимающего поршня происходит в соответствии с точными решениями, явно полученными в более ранних работах [4,16,17].
Основная цель выполнения данной работы обусловлена требованиями со стороны физиков-прикладников, которые чётко сформулировал академик Е.Н. Аврорин: «Какие конкретные рекомендации для физического эксперимента дают результаты вашего математического моделирования сильного сжатия?» Эта цель формулируется так: с использованием ранее опубликованного алгоритма построить закон движения сжимающего поршня и расчётами в прямом направлении изменения времени конечно-разностным методом «Ромб» [18] определить возможность сжатия газа с помощью найденного закона движения сжимающего поршня до нужных степеней плотности.
В работе кратко описан алгоритм расчёта в обратном направлении изменения времени безударного сильного сжатия изнутри (то есть при увеличении радиуса сжимаемого слоя) одномерных слоёв первоначально однородного и покоящегося газа в диапазоне от р0 = 1 до любой наперёд заданной конечной плотности р > 1. Представлены тестовые расчёты в обратном направлении изменения времени для плоского случая V = 0, а также некоторые результаты расчётов сжатия цилиндрических (V = 1) и сферических (V = 2) слоёв фиксированной массы до наперёд заданного конечного значения плотности. Приводятся некоторые сведения о построении конечно-разностного метода «Ромб», необходимые для понимания сконструированного алгоритма.
Основным результатом данной работы являются численные расчёты задачи сильного сжатия плоского симметричного слоя в прямом направлении изменения времени конечно-разностным методом «Ромб». Расчёты при использовании численно восстановленного закона движения сжимающего поршня.
1. Математическая постановка одномерной задачи безударного сильного сжатия газовых слоёв
Далее рассматривается изоэнтропическое сжатие одномерного газового слоя изнутри. Подобные течения являются решениями системы уравнений [3]:
а г + иаг +
(7 - 1)
и
а ( иг + V— г
0,
иг +
2
(1)
(7 - 1)
ааг + ииг = 0.
Здесь а = р(7-1)/2 - скорость звука в газе, р - плотность, 7 > 1 — константа в уравнении состояния р = р7/7,р — давление; и = (а, и) — искомые функции. Значения V = 0,1, 2 соответствуют случаям плоской, цилиндрической и сферической симметрии.
В случае V = 0 в (1) г = хь и = у1 — проекция вектора V на ось Ох1.
1^+1
В случае V = 1, 2 г = \ х2 ^ 0, а скорость газа и есть проекция вектора
V на радиус-вектор в плоскости ж1Ож2 при V =1 или на радиус-вектор в пространстве переменных х1, х2, х3 при V = 2 соответственно.
Пусть в некоторой окрестности точки (£ = ¿*, г = г*) при г* > 0 задано какое-либо фоновое течение И = И0(£,г), где компоненты вектора И0(£,г) — аналитические функции а0(¿,г), и0(£,г) — являются решением системы (1). На рис. 1 в области 0 находится фоновое течение, а в области 1 требуется построить искомое течение с наперёд заданным распределением одного из газодинамических параметров, которое будет сопряжено с фоновым течением через слабый разрыв.
Рис. 1. Области однородного покоя и ЦВ сжатия
Слабым разрывом будет являться звуковая характеристика С+, однозначно определяющаяся при решении задачи Коши:
^г
- = и0(£,г) + а0 (¿,г),г(^ = г*. (2)
Из аналитичности фонового течения следует существование и единственность решения данной задачи — аналитической функции г = г00(£), описывающей закон движения характеристики С+. Далее везде считается, что функция г = г00(^) известна и, следовательно, известны значения газодинамических параметров фонового течения на этой С+- характеристике:
и^г)1г = г00(^) = и00(;),
т.е.
= и|с+ = иоо(£). (3)
Задача (1), (3) является характеристической задачей Коши, для существования и единственности решения которой необходимо задать дополнительное краевое условие.
1.1. Задача о получении наперёд заданного непрерывного распределения плотности
Рассмотрим случай, когда звуковая характеристика
C+ : г = (t - t*) + г*
разделяет однородный покой и искомое течение. Тогда на ней определены условия:
u|с + = 0, a\C+ = 1.
Область определения искомого течения ограничена прямой t = t* и звуковой характеристикой C+ фонового течения.
Требуется получить в наперёд заданное распределение плотности p\t=tt = = р*(г). К системе (1), (3) добавляется следующее краевое условие
a\t=t* = а*(г) = piY-1)/2(r). (4)
Предполагается, что наперёд заданное распределение плотности в момент времени t = t* в точке г = г* непрерывно стыкуется с плотностью газа фонового течения.
Задача (1), (3), (4) является ХЗК стандартного вида и имеет в окрестности точки (t = t*, г = г*) единственное аналитическое решение [3,9].
Далее предположим, что в момент времени t = t* необходимо в точке г = = г* получить вертикальное распределение плотности. Тогда вместо условия (4) необходимо задать условие вертикали. Такая задача, далее называемая ХЗК1, имеет единственное аналитическое решение в окрестности точки (t = t*, г = = г*) [3].
Течение в области 1 (см. рис. 1) имеет особенность в точке (t = t*, г = г*) и описывается ЦВ Римана в случае плоской симметрии (v = 0)
2 2
и(^г) =-a--const,
7 - 1 7 - 1 (5)
г - г* Y +1 2
-const
t - t* y - 1 y - 1
или её обобщением при v = 1, 2 (см. [3]).
Далее предполагаем, что искомое течение является составным. Область определения решения состоит из двух частей: нижнего треугольника (обозначен цифрой 1 на рис. 2) — области определения обобщения ЦВ Римана и верхнего треугольника (обозначен цифрой 2 на рис. 2) — области определения течения, имеющего в момент t = t* наперёд заданное распределение a = а*(г), например, постоянное — a*(г) = a* = const.
Рис. 2. Области определения ХЗК1 и ХЗК2
Характеристика С+ определяется решением соответствующей задачи Коши, которая имеет единственное аналитическое решение, поскольку ХЗК1 имеет аналитическое решение. На этой характеристике значения параметров газа являются аналитическими функциями
a = a1(t), u = u(t, a1(t)) = u1(t),
(6)
определяемыми из решения ХЗК1. В области 2 ставится новая ХЗК, называемая ХЗК2. К системе (1) добавляются условия на характеристике С+ и краевые условия
и(^г)\г=г1№
a1(t), u(t,a1(t)), a*(г*) = a1(t*).
(7)
Аналитическая функция г^Ь) считается известной, для неё выполняется соотношение
dt
U1 (t) + a1(t).
(8)
Функция и = а*(г) считается аналитической, известной и связанной соотношением и*(г) = [р*(г)](7-1)/2 с наперёд заданным в момент Ь = Ь* требуемым распределением плотности р = р*(г). Тогда ХЗК2 (1), (8), (7) имеет в окрестности точки (Ь*,г*) единственное аналитическое решение [3,9].
1.2. Задача Коши для восстановления траектории поршня
Требуется найти закон движения поршня, безударно сжимающего однородный покоящийся газ массой m* с плотностью р =1 до плотности р = p(t*,r) = = р* = const, постоянной в сжатом слое шириной d* = m*/p* (рис. 2). Эта траектория определяется при решении в обратном направлении изменения времени следующей задачи Коши:
dr , .
dt = U(t,r), (9)
r|t=i* = r* — d*.
В (9) в правой части уравнения стоит функция u(t,r), определяющаяся при решении ХЗК2.
Теоремы, доказанные в [3], утверждают, что существует ненулевая масса газа m*, расположенная в области между траекторией движения поршня r = rp(t) и прямой r = r*, которую можно безударно сжать до любой наперёд заданной плотности p*(r). Но эти теоремы не определяют максимально возможное значение m*, а следовательно и максимальную ширину слоя d*.
2. Алгоритм численного решения одномерной задачи сильного сжатия газа методом характеристик в обратном направлении изменения времени
Приведём далее алгоритм, основанный на стандартном методе характеристик с пересчётом [19] и ранее подробно изложенный в [3,10]. Здесь приводятся лишь узловые моменты алгоритма и схема расчёта.
1. С помощью введения инвариантов Римана система (1) сводится к системе двух обыкновенных дифференциальных уравнений (ОДУ). К системе присоединяются дифференциальные уравнения для характеристик из (2), а также учитываются начальные условия (3) и краевое условие.
2. Система ОДУ решается отдельно в области 1 и в области 2.
2.1. Алгоритм построения численного решения в областях 1 и 2
Расчёт сетки в нижнем треугольнике (область 1 на рис. 2) и вычисление значений параметров газа в её узлах происходят в обратном направлении изменения времени. Схема расчёта проиллюстрирована на рис. 3.
A. В точке (г*,Ь*) интервал [1,а*] изменения а[1,а*] разбивается на N промежутков, здесь а* = р*7-1)/2. Для каждого значения а» = 1 + ¿Да, Да = = а*]у 1,^ = 0,1,...,^, используя формулу обобщённой ЦВ Римана [3] при Ь = Ь*, получаем значения скорости и».
B. Из точки (г*,Ь*) выпускаются N штук С+ -характеристик.
C. Задаём число шагов по времени М, таким образом определяется шаг по времени ДЬ = Ь*/М.
D. Делается шаг Д£ в обратном направлении изменения времени вдоль характеристики С+ (обозначена на рис. 3 как С+, у которой вторая цифра в индексе обозначает принадлежность характеристики к области 1). Из точки на этой характеристике, соответствующей моменту времени (£*—гДЬ) выпускается С--характеристика до пересечения с линией С+. Здесь и далее г = 1 ,...,М — номер шага по времени.
Рис. 3. Расчёт нижнего (слева) и верхнего (справа) треугольников
E. Используя разностные аналоги уравнений (4) и значения газодинамических параметров в начальных точках соответствующих характеристик, получаем газодинамические параметры в точке пересечения.
F. Повторяем процедуру, выпуская С--характеристику из полученной точки пересечения (см. рис. 3). После того, как найдена точка пересечения с С+ -характеристикой, переходим к пункту С, делая новый шаг по времени. Если достигнуто значение £ = 0, значит решение в области 1 построено.
Расчёт сетки в верхнем треугольнике (область 2 на рис. 2) также ведётся в обратном направлении изменения времени. Схема расчёта полностью аналогична нижнему треугольнику, с той лишь разницей, что характеристики С+, С+,... выпускаются с прямой £ = £*, на которой известно значение а = а*.
2.2. Алгоритм построения траектории сжимающего поршня
Построение траектории движения сжимающего поршня в виде г = гр(£), как и построение характеристической сетки, происходит в обратном направлении изменения времени, т. е. при £ ^ £*.
A. Траектория г = гр(Ь) движения сжимающего поршня строится исходящей из точки (£ = £*,г = г* — й*), где й* определяет ширину уже сжатого слоя газа. Из этой точки выпускаем прямую, определяемую разностным аналогом задачи Коши (9).
B. Находим точку пересечения этой прямой с одной из линий характеристической сетки, построенной в соответствии с вышеописанными алгоритмами. В найденной точке пересечения скорость газа и находится линейной интерполяцией по значениям в ближайших узлах. Пусть найденная точка пересечения
имеет координаты (Ь = ¿, г = г). Выпускаем из неё прямую, определяемую разностным аналогом задачи Коши
до пересечения этой прямой со следующей линией построенной характеристической сетки.
С Продолжаем вычисления по этому алгоритму до тех пор, пока либо траектория поршня не пересечёт характеристику С+ нижнего треугольника, либо не выйдет за пределы области, в которой существует решение.
D. В последнем случае следует уменьшить ширину сжатого слоя и вернуться к началу алгоритма.
Завершая вторую часть статьи, необходимо сделать некоторые выводы:
1. Задача в обратном направлении изменения времени с математической точки зрения поставлена правильно, в том числе существует единственное решение для двух ХЗК и существует единственная траектория сжимающего поршня, дающая в момент времени Ь = Ь* требуемое распределение плотности для нужной массы газа.
2. Закон движения сжимающего поршня численно восстановлен, т.е. возможно сформулировать конкретные рекомендации к физическому эксперименту.
3. Требуется решение другой начально-краевой задачи в прямом направлении изменения времени при заданном законе движения сжимающего поршня. Необходимость дальнейшего исследования обусловлена физической сутью изучаемой проблемы. Требуется ответить на вопрос: переводит ли, при возрастании времени, сжимающий по полученному закону поршень исходную массу газа в нужное физическое состояние, при котором достигается необходимое значение плотности.
3. Конечно-разностный метод «Ромб» для численного решения одномерных уравнений газовой динамики
Рассмотрим одномерную систему уравнений газовой динамики в лагранже-вых переменных.
дЬ + ^ дх 0
1
где ¿,ш, р, V = -,е,р — соответственно: время, массовая лагранжева коор-Р
дината, плотность, удельный объём, удельная внутренняя энергия, давление, и2
Е = е + ^--удельная полная энергия.
Система замыкается уравнением состояния идеального газа р = (7 — 1)рс^Т, е = с€Т, где с„, Т — теплоёмкость при постоянном объёме и температура газа соответственно.
К системе (10) присоединяются граничные и начальные условия в следующем виде
и|Ж=Ж; (t) Up(t)
u|x=xr (t) = 0,
u(x, 0) = 0,
(11)
где xi,xr — левая (left) и правая (right) границы отрезка, на котором решается система (10). В момент времени t = 0 xl = rp(0),xr = г*. up(t) — скорость сжимающего поршня, полученная в соответствии с алгоритмом из пункта 2.2.
Применив для (10) неявную аппроксимацию и введя шаг по времени т, получим:
т
ио+1 - -un
т
En+1 — E n
5u \ dm /
n+1
n+1
+ v
n+1
5p \ dm у 5 (pu)
(12)
т
dx
n+1
0.
Систему (12) будем решать, разбивая исходный отрезок [х , Хг] на N ячеек (пронумеруем узлы полученной сетки х0 = х, хь х2,..., х^ = хг). Граничные условия на левой и правой границах запишем так:
аоРо + вои° = г p ®NpN + вмuN =
rp(tn),
(13)
Система (12) нелинейная, для её разрешения линеаризуем давление и удельную внутреннюю энергию по Ньютону, введём итерационный процесс. После введения через соотношения связи новых неизвестных (давления Р и скорости и, определяемых в узлах пространственной сетки) и ряда преобразований получим СЛАУ относительно Р, и с четырёхдиагональной матрицей, которая решается методом, аналогичным методу потоковой прогонки [20].
Затем определяются новые координаты хм+1, плотность рм+1, ячеечные скорости им+1, удельная внутренняя энергия ем+1, из УРС восстанавливаем давление в ячейке рм+1 на новой итерации ^ + 1. Подробный вывод формул и обоснование метода изложены в [18].
Определив все величины на текущей итерации по нелинейности p, проверяем условие сходимости итераций
Ьт/2 - Рт/21 ^ ¿сои(1 + РЧ+ 1/2). (14)
Если условие (14) выполнено во всех ячейках, то итерации заканчиваются. Далее переходим на новый временной шаг.
и(*°) = 0; ь(1°) = 1; Е(1°) = 1/(7(7 - 1)). (15)
4. Результаты расчётов
В настоящей статье было описано несколько различных численных алгоритмов. Сначала будут приведены результаты расчётов в обратном направлении изменения времени. Эти результаты приведены для демонстрации работоспособности и адекватности алгоритмов. Затем будут приведены результаты расчётов в направлении возрастания времени методом «Ромб».
4.1. Расчёты в обратном направлении изменения времени
Одним из поводов для этой работы стала возросшая с момента публикации [11] вычислительная мощность компьютеров. Изложенные выше алгоритмы были реализованы в виде программы для PC на языке С++. Далее приводятся в том числе и результаты на более подробных сетках, чем использованы в работе [11].
_____________
-2 0 2 4 6 8 10
Г
Рис. 4. Пример характеристической сетки
Для проверки работоспособности программ сначала были посчитаны задачи сильного сжатия для случая плоской симметрии, когда известно точное решение [5,6]. Рассчитанные поля совпали с этим точным решением с точностью, близкой к машинному нулю. На рис. 4 приведён пример рассчитанной характеристической сетки и траектории поршня. По оси абсцисс отложена координата г траектории частицы газа, полученная в процессе решения, по оси ординат — время ¿. Выделенная прямая в верхней части рисунка — характеристика С+ верхнего треугольника, которая пересекает траекторию поршня в момент времени (¿* — 0.0001).
В следующей таблице приводятся результаты сравнения четырёх численно построенных решений с точным решением для плоской симметрии (V = 0).
Таблица 1. Результаты тестовых расчётов
y m* Pmax || cex cnum|| C || uex unum|| C
1.4 10 104 1.66355818 ■ 10-12 8.37019342 ■ 10-12
1.4 100 105 2.25801919 ■ 10-11 1.12819975 ■ 10-10
5/3 10 104 4.13891143 ■ 10-13 1.24167343 ■ 10-12
5/3 100 105 1.81159976 ■ 10-10 5.43906253 ■ 10-10
Варьировались следующие параметры: y — показатель политропы идеального газа, m*, pmax — масса сжимаемого слоя и максимальная плотность слоя газа. Были посчитаны и приведены в таблице нормы разности численно-полученных и точных значений скорости звука и скорости газа. Используются обозначения: cex,cnum,uex,unum — точное (exact) и численное (numerical) значения скорости звука и скорости газа, ||F||C = maxj(|Fj|) — норма.
Построенная численно траектория сжимающего поршня также удовлетворительно согласуется с точным решением [3]. Поточечная норма разности убывает обратно пропорционально возрастанию числа точек в расчёте. В таблице 2 приведены нормы ошибки при сгущении сетки для случая y = 1.4, d* = 10,pmax = = 104.
Таблица 2. Результаты тестовых расчётов траектории сжимающего поршня
N xM 250x250 500x500 1000x1000 2000x2000
|| rex rnum|| C 4.41842 ■ 10-3 1.09307 ■ 10-3 2.71721 ■ 10-4 6.78031 ■ 10-5
Основываясь на данных из таблиц 1 и 2, можно сделать вывод об адекватном описании предложенным алгоритмом решения задачи о сильном сжатии одномерного газового слоя в случае плоской симметрии.
Далее приводятся результаты численных экспериментов в случае V =1 и V = 2, для которых не существует точного решения в аналитическом виде.
В качестве критерия точности полученного решения задачи мы, вслед за [3, 10], выбрали следующий: масса несжатого слоя шириной равна массе
сжатого слоя шириной ¿*, т. е. проверялось равенство масс в начальный и конечный моменты сжатия газа. Кроме того, нами была посчитана работа поршня по сжатию газа и проведено сравнение с изменением полной энергии в сжимаемой системе.
Таблица 3. Результаты численных экспериментов
V 7 Ь * т* Ртах NxM Ат ,%
0 1.4 10.01166 10 104 500x500 0.117
1 1.4 1.04242 10 104 500x500 0.131
2 1.4 0.44934 10 104 500x500 0.132
0 5/3 10.04016 10 104 500x500 0.402
1 5/3 1.03799 10 104 500x500 0.502
2 5/3 0.45187 10 104 500x500 0.516
0 1.4 100.0202 100 105 2000x2000 0.02
1 1.4 4.4574 100 105 1000x1000 0.103
2 1.4 1.5411 100 105 1000x1000 0.092
0 5/3 100.1122 100 105 2000x2000 0.112
1 5/3 4.5588 100 105 1000x1000 0.844
2 5/3 1.5656 100 105 1000x1000 0.735
Такая поверка выполнялась и в промежуточные моменты сжатия. В приведённых в таблице 3 данных различие масс Ат сжатого и несжатого газа колеблется от 0.02% до 0.85%, относительная погрешность сохранения полной энергии не превышает сотых долей процента. Точность можно повышать, уменьшая АЬ и увеличивая N и М.
По представленным результатам можно сделать вывод, что предложенный алгоритм решения задачи сильного сжатия строит решение и адекватно описывает процесс сжатия с использованием найденной траектории сжимающего поршня. Времена счёта на PC не превышают нескольких минут для задач с десятками тысяч расчётных точек.
4.2. Результаты расчётов в прямом направлении изменения времени методом «Ромб»
Для реализации основной цели исследования необходимо использовать построенную траекторию движения сжимающего поршня как граничное условие при расчётах в прямую. Для этого был запрограммирован надёжный безусловно-устойчивый конечно-разностный метод «Ромб». Ниже излагаются первые результаты представляемого подхода. Для проверки работоспособности избранного метода при расчётах задачи сильного сжатия была выбрана одна из задач в плоскосимметричной постановке: 7 = 1.4, ¿* = 10, ртах = 104.
Сначала задача решается в обратом направлении изменения времени. Затем из траектории сжимающего поршня выбираются левая х = тр(£0) и правая хг = тр(г*) границы. На границах ставятся условия: щ = щр(£),щг = 0. Начальные условия щ0 = 0, V) = = 1,£0 = ^(^1—1). Далее решается начально-краевая задача, при этом на каждом временном шаге левое граничное условие интерполируется, т.к. значение скорости и координаты сжимающего поршня известны только в фиксированные моменты времени.
Таблица 4. Погрешности метода «Ромб» при расчёте плоскосимметричной задачи
tend 11 uex unum11C 11 @ex &num11C @ ex &num C 100% &num C
0.5 • tex 0.0018 3.6304 • 10-4 0.031
0.8 • tex 0.0043 8.8103 • 10-4 0.0673
0.9 • tex 0.0088 0.0018 0.12
0.99 • tex 0.0835 0.0170 0.78
На рисунках 5-6 представлены результаты расчёта вышеупомянутой плоскосимметричной задачи на различные моменты времени.
В таблице 4 представлены погрешности расчётов при варьировании момента окончания счёта tend. Используются обозначения: uex,aex — точные значения скорости и скорости звука (ЦВ Римана), unum,anum — значения скорости и скорости звука, полученные численно.
Из рисунков 5-6 и таблицы 4 можно сделать вывод, что метод «Ромб» хорошо описывает распространение центрированной волны сжатия по покоящемуся газу до момента времени tend = 0.999t *
Из таблицы 5 видно, что в вышеприведённые моменты времени сжатия существенной массы газа не происходит. Требуется считать дальше. На следующем рисунке представлены результаты расчётов до более поздних момен-
Рис. 6. Скорость и скорость звука на момент времени £ = 0.9£*
Таблица 5. Время до сжатия части массы до требуемой плотности в долях от
V Время до сжатия части массы до требуемой плотности, доля от ¿*
0.1 ■ т* 0.3 ■ т* 0.5 ■ т* 0.7 ■ т* 0.9 ■ т*
0 0.999980525 0.9999889 0.9999921 0.9999953 0.9999984257
1 0.999983714 0.9999918762 0.9999918762 0.999996381 0.9999988571
2 0.9999899 0.9999921 0.9999944 0.9999966 0.9999989
тов времени. Из рисунка 7 видно, что на левой границе отрезка формируется область газа, сжатого до требуемой плотности. Точное решение — сплошная линия, численное — пунктир.
Стоит отметить, что в районе состыковки ЦВ с однородным покоем значительным становится «размазывание» волны, что сказывается на точности результатов. Более того, численное возмущение, приходящее на правую границу, вынуждает отодвигать её, т.к. иначе погрешность расчёта в области правой границы делает решение бессмысленным. Следовательно, необходимо каким-то образом модифицировать алгоритм для более аккуратного описания течения в моменты времени, близкие к ¿*. Планируется несколько путей модификации: либо разбиение исходного пространственного отрезка на две счётные области, либо исходное неравномерное разбиение отрезка интегрирования для более детального описания интересующей области.
Выводы
Представленные результаты позволяют сделать некоторые выводы.
1. Мы получили закон сжатия исходной массы газа до заданной плотности, решая задачу в обратном направлении изменения времени. Показана работоспособность предложенного алгоритма расчёта одномерных задач
Рис. 7. Скорость и скорость звука на при £ = 0.99999£Ф (вверху) £ = 0.999995£Ф (внизу)
сильного сжатия. Точность и скорость расчётов в обратном направлении изменения времени позволяет восстанавливать искомую траекторию сжимающего поршня за 1-2 минуты. Таким образом, рекомендации по проведению физического эксперимента могут быть сформулированы. Однако решение задачи сильного сжатия в обратном направлении изменения времени с явным раскрытием особенности в точке (г*,Ь*) достаточно сложно и не является общепринятым.
2. Реализован и оттестирован безусловно устойчивый конечно-разностный метод «Ромб» для решения одномерной СУГД в лагранжевых переменных.
3. Проведены расчёты плоскосимметричной задачи сильного сжатия газа при возрастании времени методом «Ромб». Результаты показывают, что численный алгоритм с хорошей точностью восстанавливает решение вплоть до момента времени 0.999 ■ ¿*. Дальнейшие расчёты демонстрируют, что при приближении к особенности в решении, когда производные в точке г* начинают быстро расти, точность расчёта ухудшается, однако, в решении реализуется требуемое физическое состояние исходной массы газа.
Автор выражает своему научному руководителю профессору С.П. Баутину признательность за внимание, помощь и поддержку.
Литература
1. Забабахин Е.И., Забабахин И.Е. Явления неограниченной кумуляции. М. : Наука, 1988. 173 с.
2. Накколс Дж.Г. Осуществимость инерциально-термоядерного синтеза // Успехи физ. наук. 1984. Т. 143, № 3. С. 467-482.
3. Баутин С.П. Математическое моделирование сильного сжатия газа. Новосибирск : Наука, 2007. 308 с.
4. Риман Б. О распаде плоских волн конечной амплитуды // Соч. М-Л. : ОГИЗ, 1948. С. 376-395.
5. Станюкович К.П. Неустановившиеся движения сплошной среды. М. : Гос. изд-во техн.-теор. лит-ры, 1955. 804 с.
6. Мизес Р. Математическая теория течений сжимаемой жидкости. М. : Изд. иностранной литературы, 1961. 588 с.
7. Баутин С.П. Математическая теория безударного сильного сжатия идеального газа. Новосибирск : Наука, 1997. 160 с.
8. Баутин С.П. Асимптотические законы безударного сильного сжатия квазиодномерных слоев газа // Прикл. матем. и механика. 1999. Т. 63, Вып. 3. С. 415—423.
9. Баутин С.П. Характеристическая задача Коши и её приложения в газовой динамике. Новосибирск : Наука, 2009. 367 с.
10. Баутин С.П., Николаев Ю.В. Об одном методе расчёта безударного сильного сжатия одномерных слоёв газа // Вычислительные технологии. 2000. Т. 5, № 4. С 3-12.
11. Николаев Ю.В. О численном решении задачи безударного сильного сжатия одномерных слоёв газа // Вычислительные технологии. 2001. Т. 6, № 2. С. 104-108.
12. Анучин М.Г. Влияние теплопроводности на неограниченное безударное сжатие плоского газового слоя // ПМТФ. 1998. Т. 39, № 4. С. 25--32.
13. Рощупкин А.В. Исследование некоторых характеристических задач Коши, возникающих при решении неодномерных задач безударного сильного сжатия газа // Вычислительные технологии. 2002. Т. 7 (совм. вып.), Ч. 4. С. 96-103.
14. Артемьев А.Ю., Делов В.И., Дмитриева Л.В. и др. Численное моделирование безударного неограниченного сжатия газа в переменных Лагранжа по методике Д. // ВАНТ. Сер. Мат. моделирование физ. процессов. 1995. Вып. 4. С. 42-47.
15. Сараев В.А. Решение задач о неограниченном безударном сжатии газа в областях специальной формы по программам комплекса МИМОЗА // ВАНТ. Сер. Мат. моделирование физ. процессов. 1998. Вып. 4. С. 3-16.
16. Сучков В.А. Истечение газа в вакуум на косой стенке // Прикл. матем. и механика. 1963. Т. 27, Вып. 4. С. 739-740.
17. Сидоров А. Ф. Некоторые оценки степени куммуляции энергии при плоском и пространственном сжатии газа // Доклады АН. 1991. Т. 318, № 3. С. 548-552.
18. Гаджиев А. Д., Писарев В. Н. Неявный конечно-разностный метод «Ромб» для численного решения уравнений газовой динамики с теплопроводностью. // Ж. выч. матем. и матем. физ. 1979. Т. 19. N. 5. С. 1288-1303.
19. Рождественский Б. Л., Яненко Н. Н. Системы квазилинейных уравнений и их приложения к газовой динамике. М.: Наука, 1968. 529 с.
20. Дегтярев Л. М., Фаворский А. П. Потоковый вариант метода прогонки. // Ж. выч. матем. и матем. физ. 1968. Т. 8. N. 3. С. 679-684.
ONE-DIMENSIONAL MATH MODELING OF IDEAL GAS STRONG COMPESSION IN R. MISES CONFIGURATION
N.S. Novakovskiy
Postgraduate Student, e-mail: [email protected]
Ural State University of Railway Transport (USURT)
Abstract. The article describes results of one-dimensional gas layers strong compression problem's math modeling. Mathematical formulation of the problem includes two characteristic Cauchy problems (cCpl and cCp2) and one Cauchy problem for compression piston's trajectory reconstruction with time decrease. The article provided numerical algorithm based on typical characteristic method for solving these problems. Algorithm for compression piston's trajectory reconstruction is also described. Finite-difference method "ROMB" for solving system of gas dynamics equations while time increasing is described too. The main results of numerical simulations are shown in graphs and tables.
Keywords: gas strong compression, the characteristic Cauchy problem, centered wave of compression, characteristic's method, finite difference method "ROMB".
