Роуч П. Вычислительная гидродинамика (1185913), страница 43
Текст из файла (страница 43)
По сравнению с прямыми методами различные итерационные методы проше с точки зрения понимания и программирования и являются достаточно гибкими. Скорость сходимости в таких методах существенно больше скорости сходпмости в старых прямых методах (Уэстлейк [!968[), так как в них используются <разреженные» матрицы. Исторически сложилось так, что итерационные методы чаше применяются в вычислительной гидродинамике и, несомненно, не утратят своей важности. Здесь итерационные методы будут рассматриваться в хронологическом порядке.
3,2,2, Метод Ричардсона н метод Либмана Как обсуждалось ранее в 9 3.1.22, решение (стационарных) уравнений эллиптического типа аналогично получению асимптотнчески стационарного решения нестационарной задачи (Фран- ') Прп вепосредствеппом прпмеиеппв; см, равд 3.2,9, 178 а2. Методы решения уравнений для функции гока кел [!950]). Предположим, что рассматривается нестационарное уравнение диффузии для тр с «источниковым членом» и с коэффициентом диффузии, равным единице: — =ч ф — ь. дФ д! (3.366) Физический смысл нестацнонарного решения здесь не играет роли, но когда решение такого уравнения диффузии приближается к стационарному, оно стремится к интересующему нас решению уравнения Пуассона (3.363).
В некоторых случаях такая аналогия выполняется точно. Для того чтобы продемонстрировать подобную эквивалентность, выведем итерационный метод Ричардсона для эллиптического уравнения Пуассона из нестациопарной схемы с разностями вперед по времени и центральными разностями по пространственным переменным для уравнения диффузии параболического типа. Применяя к уравнению (3.366) разностную схему с разностями вперед по времени и с центральными разностями по пространственным переменным, получаем „~ее~,~ь ' бз,ра б~,~е 1,! с! б! бх' бу' — = — + — — ~, . (3.367) Покажем сначала, что ьг, ! не оказывает влияния на устойчивость уравнения (3.368).
Следуя работе Шортлн и Уэллера (1938], будем рассматривать граничные условия Дирихле и обозначим через ф точное конечно-разностное решение ') конечноразностного уравнения Пуассона (3.364). Тогда ошибка ер ! значения ф~ ! после итерации Й будет составлять а» а е !! '! $!' (3.369) Подставляя полученное отсюда выражение для зуе в уравнение (3.364), будем иметь бззр б'ее бтФ беса — — — + бх' бха бу' буз (3. 370) ') Формальное заключение, что точное решение получается при Ь- са, является скорее предположкительым, чем строгим, предел пш Фл будет а.+ отличаться от точного решения уравнения (З.збб) из-за ошибок округления. Для простоты временно ограничимся случаем, когда Лх = = Лу = Ь и уравнение (3.367) принимает следующий вид: зйал'=зре + —,]зря, .+ фа, +тру, + фа,— 4ф" . — Лз~ (3.368) 179 З.2.2.
Метод Ричардсона и метод Либмана Так как ф в точности удовлетворяет уравнению (3.364), уравнение (3.370) сводится к уравнению Лапласа 5'еа б'еа йлз Ьу» + — = О. (3.37!) Поскольку граничные условия заданы, на всех границах три = или е = О. Тогда аналогично итерационному конечноразностному уравнению (3.368) можно записать для ошибки е следующее уравнение: еа+.'=е,'1+-~-~еь, г +е,'., +ел, +ел, — 4еа 11. (3.372) И Таким образом, итерационное уравнение (3.368) для зр эквивалентно итерационному уравнению (3.372) для е, причем последнее, очевидно, не зависит от Ь.
Условие устойчивости уравнения (3.368) или (3.372) (см. равд. 3.!.5.д) имеет вид с( = аЛ!/Ла ( '/„ или при сс = ! Л! ( Лз(4. Поскольку желательно достичь асимптотического решения как можно скорее, выберем максимально возможное значение Л! = Лз/4. Подстановка этого значения в уравнение (3.368) дает ! трг,+у зри(+ 4 ~три+с 1+ зрг-ь 1+ трг, (тг + + тра 1, — 4«раг — Лзьг 11. (3.373) Приводя подобные члены, содержащие трьг Р получаем тйа ' = 1 У + зйа + тйа + зР— Лз~ 1. (3.374) '! Этот метод известен также как метод Якоби (Сальвадори и Баран [1961)), как метод «итераций полными этапами» (Креиделл [1956) ! и как метод «одиовремеииых смсшеияй» (Яиг [1954)); последнее название связано с тем, что каждое значение хр»»' вычисляется в точке, которую можно брать независимо в любом порядке из последовательности точек (1,1), и поэтому все точки можно в известном смысле рассчитать одиовремеиио.
Такое свойство алгоритма очень важно при оценке методов для вычислительных машин с параллельными процессорами. Этот алгоритм представляет собой метод Ричардсона ') при Лх = Лу. В частном случае для уравнения Лапласа (йсг= 0) он сводится просто к требованию, чтобы значение тр на новой итерации равнялось среднему арифметическому значению в четырех соседних точках. Уравнение (3.374) представляет собой тот же результат, который можно получить, разрешая стационарное уравнение (3.365) эллиптического типа относительно зри; при условии, что этот член в левой части уравнения берется на итерации й + 1, 180 2.2. Методы решения уравнений с)ля функции токо а все члены в правой части берутся на итерации /с.
Введем величину отношения размеров шагов сетки р = Лх/Лйч тогда 1 ""/,/ 2(1+ Р») Ь'Р/Е1 /+ Р/ 1 /+ г "Р/, /Ч-~ + с т/ /-~ »/ /З' (3.375) Этот алгоритм представляет собой метод Ричардсона при Лх Ф Лу. Исследование скорости сходимости методов можно провести так же, как и исследование устойчивости для уравнения переноса вихря, подставляя в уравнение для ошибки (3.372) )/ + =с/ра (3.376) или раскладывая оператор Лапласа в дискретизированной форме по его собственным функциям; результаты оказываются идентичными при условии, что граничные условия правильно учтены'), Представим здесь полученные таким способом (но в других обозначениях) результаты Франкела [1950[.
Компоненты ошибки с наибольшей и наименьшей длинами волн затухают наиболее слабо (т. е. имеют наибольшую величину ) с/(О) [). Таким образом, независимо от начального распределения ошибки по 6 данные компоненты будут доминировать при асимптотичсски больших /с. Обе эти компоненты имеют одну и ту же величину [с/(6) [, но коротковолновая компонента ошибки (Л = 2Лх) совершает знакопеременные колебания при с/ с. О, что в некотором смысле оптимально (Франкел [1950[). Уравнение (3.375) является двухслойным уравнением, и в памяти вычислительной машины необходимо хранить массивы для величин зйаы и з)ь. Если обход расчетных точек вести в направлении возрастания 1 и 1 и в уравнении (3.375) использовать, где это возможно, уже полученные новые значения зрль', то получится следующая схема: 1 тр/,/ = 2(! + 11») [.зр/+ь /+ "р/ ! /+ н тг/ /+/+ + [)хзРа/+/', — Лхзь/ /1, (3,377) которая называется методом Либл/аназ).
При программировании этого метода достаточно в памяти хранить лишь один мас- ') Более нзяшные матричные методы опвсаны в книгах Вазова и Форсайта [1960) н Эймса [1969). ') Этот метод известен также как метод Гаусса — Зейдсля (Сальвадора и Барон [1961)), как метод «птерапий неполными агапами» (Кренделя [1956)) и как метод «последовательных смешений» (Явг [1954)). Этот метод полон/ на методы саульева для нестапионарных уравнений (см.
равд.з.1.17), но пе идентичен им. 1ВГ д 22Н Метод релаксации невлзки Саусвелла сив для величины зр в уравнении (3.377), предусмотрев в программе оператор замещения. Кроме того, Франкел [1950] показал, что для наиболее «стойких» компонент ошибки с большими н малыми волновыми числами в методах Либмана и Ричардсона имеет место соотношение гл (Либмаи) = [6 (Ричардсон)]з. (3. 378) Асимптотнчески при достаточно большом числе итераций й итераций по методу Либмана эквивалентны 2Й итерациям по методу Ричардсона; кроме того в методе Либмана требуется вдвое меньший объем машинной памяти '). Более сложный вариант метода Либмана был предложен Шортли и Уэллером [1938].
3.2.3. Метод релаксации невязни Саусвелла Метод релаксации невязки Саусвелла (Саусвелл [1946] ) применялся в течение многих лет при расчетах вручную для получения численных решений важных технических и научных задач, вклгочая одно пз самых ранних решений задачи о течении при большом числе Рейнольдса (Аллен и Саусвелл [1955]). Первоначально он назывался просто «релаксациоиным методом», но здесь это название заменено на «метод релаксации невязки», чтобы было можно отличать его от метода Либмана и других итерационных методов, которые в настоящее время иногда называют релаксационными методами.
Простейшая форма метода релаксации невязки Саусвелла основана на том же уравнении (3.375), что и метод Ричардсона для вычисления новых значений трьь', Различие заключается с/' в том, что уравнение (3.375) не используется во всех без исключения узловых точках сетки. Здесь невязка гкг определяется из уравнения зРег~ г+зРг г г — 24гг г зРг г ~+зуг г ~ — 2зуг г дх' (3.379) Когда гь г = О, уравнение Пуассона (3.365) удовлетворяется, но только в точке (1,1). Таким образом, величина [гц;[ показывает, как на ошибку в точке (51) влияют все текущие прибли- '1 В деяствнтельностн для методов, подобных методу Ричардсона, потребность в вслнчнне объема памятн ЭВМ можно сократить за счет нехнт.
рых программистских прнемов н использования несколько большего числа операторов замещения. !82 З,2, Методы решения уравнений для функции тока женные значения >[>т,ь Затем перебором находят ту узловую точку, где величина [гь >[ является наибольшей, полагают в этой точке гь ! — — 0 и при помощи уравнения (3.375) находят отсюда новое значение фн ь Это в свою очередь приводит к изменению невязок г во всех соседних точках, и процесс повторяется. Метод Саусвелла не применяется на современных электронных вычислительных машинах, так как время, нужное для нахождения наибольшей величины [ть>[ и для пересчета невязок т в соседних точках, не отличается существенно от времени, необходимого для непосредственного применения схемы (3.375).