v = max

[0. с, (i-ir);

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

Получающееся разностное уравнение можно решать различными способами. Холст и Боллхауз делали это как обычным методом последовательной верхней релаксации по линиям, так и с использованием схемы приближенной факторизации. Обсудим схему приближенной факторизации Холста, называемую ПФ2-схемой, поскольку методы последовательной верхней релаксации были описаны при рассмотрении других систем уравнений.

Схема приближенной факторизации (и многие другие) в случае релаксационной задачи, описываемой уравнением вида

=0, где L - некоторый дифференциальный оператор, для уравнения полного потенциала может быть записана в виде

МС + а1ф = 0, (6.172)

В этом выражении контравариантные компоненты скорости имеют вид

/ + 1/2. / = (l)/ + l/2. / {Фш, / /) +

/+1/2 T ("2)/, / + 1/2 (/ + 1, / + /""/-1, /) +

+ (зк/+1/2(/./ + 1""/,/)

где значения в полуцелых узлах получают усреднением, а величины А определяют по уравнению (6.152). Такое представление справедливо для дозвуковых областей течения. Теперь введем смещение назад при аппроксимации плотности и заменим уравнение (6.167) следующим:

Jp(t) 1 + (рт) =0, (6.169)

Р/ = (1 - i+k. /) Р/+1/2. / + "i+k, i9i+2k-U2, г (6.170)

и, С + 1/2./<0.

Для искусственной вязкости выбираем выражение

где о - параметр релаксации, С - поправка (*+-ф), 1ф - невязка (приближенное решение не удовлетворяет дифференциальному уравнению в частных производных) и Л -оператор, определяющий тот или иной итерационный метод. В схемах приближенной факторизации представляют в виде произведения двух или более операторов: N = N\N2 Операторы N\ и N2 следует выбирать так, чтобы их произведение аппроксимировало L. При этом годятся только простые матричные операции и схема в целом будет устойчивой. В ПФ2-схеме, которую использовал Холст, оператор представляется в виде

aNCS., = - [а - 4, (if) J [aV, - V.p, {*f) aJc," ,.

(6.173)

где a -свободный параметр, который можно интерпретировать как (АО"". Некоторый набор а используется во время вычислений для подавления высоко- и низкочастотных ошибок в решении. Холст и Боллхауз [Hoist, Ballhaus, 1979] привели один такой набор. Его можно использовать и для других расчетов обтекания трансзвуковых профилей. При этом ПФ2-схема работает почти оптимально.

В ПФ2-схеме уравнения (6.173) решение определяется при помощи двухшаговой процедуры:

(6.174)

где f - результат на промежуточном шаге. На первом шаге fi,f получают из решения двухдиагональной системы, аналогичной той, которая встречалась в факторизованной схеме с расщепленными потоками в § 6.3. На втором шаге требуется решать уже трехдиагональную систему. Направление расчета на шаге 1 - от профиля, на шаге 2 - к профилю. В этой схеме нет ограничения на направление расчета в зависимости от направления потока. Однако при использовании других методов, например метода последовательной верхней релаксации, необходимо следить, чтобы в сверхзвуковых зонах направление расчета совпадало с направлением потока. Это как бы соответствует введению стабилизирующего члена Фlt,

Джеймсон [Jameson, 1974] подчеркивал, что при решении уравнения потенциала вблизи звуковой линии возникают трудности с устойчивостью. Чтобы обойти их, в разностную схему добавляют времениподобные члены. Последние имеют вид ф

И и обычно включаются в операторы в релаксационных схемах. В нашем случае к операторам на втором шаге уравнения (6.174) добавляются члены типа a/(iV, aKiA, Используются только разности вверх по потоку и обычно только в сверхзвуковых областях. Знак этих членов выбирается так, чтобы величина диагонального члена росла на втором шаге уравнения (6.174), что обеспечивает диагональное преобладание. Было показано [Hafez et al., 1979], что добавлением этих времениподоб-ных членов можно произвести модификацию зависимости для плотности. Она осуществляется добавлением члена фt к квадрату скорости в уравнении для плотности (6.153). Такой подход также обеспечивает стабилизирующую добавку, вид которой близок к той, что предложил Джеймсон.

При постановке граничных условий на профиле используются фиктивные точки внутри стенки, чтобы можно было задать условие отражения. Условие V = О на профиле задается как

/РКЧ (6.175)

в фиктивных точках, где Л/ -точка на поверхности профиля. В этих расчетах граничное условие задается явным образом.

Если подъемная сила профиля равна нулю, значения потенциала скорости и плотности на внешней границе сохраняются постоянными и равными значениям в набегающем потоке. Для ненулевой подъемной силы эти величины на внешней границе должны давать циркуляцию, совместимую (совпадающую) с завихренностью, найденной из решения, и обновляемую в конце каждой итерации. В конце каждой итерации циркуляция рассчитывается по скачку потенциала скорости на задней кромке:

Г = фи,-ф1,. (6.176)

В начале следующей итерации скачок потенциала скорости задается вдоль всей вихревой пелены. Поправка определяется как разность при переходе через вихревую пелену

Г~Г = С2~С?, (6.177)

= 3 (Г ~ Г "О + Г"-. (6.178)

На рис. 6.28 показаны результаты расчета трансзвукового обтекания профиля NACA0012. На нем приведено сравнение коэффициентов давления, вычисленных по ПФ2-схеме, с коэффициентами, взятыми из работы Лока [Lock, 1970]. ПФ2-схема работает хорошо, о чем свидетельствуют полученные результаты.