ABOUT PERMISSIBLE LIMITS AT ASSIGNMENT OF BOUNDARY EFFECT FOLLOWING THE SHOCK WAVE

Abstract

The author considers influence of boundary effect following the shock wave on evolution of the shock wave itself, when it penetrates the shift layer. Permissible limits at assignment of boundary effect are revealed. Concepts of Uskov’s invariant and angle are introduced.

Full Text

Рассмотрим задачу о проникновении скачка уплотнения (СУ) в сдвиговый слой (СС) невязкого [1-3] или вязкого теплопроводного газа при больших числах Рейнольдса [3-5]. Криволинейный СУ схематизируем поверхностью сильного газодинамического разрыва, на которой выполняются либо обычные, либо обобщенные соотношения нулевого и первого порядков, асимптотически учитывающие действие фактора вязкости-теплопроводности (ВТ) в приближении СС; в газодинамическом отношении эти соотношения являются точными. Приведем компактную матричную запись обобщенных дифференциальных соотношений (ОДСС) [3-5], выполняющихся на таком СУ (рис. 1): ai (ФJ)» + b,Kw+ Ci (Ф) + ^9 + Х V + f,Kw' + Si] (фj )m + Ni J nn nn eij \ф j = 0, (1) где Ф = (W, e, p, h)T, здесь p - давление, h - энтальл пия, e и e - углы наклона линии тока к оси OX декартовой (цилиндрической) системы координат XOY до и за СУ соответственно; 5 = 0(5 = 1); е2 = ( Rem) 1; N=N< |Гфj 1,К), фj, Ф, Kw]; 1,] =1,-,4. Рис. 1. Системы собственных координат, связанные со скачком уплотнения и линией тока В (1) применено суммирование по повторяющемуr' d 2Q ся индексу j; величины Kw и Kw = ——, где w = т dw координата в продольном к СУ направлении (в отлил л чие от модуля скорости W); ^=е+ст=е+ст, представляют собой продольную кривизну СУ (в отличие от радиальной кривизны 1/y в осесимметричном случае) и производную от нее; вектор N явно объединяет нелинейные слагаемые; все производные представлены в локальной системе естественных (собственных) л л координат (s, n) - до СУ, (s, n) - за СУ, причем первые считаются известными. Предполагается, что все величины являются безразмерными. Острый угол ст между вектором скорости W и СУ, т. е. между s и Ts (см. рис. 1), показывает интенсивность СУ: |ст| > ам , где ам = arcsin (1/M) - угол Маха, а знак x этого угла определяет семейство, к которому принадлежит СУ (характеристика при ст =хам) в предельном невязком случае (е 2 = 0 в (1)). Функциональные матрицы (векторы) коэффициентов A-G сложным образом зависят от газодинамических величин по обе стороны СУ. В невязком случае A-D в (1) после соответствующих допустимых преобразований в точности приводятся к аналогичным коэффициентам, полученным в [1], где использована другая группа зависимых переменных. Практическое применение ОДСС (1) предполагает их замыкание, для чего требуется задать какую-либо из производных за СУ или их комбинацию, т. е. определить краевой эффект (КЭ) за СУ. В качестве такой комбинации может быть использована дополнительная дифференциальная связь, в частности, исключающая КЭ [2-5]. Необходимость введения термина «краевой эффект» вызвана прежде всего нестандартной (не краевой) постановкой задачи для исследуемого стационарного ударного течения и лишь частичным смысловым совпадением с общепринятым термином «краевое (граничное) условие», обычно используемым при решении краевых (начальнокраевых) задач. Краевой эффект за СУ имеет производный смысл от термина «краевое условие» приме- 6 Вестник Сибирского государственного аэрокосмического университета имени академика М. Ф. Решетнева нительно к задней поверхности скачка; КЭ имеет отношение к продолженной системе соответствующих законов сохранения и необходим для замыкания соотношений первого порядка на СУ; КЭ может совпадать или не совпадать с реальным краевым условием того же порядка, например когда СУ оказывается присоединенным к границе области; КЭ имеет смысл локального источникового члена в каждой расчетной точке задней поверхности СУ. Такая постановка задачи позволяет рассматривать фронт СУ, проникающий в слой идеального или вязкого теплопроводного газа, как траекторию условной материальной частицы в подходе Лагранжа и произвести редукцию исходной начально-краевой задачи для уравнений в частных производных Эйлера или Навье-Стокса к задаче Коши для системы обыкновенных дифференциальных уравнений [3-5]. Система (1) нелинейных (относительно первых производных за СУ) уравнений после ее замыкания может быть проинтегрирована без приведения к нормальному виду по причине отсутствия старших производных в первом уравнении системы - законе сохранения массы, для чего разработаны эффективные итерационные методы [5]. Существенно, что в процессе интегрирования указанной системы удается обходиться без поточечной постановки обобщенных или обычных условий на скачке: эти условия с высокой точностью будут выполняться уже автоматически, что существенно отличает данный подход от других. Постановка этих условий требуется только в начальной точке СУ, откуда начинается его эволюция. Используя математический язык сплайн-функций, можно сказать, что локально применяемый аналитический аппарат в виде ОДСС (1) в гладком представлении фронта СУ позволяет перейти к сплайну более высокого - третьего порядка, т. е. к кусочно-кубическому представлению скачка. В предельном невязком случае (при е 2 = 0 в (1)) ОДСС переходят в невязкие дифференциальные соотношения (ДС), но проблема, связанная с их замыканием, сохраняется, поэтому вполне разумно поиск недостающего уравнения, определяющего КЭ, осуществить вначале в более простой постановке - при отсутствии фактора ВТ. В этом случае после замыкания ДС полученная система линейных уравнений становится однозначно разрешимой относительно производных за СУ и его кривизны уже в аналитическом виде, при этом итерации не требуются. Обобщив возможные КЭ за СУ, представим расширенную дифференциальную связь, с помощью ко- А торой p n исключается из (1), в виде m X 1~2 Ш -1 • p л = c , Г у, = -X tg(a M> (2) При подстановке (2) в (1) подразумевается, что л тангенциальная производная от давления p's в (2) заменена через уравнения Эйлера (УЭ) нормальной А производной от угла наклона ^ . Таким образом, количество компонент неизвестного вектора производных X = 1 W л, 0л, Kw , hn7 n7 w 7 ) в (1) становится равным количеству уравнений и при невырожденности матсистема будет однорщы А = (ai j =1,2, b , ai4 )i=L_ значно разрешимой, при этом предполагается, что матрица уже преобразована подстановкой (2), исклюл чающей p n как неизвестную величину. Коэффициент m в (2) выбран не случайно, поскольку в работе [2] строго доказано, что случай b = 1, c = 0 (это случай изобарической связи) точно соответствует рефракционной модели и догоняющие (СУ) возмущения здесь отсутствуют. Этот коэффициент фиксирует за СУ наклон характеристики противоположного ему семейства и, таким образом, (2) означает равенство нулю полной производной от давления в указанном направлении или локальную для точек на задней поверхности СУ изобаричность в данном направлении. Здесь важно отметить, что малые отклонения значений параметра b от единицы также приводят к малым отклонениям возмущенного таким образом решения. В этом случае изобарическая форма задания КЭ за СУ удобна по причине явного задания наклона изобары за скачком, а не какой-либо из производных за ним. Очевидно, что с помощью управления параметрами b и c (b Ф 0 по техническим соображениям) могут быть определены все возможные КЭ за СУ, включая реальные, когда известны кривизна стенки, от которой отражается СУ, или закон изменения давления вдоль границы перерасширенной струи [2]. Случай с = 0, b и -1 представляется важным и требует дополнительного рассмотрения. Как показано В. Н. Усковым [1], на криволинейном СУ конечной интенсивности выполняется соотношение N — CA 115 C2 A25 4 ( Л л,- A5 А V A25 N, (3) где b и с - параметры связи; ам - угол Маха за СУ. где N,, N j - неравномерности течения в текущей точке до и за СУ соответственно, в частности: N1 = д ln p/ds, N2 = c0/ds, N3 = dln p0/dn, N4 = 8/y - поперечная кривизна СУ, здесь 5 = 0 соответствует плоскому, 5 = 1 - осесимметричному случаю течения. Продольная кривизна СУ Kw = N5 исключена из (3), а коэффициенты Cj, A, зависят как от числа Маха перед СУ, так и от его текущей интен- Л / сивности J = p/p. Обратим внимание на то, что -х Л Л 7 Математика, механика, информатика обозначения коэффициентов Aji и величин Ni в (3) никак не связаны с аналогичными обозначениями матрицы А и Ni в (1). Теперь рассмотрим, как зависят друг от друга неравномерности на слабом разрыве (СЛР), в который вырождается СУ при J ^ 1. Как следует из решения соответствующей задачи линейной алгебры [1], при СЛР сохраняется простой комплекс: А А А N j -хГ(м) N 2 = +-хГ(м) N2), J = 1 Г(м) = Г(м) = у мУх/м2-!, т. е. = 0, N2 = -1 д ln p ум2 д n л -1 д ln p N 2 =—г— ум д n = + + m !• pn), I 2 м -1 4м2-1 Пользуясь (5), приведем левую часть (3) для СУ произвольной интенсивности к форме левой части (6) для СЛР. Тогда вместо ^ получим лх 1У =+C1A15 )/ C2A25 Ум 2 Л где показатель направ- (4) (n1 - X Г(м) N2 )=1 или [N1 ] J=j- хГ(м)[N2 ] J= = 0, где квадратные скобки означают разрыв соответствующей величины. Другими словами, этот левый комплекс (в силу знака перед х) является инвариантным (как инвариант Ускова) на СЛР [2] и не означает тривиального равенства N1 = N1, N2 = N2 если это имеет место, то СЛР представляет собой обычную характеристику, на которой производные не рвутся. Выразив неравномерности через нормальные производные ления х входит в коэффициенты А, Сравним между собой пары коэффициентов л хГ(м) и (C^_A15 )/(C2А25) в соотношениях (4) и (3) х л X для СЛР и СУ соответственно и 'm и m У , а также л х х соответствующие им угол Маха ам = arctg(m ) (в этой записи уже учитывается знак угла) и угол лх лх а У = arctg(m У). Последний, согласно [1], определим как угол Ускова. Заметим, что сравниваемые пары коэффициентов или углов зависят как от числа Маха перед СУ, так и от его текущей интенсивности. Необходимость в таком сравнении вызвана прежде всего тем, что строгое обнуление определителя матрицы А = (a, j =12,b,,a,4) ^ 4 (см. ОДСС (1) при е2 = 0 [3-5]) при линейной части неизвестного вектора производных X = 1 W ,K„,, h' происходит при исключении величины p n с помощью дифференциальной связи, основанной на левом комплексе, но лх (5) с направляющим коэффициентом m У а не л X из соответствующих (перед и за СУ) уравнений импульса [2; 5] и подставив (5) в (4), перепишем (4) в терминах полных производных от давления, учитывая равенство величин по обе стороны СЛР: (6) m = tg(a^, как можно было бы ошибочно предположить, проводя аналогию с соотношениями известного метода характеристик [1]. Однако, как следует из несложных расчетов, в широком диапазоне интенсивностей СУ и чисел Маха перед ним указанные направляющие коэффициенты (равно как и значения определителя матрицы А) имеют близкие, иногда даже совпадающие значения. В связи с вышеизложенным закономерны следующие вопросы: - как соотносятся направляющие коэффициенты лх х лх лх А л ч_/ mv и или углы а У и ам при конечной интенm Таким образом, имеет место равенство полных производных от давления в указанном направлении по обе стороны поверхности СЛР, при этом сами однотипные производные, как и в (4), могут отличаться. Существенно, что направляющий коэффициент mX = X • tg^M), где аМ - угол Маха, в каждой точке СЛР совпадает с наклоном характеристики того же семейства. Анализируя левые части (6) и (2) при b = 1, можно сделать вывод, что ^ = - m %, поэтому, в противовес (6) и согласно терминологии [2], левую часть (2) при b = 1 определим как правый комплекс, соответствующий рефракционной модели. сивности СУ (J Ф 1); имеет ли место предельный переход m У ^ m л X л X (а У ^ ам ) при постепенном вырождении СУ в СЛР (J ^ 1) ? Для ответа на эти вопросы приведем расчетные графики (рис. 2), характеризующие различие указанных углов, при числе Маха невозмущенного потока м = 2 и изменении интенсивности СУ в диапазоне значений J е [1, Js], где Js - звуковая интенсивность СУ (м| j =1) [1]. х 1 8 Вестник Сибирского государственного аэрокосмического университета имени академика М. Ф. Решетнева Рис. 2. Углы Маха и Ускова, град, за СУ (слева: 1 - угол Маха; 2 - угол Ускова) и их разность, град (справа: 1 - разностный угол; 2 - нуль), для м = 2, у = 1,4, х = -1 Тогда при вырождении СУ в СЛР (J ^ 1) в пределх X ( 1 Л ле будем иметь ам = ам = х • arcsinI ^ I = -30° для X = -1. Из расчетов следует, что в широком диапазоне интенсивностей СУ (исключая окрестность звуковой интенсивности) и чисел Маха перед скачком указанные углы и соответствующие им коэффициенты имеют близкие значения, а в пределе (J ^ 1) совпадают. Близость этих углов позволяет рассматривать угол Ускова как двойник хорошо известного угла Маха, однако угол Ускова в отличие от угла Маха имеет отношение исключительно к задаче первого порядка (дифференциальной задаче) для криволинейного СУ, хотя формально может быть рассчитан (но не востребован) и за плоским СУ. Такая ситуация возникает, например, в методе характеристик [1], когда криволинейный СУ поточечно выстраивается кусочно-плоскими элементами. В качестве кривизны такого СУ может выступать лишь ее приближенный разностный (например, трехточечный) аналог, но не сама кривизна скачка, как в предлагаемом подходе. Вернемся к заданию КЭ за СУ. С учетом сказаннолх го выше, наклон изобары m У , определяемый углом Ускова, является запрещенным при задании КЭ, однако имеется возможность несколько отступить от этого направления. Продемонстрируем на конкретном примере принципиальное отличие правого и левого комплексов при задании с их помощью КЭ за СУ в рамках изобарической дифференциальной связи (с = 0 в (2)). Как уже было замечено, при использовании правого комплекса даже значительные вариации коэффициента m % приводят к малому изменению соответствующего решения. В случае же использования левого комплекса отмеченной близости соседних симметричных (±Д) решений при вариациях направляющего коэффициенл х та m У +±Д); УД >0, Д ~ 0 не существует и эти решения ведут себя по-разному. Так, ветка решений, л х определяемых коэффициентом mУ (1 -Д), соответствует резкому принудительному усилению СУ в случае, если число Маха в неравномерном невозмущенном СС убывает (как в пограничном слое) в процессе проникновения в него СУ, а ветка решений, определ х ляемых коэффициентом m У (1+ Д), - принудительному вырождению СУ до СЛР. Если невозмущенный СС имеет возрастающий профиль числа Маха, то ситуация с этими решениями меняется на противоположную. Более того, чем меньше величина параметра Д , тем более интенсивно происходит любой из указанных процессов, что связано с большими градиентами газодинамических величин, т. е. разворотом СУ в ту или другую сторону в меньшей окрестности, дающим в пределе угловую точку на скачке. Представим результаты вычислительного моделирования проникновения СУ в плоский вихревой слой для значения малого параметра Д = 0,0005 в виде расчетных графиков (рис. 3, 4), где по вертикальным осям отложены значения безразмерной нормированной на толщину слоя координаты y, по горизонтальным осям - функции: на рис. 3 - скорость и энтальпия, на рис. 4 - характеристический коридор. Конкретные значения безразмерных газодинамических параметров во внешнем потоке и на плоскости симметрии, совместно определяющие невозмущенный, рассчитанный разностным методом изоэнталь-пийный сдвиговый (вихревой) слой, были следующими: у = 1,4; м* = 2,28; р* = 1,18; W* = 0,95 при у = 1 и W0 = 00 = 0 при у = 0. Начальная интенсивность J *+% ) = j/%/W0% падающего СУ составляла t% = 20 % от логарифма Л* = ln(Js) звуковой интенсивности скачка. Характеристический коридор представляет собой «распечатку» в каждой расчетной точке СУ величин |хам|, |ст|. Р + Хам (см. рис. 1) модулей значений углов наклона СУ и характеристик одного с ним семейства по обе его стороны к текущей линии тока в невозмущенном течении (см. рис. 4). л 9 Математика, механика, информатика hji Рис. 3. Распределения скорости (слева) и энтальпии (справа) вдоль передней и задней поверхности СУ для различных невязких моделей: 1 - рефракционная модель; 2 - модель с принудительным вырождением СУ; 3 - модель с принудительным усилением СУ; 4 - невозмущенное решение , • Г."" ,4 t **.**' 1 1 ....... •• 1 % . *. 1 i 1 чч 1 1 Ч. С 1 \ 1 '( \ \ \ , \ \ Ч \ S ч \ \ \ % \ ч N. N ч V V Ч V ч v Ч \ ' % ' V ' ч 1: С¥; КЭ:0 . , £\ЕГ+ — 4: до С?; КЭ:0 I I ш.ш 5 : за СЭ; КЭ: О - 6: до С¥; КЭ: - • •• • 7: за С¥; КЭ:-8 : до С¥; КЭ: + --- 9: за СП; КЭ: + Рис. 4. Характеристический коридор, град: 1, 4, 5 (эквидистантные кривые) - рефракционная модель; 2, 6, 7 - модель с принудительным вырождением СУ; 3, 8, 9 - модель с принудительным усилением СУ; 4, 6, 8 - невозмущенное решение, соответствующее различным моделям Кривые, соответствующие КЭ, усиливающему СУ (mу (1 -А)), помечены как КЭ: + ; вырождающемуся л х СУ (mу (1+ А) ) - как КЭ: - ; рефракционной модели, исключающей КЭ за СУ (m * (1 ± 0) ) - как КЭ: 0. Крайнее левое семейство слившихся кривых соответствует характеристикам в невозмущенном потоке, крайнее правое (веер) - в возмущенном потоке (за СУ), среднее семейство - самому СУ. Очевидно, что при полном вырождении СУ все три однотипные кривые из различных семейств, соответствующие этому решению, пересекутся в одной точке внизу и тогда характеристический коридор «схлопнется», как на рис. 4. В случае КЭ, усиливающего СУ, этот коридор, наоборот, расширится, и только в нейтральном случае останется практически без изменений (кривые эквидистантны). Таким образом, при А Ф 0 базовый левый комплекс может эффективно использоваться в качестве математического инструмента управления поведением СУ: либо для его усиления, либо для перевода в разряд СЛР. Следует отметить, что при использовании модели вязкого теплопроводного газа величина малого параметра А не может быть слишком малой, посколь 10 Вестник Сибирского государственного аэрокосмического университета имени академика М. Ф. Решетнева ку градиентным процессам разворота СУ на малом участке физического пространства соответствуют малые радиусы кривизны СУ, т. е. большие кривизны разных знаков, что неизбежно входит в противоречие с асимптотическим учетом фактора ВТ и самой схематизацией скачка [4]. Проведем аналогичные исследования с учетом фактора ВТ. В ходе расчетов, где использовалась стационарная модель на основе ОДСС (1), было замечено, что при КЭ, ослабляющем или усиливающем СУ, влияние фактора ВТ на процесс проникновения скачка в слой больше и кривые, соответствующие вязким и невязким решениям, расходятся сильнее. Это объясняется тем, что исключение КЭ (фактически - разглаживание течения за СУ) уменьшает влияние фактора ВТ: градиенты газодинамических параметров входят в величины вязких напряжений и теплового потока. Для иллюстрации данного факта на рис. 5 показан характеристический коридор для вязких и невязких решений с КЭ и без него для Jm (t% = 63,7 %), Pr = 1, Rex = 1,3-103 в начальной точке падения СУ на СС со следующими параметрами: у = 1,4; Мм = 2,5; рм = 1; Wm = 1 при у = 1 и W0 = 90 = 0 при у = 0. Как и на рис. 3, 4, по вертикальной оси графиков на рис. 5 отложены значения безразмерной нормированной на толщину слоя координаты у. Анализ полученных графиков показывает, что вне зависимости от знака КЭ всегда происходит частичное, но более значительное, чем при исключенном КЭ, сужение характеристического коридора за счет ослабляющего действия фак тора ВТ. Таким образом, если имеет место КЭ, ослабляющий СУ, то в результате суммарного однонаправленного действия обоих факторов сужение коридора оказывается более сильным. При произвольном числе Прандтля Pr Ф 1 действие фактора ВТ может носить более сложный характер, поскольку вязкость и теплопроводность влияют на ударный процесс по-разному.