Роуч П. Вычислительная гидродинамика (1185913), страница 39
Текст из файла (страница 39)
Хотя сама идея иестацпонарного подхода привлекательна (Харлоу и Фромм [1965[, Макано [19651), естественно возникает следующий вопрос; зачем возиться с нестационарным решением, если интерес представляет только возможное стационарное решение? Почему бы не положить дьй/д1 и не работать с уравнениями, описывающими стационарное течение? Такой подход на самом деле с успехом использовался многими авторами, Однако в обшем случае пока рекомендуется нестационарный подход. В качестве первого решающего аргумента в его пользу продемонстрируем эквивалентность простейшей схемы для стационарного уравнения и одной разностной схемы для нестационарного уравнения.
Одномерное стационарное модельное уравнение переноса получается из модельного уравнения переноса параболического типа (2.18) и имеет вид дч деч и — =а —. дх дх' ' (3.354) Это уравнение можно представить в конечно-разностной форме, используя центральные разности по пространственной переменной: — 2г,.
+ й. (3. 355) 2 Дх Лхе Для решения этого уравнения эллиптического типа можно воспользоваться какой-либо из итерационных схем, которые будут обсуждаться в равд. 3.2. В простейшей итерационной схеме (Ричардсона или Якоби) уравнение (3.355) разрешается относительно ~,. Считая известным некоторое начальное приближение Ц для всех с, новые значения се+' на (й+1)-й итерации определяются по стояшим в правой части старым значениям Эд, Методы решения уравнения переноса вихря 162 после н.й итерации: ~~" = — 4 (С'.1 — Ь*;-~)+ 2 (~,'+,+~,',). (3.356) Итерации продолжаются до тех пор, пока не будет выполнен некоторый критерий сходимости (разд. 3.4). Если из обеих частей уравнения (3.356) вычесть ЬЯ, то итерационная вычислительная схема не изменится или '3 3 8) ! е 2а 2Ьх 2 Ьх Заметим теперь, что уравнение (3.358) будет эквивалентно уравнению (3.18) (неста14ионарному конечно-разностному уравнению с разностями вперед по времени и с центральными разностями по пространственной переменной), если определить шаг по времени в уравнении (3.18) как Мн = 1, итерационную скорость конвекции как и" = иЛха/(2а) и итерационный коэффициент диффузии как схн = Лхв/2.
Из анализа уравнения (3.18) известно, что для сходнмости требуется, чтобы итерационное число Куранта С" не превышало единицы, т. е. чтобы или (3.359 б) Ке,<2. Значит, при сеточном числе Рейнольдса ибх/сс ) 2 данная итерационная схема будет расходиться, Этот пример показывает, что итерационная схема Ричардсона в точности эквивалентна нестационарной схеме и является ограниченной. Другие итерационные схемы для уравнений эллиптического типа эквивалентны или по меньшей мере аналогичны нестационарным схемам для уравнений параболического типа.
Впервые на такую аналогию указал Франкел [1950). Впоследствии Ходжкннс [1966[ установил соответствие между полуаналитическим методом Чебышева и решением не- стационарного уравнения гиперболического типа. Хейвуд [1970[ исследовал связь между решением уравнений для стационарного течения и пределом решения нестационарного уравнения. В схемах нижней и верхней релаксации существенную роль играет умножение членов, стоящих в правой части уравнения (3.357), па множитель г, причем г ( 1 соответствует нижней релаксации, а г > 1 — верхней релаксации.
Например, для схемы (3.358) уменьшение г, очевидно, эквивалентно уменьше- 3Л.22. Слежы олл расчета стационарньгл течеяиа !83 нню величины шага по времени. Как следует ожидать из анализа для нестационарных уравнений, прн увеличении числа Ке требуется нижняя релаксация. Текстор [1968) н Техейра [1966] экспериментально установили, что нужная величина г — 1Яе в соответствии с полученным ранее ограничением на число Куранта. Но начиная с некоторого Ке н выше, сходимость не достигается даже при произвольно малых г').
Это подтверждает опыт многих исследователей, которые (успешно) использовали схему такого типа; см., например, Том и Орр [1931], Том [1933], Том и Апельт [196!], Кавагути [1953, !961, 1965], Бургграф [1966), Майкл [1966], Гамнлец с соавторами [1967а, 1967б], Фридман с соавторами [1968], Деннис с соавторами [1968], Фридман [!970], Ли и Фын [!970). Ограничение на Ке, сохраняется даже при использовании некоторых (но не всех) приближений пограничного слоя (Плоткин [1968]). При итерационном решении стационарного уравнения это ограничение удается устранить, вернувшись для аппроксимации конвективных членов к разностям против потока, как это сделали Ранчел и Вольфштейн [1969], Гринепэн [1969а, 1969б], Госмэн и Сполдинг [1971).
Возможны и другие итерационные схемы, которые могут оказаться вполне эффективными; в дополнение к отмеченным работам укажем здесь ряд работ. Аллен и Саусвелл [1955] использовали релаксацнонную разностную схему Саусвелла (равд. 3.2.3) для решенця задачи об обтекании цилиндра при Ке = 1000; Гриффнтс с соавторами [1969) применял линейный метод последовательной верхней релаксации (равд. 3.2.4) в цилиндрической системе координат; Катсанис [1967) и Бреди [1967] решали итерационным методом стационарные уравнения, описывающие потенциальное течение. Отметим интересный исторический факт, заключающийся в том, что большинство исследователей, применявших итерационные схемы для решения стационарных задач, пе занималось анализом устойчивости н скорости сходимости своих схем, а определяло характеристики эмпирически, хотя уже в то время из ранней работы фон Неймана был известен метод исследования устойчивости для уравнений, описывающих нестационарное течение.
Возможное объяснение этого факта заключается в том, что методы расчета стационарных течений развивались из раздела численного анализа, относящегося к решению уравнения Пуассона, для которого простейшие итерационные методы не имеют ограничений, связанных с устойчивостью. ') Несмотря на то что в лействнтельиости зто ограничение, наложенное на Пс, может быть связано с линейной неустобчнвостшо в итерационном процессе, аналогичном процессу нестациопарного решения, оио может быть связано н с влиянием границ, которое рассматривается в равд, 3.3.8. 8.Д Методы регивния уравнения яервноса вихря 164 Рассмотрение нелинейных конвектпвных членов может изменить строгую эквивалентность между итерационной схемой Ричардсона и нестационарным подходом.
При пгстациопарном подходе на каждом шаге по времени решается уравнение переноса вихря и обычно итерируется до сходпмости уравнение Пуассона эллиптического типа, Чача = ь. При стационарном же подходе каждое из этих двух уравнений можно (хотя и не обязательно) итерировать последовательно. При таком подходе с «комбинированным итерированием» сходимость итераций, если опа имеет место, может быть достигнута за меньшее число шагов; Госмэн (личное сообщение; см. также Госмэн и Сполдинг [1971]) указывал, что затраты машинного времени при этом сокращаются на 40%. Однако, поскольку для уравнения Лтазь = Ц не обеспечивается сходимость на каждой итерации для уравнения вихря, плохой выбор исходных значений при стационарном подходе может привести к неустойчивости, обусловленной нелинейностью (Техейра [!966], Лил [1969], Ричарде [1970] ) из-за конвективных членов и~или из-ча граничных условий.
Опыт показывает, что если граничные условия итерируются с достаточной степенью точности, то такая трудность не возникает ни при использовании явных нестационарных схем, ни прв использовании неявных схем. Более того, программа на Фортране может быть написана для нестационарных уравнений, а затем приспособлена для стационарного подхода. При «комбинированном итерировании» уравнения Пуассона и уравнения переноса вихря можно пользоваться простым критерием сходимости для уравнения Пуассона.
(Эту процедуру действительно можно рекомендовать для расчетов; см. равд. 3.4), Преимущество, присущее итерационному методу Либмана (методу Гаусса — Зейделя) или итерационному методу последовательной верхней релаксации (будут рассмотрены в равд. 3.2), которые аналогичны нестационарным явным схемам метода чередующихся направлений (равд. 3.!.17), можно обеспечить простым добавлением в программу оператора Е4,1Ш'тгАЕЕг(СЕ для массивов ~"+' и ",". На практике использование меньших значений параметра нижней релаксации вблизи грашец (Фридман [1970] для расчетов в граничных точках брал параметр г приблизительно равным одной трети от его значения, принятого для внутренних точек) может быть реализовано введением переменного в пространстве ') шага сзй '] Таким образом, локальное значение Л!(х, у) ьюжно определять в зависимости от локального критического значения ЛД а затем брать его локально дли продвижения с посредствам д;/д1 Нестаггнанарное решение нри атом, очевидно, не будет иметь смысла, йа скорость сходимости итераций к стационарному состоянию может быть увеличена.