dt дх и u{p) = df{p)/dp.

TVD-схема (14.88), (14.89) может быть представлена в виде

/,+1/2 = 0.5 (fj + f,+i) - 0.5 [ф (г) C]/+,/2Af/+i/2 -

- 0.5а [1 - (r)]/+,/2Af/+,/2, (16,160)

Можно напомнить (п. 9.1.2), что производные четвертого порядка в правой части уравнения (18.6) приводят к положительной диссипации, если коэффициент при них меньше нуля.

Диссипация четвертого порядка (18.155) использовалась также в схеме искусственной сжимаемости (17.51) для исследования несжимаемых течений, в маршевой по пространству схеме (16.145) - для расчета сверхзвуковых вязких течений и для нахождения стационарных решений уравнений Эйлера (п. 14.2.8).

18.5.2. Ударные волны

Если в расчетной области имеются ударные волны, численной диссипации, описанной в п. 18.5.1, недостаточно для предо-тврапхения дисперсионных осцилляции вблизи скачка - областях преимупхественно невязкого течения. Вблизи твердой поверхности физические диссипативные механизмы уменьшают большие градиенты, связанные со скачком, и, если только скачок не слишком сильный, введения дополнительной численной диссипации не требуется.

Явное введение численной диссипации (искусственной вязкости) для сглаживания профилей скачков описано в п. 14.2.3. Это может привести к существенному размазыванию скачка (рис. 14.18). Алгоритмы FCT и схемы TVD (п. 14.2.6) можно рассматривать как схемы, вводящие искусственную вязкость до тех пор, пока профиль скачка не станет монотонным, а затем выборочно убирающие численную диссипацию для получения более резкого профиля скачка. В данном разделе типичная схема TVD будет рассмотрена как трехточечная центрально-разностная схема (второго порядка) с добавлением численной диссипации. Такое рассмотрение позволяет построить более эффективные неявные алгоритмы расчета стационарных сжимаемых турбулентных течений со скачками.

Нелинейный скалярный закон сохранения можно записать в форме

Р+%>- = 0, (18.158)

Cj+u2 = UfjU2At/Ax, a = signC/+i/2 и (г) - ограничитель. Возможны различные выборы ф{г) \ один из них определяется формулой (14.86).

Если (18.160) и эквивалентное выражение для -1/2 подставить в уравнение (18.159), то окажется, что члены вида 0.5 ( + +i) образуют трехточечную центрально-разностную аппроксимацию производной df/dx. Остальные члены (18.160) лри подстановке в (18.159) образуют численную диссипацию.

Уравнение (18.159) является явным алгоритмом, пригодным для расчета нестационарных задач. Для решения задач, описываемых стационарными уравнениями Эйлера или Навье- Стокса, желательно использовать неявные алгоритмы, а для численного представления потоков целесообразно применять формулы, не зависящие от шага по времени. В этом случае стационарное решение не будет зависеть от шага по времени.

Это осуществлено в работе [Yee, 1987], в которой соотношение (18.159) заменено уравнением

/+ Р лГ(ж/2 ~ - 2) = р - (J ~ Р) КГ(f/+i/2 -

(18.161)

где р играет ту же роль, что и в уравнении (18.73). Формула (18.160) заменяется следующим выражением:

+ U2 = 0-5(f/ + f/ + ,)~0.5U;+i/2[I-(r)],/2Ap/H-l/2- (18.162)

Оказалось, что для того чтобы схема (18.161) была схемой TVD (т. е. удовлетворяла условию (14.81)), должно выполняться условие типа КФЛ

При расчете стационарных задач рекомендуется всхеме (18.161) г.спользовать значение Р = 1. При таком выборе р условие (18.163) выполняется при любом А.

Чтобы иметь возможность использовать метод решения систем уравнений с трехдиагональной матрицей (п. 6.2.3), необходимо с целью получения линейной относительно Ар + системы линеаризовать (18.161) относительно временного слоя п. В результате линеаризации получается трехдиагональная система уравнений

/ = 0-5р-(-/-1-/-1/2), В/=1+0.5р~(£/-,/2 + £/+1/2)>

(18.165)

Из-за конкретного вида коэффициентов в], В), В/схема (18.164) является пятиточечной, хотя и трехдиагональной относительно Др. Для решения методом установления чисто стационарных задач более экономно [Yee, 1987] опустить ограничитель в выражении для Е в (18.165). Таким образом, £/4.1/2 = t/y+i/2. Однако ограничитель сохраняется в правой части уравнения (18.164), в результате чего схема имеет второй порядок точности по пространству.

Обобщение схемы TVD на системы уравнений, подобные уравнениям Эйлера, осуществляется путем замены уравнений для отдельных компонент системой характеристических соотношений (см. формулы (14.91), (14.92) (п. 14.2.6)). Описанный выше алгоритм применяется к каждой характеристической переменной. После суммирования вкладов отдельных компонент вместо (18.164) получается блочно-трехдиагональная система.

Все сказанное выше можно проиллюстрировать на примере суммарной системы уравнений Эйлера (14.43)

l + -L=.o. (18.166)

Характеристические переменные получаются в результате замены (18.166) уравнением

1 + А-?- = 0. (18.167)

где А = dF/dq. Элементы матрицы Якоби А для эквивалентной двумерной задачи определяются формулой (14.99). Матрица А разлагается на множители так же, как (18.60):

А = ТлАдТл. (18.168)

Диагональная матрица Ал составлена из собственных чисел матрицы А, т. е. diag Ал = {и, и+ а, и - а}.

Эквивалентный (18.164), (18.165) неявный алгоритм для решения уравнения (18.166) имеет вид

B}Aq?+l + BfAqr + B?Aq?tl = -(¥%г - PUn), (18.169)