CC BY
102
Вычислительные технологии
Том 17, № 3, 2012
Представление приближенных решений
полной системы уравнений Навье — Стокса
>к
в одномерном случае*
С. П. Бдутин, В.Е. Замыслов Уральский госудаpственный университет путей сообщения, Екатеринбург, Россия e-mail: SBautin@math.usurt.ru, VZamislov@math.usurt.ru
Рассматривается полная система уравнений Навье — Стокса, решения которой описывают течения сжимаемого вязкого теплопроводного идеального газа при постоянных значениях коэффициентов вязкости и теплопроводности. От традиционного использования плотности и температуры в качестве независимых термодинамических переменных осуществляется переход к удельному объему и давлению, через которые записываются уравнения состояния и полная система уравнений Навье — Стокса. Этот переход позволяет записать систему с частными производными в нормальном виде относительно производных по времени и с полиномиальной правой частью. В случае одномерных плоскосимметричных течений решения выписанной системы уравнений строятся с использованием идей метода Галёрки-на в виде бесконечных сумм от гармоник по пространственной переменной. Для данного представления на границах отрезка пространственной переменной выполняются условия прилипания и теплоизоляции. Коэффициенты бесконечных сумм есть искомые функции, зависящие от времени. Для них выписана бесконечная система обыкновенных дифференциальных уравнений. При учете конечного числа гармоник численно построены решения соответствующей конечной системы обыкновенных дифференциальных уравнений. В частности, описано движение бегущего ударного перехода, не имеющего существенных осцилляций в районе его фронта.
Ключевые слова: сжимаемый вязкий теплопроводный газ, одномерные течения, полная система уравнений Навье — Стокса, приближенные решения, метод Галёр-кина.
Введение
Использование для описания движения несжимаемой вязкой теплопроводной среды традиционной системы уравнений Навье — Стокса [1] (при постоянной плотности р = const) приводит к тому, что для построенных решений этой системы не выполняются (см. [2]) либо закон сохранения энергии при отсутствии в системе уравнения для температуры, либо основное термодинамическое тождество при включении в систему уравнения для температуры.
Для описания течений сжимаемого вязкого теплопроводного идеального газа [3] необходимо исследование решений полной системы уравнений Навье — Стокса (ПСУНС) [1-4] и в качестве независимых термодинамических переменных можно использовать
* Работа выполнена при финансовой поддержке РФФИ (грант № 11-01-00198).
плотность р и температуру Т [2]. В этом случае запись системы уравнений с частными производными в нормальном виде относительно производных по времени приводит к тому, что в правых частях уравнений движения и уравнения энергии плотность р стоит в знаменателе, что, естественно, затрудняет применение аналитических подходов для представления решений ПСУНС.
В [5-7] в качестве независимых термодинамических переменных берутся удельный объем 8 = 1/р и давление р. Тогда при записи ПСУНС в нормальном виде относительно производных по времени правые части всех уравнений имеют полиномиальный вид с квадратичными нелинейностями. Это дает возможность строить приближенные решения ПСУНС с использованием идей метода Галёркина [8].
В настоящей работе предлагаемый подход к построению приближенных решений ПСУНС демонстрируется в случае одномерных плоскосимметричных нестационарных течений.
1. Полная система уравнений Навье — Стокса
В случае постоянных значений у, к — коэффициентов вязкости и теплопроводности, а также равенства нулю второго (объемного) коэффициента вязкости: у' = 0 — для идеального газа с уравнениями состояния
p = RpT, e = cvoT, R, cvo — const > 0,
где e — внутренняя энергия, при одномерных плоскосимметричных течениях ПСУНС в безразмерных переменных имеет вид [5-7]
8t = 8ux — u8x,
ut = —uux — 18px + уо 8uxx, (1)
pt = —upx — Ypux + ко (8p)xx + yoY(Y — l)u2
x
(здесь t — время, x — пространственная переменная, u — скорость газа), а уравнения состояния становятся следующими:
T = 8p, e = T.
В системе (1) и в уравнениях состояния при введении безразмерных переменных масштабные значения р00, T00, x00, u00, t00, p00, e00 связаны соотношениями [2]
uoo = \J R7T00, t00 =-, p00 = RP00T00, e00 = CvoT»,
u00
где 7 = 1 + R/cvo. Из приведенных равенств следует, что за масштабные значения скорости, давления и внутренней энергии выбраны соответственно скорость звука, давление и внутренняя энергия газа с параметрами р00, T00. При этом значения положительных констант у0 и к0
4 у к
У0 =о-, К0 =-.
3 P00u00x00 CvoP00u00x00
При исследовании течений вязкого теплопроводного газа часто вводят числа Рей-нольдса и Прандтля
P00u00x00 ^ CvoYy
Re =-, Pr =-.
ук
Тогда ^о = 4/(3И,е), Рг = 37^0/(4к0). Для воздуха обычно полагают Рг = 0.72, 7 = 1.4, следовательно, ко ~ 1.458333д0.
Система (1) смешанного типа и состоит из трех уравнений: первое — уравнение неразрывности — образует гиперболическую часть системы, так как определяет в течениях сжимаемого теплопроводного вязкого газа наличие слабого разрыва на контактной поверхности [9]; второе и третье — уравнения движения и энергии — составляют параболическую часть системы, поскольку содержат вторые производные скорости и температуры по пространственной переменной.
2. Представление решения полной системы уравнений Навье — Стокса
Используя идеи метода Галёркина [8], решение системы (1) строится в виде
те
8(t, x) = 1 + ^^ 8k (t) cos kx, k=l
X
u(t, x) = uk(t) sin kx,
k=l
(
p(t, x) = 1 + ^^ pk (t) cos kx.
(
,x) = 1 + ^ Pk (t) cos I k=0
Для приведенных представлений (2) при x = 0 и x = п автоматически выполняются условия прилипания и теплоизоляции
u|x=0, x=n = 0, Tx|x=0,x=n = °. (3)
Чтобы получить уравнения для коэффициентов 8k (t), Uk(t), pk(t), соотношения (2) подставляются в (1) и каждое уравнение полученной системы проецируется на свою систему базисных гармоник, а именно, умножается соответственно на cos tx, sintx и cos tx (t = 1, 2,... для всех трех уравнений и, кроме того, t = 0 для третьего уравнения) и интегрируется на отрезке [0, п]. В результате получается следующая бесконечная система обыкновенных дифференциальных уравнений (СОДУ) для бесконечного числа
искомых коэффициентов 8¿(t), ui(t), p0(t), pi(t), t = 1, 2,...:
2(
8'f (t) = tui(t) + - У^ (makme + kbkmi) 8k(t)um(t), (4)
п
k,m=l
2 ( 1 2 (
ui(t) =--mbkimuk(t)um(t) + -tpi(t) +--mbmik8k(t)pm(t) - ^0t24f(t)-
п ^ Y Yn
k,m=l k,m=l
2(
-^0- У^ m2bmik 8k (t)um (t), п
k,m=l
1 ( 1 (
p0 (t) = 1(1 - Y) ^ kuk(t)pk(t) + -m(Y - 1) ^ k2uk(t),
k=l k=l
2 те
p'e(t) = 2 У^ (mbkme - lkakme) uk(t)pm(t) - y[1 + po(t)]£ue(t)-
k,m=l
2 те
-Ko£2 {[1 + po(t)]5¿(t) + pe(t)} - Ко - [(m2 + k2)akme - 2kmbkme] 5k (t)pm(t) +
k,m=l 2 ^
+Уо Y (y - 1)- > kmakmeUk (t)um(t).
IT ' J
-
k,m=l
Здесь
n n
akm£ = J cos kx cos mx cos ixdx, bkm£ = J sin kx sin mx cos ixdx, oo а константы k, m, € независимо друг от друга принимают целые положительные значения: 1, 2, ... Начальные данные для системы (1) также представлены в виде, аналогичном виду (2):
те
£0 ,
S(t , x) |t=o — 1 + ^ ] Sk cos kx k=l
те
u(t, x) |t=0 = uk sin kx,
k=l
те
p(t, x) |t=o = 1 + Y1 pk cos kx,
k=o
где 500, uk, pk — заданные константы:
4 (t)|t=0 = Uk (t)|t=0 = uk, Pk (t) |t=0 = pk, к = 1, 2,..., po(t)lt=o = p0. (9)
Решение начально-краевой задачи (1), (3), (8) описывает процесс стабилизации одномерного течения при t ^ в области 0 < x < п от начального неоднородного состояния (8) к состоянию однородного покоя.
В работах [10, 11] доказано, что при определенных условиях на начальные данные поставленная начально-краевая задача для ПСУНС имеет в соответствующем функциональном пространстве единственное решение. В [5-7] при использовании на отрезке [—п, +п] полных рядов Фурье для представления периодических по пространственным переменным решений задачи Коши в общем трехмерном нестационарном случае доказано существование решения в малом по времени при определенной малости начальных условий, а в одномерном плоскосимметричном случае при условии, что с ростом времени u ^ u* = const = 0, доказано существование решения при всех t > 0. Это решение с ростом времени стабилизируется к однородному потоку.
В [12, 13] (см. также [2]) с использованием методики представления решений характеристической задачи Коши в виде рядов по степеням формальным образом введенного малого параметра построены приближенные решения ПСУНС, в частности, периодические по пространственной переменной.
Теоремы о существовании решений у начально-краевой задачи (1), (3), (8) в случае использования представлений (2) при условии сходимости представлений (8) авторам настоящей работы не известны.
Для построения приближенного решения рассматриваемой задачи можно ограничить в (2) и (8) число гармоник, например, положив во всех суммах в качестве верхнего предела суммирования конечное число К. В этом случае СОДУ (4)-(7) станет конечной, состоящей из (3 х К + 1) уравнений для такого же числа искомых функций 8 к(¿), пк (¿), Ро(^), Рк (¿), к = 1, 2,..., К. Естественно, что решать получившуюся задачу Коши (4)-(7), (9) можно только численно.
3. Результаты расчетов
Проверка точности получаемых приближенных решений проводилась стандартным образом по следующим направлениям.
Во-первых, при изменении значений начальных данных с помощью констант (9) выяснилось, что при достаточно больших по модулю значениях этих констант решение строится только на достаточно ограниченном промежутке времени. В качестве примера подобной ситуации на рис. 1, а показаны графики коэффициентов pk(t), 1 < k < 5, при K = 50, y = 1.4, = 0.001, ко = 1.4583^0, ¿|t=o = 1, u|t=o = 0, p|i=o = 1 + 0.7cosx. Пример поведения коэффициентов pk(t), 1 < k < 5, в случае, когда при t ^ течение стабилизируется к состоянию однородного покоя, представлен на рис. 1, б. Здесь все входные данные при решении СОДУ были такими же, как в предыдущем примере, за исключением начального распределения давления — в данном случае p|t=0 = 1 + 0.1 cos x.
Pk(f) 0.5 0
-0.5 -1 -1.5
/5 .......
2
0
0.5
1.5
2.5
10 20 30 40 50 60 70 80 90 t
100 200 300 400 500 600 700 800 900 t
a
Рис. 1. Зависимость первых пяти коэффициентов разложения давления от времени при их достаточно больших (а) и достаточно малых (б) по модулю начальных значениях
Во-вторых, для тех значений начальных констант (9), при которых решение восстанавливалось фактически до однородного покоя, осуществлялась проверка закона сохранения массы: значение интеграла
должно быть постоянным. Для варианта, приведенного на рис. 1, б, график функции относительной погрешности
при вычислении массы в зависимости от времени, представлен на рис. 2, из которого следует, что погрешность выполнения закона сохранения массы в данном варианте составляет величину < 1 %.
В-третьих, для проверки точности получаемых решений проводились расчеты с разным числом гармоник, используемых в представлении (2): К = 10, 20, 30, 40, 50. При увеличении числа слагаемых в (2), т. е. при увеличении числа К, решения равномерно и по х, и по Ь сходятся к некоторому конкретному решению. Этот вывод иллюстрируют на рис. 3 графики р(Ь,х) в момент времени Ь = 33 при К = 10, 20, К = 20, 30, К = 30, 40, К = 40, 50. Сравнение проводилось для входных данных варианта, представленного на рис. 1, б, на временном промежутке 0 < Ь < 100 и для всех х Е [0,п]. При этом максимальные величины относительных различий значений давления для указанных парных значений К получились соответственно следующими: 4.33, 2.05, 1.12 и 0.56%.
Выявленное различие значений газодинамических параметров, рассчитанных при разных К, связано с принципиальным, на наш взгляд, свойством решений СОДУ (4)-(7) (естественно, при суммировании не до бесконечности, а до конечного числа К). Так, пусть в начальный момент времени модули значений констант (9) при гармониках с малыми номерами больше, чем при гармониках с большими номерами. Например, при 1 < к < П\ начальные условия для 8 к, Щ, Рк — не нули, а при П\ < к < К — нули. Тогда
о
а, %
0.9
0.8
0.7 0.6 0.5 0.4 0.3 0.2 0.1
HH.HWUUU.UVU.VVV.VV.V1
О 10 20 30 40 50 60 70 80 90 г
Рис. 2. Относительная погрешность выполнения закона сохранения массы для одного конкретного решения
1.2 1.1 1
0.9
0.
1.2 1.1 1
0.9 0.
К = 10, 2С
0 0.5 1 1.5 2 2.5 3 V
К = 30, 4( )
f
/
О 0,5 1 1.5 2 2.5 3
1.2 1.1
1
0.9
X 0.:
1.2 1.1 1
0.9 .г 0.8
К = 20, 30
О 0.5 1 1.5 2 2.5 3 Р
II Й 40, 5 Э
/
О 0.5 1 1.5 2 2.5 3
х
Рис. 3. Распределения давления в конкретный момент времени при разном числе К
с ростом времени происходит передача начальных возмущений от начальных гармоник на гармоники с большими частотами: модули их амплитуд начинают расти. При этом данный процесс происходит достаточно равномерно по времени. Пусть общая стабилизация потока к однородному покою идет относительно долго (время порядка 1/у0). Тогда амплитуды гармоник с последними номерами (близкими к K) становятся соизмеримыми с амплитудами низкочастотных гармоник. Это является результатом того, что последним гармоникам некуда передавать возмущения: в используемом представлении (2) нет гармоник с частотами, большими числа K. В таком случае используемые представления хуже передают искомые течения. Данный факт иллюстрируют результаты расчета при K = 50, y = 1.4, у0 = 0.0001, к0 = 1.4583у0, í|t=0 = 1, u|t=0 = 0.1 sinx, p|t=o = 1. На рис. 4,1 для момента времени t = 29 приведены распределение скорости газа как функция от x и значения модулей коэффициентов для скорости газа.
В итоге получается, что для конкретного значения у0 найдется такое число K0 — частота последней используемой гармоники, что общий процесс стабилизации потока произойдет раньше, чем значительные возмущения достигнут гармоник с частотами, близкими к K0. Следовательно, невозможность передачи возмущений на гармоники с отсутствующими частотами, большими K0, не приведет к искажению решения. Этот факт иллюстрируют результаты расчета варианта, приведенные на рис. 4, II, отличающегося от предыдущего только значением у0 = 0.001. А если для заданного у0 число K используемых гармоник существенно меньше требуемого K0 (например, K = K0/2), то невозможность процесса передачи возмущений с последних используемых гармоник на отсутствующие высокочастотные гармоники приведет к искажению решения (см. рис. 4, I).
Заметим, что если ненулевые начальные данные будут только у гармоник с номерами U\ < k < n2, U\ > 1, то распространение возмущений и рост модулей амплитуд наблюдаются не только у гармоник с номерами k > n2, но и у гармоник с уменьшающимися частотами от k = n — 1 до k = 1. Такое перераспределение начальных возмущений от одних гармоник на другие в случае многомерных течений будет приводить к достаточно сложным картинам течения даже при соответствии начальных данных простым течениям.
0.2 0.1 0
-0.1 0.2
0.024 0.016 0.008
О 0.5 1 1.5 2 2.5 3
х
»1-1
I
In iiiitiii mnt^jvli.jr п.. JItfL ll
О 10 20 30 40 50
0.024 0.016 0.008
|
i ll г
Шкшл fB 11 ..
О 0.5 1 1.5 2 2.5 3
О 10 20 30 40 50'
б
Рис. 4. График скорости (а) и модули коэффициентов её разложения (б) в момент времени г = 29 при ^о = 0.0001 (I) и ^о = 0.001 (II)
Рис. 5. Распределения давления в разные моменты времени для одного конкретного решения
I
а
Наконец, еще один важный вывод, следующий из результатов расчетов. Некоторые начальные возмущения с течением времени трансформируются в один или в несколько бегущих ударных переходов, постепенно затухающих во времени. При этом имеются варианты, когда несколько ударных переходов почти не догоняют друг друга, и варианты, когда они догоняют друг друга, сливаясь в один бегущий ударный переход. На рис. 5 приведены распределения давления для течения с входными параметрами K = 50, Y =1.4, ^о = 0.001, к0 = 1.4583^о, ^|t=o = 1, u|t=0 = 0, p|t=0 = 1 + 0.1 cos x — 0.1 cos 3x в моменты времени t = 0, 10, 25, 75.
Выводы
1. В одномерном плоскосимметричном случае построены приближенные нестационарные решения начально-краевой задачи для ПСУНС. Соответствующими вычислительными приемами и сравнением с другими известными приближенными решениями [2, 12, 13] подтверждено, что решения, построенные в случае относительно небольших начальных данных, с достаточной точностью передают одномерные течения сжимаемого вязкого теплопроводного газа.
2. Предложенный подход позволяет получать решения на достаточно больших промежутках времени, в том числе, когда течение практически стабилизируется — в данном случае к однородному покою. Расчеты показали, что время стабилизации потока составляет ~ 1/у0.
3. Приведены примеры достаточно нетривиальных решений.
4. Установлен факт перераспределения начальных возмущений с одних гармоник на другие и сделаны предположения о механизме такого перераспределения через кратные частоты.
5. Очевидно, что применение данной методики возможно и в неодномерном случае, и при других финальных (т. е. при t ^ распределениях газодинамических параметров. Однако в общем трехмерном нестационарном случае соответствующие СОДУ будут содержать принципиально большее число уравнений: ~ 5 х K3.
Список литературы
[1] Лойцянский Л.Г. Механика жидкости и газа. М.: Наука, 1970.
[2] Баутин С.П. Характеристическая задача Коши и ее приложения в газовой динамике. Новосибирск: Наука, 2009.
[3] Рождественский Б.Л., Яненко Н.Н. Системы квазилинейных уравнений и их приложения к газовой динамике. М.: Наука, 1968.
[4] Ковеня В.М., Яненко Н.Н. Метод расщепления в задачах газовой динамики. Новосибирск: Наука, 1983.
[5] Титов С.С. Пространственно-периодические решения полной системы Навье — Сток-са // Докл. АН. 1999. Т. 365, № 6. С. 761-763.
[6] Титов С.С. Решение нелинейных уравнений в аналитических полиалгебрах. I // Изв. вузов. Математика. 2000. Т. 1(452). С. 66-77.
[7] Титов С.С. Решение нелинейных уравнений в аналитических полиалгебрах. II // Там же. 2000. Т. 6(457). С. 45-52.
[8] Флетчер К. Численные методы на основе метода Галёркина. М.: Мир, 1988.
[9] Баутин С.П. Представление решений системы уравнений Навье — Стокса в окрестности контактной характеристики // Прикл. математика и механика. 1987. Т. 51, вып. 4. С. 574-584.
[10] Антонцев С.Н., Клжихов А.В., Монлхов В.Н. Краевые задачи механики неоднородных жидкостей. Новосибирск: Наука, 1983.
[11] КАжихов А.В. Избранные труды. Математическая гидродинамика. Новосибирск: Ин-т гидродинамики СО РАН, 2008.
[12] Блутин С.П. Об одном классе решений системы уравнений Навье — Стокса, периодических по пространственным переменным // Моделирование в механике. 1987. Т. 1(18), № 1. С. 3-15.
[13] Блутин С.П. Аналитическое построение течений вязкого газа с помощью последовательности линеаризованных систем Навье — Стокса // Прикл. математика и механика. 1988. Т. 52, вып. 4. С. 579-589.
Поступила в 'редакцию 6 февраля 2012 г., с доработки —12 марта 2012 г.
