Роуч П. Вычислительная гидродинамика (1185913), страница 83
Текст из файла (страница 83)
Метод частиц в ячейках и метод жидкости в ячейках 361 ные осцилляцин величин в ячейках. Эти осцилляции моделируют поведение молекул газа, но при очень малом количестве «вычислителыгых молекула. В обоих рассмотренных методах для аппроксимации конвективных членов используется схема с допорными ячейками (вторая схема с разностями против потока) и, следовательно, в обоих методах имеется схемная искусственная вязкость (см.
равд, 5.5.1, 5.5.2). Джептрп, Мартин н Дали [1966] указали, что наличие в обоих методах искусственной вязкости г) — [и] означает пеипварнантность искусственной вязкости относительно преобразования Галилея, т, е. невозможность использования в этих методах преобразования, состоящего в обращении потока '). Кроме того, как отметили Эванс и Харлоу [1958, 1959], а также Лонгли [1960], без введения явной искусственной вязкости метод будет локально неустойчив в точках торможения потока, так как здесь схемная вязкость г) — [и] стремится к нулю; см. также формулу (5.25) н далее.
В исходных работах оба метода были записаны как в декартовых, так и в цилиндрических координатах. Метод частиц в ячейках становится наиболее эффективным прн решении задач с поверхностями раздела (свободные поверхности и многокомпонентныс среды), поскольку отдельным частгщам можно приписать различные массы, удельные теплоемкости и т.д. в целях моделирования двухжидкостной среды, свободной поверхности жидкости и даже границы жидкости с деформируемым телом. За годы успешных приложений этого метода постепенно разрешались вопросы, связанные с пустыми ячейками, граничными условиями, деталями процедуры осреднения параметров частиц (Эванс и Харлоу [1957, 1958, 1959], Эванс с соавторами [1962], Харлоу [1963, !964]). Обзор этих методов был сделан Амсденом [1966]. Мадер [1964] в разработанном им методе взрыва в ячейках (метод Е1С) обобщил данный подход, рассмотрев течение с химическими реакциями, Херт [1965] также по методу частиц в ячейках проводил расчеты детонационной ударной волны при взрыве.
Дикмен с соавторами [1969], а также Морс и Нилсон [1971[ при помощи метода частиц в ячейках рассчитывали неустойчивость плазмы. Армстронг н Нилсон [1970], рассматривавшие нелинейное развитие сильной двухпотоковой неустойчивости в плазме, продемопстрировалп хорошее согласование результатов нестационарного расчета по методу частип в ячейках и расчета по методу интегральных преобразований. Точность рассматриваемого подхода была также продемонстрирована с ') Это справедлнво также длв любого метода, в котором применяЕтся гхема с коне шымп разностямн против потока. 362 б.а Схемы с каленой искусственной энэкостыо помощью некоторых алгоритмов типа метода частиц в ячейках для мпогокомпонеитных сред, разработанных в Р)туз)сз 1п1егпаБопа1 Согпрапу (Бакипгем и др.
[1970], Уотсон и Годфри [!967], Уотсон [!969] ). Амсден и Харлоу [!965] рассчитали крупномасштабпые характеристики сверхзвукового турбулентного течения в донной области. Крейн [1968] пытался точно рассчитать гиперзвуковое течение о ближнем следе за телом, примепяя метод частиц в ячейках к уравнениям для певязкого газа; однако этот метод не подходит к подобной задаче, и поэтому расчет окончился неудачей. Точность метода жидкости в ячейках независимо подтвердили Гурурая и Деккер [1970] на расчетах сложных двумерных задач о распространении ударных волн и Сатофука [1971] на расчетах двумерных плоских и осесимметричных задач о течении в ударной трубе.
Алгоритм типа метода Р).1С под названием Т011 был разработан Джонсоном [!967]; см. также работы Хилча и Ларсена [1970] и Рейнольдса [1970]. Ссылки на другие работы с расчетами по методам Р!С и Р(.1С, выполненными в Лос-Аламосской лаборатории, можно найти в работе Харлоу и Амсдена [!970а]. Батлер [1967] включил в алгоритмы методов Р1С и Р1.1С учет вязкости и теплопроводпости и обнаружил, что при этом оба метода дают сопоставимые результаты.
б,б,4, Схема Лампа 17т = 1)с — б! — ] н+з о бт 1 » бх ]з (5.31) он заменил (в правой части этого уравнения) величину ()т ее средним по пространственной переменной значением на п-м слое по времени: э+1 1 л э др !л (7, = — ((Л„+ (7,,) — Л1 — ]! . 2 бх ]с' (5.32) ') Обычно эта схема называется схемой Лакса. Впервые она была опуб. лнкована в открытой лнтсратуре в работе Куранта с соавторами [19621 (в подстрочном примечании) н названа там схемой Келлера н Лакса, Рнхтмайер 119631 в связи с этой скемой упоминает также Фридрихса.
Эта схема ') была опубликована Лаксом [!954] в его фундаментальной работе, посвященной консервативным уравнениям, Ланс в первую очередь интересовался законами сохранения и лишь во вторую очередь — конечно-разцостными схемами. Для устойчивости расчета одномерного течения невязкого газа по уравнениям (4.66) при помощи схемы с разностями вперед по времени и центральными разностями по пространственным переменным типа 5.5.4. Схема Лакеи 363 Эта простая и исторически важная схема обладает рядом поучительных свойств. Производные по пространственным переменным аппрокспмируются здесь центральнымп разностями, что обеспечивает второй порядок точности, однако схема обладает диффузней (Рихтмайер [!963) называет ее диффузионной).
Рассмотрим модельное уравнение (5.1) прп а = О; для нсго схема .г!акса дает дх и (5.33) Исследуем схему на устойчивость по методу Хсрта (см. равд. 3.1.5.в), проводя разложения в ряды Тейлора: ~+ ~~ Л1+ —,' ",",,' Л1а+ 0 (Л!') = деи 2Е' дх 2 дх' ' О дх' — й — ~ + 0 [Лт') (5.34) или С учетом соотношения деи д ди - д ди д'и получаем ди ди г ахе Ы ,~ деи — = — й — + ( —, — — й-х! — + 0(Лх-) .
(5.37) дГ дх (,2М 2 ) дхе Из этого нсстационарного анализа видно, что в схеме Лаков вводится эффективный коэффициент искусственной диффузии или Лхе а = — [1 — Са] 2аг (5.39) Как обычно, для устойчивости расчета по модельному уравнению требуется, чтобы а. ) О или С ( 1.
При С = 1 получается точное решение молельного уравнения. Поскольку схема применяется ко всем переменным У=(р,ри,Е,), искусственная 35. Следы с неявкой искусственной вязкостью диффузия вводится не только в виде искусственной вязкости, но и в виде искусственных диффузии массы и теплопроводности '). Из равенства (5.3?) следует, что ошибка аппроксимапни в этой схеме является величиной (5.40) Е = О (Л., ЛЛ ЛхзгЛ1), Такое выражение показ1лвает, что ошибка аппроксимации неограниченно возрастает прп фиксированном Лх п Л1 — »О. Этот факт имеет большое значение. Именно он привел в сильнейшее замешательство автора настоящей книги, когда он случайно запустил программу расчета распространения ударной волны с шагом Л1 = 0 и прп этом обнаружил, что ударная волна все же распространяется (рассмотритс выражение (5.32) прп Л1= 0).
В действительности н этом случае возмущение пе распространяется вместе с фронтом волны прн его движении, а диффундирует от начального разрыва при г = О. Очевидно, что при достаточно малых Лг в схеме Лакса получается большой коэффициент а„ достаточный для стабилизации расчетов течений с сильными ударными волнами. ПриС=( демпфирование исчезает и схему нельзя применять для расчета течений с ударными волнами. Схема Лакса очень просто обобщается на случай двух или трех пространственных измерений, а именно оул и",,,','=ф?~,ь,+(?", и,+и,,„,+(?сг Д+Л~ — ',",, (5,4П л+1 1 Г л л л л л иид.--,"(?„ьг,.+и.
сл,+ис1+с,+ист с.+исг„,+ кг?л + (/1, л а 1) + М вЂ” . (5.42) При этом соответствующие условия устойчивости имеют вид (5.43) С -(р+ )Л ъ'дхау +ауде+халк <) зв а (5.44) Таким образом, в трехмерном случае при Лх = Лд = Лг макси- мально возможный шаг по времени умепыпается в хгЗ раз. Упражнение. Получить выражения дхя а. в случае двух илн грех про. странственных переменных.
') Диффузионная схема дважды нарушает свойство траиспортнвности. Если, например, схема «чехарда» (см. разд, Зд.б) переносит возмущения вверх по потоку против направления скорости, то диффузионная схема пере. носит их, кроме того, и в направлениях, перпендикулярных скорости. 5.5.5. Схелгп .токга — Вендроффа йядажяепис Определить условия, при когорьы схема Русакова сводится к схеме Лакса. Моретти и Аббетт [1966а] применяли двумерный вариант схемы Лакса в сочетании с методом характеристик и выделеиием скачков для попытки расчета течения в доипой области.
Оии обнаружили явлеиие, названное имп «з!а11!пд» и состоящее в следующем. При наличии градпеита по пространственным перемеииым, такого, что (г7 Ф (г7 = ~(4 ((г7тн г + (7'-н г+ (тг г ь1+ (I";, г 1), (545) зависящее от времени решеиие определяется уравнением (7' — (("; = (б(l";(б() ЛБ (5.46) я-~! я так что ((г ~(l, для всех л.
Положение дел может измениться при измеиеиии Ай Конечно, схема Лакса пе прсдвазиачеиа для расчета дозвуковых безударных течеппй, однако этот пример демонстрирует еще один иедостаток схемы. Несмотря иа свои недостатки, схема Лакса имеет важное преимушество — простоту, Оиа также легко приспосабливается к цилипдрическим и сферическим координатам, а также к трехмерным задачам.
По-видимолчу, эти обстоятельства и явились основной причиной того, что ее применяли Богачсвский и Рубии [1966], Богачевский и Мейгс [1966], Богачевский и Костофф [1971], Бариуэлл [!967], Ксерикос [1968], Эмери и Ашерст [1971]. Кеицер [1970б] исследовал схему Лакса и схему «чехарда» (см. равд, 3.!.6) в различных комбинациях и при различных шагах по времени применительно к двумерной задаче, в которой ударная волна принималась за границу расчетной области. Из-за простоты программирования и надежности схема Лакса может применяться иа ранних стадиях разработки алгоритма решеиия задачи с тем, чтобы впоследствии заменить ее более сложной схемой.
Упражнение. Показать, что в модельиом уравнении представление коивектпвиого члена по схеме Лакса, а диффузиоииого члена центРальными разиостями приводит к безусловио иеустойчпвой схеме. Указание. Исполь. зовать иссчедоваипе схемы с разностями вперед по времеви п центральными разиостями по пространственной переменной, заменив а иа и + а,. Упражнение. Показать, что в схеме Лакса для стационарного уравиеиия получается величива и, = Ьхт/(2Ы). б.б.б. Схема Лаиса — Веидроффа Лакс и Веидрофф [1960, 1964] разработали класс схем, который позволил достигиуть значительного прогресса в теоретинском изучеиии разиостпых методов и привел к разработке ББ. Схемы с неявной искусственной вязкостью звб класса двухшаговых схем (см.