Роуч П. Вычислительная гидродинамика (1185913), страница 26
Текст из файла (страница 26)
Однако анализ показывает, что схема Дюфорта — Франкела устойчива, в то время как схема Ричардсона безусловно неустойчива. Разрешение этого парадокса связано с идеей корректной постановки задачи или чувствительности решения к начальным данным. Поскольку некоторые специальные случаи несущественны для практических расчетов, мы ие будем рассматривать уст >йчивость в этих случаях. Тогда, если !о! 3.!Я.
Первая схема с разностями против потока странственной переменной. Но при с(че '/я эти две схемы, полученные по ним результаты и свойства устойчивости оказываются различными. 3.1.8, Первая схема с разностями против потока. Ошибки, обусловленные схемной искусственной вязкостью Одношаговая явная двухслойная по времени схема, обеспечивающая статическую устойчивость для конвектпвных членов, основана на использовании односторонних, а не центральных разностей по пространственным переменным. Когда скорости положительны, то используются разности назад, н наоборот').
Таким образом, односторонняя разность всегда берется против потока, т. е. в направлении вверх по течению от точки, в которой вычисляется бс/61з). Данная схема имеет ошибку аппроксимации Е = 0(Ж, Лх) и записывается так: — при и>0, Ьх ~и+1 йв (3.176) иь'!' ! — иь" при и<0. ох Данная схема часто применяется в литературе под различными названиями и с разными объяснениями. Метеорологам давно известно стабилизирующее влияние разностей против потока (Лилли [1965], Вазов и Форсайт [1960]) или «наветренных (Франкел [1956]) разностей '), и они используют их для решения задач о течении несжимаемой жидкости и о течении в приближении Буссинеска.
Математики относят эту схему к разностным уравнениям «с положительными коэффициентами» (Вазов и Форсайт [1960], Моцкин и Вазов [1953]); употребляется также термин «несимметричные разности» (Ломекс с соавторами [1970]). Эту схему Рихтмайер [1957] сначала приписывал Лелевье, впервые применившему ее для исследования течения невязкой сжимаемой жидкости при условии симметрии слоя. Робертс и Вейс [1966], Курцрок [1966] и Крокко [1965], по-видимому следуя Рихтмайеру, тоже связывали ее с именем Лелевье. В более поздней работе того же Рихтмайера (Рихтмайер [1963]) эта схема была названа «схемой с разностями против потока» ) Термины «назад» нлн «вперед», очевидно, имеют смысл только по отношению к скорости. ') С учетом этого прн переводе принято название «схема с разностями против потокам †Пр. ред.
') «Наветренньгм разностям» противоположны безу«лов!!о неустойчнныг «подветренные разности» (Франкел [!9бб1) 102 3.1. Методы решения уравнения переноса вихря и указано (в этом к Рихтмайеру присоединились Стоун и Брайен [1963]), что она восходит к статье Куранта, Изаксона и Риса [1952). В этой статье была впервые продемонстрирована тесная связь данной схемы с теорией характеристик, а сама схема применена для плоских течений невязкой сжимаемой жидкости. Схема «типа П» Лонглн [1960] и первая схема Филлера и Ладлоффа [1961) представляют собой применение обсуждаемой схемы для исследования одномерных течений сжимаемой жидкости с учетом вязкости. В методе Г) 1С (метод жидкости в ячейках) Джентрп, Мартина н Дали [1966] данная схема называется «разностной схемой с донорнымн ячейками» '). Курцрок [1966] применил эту схему для исследования плоского течения вязкой сжимаемой жидкости и нашел критическое условие для величины шага по времени в этом случае.
В литературе можно найти и много других приложений этой схемы. Переписывая уравнения (3.176) и вводя число Куранта С = = исх1/свх с положительной постоянной скоростью и, получаем (3. 177) втпя-! — втп С (~п втп ) При С=! эта схема дает йп,+'= Ц' ц что соответствует точному решению (см.
равд. 3.!.6). Условие С=! является также предельным условием устойчивости (см. предыдущее упражнение). При С (1 схема вносит искусственное затухание; при этом исследование устойчивостц по методу фон Неймана показывает, что матрица перехода имеет собственные значения )х ( 1. Любая схема для уравнения с одним только конвективным членом в невязком случае при д ( 1 обладает таким схемным искусственным затуханием, а разложение в ряд Тейлора (так же, как и применение метода Херта исследования устойчивости) показывает, что уравнение (3.176) эквивалентно уравнению — = — — + а,—, + ЧВП + ПВП, дв д (ив) дяв (3.
178) где ПВП вЂ” производные высшего порядка, а ссе /аи '-'-т (! С) (3. 179) Поскольку здесь появляется нефизический коэффициент сс, при производной дяь/дхя, это объясняет не только искусственное затухание, но, говоря конкретнее, и то, что схема с разностями ') Схему с донорнымн яясйками Джентри, Мартина и Дали мы будем рассматриват~ ниже, называя ее второй схемой с рааностями против потока (см.
равд. 3.1.11). 8д.а лервая схема с разностями против потока 108 против потока обладает схемной искусственной диффузией (Нох и Проттер [1963!), или схемной искусственной вязкостью '). Интерпретация коэффициента ал в случае многомерных и вязких течений не столь очевидна, как это могло бы показаться. Рассмотрим, например, случай, когда достигается стационарное состояние. Тогда левая часть уравнения (3.176) обращается в нуль и можно уменьшать Лс, не меняя при этом решения конечно-разностного уравнения. Уравнение же (3.179) наказывает, что уменьшение Л( приводит к увеличению а, (через С).
Если понятие схемной вязкости а, имеет какой-либо смысл, то решепне конечно-разностиого уравнения, казалось бы, должно было зависеть от величины а,. Однако если вместо исследования нестацноиарного уравнения положить дь/д( = 0 в уравнении (3.176) и воспользоваться разложением в ряды Тейлора, то получится а, = — '/зи сзх. (3.180) При этом а, не является функцией от Л( и стационарное решение не зависит от М. Противоречие между выражениями (3.179) и (3.180) для а, можно объяснить, вспомнив, что для модельного уравнения, содержащего только конвективный член, единственным возможным стационарным решением при и = сопз1 является тривиальное решение,",е' = !.лс = сопз(, В этом случае бтс)без = дтй/дхт = О, что допускает произвольный вид коэффициента а,.
Рассмотрим теперь применение схемы с разностями против потока к уравнению плоского течения, учитывающему как конвекцию, так и (физическую) диффузию. При постоянных положительных ссь о; получаем злят! ьлл ~л Зл Зл,ел ( ~п 1 ~л Зйл ~л 1 ~л З~лс охз Луз причем условие устойчивости накладывает следующее ограничение на сз(: б(~~ [2а( + з ) + д + л 1 ' (3'182) Разложения в ряды Тейлора дают д~ д (и~) д (о~) дзь — = — — — — +(а+а ) — + де дк ду ех дх' +(а+а,„) д, + ЧВП+ ПВП, (3.183) ') Кан мы увидим позже, представление коэффициента схемной вязкости а, для схемы с разностями протнв потока не едияственно.
104 8.1. Методы рви«ения уравнения переноси вихря где в случае исследования нестационарного уравнения а,„= '/ии Лх (1 — С,), а,„= '/ип Лй (1 — Се), (3.184) а в случае исследования стационарного уравнения ае„= '/еи Лх, а„= — '/еи Лу, (3.185) Автор настоящей книги (см. приложение Б) показал, что результат (3.185) действительно соответствует стационарным решениям. Для нестационарных решений формулы (3.184) показывают, что влияние схемной вязкости будет минимальным, если С н С, по возможности близки к единице.
Однако на практике невозможно добиться, чтобы эти две величины одновременно были близки к единице во всех частях области течения, поэтому схемная вязкость обязательно будет входить в расчеты. Заметим также, что коэффициенты схемной вязкости зависят от составляющих скорости и и о, которые рассматриваются относительно неподвижной эйлеровой системы координат. Это приводит к нарушению принципа инвариантности Галилея, т.
е. преобразование, связанное с обращением скорости иевозмущенного потока и допустимое для дифференциальных уравнений, неприменимо к этим консчно-разностным уравнениям, за исключением случая, когда Лх-~ О, Лу — н О. Термины «искусственная вязкость» и «схема первого порядка» часто используются как синонимы, но в действительности это не так. Например, можно просто добавить в схему второго порядка дополнительный член ав„д'ь/дх' с явной искусственной вязкостью авх — Лх'. Такая схема с явной искусственной вязкостью имеет второй порядок; на ней, в частности, основывается метод фон Неймана — Рнхтмайера для расчета ударных волн (см.
равд, 5.4.! ). Обсудим некоторые соображения относительно того, что точное решение невозможно до тех пор, пока не выполнено условие а, « а. Из соотношений (3.185), полученных прн исследовании стационарного уравнения, видно, что для выполнения этого условия должно быть иЛх/а « 2 н оЛу/а « 2, т. е. сеточные числа Рейнольдса по различным направлениям должны быть много меньше 2. Эти условия являются требованием формальной точности, но па практике положение оказывается не столь уж плохим. Рассмотрим некоторую область, где применимо приближение пограничного слоя (см.
Шлихтинг [1968] ). Тогда производная де~/дх' будет мала н вклад члена с коэффициентом а+ ае. в уравнение (3.183) будет мал. Кроме того, величина о также мала, поэтому а,„в (3.184) н (3.185) может быть меньше, чем а (см. задачу 3.9). Исследования Ранчела и Вольфштейпа (1969] (см. также Вольфштсйн !1969] ) показывают, что для плоского течения 105 ЗЛ.В.
Первая схема с разностями против потока т, ж >)зиЛх ейп 26, где 6 — угол, который линия тока образует : осью абсцисс. Эти авторы использовали вторую схему с раз>остями против потока '), которая будет рассмотрена ниже, для ласчета течения внутри замкнутой прямоугольной области с >дной подвижной границей. При Ке= 100 на неравномерной :етке 13;к', 13 максимальное сеточное число Рейнольдса было >коло 20. Тем не менее полученные здесь результаты достаточно хорошо согласуются с решением этой задачи, полученным при помощи схемы второго порядка на сетке 5!;>с', 5!. Аналогично, результаты расчетов естественной конвекции, выполненных Торрансом [1968) также с использованием второй схемы с разностями против потока при большом числе Грасгофа [эквивалентном Ке = 300), отличаются от решения, полученного при помощи схемы второго порядка, менее чем на 5о7о.