Anderson-et-al-1 (1185923), страница 57
Текст из файла (страница 57)
9 6.3. Методы сквозного счета Методы сквозного счета широко применяются для расчета течений невязкого газа с ударными волнами. В этом случае уравнения Эйлера записываются в дивергентной форме и для расчета скачков уплотнения или любых других разрывов не принимается никаких специальных мер, они получаются как часть решения всей задачи.
Рассчитанные этими методами скачки <размазываются» на несколько ячеек сетки, но простота этого подхода может оказаться решающим фактором при сравнении со ехемал1и с выделением скачков, несмотря на то что последние дают несколько более точные результаты. При расчетах течений методами сквозного счета часто пользуются выделением скачка на границе. Поскольку на границах скачки могут быть выделены по любой из обычных схем, обсуждаемых в настоящем параграфе и в 9 6.4, то на самом деле преимущество схем сквозного счета становится еще более ощутимым, когда имеется сложная структура'скачков. В этом случае удается уловить структуру скачков и отпадает необходимость в специальной обработке каждого скачка. Лаке [(.ах, 1954) показал, что скорость и интенсивность скачка рассчитываются верно, если уравнения Эйлера записаны в дивергентной форме.
Это означает, что физически корректное слабое решение, соответствующее соотношениям Гюгонио — Рэнкина на скачках, получается в случае, когда используется дивергентная форма записи уравнений Эйлера и последние дискретизируются консервативным образом. В 9 6.4 мы получим сла- э 6.3. Методы сквозного счета 317 (6.23) Для стационарного изоэнергетического течения полная энтальпия сохраняется вдоль линий тока.
В этом случае уравнение энергии не интегрируется, а используется аналитическая интегральная форма иа+ оа Н= ~ ~+ = сопзй у — 1 р 2 (6.25) Уравнения (6.22), (6.24) при постоянной полной энтальпии образуют в случае сверхзвукового течения гиперболическую систему.
Ее можно решать маршевым методом по координате х, т. е. последовательно интегрируя уравнения вдоль направления х, начиная с поверхности задания начальных данных. Геометрия такой задачи показана на рис. 6.6. Начальные данные задаются на линии х = 0 и решение находится шаг за шагом при продвижении по координате х с учетом граничных условий на поверхности тела и соответствующих условий при у,„. Уравнения (6.22) — (6.24) можно записать в виде дЕ дг — + — =О, дх ду (6.26) где Е = р+ ри' Г = рио Уравнение (6.26) можно решать любым методом, пригодным для гиперболических уравнений в частных производных из числа бое решение уравнений Эйлера, используя их недивергентную форму записи, и покажем, что в этом случае разрыв решения с физической точки зрения уже не интерпретируется как ударная волна.
Таким образом, вид получаемого решения зависит от формы представления решаемых уравнений. Дивергентную форму записи проиллюстрируем на примере сверхзвукового обтекания совершенным газом некоторой двумерной поверхности. Если считать, что ось х совпадает с поверхностью тела и эта координата является маршевой, то уравнения, описывающие это течение, являются двумерной версией уравнений (5.192) и их можно записать в виде — + — =О, дри дро дх ду д„(Р+Рй)+ д (Рио)=О, д е д — (рио) + — (р+ ро') = О. (6.24) 318 Гл.
6. Численные методы решения уравнений течения представленных в гл. 4. Схема Мак-Кормака будет хорошим выбором. Вариант схемы Мак-Кормака с разностями вперед на шаге предиктор и с разностями назад на шаге корректор в применении к уравнению (6.26) может быть записан в виде Е1+' = Е1 — — (Р1+! — Р1) (6.27) Е1+ = 2 ~Е1 + Ет+ — ь (Р1 +' — Р1+!')(. ар По завершению шагов предиктор и корректор от вектора Е следует перейти к примитивным переменным р, р, и, п. Только после шлрааланпе Поаерхлосшь тела Рис. 6.6.
Система координат для маршевой задачи. этого можно составить новый вектор потока для следующего шага интегрирования. После получения решения на каждом шаге по маршевой координате сразу находим компоненту у скорости: о = Ез/Еь где нижние индексы обозначают элементы вектора Е. Для отыскания компоненты х скорости необходимо решить квадратное уравнение.
Для исключения р можно скомбинировать Е, с уравнением энергии; тогда имеем Ез ит + [(у — 1)/2у) (2Н вЂ” иа — с') Исключим р, используя равенство р = Е!/и. (6.28) Это приводит к квадратному уравнению для и с корнями и = —, — +. т т ~ — — ~ — — (20 — р ). (6.29) у Ея Iт' у Еа Хя у — 1 2 у+! Еь Ъ'чу+1 Е1 ) у+1 3!9 % 6 3. Методы сквозного счета Обычно радикал следует брать со знаком плюс.
Плотность теперь рассчитывается через Е, по формуле (6.28), а давление выражается через Ет в виде Р = Е, — риа. (6.30) Проделав эти вычисления, можно вычислить новое значение Г и перейти к следующему шагу интегрирования. Оно приводит к расчетной сетке, изображенной на рис. 6.8. Теперь задача обтекания клина решается без особого труда, поскольку линии постоянства т1 растут линейно с х. Так как основ- Пример 6.2. 1зассчитайте обтекание двумерного клина (полу- угол при вершине равен 15'), движущегося со скоростью М = 2. 1-ассмотрите течение невязкого совершенного газа. Решение. В задаче требуется определить положение скачка уплотнения и его интенсивность, а также картину течения вблизи клина.
Клин и картина его обтекании показаны на рис. 6.7. Обтекание двумерного клина с присоединенной ударной волной является коническим течением. Это значит, что парамет- Эй =2 15' ры потока вдоль лучей, исходящих из вершины клина, постоянны (Апбегзоп, 19821. Это существенно облегчает задачу. Определяющими уравнениями задачи являются уравнения (6.22) — (6.24) и уравнениеэнер- рис. вд, Оотекание клина с нрисогии (6.25), граничными уело- единенной ударной волной. виями — условие скольжения на поверхности клина и условия в свободном потоке вне ударной волны.
Конечно, можно направить ось х вдоль поверхности клина и интегрировать уравнения шаг за шагом вдоль этой оси до тех пор, пока число Маха в ударном слое остается больше единицы. К сожалению, по мере продвижения вниз по потоку ударный слой утолщается и в конце концов это приводит к взаимодействию внешней границы нашей расчетной области (при у = = у „) с ударной волной. Эти трудности можно легко преодолеть, если учесть, что ударная волна прямая, а толщина ударного слоя растет линейно с увеличением х.
Выполним преобразование переменных й=х, т(=у/х. (6.31) 320 Гл 6 Численные методы решения уравненнб течения ныс уравнения являются гиперболическими относительно координаты ~, то начальные данные необходимо задавать на некоторой линии, не являющейся характеристической. В качестве последней удобно выбрать линию $ = !. Уравнения в частных производных интегрируются в направлении $ с произвольно заданными начальными условиями. Так как течение вокруг двумерного клина является коническим, то решение для него будет достигаться при больших ~(асимптотически). вопя! вб =2 Рис. 6.8. Клин с ударным слоем после преобрааования.
После преобразования координат от (х, у) к ($, а)) уравнения принимают вид — + — =О, дЕ дР (6.32) д$ дч где Е = $Е, Р = Р— т!Е. Можно избежать дополнительных трудностей, если в этой задаче использовать свойства конического течения. Устойчивость схемы, применяемой для решения уравнений (6.32), зависит от собственных значений матрицы [А) расширенной системы, записанной в координатах Я,т)): — + [А] — + Н =О. дв дп (6.33) В этом уравнении вв есть вектор примитивных переменных, а Н вЂ” источниковый член в расширенной системе.
Оказывается, что собственные значения [А) явно зависят от ~, т. е. ~ входит в выражения для собственных значений. Когда при решении задачи мы продвигаемся вниз по потоку вдоль маршевой координаты $, допустимый размер шага должен изменяться, если применяется явная схема, например схема Мак-Кормака. Если при увеличении 8 этого не происходит, могут возникнуть трудности с устойчивостью разностной схемы.
Их можно избежать, 321 3 6.3. Методы сквозного счета если интегрировать уравнения от й = 1 до $ = 1 + Л$ путем итераций до достижения сходимости. Тщательного рассмотрения требуют граничные условия. Количество точек в направлении т1 должно быть достаточным, чтобы ударная волна могла формироваться естественным образом и не взаимодействовала с граничными условиями в свободном потоке, которые ставятся при т1 = т1,„.
Например, если в нашей задаче ударная волна образует с поверхностью клина угол в 20' и в ударном слое мы размещаем 10 точек, то т1,,„= 1п (20') = 0.3640, Ьт1 = 0.3640/(10 — 1) = 0.0404. Пусть при том же шаге мы добавляем еще пять точек; тогда т1 „= О. 0404 (15 — 1) = 0.5662 Прежде чем мы закончим с задачей обтекания клина, стоит заметить, что ее решение можно было бы получить, считая его зависящим от времени.