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

Вычислительные технологии

Том 5, № 4, 2000

ОБ ОДНОМ МЕТОДЕ РАСЧЕТА БЕЗУДАРНОГО СИЛЬНОГО СЖАТИЯ ОДНОМЕРНЫХ СЛОЕВ ГАЗА

С. П. Блутин, Ю. В. Николаев Уpальский государственный университет путей сообщения

Екатеринбург, Россия e-mail: bautin@ac.usart.ru

A mathematical description of a shockless powerful compression of an ideal gas to any preassigned density (including the infinite value) is important in connection with the problem of laser thermonuclear synthesis. The description of the algorithm for calculating the same shockless powerful compression of unidimensional gas layers is presented here. The algorithm efficiency has been tested on specific calculation examples.

Математическое описание процесса безударного изэнтропического сжатия идеального газа до любого наперед заданного значения плотности, в том числе до бесконечной плотности (подробную библиографию см. в [2]), представляет интерес в связи с проблемой лазерного термоядерного синтеза [4, 9]. В случае плоскосимметричных течений (v = 0) простая центрированная волна Римана описывает сжатие плоского слоя газа в конечный момент времени t = t* до бесконечной плотности [13]. Состыковка центрированной волны Римана с однородным потоком газа дает решение задачи о получении в сжатом плоском слое любого конечного значения плотности [8]. На рис. 1 изображена область определения этого кусочно-составного решения. На рис. 2 — поверхность c = c(t,xi) скорости звука газа у этого кусочно-составного решения.

В случае цилиндрически (v =1) и сферически (v = 2) симметричных течений автомодельные решения Л. И. Седова [11] описывают безударное сильное сжатие идеального газа, первоначально однородного и покоящегося внутри цилиндра или сферы (см., например, [5, 6]). Задача о выборе оптимальных законов движения непроницаемых поршней, безударно сжимающих одномерные слои газа, рассмотрена в [7, 12]. В отличие от предложенной в [12] схемы движения поршня, которая совпадает с приведенными ранее [8, 13] (случаи сжатия до бесконечной или до конечной плотностей соответственно), в работе [7] предложена принципиально другая схема движения поршня: точка, в которой центрирована волна сжатия, лежит на поршне. В работах [7, 12] отсутствует строгое математическое обоснование существования решений.

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

© С. П. Баутин, Ю.В. Николаев, 2000.

Рис. 1. Область определения решения: 0 — однородный покоящийся газ, 1 — центрированная волна Римана, 2 — однородный движущийся газ.

V = 0,1 < I

Рис. 2. Поверхность скорости звука составного течения.

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

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

все течение до момента сжатия определяются однозначно как решение характеристической задачи Коши стандартного вида [1, 2].

Доказанные в [2] теоремы утверждают, что существуют цилиндрические и сферические слои с ненулевой массой газа, которые можно безударно сжать до любой плотности. Однако эти теоремы не позволяют определить — предельную ширину исходных слоев, которые при фиксированных V, 7 можно безударно сжать до заданной плотности р*. Такое возможно в случае, если течение строится не только в некоторой окрестности точки (£ = ¿*, г = г*), но и в существенно большей области.

В работе описан алгоритм расчета безударного сильного сжатия одномерных слоев первоначально однородного и покоящегося газа в диапазоне от р0 = 1 до любой наперед заданной конечной плотности р* > 1, а также иллюстрация его "работоспособности" на примере конкретных расчетов. В основе численного метода лежит известный [10] метод характеристик, позволяющий достаточно точно рассчитывать одномерные нестационарные течения газа.

1. О математическом pешении задачи безударного сильного сжатия одномерных газовых слоев

Одномерные симметричные течения для политропного газа являются решениями системы уравнений:

(y — 1) ( и\

at + uor +--^-° l^U- + v-J = 0,

2 2 2 2 (1-1)

ut + uur + ---s aar +— a ssr = 0,

(Y — 1) Y

k St + usr = 0.

Здесь r = x2 > 0, a = p(Y-1)/2, p — плотность, y — константа в уравнении состо-

яния p = A2(S)pY/y, Y > 1, P — давление, S — энтропия, s = A(S), c = sa — скорость звука, U = (a, u,s) — искомые функции. При v =2 (случай сферической симметрии) u есть проекция V вектора скорости газа в точке M G R3 на радиус-вектор этой точки в R3, и предполагается, что V параллелен этому радиус-вектору. При v =1 (случай цилиндрической симметрии) u есть проекция V на радиус-вектор этой точки в R2, и при этом проекция V на ось Ox3 считается равной нулю. При v = 0 (случай плоской симметрии) u есть проекция V на ось Ox1 и проекции V на другие декартовы оси равны нулю.

Далее будет рассмотрено сжатие слоя газа изнутри. Пусть в некоторой окрестности заданной точки (t = t*,r = r*), r* > 0, известно аналитическое решение системы (1.1) U = Uo(t,r) — фоновое течение, у которого Co(t*,r*) = so(t*, r*)ao(t*, r*) > 0. В расчетах в качестве фонового течения будет браться однородный покоящийся газ и поэтому далее s = 1 , а a — скорость звука в газе.

С фоновым течением U0(t,r) будет сопряжено через слабый разрыв искомое решение, являющееся волной сжатия и распространяющееся по фоновому течению.

Слабым разрывом будет звуковая характеристика семейства C + фонового течения, проходящая через точку (t*,r*). Характеристика C + однозначно определяется при реше-

нии задачи Коши

dr

— = u0(t,r) + so(i,r)ao(i,r) , r(t)|i=it = r*.

Из аналитичности фонового течения следует существование и единственность решения данной задачи — аналитической функции г = Го (¿), описывающей С+ -характеристику. Далее везде считается, что функция г = г0(£) известна и, следовательно, известны значения газодинамических параметров фонового течения на этой С+-характеристике:

ё/со(*)= ад, г) цго(4).

С учетом исходного предположения, эти параметры также являются аналитическими функциями в некоторой окрестности точки £ =

Поскольку здесь решается задача о безударном сжатии газа до бесконечной плотности, то у искомого течения при £ ^ — 0 производная аг должна обращаться в бесконечность.

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

г = г(£', а) ,

£ = И .

Теперь ¿',а — новые независимые переменные, а г, и, 5 — искомые функции новых независимых переменных.

С учетом замены система (1.1) примет вид

(т — 1)

г (и — г4) +---— а(гмст + ) = 0,

raut + (u - rt)uCT + 2a2ssCT +

Y

rCTSt + (u - rt)sCT = 0.

2

(Y -1)

as2 = 0,

(1.2)

В системе (1.2) и ниже штрих у £ опущен.

Искомое решение системы (1.2) должно примыкать к фоновому течению через С+ -характеристику, которая в пространстве (£,а) задается соотношением а = а00(£). Поэтому нужное нам решение должно также удовлетворять условиям

r(t,a)Uaoo(t) = ro(t)>

u(t,a)Uaoo(t) = u00 (t), s(t,a)|a=aoo(t) = S00(t).

Кроме (1.3), искомое решение должно удовлетворять "условию вертикали"

r(i,a)|t=t+ = r* , r* = const > 0 .

(1.3)

(1.4)

"Условие вертикали" передает следующее свойство решения: в физическом пространстве кривые а(£, г)|4=сопв4<^ = а(г) при £ ^ — 0 переходят в вертикальную прямую г = г*. Только при таком поведении а(£, г) при £ ^ — 0 может получиться сжатие до бесконечной плотности.

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

Теорема 1. Если компоненты вектора и0(г,г) являются аналитическими функциями в некоторой окрестности точки (г*,г*), то задача (1.2) - (1.4) имеет единственное решение, аналитическое в некоторой окрестности точки (г = г*, а = а0), где а0 = а00(г*). Доказательство теоремы приведено в [2]. Решение задачи (1.2)-(1.4) имеет вид:

и(г, а)

2

_7 - 1 г(г, а) = г* +

з°а + и0

+ (г - г*)и(г,а)

) 50а + и0

7 - 1

в(г,а) = а0 + (г - г*)г(г,а)

(г - г*) + (г - г*)2г(г,а)

(1.5)

где и(г,а), г (г, а), §(£, а) — некоторые функции, аналитические в окрестности точки (г = г*, а = а0).

Если известную центрированную волну Римана

представить в виде

0 + ^ - 1 а = аи +--2— и ,

и =

2 (х - х*) 2

7 +1 (г - г*) т +1

а

и

х

а-

7 - 1 7 - 1

а

х* +

7 + 1 7 - 1

а-

7 - 1

а0

(г - г*)

5 = 1

то можно заметить, что:

1. Вне зависимости от V главная часть решения задачи (1.2)-(1.4) в виде (1.5) при г ^ г* полностью совпала с центрированной волной Римана.

2. Построенное решение (1.5) так же, как и центрированная волна Римана, при г < г* описывает волну сжатия, а при г > г* — волну разрежения.

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

С помощью решения (1.5) поставим задачу Коши для обыкновенного дифференциального уравнения

п(г,а1(г)) + гст (г,ах(г))

и(г,а1(г)) + а1(г) , а1(г)|4=, = а

где а* = р*7 1)/2, а константа р* равна значению р*(г*). Функция р* (г) и является тем наперед заданным распределением плотности, до которого нужно в момент г = г* сжать исходный слой газа. Функция а = а1 (г), являющаяся решением этой задачи, задает в

1

5

2

2

2

*

пространстве независимых переменных а характеристику С+ волны сжатия (обобщения центрированной волны Римана), выходящую из точки (£ = £*, а = а*). Эту кривую а = а^) необходимо строить при £ < £*. Функция г = г(£, а^)) определяет эту же характеристику С+ в пространстве независимых переменных г — здесь она выходит из точки (£ = £ *,г = г *) в обратном направлении изменения времени. На этой характеристике С+ значения параметров газа являются аналитическими функциями

а = а1(£) , и = и(£,а1(£)) = и1(£) , 5 = з(£,а(£)) = 51 (£) ,

определяемыми из решения (1.5). Еще раз подчеркнем, что если течение в нижнем треугольнике строится в пространстве независимых переменных а, то оно особенностей не имеет. В пространстве переменных г это течение является обобщением центрированной волны Римана и имеет в точке (£ = £ *, г = г *) особенность — предельные значения газодинамических параметров в этой точке различны на разных прямых, входящих в эту точку.

Теорема 2. Задача

(7 — 1) ^ и\ а4 + иаг +--^-а ^иг + V-^ = 0 ,

и + ииг +

52ааг +— а255г = 0

(7 — 1)

7

5 + и5г = 0

а(£,г)|с+ = а1(£) , и(£,г)|с+ = и1(£) , 5(£,г)|

, а(£,г)|4=4, = а* (г) , а* (г*) = а1(£ * )

С++

51(£)

является характеристической задачей Коши стандартного вида, у которой в некоторой окрестности точки (£ = £ *,г = г *) существует единственное аналитическое решение, если функция а * (г) является аналитической в некоторой окрестности точки г = г *.

Именно так формулируется задача о получении в верхнем треугольнике наперед заданного распределения функции а|4=4н, = а * (г) (что эквивалентно заданию плотности газа р|4=4н_ = р * (г)), непрерывно примыкающего к заданному фоновому течению — в данном случае — к обобщению центрированной волны Римана.

Решение последней задачи и в самой рассматриваемой точке (£ = £ *, г = г *), и в некоторой окрестности этой точки особенностей не имеет, в том числе и на характеристике С+ : г = г(£,а1(£)). Через эту характеристику состыкованы решения в нижнем и верхнем треугольниках.

2. Алгоритм расчета течения методом характеристик

Будем считать, что фоновое течение С/0 (£, г) является однородным покоем, т.е. задается соотношением

^ Г и0(£,г) = 0 ,

ад, г) = I а0(£,г) = 1 ,

[ ад, г) = 1.

Тогда все рассматриваемые течения в газе будут изэнтропическими, т. е. s(t, r) = const = 1.

С учетом этого система (1.1) примет вид:

Y — 1 y — 1 au at + uar +---— aur = — v—---,

2 2 2 Г

ut +--aar + uur = 0 .

Y — 1

Рассмотрим случай сжатия слоя газа изнутри, когда звуковая характеристика, разделяющая фон и искомое течение, является C+- характеристикой r = (t — t*) + r*, на которой: u|C+ = 0 и a|c + — 1. Область определения искомого течения ограничена прямой t — t* и звуковой характеристикой C+ фонового течения. Эта область состоит из двух частей: нижнего треугольника — области определения обобщения центрированной волны Рима-на и верхнего треугольника — области определения течения, имеющего в момент t = t* наперед заданное распределение a = a*(r) (например, постоянное — a*(r) = a* = const).

Расчет течения производится стандартным методом характеристик (с пересчетом) [10], с построением характеристической сетки по слоям. Сначала построим характеристическую сетку в нижнем треугольнике. Для этого разобьем отрезок [1, a*] на n равных частей. Пусть ai задают значения a в точках разбиения: ai = 1 + i • Aa (i = 0,1, ... , n), где Aa = (a* — 1)/n. Искомое течение при t ^ t* — 0 описывается формулами (1.5), поэтому при t = t*

o a — 1 u¿ = 2-- .

Y — 1

По ui и ai вычисляются значения инвариантов Ri,L¿. Таким образом, точка (t*,r*) определяет нулевой слой, в котором известно n +1 значений Ri и L¿. По этим значениям однозначно находятся в плоскости переменных t,r звуковые характеристики C+:

dr

— = ui + ai i = 0, ... , n, dt

выходящие из точки (t*, r*). На рис. 1 последняя используемая С+-характеристика этого течения обозначена как C+.

Из точки, лежащей на С+-характеристике и соответствующей моменту времени tj = t* — At (At = const > 0 — шаг по времени), выпускаем C- до пересечения с C+. Используя формулы метода характеристик, находим как точку пересечения этих характеристик, так и значения инвариантов R, L в этой точке. Из найденной точки продолжаем C- до пересечения с C+. Определяем точку пересечения C- и C+, а также значения инвариантов в этой точке и так далее, до C+. Построен слой с номером 1. Теперь, используя точки слоя 1 и точку на C+, соответствующую моменту времени t2> = t* — 2At, строим слой с номером 2 и т. д.

Таким образом, расчет обобщения центрированной волны Римана — течения, определенного в нижнем треугольнике, — происходит в направлении убывания переменных t, r от точки (t*, r*).

Расчет верхнего треугольника производится аналогичным образом. В качестве нулевого слоя берется точка (t*,r*). В ней заданы a = a* и u = u* = 2(a* — 1)/(y — 1), а на линии t = t* полагаем a = const = a*. Из точки (t^,rj) (из сетки в нижнем треугольнике) выпускаем C- до пересечения с прямой t = t*. В этой точке значения инвариантов считаются следующим образом. Вдоль C- можно перенести только инвариант L, но на линии

£ = £ * известно значение а = а *, тогда Я в точке пересечения вычисляется через Ь и а * по известной формуле

Я = Ь + (1.6)

7 — 1

Таким образом строится первый слой в верхнем треугольнике, который является продолжением первого слоя из нижнего треугольника. Второй слой в верхнем треугольнике содержит на одну точку больше и строится следующим образом: из точки пересечения прямой £ = £ * и С-, проходящей через первый слой, выпускаем С+-характеристику. Из (£п, гП) выпускаем С- до пересечения с С +. Затем С— продолжается до пересечения с прямой £ = £ *. Значения инвариантов в точке пересечения считаются с помощью формулы (1.6). Слой с номером 2 построен. Слой 3 содержит уже четыре точки и строится аналогично слою 2. Таким способом последовательно определяются слои расчетной сетки в верхнем треугольнике, являющиеся продолжением слоев с теми же номерами расчетной сетки из нижнего треугольника.

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

Одной из целей исследований безударного сильного сжатия является определение максимальной массы газа, которая сжимается безударно до заданного значения а *. Можно ожидать, что при больших а в течении газа будут возникать особенности и это приведет к тому, что с ростом значений а *, 7 и V можно ожидать уменьшения

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

Кроме этого, необходимо следить за размерами сторон каждой из получаемых ячеек. Максимальная длина стороны ячейки устанавливается от 1.5Д£ до 2.0Д£ (возможно и меньше). Если возникнет ячейка больших размеров, то далее, начиная с такой ячейки, шаг Д£ уменьшается, и вводятся дополнительные слои и лучи.

Построение траектории движения поршня, так же, как и построение характеристической сетки, происходит в обратном направлении изменения времени, т. е. при £ < £ *.

Траектория движения строится исходящей из точки (£ = £ *, г = г * — d*), где d* определяет ширину уже сжатого слоя газа. Из этой точки выпускаем прямую, определяемую разностным аналогом следующей задачи Коши:

dг

И

г|4=4* = г * — ^ .

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

dг

^ = Ч=>= , г|4=* = г

до пересечения этой прямой с одной из линий построенной характеристической сетки и т.д.

Построенная таким образом траектория движения непроницаемого поршня может достичь характеристики фонового течения — С+-характеристики. На этом завершается построение траектории и решения всей задачи о безударном сильном сжатии в целом, поскольку определится точка (£ = £0,г = г0), которая лежит на С+-характеристике. Из нее в момент £ = ¿0 стартует непроницаемый поршень, безударно сжимающий однородный покоящийся в момент £ = £0 слой газа шириной d0 = г * — г0 с р = 1 до слоя шириной d* в момент £ = £ * с плотностью газа р = р *.

Однако, возможна ситуация, когда траектория движения поршня войдет в область, где характеристическая сетка перестает быть регулярной. В этом случае дальнейшее построение траектории движения поршня из точки (£ = £ *,г = г * — d*) невозможно, и при фиксированных V, 7, а * следует уменьшать d *.

3. Пpимеpы pасчетов

Описанный выше алгоритм был реализован в виде программы для персональной ЭВМ. Расчеты проводились при Д£ = 0.01 ^ 0.0001 и при числе разбиений отрезка [1, а * ] п = 100 ^ 1000. В качестве критерия точности полученного решения задачи был принят следующий: масса несжатого слоя шириной d0 равна массе сжатого слоя шириной d*, т. е. про-

Рис. 3. Характеристическая сетка точного реше- Рис. 4. Характеристическая сетка одного ва-

ния в случае V = 0.

рианта в случае V = 2.

верялось равенство масс в начальный и конечный моменты сжатия газа. Такая проверка выполнялась и на промежуточных стадиях сжатия. В приведенном ниже втором варианте различие масс сжатого и несжатого газа получилось меньше 0.1%. Точность можно повышать, уменьшая Д£ и увеличивая п. При этом, естественно, возрастет время счета.

В качестве примеров были рассчитаны следущие варианты.

Первый вариант: v = 0 (плоская симметрия); 7 =1.4; а* = 50. Рассчитанные поля течений до d0 = 1 совпали с известным точным решением [8]. Следует отметить, что в этом варианте происходит сжатие до р* ^ 3 • 108. На рис. 3 приведены часть характеристической сетки и траектория поршня.

Второй вариант: v = 2 (сферическая симметрия); 7 = 1.4; а* = 25. В нем удалось сосчитать до d0 = 0.1. Отметим, что газ в этом варианте сжался до р* ^ 107. На рис. 4 приведена часть характеристической сетки для этого варианта.

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

[1] Баутин С. П. Характеристическая задача Коши для квазилинейной аналитической системы. Дифференциальные уравнения, 12, №11, 1976, 2052-2063.

[2] Баутин С. П. Математическая теория безударного сильного сжатия идеального газа. Наука, Новосибирск, 1997.

[3] Баутин С. П. Асимптотические законы безударного сильного сжатия квазиодномерных слоев газа. Прикладная математика и механика, 63, вып. 3, 1999, 415-423.

[4] Забабахин Е. И., Забабахин И. е. Явления неограниченной кумуляции. Наука, М., 1988.

[5] Забабахин И.Е., СимонЕнко В. А. Сферическая центрированная волна сжатия. Прикладная математика и механика, 42, вып. 3, 1978, 373-576.

[6] Каждан Я. М. К вопросу об адиабатическом сжатии газа под действием сферического поршня. Журнал прикладной механики и технической физики, №1, 1977, 23-30.

[7] Крайко А. Н. О неограниченной кумуляции при одномерном нестационарном сжатии идеального газа.Прикладная математика и механика, 60, вып. 6, 1996, 1000-1007.

[8] Мизес Р. Математическая теория течений сжимаемой жидкости. Изд-во иностр. лит-ры, М., 1961.

[9] НАкколс Дж. Г. Осуществимость инерциально-термоядерного синтеза. Успехи фи-зизических наук, 143, №3, 1984, 467-482.

[10] Рождественский Б. Л., ЯнЕнко Н. Н. Системы квазилинейных уравнений и их приложения к газовой динамике. Наука, М., 1968.

[11] Седов Л. И. Методы подобия и размерностиi в механике. Наука, М., 1981.

[12] Сидоров А. Ф. Безударное сжатие баротропного газа. Прикладная математика и механика, 55, вып. 5, 1991, 769-779. Письмо в редакцию. Прикладная математика и механика, 56, вып. 4, 1992, 698.

[13] СтАнюкович К. П. Неустановившиеся движения сплошной среды. Гос. изд-во техн.-теор. лит-ры, М., 1955.

Поступила в редакцию 23 августа 1999 г.