Роуч П. Вычислительная гидродинамика (1185913), страница 23
Текст из файла (страница 23)
Но при применении только к конвективным членам (т. е. прн а = О) эта схема, называемая схемой со средней точкой (Лилля (!965)) или схемой «чехарда со средней точкой» нли — чаще всего — просто «чехарда> (Курант, Фридрихс, Леви [1928) ), обладает нужными свойствами устойчивости. Таким образом, уравнение дй д (ий) д( дх (3.144) ') Схемы, в которых в правую часть уравнения входят значения Цг+г и ь, и называются неявными, и при этом для вычисления значений лч! лгг на новом временнбм слое, вообще говоря, требуется обращение матрицы. г) Схемы, которые в этой книге именуются однослойными, некоторые авторы называют одношаговыми.
по этой схеме представляется в виде лег л-Г л л (3.145) 2о( 2 ах Данная схема имеет второй порядок точности по пространству и по времени; зто одношаговая явная трехслойная по вре- 8.д Методы решения уравнение переноса вихри мени схема. Значит, для вычисления новых значений на слое п+ 1 в этой схеме необходимы значения на слоях и и и — 1. Заметим, что новые ь на четном временном слое вычисляются как значения ~ на предыдущем четном временном слое плюс некоторое приращение, а предыдущий нечетный временной слой при этом «перепрыгивается» (отсюда и название схемы — «че.
харда»). »та+1 т»п-1 С (ьт» ьтп ) (3.145) где С = пб1/Лх — число Куранта. Подставляя сюда фурье-ком- поненты, получаем связь для амплитуд (ти~-1 ) п+ 1,н-1 (3.147) где а = — 27С з(п О. (3.148) Чтобы получить матричное уравнение, добавим к (3.147) тождество 'е'"=1 ° Рн+ 0 ° )ти '. (3.149) Рассматривая это уравнение совместно с уравнением (3.147), получаем (3.150) где множитель перехода 6 теперь представляет собой матрицу -ИЛ (3.151) Для изученной ранее однослойной схемы множитель перехода 6 был просто числом, а условие устойчивости имело вид ~ 6( ( 1.
В данном случае, когда 6 представляет собой матрицу, условие устойчивости (согласно фон Нейману) имеет следующий вид: (3.152) Метод фон Неймана исследования устойчивости для этой и других многослойных схем применяется следующим образом. Используя те же определения и предположения, что и в предыдуших примерах, уравнение (3.145) можно записать в следующем виде: Зд 6. Одношаговые явнвге схемы вт де Л вЂ” все возможные собственные значения матрицы 6 '). ,Собственное значение Л матрицы определяется как корень характеристического уравнения, которое получается приравниванием нулю определителя матрицы, у которого из каждого диагонального элемента вычитается Л.) Таким образом, характеристическое уравнение для определения матрицы 6 записывается так: — в (3.153) Когда 6 представляет собой просто число, как в предыдущих римерах, оно рассматривается как одномерная матрица.
Тогда равнение для определения собственного значения принимает нд 6 — Л =0 илп Л =- 6, поэтому условие (3.152) сведется , предыдущему условию (3.102), а именно ~ 6! =!.) Раскрывая пределнтель и решая полученное квадратное уравнение для Л, находим два решения: Л = й).а~ 'у'а'+4). (3.154) Учитывая, что а = — 2г'С э!п 0 и ат = — 4С з!па О, имеем в,= — гс ь вь вз."с"~в. (3.155) В тех случаях, когда Ст э!пз 0 ) 1, квадратное подкоренное выражение будет отрицательным, и тогда Ле — — !1 — С з!и 0 ~ тУСз вбп'0 — 1]. (3,156) При этом абсолютная величина )Л!) 1, что означает неустойчивость.
В том же случае, когда Са з)пз 0 = 1, для чего, вообще говоря, требуется выполнение условия С ~~1, (3.157) для абсолютной величины Л получаем )Л !'=С'э)п'О+(1 — С'з!п'О), (3.158) 1 Лм ! = ! при С ~(1. (3. 159) Это, очевидно, удовлетворяет условию устойчивости (3.152) в предельном случае равенства. Аналогичный результат получается также при двух пространственных переменных, но здесь для устойчивости требуется выполнение неравенства С„+ С„= 1, ') Часто для удобства это условие формулируется так: для устойчивости схемы спектральный радиус матрипы 6 не должен превышать единицы, т, е.
р(б) ( К где р(б) = шах)Л,), а Лр есть р-е собственное значение 6. Спектральный радиус, очевидно, есть радиус круга в комплексной плоскости, центр которого находится в точке (ц О) и внутри которого лежат все собстненпые значения. 88 Зд. Методы решения рраенения лереноеа вихря Можно думать, что даваемое (3.159) значение !Х )= 1, отвечающее границе устойчивости, является приемлемым лишь в крайнем случае, но в действительности это значение даже весьма желательно, если решение исходного дифференциального уравнения не является затухающим.
В самом деле, уравнение конвекции (3.144) при отсутствия вязкости и постоянном и выражает тот факт, что произвольное начальное распределение функции ~(х,0) просто сдвигается со скоростью конвекции и; значит, для любого сдвига т по времени решение этого уравнения имеет вид ~ (х, 1+ т) = !", (х — ит, 1). (3.160) Таким образом, метод фон Неймана показывает, что схема «чехарда» правильно моделирует одно пз свойств, присущих решению исходвого дифференциального уравнения при и = = сопз(, а именно отсутствие затухания.
Любая разностная схема для решения уравнения для невязкой жидкости, такая, что ~0)(1, обладает ошибкой, обусловленной искусственным или численным затуханием, В любой сходяшейся разностной схеме эта ошибка должна, конечно, стремиться к нулю при Лх-+-0, Л1-ы0, но применение метода фон Неймана показывает, что схема «чехарда» при и = сопз! и С ~ 1 обладает нулевой ошибкой, обусловленной затуханием, даже при конечных Лх и Лйь Действительно, в частном случае С = 1 рассматриваемая схема дает точные результаты.
Полагая т = Лй имеем х — ит = х — СЛх, поэтому при С = 1 решение (3.160) можно переписать так: »та+1 т»л (3.! 61) Следовательно, точное решение дифференциального уравнения, если его рассматривать в узловых точках конечно-разностной сетки, выражает перенос величины ь нз точки с' — ! на слое и в точку ! на слое п+!. Величина ь за один шаг по времени переносится с конвективной скоростью и на расстояние иЛ(, а при С = 1 расстояние иЛ! равно Лх. Через два временных шага точное решение будет »та+1»тл-! (3.! 62) В результате применения схемы «чехарда» (3.!46) при С=1 получаем ~л~-1 л-1 ~п + ~л (3.! 63) Задав правильные начальные значения согласно уравнению (3.161), а именно ~",.+,— — й",-' и Ц,=~",—,', получаем, что разностное решение (3.163) совпадает с точным решением (3.162).
Таким образом, при и = сопз1 и С = 1 схема «чехарда дает точное решение 89 З.пб. Одношаговые явные схемы Однако при практических гидродинамических расчетах, когда скорость меняется в пространстве, ограничение на шаг Л1 будет определяться (если не учитывать дополнительные усложнения, обусловленные нелинейностью) наибольшим значением и в узловых точках сетки. Значит, вообще говоря, нельзя получить С = 1 во всех точках. Но при С ( 1 разностная схема «чехарда» уже не будет давать точного решения.
Прежде всего вопреки результатам метода фон Неймана рассматриваемой схеме может быть присуще некоторого рода численное затухание, хотя это обычно не допускается. На рис. 3.10 представлены трехмерные графики величины 9(х,(), рассчитанной по схеме «чехарда» при синусондально меняющейся на входной границе потока величине ~(0,1) = з(п й При С = 1, как видно на рис. 3.10,а, получается точное решение с синусоидальным законом на входной границе, которое переносится за счет конвекцни без затухания. На рис.
3.10, 6 построена диаграмма, рассчитанная с вдвое меньшим шагом по времени, т. е. при С = '/а. Ясно видно, что в этом случае максимум амплитуды первого горба уменьшается по мере движения вниз по потоку. Рис. 3.10, и снова соответствует С = '/з, но период изменения ~ на входной границе и величина шага по времени выбраны таким образом, чтобы их отношение (и шаг Лх) были такими же, как и для рис. 10,и; в этом случае затухание очень сильное. Здесь возникает вопрос относительно применения термина «затуханием Метод фон Неймана стал настолько широко известным, что «затухание» обычно понимается в смысле поведения фурье-компонент как случай, когда ~Х~ ( 1, Очевидно, если термин «затухание» означает по определению, что 1Х1 ( 1, то схема «чехарда» по определению не обладает затуханием.
Уменыпение же экстремальных величин амплитуды, которое видно на рис. 3.10, 6 и 3.10, в, правильно связывать с дисперсионной ошибкой '), которая проявляется в методе фон Неймана и будет кратко обсуждаться в дальнейшем. Такая терминология целесообразна, н ее можно даже рекомендовать в тех случаях, когда имеется в виду, что схема в самом деле может вызывать уменьшение экстремальных величин амплитуд волн.
В обычной речи такое свойство называется «затуханием»; что же касается рис. 3,10,в, то лучше было бы говорить, что волна ') Днсперсиовную ошибку мы будем обсуждать на стр. !23 — 124, а пока ограничимся кратким замечанием. Ошибка в волновой скорости различна для различных фурье-компонент, но»тому каждая из них имеет рззличную скорость конвекции. Таким образом, фурье-компоненты исходного распреде. лепна имеют тендениии «размазываться» нли днспергировать в пропсссе решения. Эта ошибка обычно больше для компонент с меньшей длиной волны и, конечно, зависит от чвсла Куравта, причем при С = 1 она не возникает.
90 3.!. Методы решения уравнения переноса вихря Рис. 3.10. Решения уравнения дь/д/ = — идь/дх, полученные при помоши схемы «чехараа». Здесь С вЂ” число Куранта, Ж вЂ” сеточная частота. Диаграммы любезно предоставлены У. Сандбергом из лаборатории Сандиа. а:С = 1, /У = 8; б; С = 0.5, /хх = 16. ДД Гнетодм решения уравнения переноса вихря не затухает, а уменьшает свою амплитуду по мере ее переноса за счет конвекции. Основываясь на рис. 3.!О, можно сделать и другое очень интересное заключение. Обычно не обращают внимания на то, что при решении конечно-разностных уравнений для задач, аналогичных представленной на рис. 3.10, существует два характеристических параметра.
Первый параметр представляет собой число Куранта, которое является единственным параметром при решении конечно-разностного уравнения во внутренних точках. Вторым параметром является сеточная частота М = = 2п/Л(, т, е. число временных слоев за период изменения функции на входной границе потока. Сравнение рис. 3.10,б и 3.10,в приводит к выводу, что для фиксированного С ( 1 затухание (уменьшение экстремальных амплитуд) ослабляется с увеличением сеточной частоты в том случае, когда сеточная частота М представляет собой целое число.