CC BY
54
________УЧЕНЫЕ ЗАПИСКИ Ц А Г И
Том VIII 1977
№ 4
УДК 533,6,011.32-/5+518.61
ЧИСЛЕННОЕ МОДЕЛИРОВАНИЕ СЛОЖНЫХ ЗАДАЧ АЭРОГАЗОДИНАМИКИ МЕТОДОМ „КРУПНЫХ ЧАСТИЦ"
II. Асимптотика звуковых течений. Методические расчеты
О. М. Белоцерковский, Ю. М. Давыдов
Рассматриваются случаи звукового обтекания плоского клина и цилиндрического торца. Приводится сравнение результатов численных расчетов с асимптотиками Овсянникова — Франкля и Гудер-лея — Фальковича. Производится тестирование для контроля правильности постановки граничных условий: расчет на различных сетках, „отодвигание* границ расчетной области и др.
Данная работа является продолжением статьи [1]. Сйисок цитируемой литературы составлен отдельно.
Прежде чем привести результаты численных исследований различных задач газовой динамики, выпишем асимптотику для внешних звуковых течений в плоском и осесимметричном случаях [2], которая также используется в дальнейшем в качестве одного из контрольных тестов для оценки надежности полученных результатов на трансзвуковых режимах.
2.1 Рассмотрим для плоских потоков асимптотику типа решения Овсянникова — Франкля в случае обтекания клиновидного профиля. Аналитические выражения здесь можно получить лишь для одного (звукового) режима из трансзвукового интервала скоростей. Все нижеприведенные рассуждения относятся к областям, достаточно удаленным от тела, поэтому перейдем к уравнению для потенциала скорости ср. С. А. Чаплыгин показал [3], что потенциал скорости f и функция тока ф на плоскости годографа удовлетворяют системе линейных однородных уравнений с частными производными первого порядка. Такое же утверждение справедливо и для плоскости канонических переменных 6, г\ [4];
<?п — %В (у) фв, <Ръ = — В (т)) (2.1)
где В (^) — некоторая функция, удовлетворяющая следующим требованиям: В(у\)^> 0; В ("У})-»- 0 при •»] -* оо; В(щ) голоморфна при
7) = 0.
Исключая из (2.1) потенциал скорости ер, придем к уравнению Чаплыгина—Франкля для одной функции тока ф
+ = а (2-2) В первом приближении В (п) = В (0) = В0. Тогда (2.1) примет вид
9ч = В0'Чфв> сРе = -5о^. (2. Г)
а уравнение (2.2) перейдет в уравнение Дарбу — Трикоми
^еа + ■= °- (2.2)
Определяя главный член решения по Ф. И. Франклю [5] и опуская промежуточные выкладки [4], получим законы убывания скорости возмущения № — С, где Ш = [и, V), С —скорость звука.
Вдоль оси симметрии (6 = 0) и — с— |х|-1/20(1) и вдоль звуковой линии (т) = 0) и — с = |х|-3/40(1), V — с = ||—3/50(1) асимптотическое поведение звуковой линии вдали от профиля имеет вид
\у | = | х |5/40 (1). Аналогично ведет себя и предельная характери-
стика
На фиг. 2.1 результаты асимптотических оценок (штриховые линии) сравниваются для клиновидного профиля с расчетами при Моо= 1 методом „крупных частиц11 (сплошные линии). На фиг. 2.1, а приводятся формы звуковой линии с и предельной характеристики Ь, а на фиг. 2.1, б — 2.1, г показаны закономерности убывания компонент скоростей (фиг. 2.1, б—вдоль оси симметрии; фиг. 2.1, в, г — вдоль звуковой линии). Отсюда следует, что уже на расстояниях двух-трех радиусов от тела имеет место весьма удовлетворительное согласование асимптотики с расчетными данными.
В осесимметричном случае получим асимптотику типа решения Гудерлея — Фальковича для звукового обтекания полубесконеч-лого цилиндрического торца. Здесь имеем
(С2 д2 ? лич) д2 У I (С1 ,т,2\ д2у . д-р „
ус и,) дг2 гиъ дгдг ^ V) дг1 ^ г дг
~2~ (и2 + V2) + Ь = -у С1~\- /г*, и==~^~>
(2.3)
где /г — удельная энтальпия; звездочка обозначает критическое состояние.
Для уравнения (2.3) ставится задача Коши с данными на оси симметрии. Решение ищется в виде ряда
? = с*
I
(2.4)
5=--------к—- п = - • (2-5)
(2т*)1'3 гп 2Р3 с2 [ ди2 3 4 '
Главный член асимптотического решения поведения /Ш|) (?) (нулевое приближение), впервые полученный К. Г. Гудерлеем, С. В. Фальковичем и И. А. Черновым [6, 7], находится путем подстановки (2.4) в (2.3) и пренебрежением всеми членами с /ш<(£), где г> 1. В результате получается нелинейное дифференциальное уравнение для /Шо(1).
йГ~^ - п2 ^ + п (5я - 4) ? -(Зя - 2)2/»о = 0.
Фиг. 2.1
Из асимптотического разложения члена /ш„(£) при £ -> + со получим, что ю0 = (2/г — 1)/(2 —А); п— 1/(2 —&), где значение А зависит от формы образующей тела.
Подробные асимптотические выражения приводятся в [2]. Здесь следует лишь отметить, что для корректного решения нашей задачи надо учитывать и второй член разложения более высокого (второго) порядка малости /Ш1(£), который „подправляет“ решение для нулевого приближения. Следует отметить, что функция /Ш1(?) не удовлетворяет условиям аналитичности на предельной характеристике, поэтому берется ф 0 лишь в области, расположенной вниз по течению от фронта ударной волны. Рассматривая область течения вверх по потоку, мы можем, таким образом, сравнение численных результатов проводить по нулевому приближению /ш0(Е), где ш0 = — 2/7, «==4/7.
В трансзвуковом потоке можно выделить несколько особых линий: 1 — звуковая линия, 2—предельная характеристика, 3—линия, на которой вектор скорости направлен горизонтально, 4 — ударная волна [6]. Вдоль каждой из этих линий переменная к постоянна, поэтому, согласно (2.5), все они имеют вид г=сопз! г41. На фиг. 2.2, а для случая обтекания осесимметричного торца асимптотическая форма этих особых линий (штриховые линии) сравнивается с результатами численных расчетов при Моо= 1 (сплошные линии), а на фиг. 2.2, б приводятся профили изменяющихся величин—-скорости потока вдоль звуковой линии (асимптотика здесь такова: ъс = д<?1дг = сг~917, и0=^1, так как в безразмерном виде и2с-\-ъ2с= 1, а г>е — мало). Видно, что хорошее совпадение наблюдается, как и в плоском случае, уже на расстояниях трех-четырех радиусов от тела.
Результаты сравнений с асимптотикой для плоского и осесимметричного случаев говорят о надежности полученных численных результатов и о правильной трактовке граничных условий.
2.2. Приведем теперь некоторые результаты методических расчетов сложных течений газовой динамики, полученные с помощью метода „крупных частиц" (численный эксперимент). Вычисления проводились на машине БЭСМ-6, при этом область интегрирования разбивалась на 40 (по вертикали) и от 20 до 60 (по горизонтали) расчетных ячеек. В большей части расчетов использовалась схема, где для вычисления АМп применялись формулы первого порядка точности (1.5) из I этой работы [1]. Машинное время расчета варианта не превышало, как правило, одного часа.
Методический счет по разработанному методу проводился на примере задачи о внешнем обтекании полубесконечного осесимметричного торца. Эта задача достаточно сложна для построения расчетного алгоритма (бесконечный радиус затупления носовой части, наличие точки излома контура, образование ударных волн у тела и в поле течения и т. д.). Особый интерес вызывает расчет обтекания этого тела на трансзвуковых режимах и при малых сверхзвуковых скоростях, когда область возмущения за отошедшей ударной волной становится значительной по своим размерам [8].
На фиг. 2.3, а приводится положение головной и внутренней ударных волн и звуковых линий для разных чисел Маха сверхзвукового набегающего потока (1,1 <^Моо-< 14,5). Ударные волны здесь определялись как линии, на которых производная плотности по одному из пространственных направлений имеет максимум. Эти
а)
М=10
— асимптотика
Фиг. 2.2
ъ 11 2/ з/ /5/6/Х' с~
1 / X
1% Е 17
А II 1 |1Ш 1 1 1
г/Л
3
2 1 0 1 2 3 г/Л
2~М^=1,2; 3~М^=1,5-, Ч-М^2,0; 5-11^3,0-, 6-М^4,1; 7-М 5~
6)
расчет методом,, нрипныи частиц ”
• энсперимент А. Д. Губ'чика
Фиг. 2.3
сл
результаты хорошо согласуются с экспериментом (фиг. 2.3, бг Моо = 4,1, сплошная линия — расчет методом „крупных частиц**, пунктир — эксперимент А. А. Губчика), а также с расчетами в случае плоского течения по схемам Русанова, Лакса и Лакса — Венд-рофа [9, 10]. Следует заметить, что метод „крупных частиц** дает меньшее „размытие** ударной волны, а также достаточно точно описывает течение и в непосредственной окрестности угловой точки [10].
Как уже отмечалось, для трансзвуковых задач возникает вопрос о корректной постановке краевых условий на границах области интегрирования (этот вопрос детально исследовался в [11]). Наибольшую погрешность, как оказалось, вносит правая „открытая** граница СО (фиг. 2.3, а). Поэтому в программе был предусмотрен ряд внутренних контрольных тестов, обеспечивающих надежность получаемых результатов.
Так, например, практиковались расчеты на разных сетках аппроксимации, а также „склейка** полей течений (какой-либо внутренний столбец ранее рассчитанного случая использовался в качестве начальных условий для нового поля и расчеты продолжались вниз по потоку вправо). По согласованию результатов в зоне „перекрытия** нового и старого полей можно судить о влиянии граничных условий (фиг. 2.4).
Далее, при расчете обтекания полубесконечного цилиндра (фиг. 2.3, а) вычисления проводились для разных видов граничных условий на СИ и при различных длинах цилиндра (тело как бы „вдвигалось** в область интегрирования) до полного установления параметров потока в интересующей нас зоне. Как правило, небольшие возмущения на „правой** границе области СЕ) затухают уже на расстоянии нескольких (двух-трех) счетных слоев.
При этом для случая обтекания полубесконечного тела условия невозмущенного течения выполняются с достаточной точностью на левой (—со по г) и верхней (+оо по г) границах области интегрирования, а при г + оо (правая граница) в окрестности тела реализуются условия равномерности потока.
На фиг. 2.5 показано изменение профилей давления по г (поперек области) для равноотстоящих сечений АА, ВВ и СС, полученных при Моо = 1 для очень большого поля на сравнительно грубой расчетной сетке (расстояние по оси абсцисс отнесено к диаметру тела). Видно, что при г -> + оо, др/дг -+0 и др/дг -> 0 (включая и пристеночный слой), что и обеспечивает параллельность потока поверхности тела на +оо по г.
Как говорилось в I ч. этой статьи [1], для единообразия вычислений вдоль всех границ области интегрирования вводились слои „фиктивных** ячеек. При условии „непротекания** на поверхности
Р
СС
о
АА ЛВС — М^1
о
1
2.
10 15 20 г/В
Фиг. 2.5
Фиг. 2.4
тела в нашей схеме нормальная компонента скорости меняет знак (в терминах численной схемы это запишется так: ип+\ = ик, г>лг+1 = = — Ъы, где N—граничная к телу ячейка, а N + 1—соответствующая ей фиктивная ячейка внутри тела), а при условии „прилипания" обе компоненты скорости изменяют знак (Ндг+1 = — 'УЛГ+1 = — 1)ц). Отсюда следует, что на границе этих двух ячеек (поверхности тела) имеет место разрыв значений скорости. В обоих случаях аппроксимационная вязкость „сглаживает" разрывы в граничных условиях. Если теперь построить по вертикали к телу профили скоростей (помня, что значения всех параметров в ячейке относятся к ее центру), то нетрудно установить, что для обоих типов граничных условий на поверхности тела реализуются значения скоростей, близкие к нулю (условия „прилипания") [12]. Этот факт, наряду с наличием в разностных уравнениях диссипации и объясняют, по-видимому, причины образования у поверхности тела развитого пограничного слоя. По существу, появление эффектов трения на поверхности тела (а следовательно, и „пограничного" слоя) связано здесь с наличием вязкости в в потоке и существованием градиента скоростей, вызванного неравномерностью поля в окрестности тела*.
На фиг. 2.6 приведены вдоль сечений профили плотности для случаев так определяемых здесь условий „непротекания“ (сплошная линия) и „прилипания11 (пунктир). Видно, что даже вблизи тела различие между этими двумя случаями несущественно; по мере же удаления от тела оно вообще исчезает.
Для оценки влияния „откры-тых“ границ области интегрирования проводились методические расчеты на трансзвуковых режимах [13]. На фиг. 2.7 показаны результаты „закритического" обтекания полубесконечного торца при Мс;о = 0,9 для различных отношений его длины I к радиусу И.
Если поле течения перед телом устанавливается достаточно быстро и в рассматриваемых случаях на расстояниях до ^-1,5/? слева от среза (г<0) оно практически не меняется, то течение правее угловой точки (г->0) стабилизируется лишь при///?~2 — 4. При дальнейшем „вдвигании11 торца начинают сказываться краевые условия левой границы.
На фиг. 2.7, г проводятся также данные вычисления на боль-
1 ~
* При отмеченной трактовке краевых условий здесь наблюдается аналогия со струей в образовании пограничного слоя.
М=2
Фиг. 2.6
М^=0,9, а)1/Л = 0,5'е, (Г)1/Я=2,0, в)1/Л = 2,72; г) 1/11=7 Фиг. 2.7
шей расчетной сетке, где пунктиром отмечена обычно используемая область расчета (60 ячеек по горизонтали и 40 по вертикали), которая здесь является уже внутренней. Сравнение обоих решений в зоне угловой точки (зона „перекрытия" расчетных областей) указывает на их хорошее совпадение на всех режимах обтекания. В частности, давление торможения в критической точке 0 (область „медленного" установления), полученное из вычислений, отличается от точного значения при 0,2<Мсо<0,9 меньше, чем на один процент.
Аналогичное тестирование производилось и для нестационарных задач газовой динамики (диффракционных [14], взрывных [15]).
Эти методические расчеты позволили определить оптимальный размер расчетной сетки и число узлов аппроксимации (обычно используется не более 2,5 тыс. расчетных ячеек).
В следующей статье этого цикла будут приведены результаты численных исследований различных задач газовой динамики (численный эксперимент).
ЛИТЕРАТУРА
1. Белоцерковский О. М., Давыдов О. М. Численное моделирование сложных задач аэрогазодинамики методом „крупных частиц". Ч. I. Метод. Исследование схем. „Ученые записки ЦАГИ", т. 8, № 3, 1977.
2. Давыдов Ю. М. Метод „крупных частиц" для задач газовой динамики. Канд. диссертация, М., МФТИ, 1970.
3. Ч а п л ы г и н С. А. О газовых струях. Диссертация, 1902, изд. ГИТТЛ, 1949.
4. Овсянников Л. В. О движении клиновидного профиля со скоростью звука. Труды ЛК ВВИА, 33, 1950.
5. Ф р а н к л ь Ф. И. Исследования по теории крыла бесконечного размаха, движущегося со скоростью звука, ДАН, 1947, 57, № 7.
6. Гу дер лей К. Г. Теория околозвуковых течений. М., Изд. иностр . лит., 1960.
7. Фалькович С. В., Чернов И. А. Обтекание тела вращения звуковым потоком газа. Прйкл. матем. и механ., т. 28, вып. 2, 1964.
8. Давыдов Ю. М. Разработка нестационарного метода .крупных частиц* и расчет обтекания цилиндрического торца на трансзвуковых и сверхзвуковых режимах. Труды Казанского авиационного института. 1971.
9. Emery A. F. An evaluation of several differencing methods for inviscid fluid flow problems. J. Comput. Phys., vol. 2, N 3, 1968.
10. Белоцерковский О. М., Давыдов Ю. М. Расчет трансзвуковых течений методом «крупных частиц*. В сб. «Численные методы механики сплошной среды, т. 1, № 6, Новосибирск, 1970.
11. Belotserkovskii О. М., Davydov Ju. М. Numerical approach for some transsonic flows. Lecture Notes in Physics, vol. 19, 1973.
12. Давыдов Ю. М. Численное исследование некоторых переходных явлений. В сб. .Асимптотические методы в теории систем*, вып. 6, Иркутск, 1974.
13. Белоцерковский О. М., Давыдов Ю. М. Расчет и исследование характеристик сложных задач газовой динамики методом .крупных частиц". Препринт ВЦ АН СССР и МФТИ, 1972.
14. Давыдов Ю. М. К расчету нерегулярного отражения ударных волн методом .крупных частиц*. Труды МФТИ, серия .Аэромеханика, процессы управления*, 1973.
15. Давыдов Ю. М., Шевырев С. П. Расчет некоторых взрывных задач методом .крупных частиц". Аэродинамика, вып. 4 (7), изд. Саратовского университета, 1975.
Рукопись поступила 13/Х 1975
