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

УДК 532.59; 532.595

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

© 2009 В.В. Никонов, В.Г. Шахов

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

Поступила в редакцию 11.09.2008

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

Ключевые слова: идеальный газ, одномерная волна, численный метод, задача акустики, схема вверх по потоку, метод Годунова, подход Лагранжа, шаг по времени, точное решение, начальные условия

1. МАТЕМАТИЧЕСКАЯ ФОРМУЛИРОВКА ЗАДАЧИ

Уравнения, описывающие одномерную задачу линейной акустики, получаются из уравнений Навье-Стокса при отбрасывании конвекционных и вязких членов. Закон изменения давления считается адиабатическим, в результате в размерных переменных будем иметь [1]

др* * ди*

—Чт + р —V = 0

д 1 дх

du c Эр —г +——— = 0

д t р д x

(1)

Здесь и * - скорость, р * - плотность, х* -координата, 1* - время, с*2 = кр* /р* - квадрат скорости звука, р - давление, к - показатель адиабаты.

Введя следующие безразмерные переменные: и = и*/и0, р = р*/р0, х = х*/Ь , г = 1*Ь/и0 , с = с /и0, где Ь - характерный размер, и0 -характерная скорость, р0 - характерная плотность, система уравнений (1) примет безразмерную форму[1]

Эр ди —+ р-= 0

д г д х '

дu + c др = о д t р д x

(2)

Сделаем следующую замену переменных

w, = u + P ,

w2 = u - P, (3)

Никонов Валерий Владимирович, кандидат технических наук, инженер НТП "Авиатехнокон". E-mail: [email protected]

Шахов Валентин Гаврилович, кандидат технических наук, профессор, заведующий кафедрой аэрогидродинамики. E-mail: [email protected]

где P (смотри [2]) определяется следующим образом

ёр

(4)

P =f c^ р

тогда систему (2) можно привести к системе уравнений, аналогичных уравнениям конвекции (переноса) [2]

д w1 д w1

д t

д x

д w2 д w2 ■ - c-

д t

д x

= 0

= 0

(5)

Здесь мы будем считать, что плотность в волне сжатия изменяется по закону изоэнтропичес-кой адиабаты, а не ударной адиабаты. Данный подход справедлив [2], если плотность газа в волне сжатия возрастает не более чем в два раза по сравнению с невозмущенным газом. В случае изоэнтропической адиабаты равенство (4) приобретает вид

2c —

p = ^Чр 2 -1). k -1

(6)

Система (5) не имеет аналитического решения, но ее можно решить [2] методом характеристик. Если задано начальное распределение значений и и Р, то на следующем слое эти величины находятся как

и(1,х) = 2К(1о,х1) + +1) ,

Р(Х х) = 2- w2(to,x1+1 Я , (7)

где

t = t0 +-

c + c

1+1

x = xi + ci(t -10). (8)

h

Здесь Ь - шаг однородной сетки, ^ - начальный момент времени. Решение, полученное с помощью метода характеристик, мы будем использовать в качестве тестового.

2. ЧИСЛЕННЫЕ СХЕМЫ ПРИМЕНЯЕМЫХ МЕТОДОВ

2.1. Применение метода Годунова для решения уравнений нелинейной акустики

Хотя в книге [3] метод Годунова для решения задач нелинейной акустики не описан, здесь мы рассмотрим применение его подхода для линейных задач, адаптированного под нелинейную задачу.

В методе Годунова [3] используется "шахматная" сетка, когда в центрах ячеек определяются "малые" переменные, а на границах ячеек -"большие" переменные. Значения величин плотности и скорости в следующий момент времени при применении метода Годунова находятся как

Pi = Pi-1

1 (U-1 - Uj-1 )

1 , V^i+0.5 ui-0.5/ h

ui = ui-1

tc2

hp

(Ri+0.5 _ Ri-0.5 )

(9)

где д t - шаг по времени,

U

j-1 _

1

J-1 , J-1 w1i + w2i+1

R

j-1 =

i+0.5

PJ-1

Li+0.5

2

k _1 ^k-1 k 1+1

2cn

PJ-1 =-Pi+0.5 2

w

j-1

_ w

У

J-1 2i+1

(10)

а w1 и w2

находятся с помощью (3). Если центр ячейки 1+1 находится "внутри" тела, то граничные условия (ГУ) на стенке определяются следующим образом:

ui+_;=-ui-1

pi+1=pi-1 •

(ii)

J J-1 дtc , w1i = w1i

(wj-1 _ whi:;)

j-ь

1 1-1 А tC . 1-1 1-1 .

№2,= + - W2i ) . (12)

ГУ на стенке накладываются согласно (11), после чего они подставляются в (3).

2.3. Применение подхода Лагранжа для решения уравнений нелинейной акустики в форме конвективных уравнений

Систему уравнений (5) можно также решить с применением подхода Лагранжа к рассмотрению движения среды. Здесь предлагается использовать метод, применяемый для моделирования конвекции в методе "вихрь в ячейке" (ВЯ) [4]. При этом моделируется движение двух волн. Одна распространяется вправо и переносит величину w1. При этом в следующий момент времени

J-1

xj = xj 1 + c д t.

(13)

Другая волна распространяется влево и переносит величину w2,причем

*j

W2i = W2i

J_1

xi = xJ 1 - c д t .

(14)

Здесь w1*i и w2i обозначают величины w1 и w2 на искаженной после перемещения сетке. Для того, чтобы получить величины w1 и w2 на основной однородной сетке, аналогично методу ВЯ производится процедура перераспределения новых значений в ячейки сетки

wuk = wuл(xk_ _ xJ),

k=w^xr1 _ xi),

w

(15)

При этом ГУ (11) подставляются в (3), а потом (3) подставляются (10).

2.2. Применение схемы вперед по потоку для решения уравнений нелинейной акустики в форме конвективных уравнений

К системе уравнений (5) удобно применить конечно-разностную схему вперед по потоку (ВП), в результате будем иметь

где в качестве интерполяционной функции Л предлагается использовать кусочно-постоянное

, (16)

или кусочно-линейное распределение [1-|z/h|, |z/h| < 1 [ 0, |z/h| > 1.

Рассмотрим ГУ на стенке для правой волны. Если ячейка при своем движении (13) оказывается "внутри" тела, то она "отражается" от его поверхности, и ее координата определяется как

Mz) =

(17)

xJ = 2xb _ x; _ c д t

(18)

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

=-чн-Р-1.

(19)

Равенство (19) получается при подстановке ГУ (10) в выражения (3). После чего полученная величина (19) перераспределяется в ячейки расчетной сетки согласно (15), (16) или (15), (17). ГУ на стенке для левой волны получается аналогично.

3. ТЕСТИРОВАНИЕ МЕТОДОВ

Описанные выше четыре численных схемы тестировались на трех задачах, отличающихся начальными условиями (НУ) и ГУ.

В первой тестовой задаче НУ соответствуют однородному потоку

и0 = 1, Р0 = 1, (20)

а в момент времени Ь > 0 в точке х = 0 ставится непроницаемая перегородка.

НУ второй тестовой задачи задавались следующим образом:

и0 =

0,

0.4 < x < 0.6 x < 0.4 и x > 0.6,

р0 =1. (22)

В тестовых расчетах использовались конечные области. В первой задаче область принималась равной x е [-1,1 , а во второй и третьей задачах области имели вид x е [0,1 .

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

Д = ka

h

(23)

Р0 =

[ 1, x < 0.5 1.2, x > 0.5 ,

(21)

твердые границы отсутствовали.

В третьей тестовой задаче твердые границы также отсутствовали, а НУ задавались как

называемому критерием Куранта-Фридриха-Леви [5], с величиной коэффициента ка = 1. Заметим, что для предложенной схемы, использующей подход Лагранжа, можно варьировать шагом по времени, и наилучшие результаты получаются, когда ка - ограниченное натуральное число.

Результаты, показанные на рис. 1-6, получены после 20 шагов по времени

(t = 5.76195268-10-4 ) для с = 347Л°5 , h = 0.01. В некоторых случаях графики результатов, получаемые с помощью разных методов,

Рис. 1. Распределение плотности и скорости в задаче (20) для всех рассмотренных методов:

численное решение,

1.2 1

0.9

о.е

П.7 0.6

- аналитическое решение

и

о -10 ■20 -30 ■40 ■50

0.2 04 0.6 0.В 1

0.2 0.4 0.6 0 3

Рис. 2. Распределение плотности и скорости в задаче (21) для метода, использующего схему Годунова: О - численное решение, - аналитическое решение

1

тах

Рис. 3. Распределение плотности и скорости в задаче (21) для метода, использующего ВП схему: О - численное решение,--аналитическое решение

Рис. 4. Распределение плотности и скорости в задаче (21) для метода, использующего подход Лагранжа с распределением (16):

численное решение,

аналитическое решение

1.2 1.1 1

0.9 О.В 0.7 0.6

Ï

-

и

с -10 ■20 -30

^о

■50

„............. ri

0.2 о* 0.6 о.е

0.2 0.4 О.Б О.е

Рис. 5. Распределение плотности и скорости в задаче (21) для метода, использующего подход Лагранжа с распределением (17):

численное решение,

Ô = max

char

u

char

•100%,

(24)

аналитическое решение

визуально не отличаются, поэтому они приведены один раз. Заметим, что на графике скорости на рис. 2 наблюдается два ряда чередующихся точек. Это указывает на "пилообразный" характер решения.

Погрешность решения определялась следующим образом:

и - ис

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

Ошибка численного решения (24) после 20 шагов по времени ( ка =1) представлена в табл. Из анализа табл. видно, что наибольшую точность имеет метод, использующий подход Лаг-ранжа с распределением (16). Худшие результаты метода, использующего подход Лагранжа с распределением (17), объясняется большой численной диффузией данного перераспределения. Значительный рост погрешности в задаче (21) объясняется тем, что при таких начальных данных начинает сказываться нелинейность задачи (изменения плотности и скорости).

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

Р 1.00Ъ 1.034 1.003 1.0П2 1.001 1

G.099 □ .098 Q.QQ7 Ü.99G й.ЭВГ,

LI

о г 0.4 0-6 0.3

1

Рис. 6. Распределение плотности и скорости в задаче (22) для всех рассмотренных методов: О - численное решение, - аналитическое решение

Таблица. Максимальная погрешность решения рассмотренных задач, полученного различными методами

su, %

Схема Схема ВП Подход Лагранжа с Подход Лагранжа с

Годунова распределением (16) распределением (17)

Задача(20) 0.12 0.0046 0.0010 1.1

Задача(21) 57 59 0.0000 65

Задача(22) 0.060 0.0060 0.0028 0.55

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

В заключение можно сделать следующие выводы:

1) Метод, использующий подход Лагранжа с распределением (16), адекватно моделируют распространение волн в одномерных задачах нелинейной акустики в рамках применимости адиабатического закона.

2) Предложенный метод позволяет варьировать в широких пределах шагом по времени.

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

1. Годунов, С.К. Уравнения математической физики [Текст] / С.К. Годунов. - М.: Наука. - 1971. - 416 с. : ил. 2. Лощянский, Л.Г. Механика жидкости и газа [Текст] / Л.Г. Лойцянский. - М.: Наука. - 1978. - 736 с. : ил.

3. Годунов, С.К. Численное решение многомерных задач газовой динамики [Текст] / С.К. Годунов, А.В. Забродин, М.Я. Иванов, А.Н. Крайко, Г.П. Прокопов. - М.: Наука. - 1976. - 400 с. : ил.

4. Никонов, В.В. Моделирование двумерного ламинарного пограничного слоя с помощью метода "вихрь в ячейке" [Текст] / В.В. Никонов, В.Г. Шахов // Вестник СГАУ. - Самара. - 2006, № 3 (11). - С. 25 - 30

5. Ferziger, J. Computational methods for fluid dynamics [Text] / J. Ferziger, M. Peric. - 3 rev. ed. - SpringerVerlag. - 2002. - 423 p.

LAGRANGE APPROACH APPLYING TO THE SOLUTION OF ONE-DIMESIONAL PROBLEM OF WAVE PROPAGATION IN AIR UNDER CONDITION OF ADIABATIC LAW

© 2009 V.V. Nikonov, V.G. Shakhov

Samara State Aerospace University

Simulation of one-dimensional waves in air using of several numerical schemes is considered. In addition to known schemes for acoustics problem solution a scheme using Lagrange approach to flow modeling is proposed. The numerical solution results for three kinds of initial conditions are compared with solution of a characteristics method. Conclusion about applying of considered schemes is made. It is shown, that time step can be varied in the proposed scheme without decreasing of solution accuracy.

Keywords: Ideal gas, one dimensional wave, numerical method, acoustics problem, upwind discrete scheme, Godunov method, Lagrange approach, time step, exact solution, initial conditions

Valery Nikonov, Candidate of Technics, engineer NTP "Aviatechnokon". E-mail: [email protected] Valentin Shakhov, Candidate of Technics, professor, Head of Aerohydrodynamic Department. E-mail: [email protected]