Федоренко Р.П. Введение в вычислительную физику (1185915), страница 72
Текст из файла (страница 72)
Лагранжевы координаты здесь неудобны: отображение (г, х, у) «е(г, $, г)) разрывно. 350 пгивлижвиныв мвгоды вычислительной +изнки [ч. и ем времени теряют точность до такой степени, что расчет не заслуживает никакого доверия. Первоначально прямоугольные ячейки постепенно деформируются, теряют даже форму параллелограммов, а в особо сложных случаях часто наблюдается «выворачиваниеь лагранжевой ячейки. Р1С-метод (метод частиц в ячейке), Р1С-метод (Рагйс!е-1пСей) был предложен и разработан Ф. Харлоу (Лос-Аламос) в 1955 г.
и является типично американским. Этот метод требует большого обьема оперативной памяти и большою числа арифметических операций, поэтому в нашей стране он нашел применение лишь в последние годы — после ввода в эксплуатацию мощных машин серии ЕС. Теоретическое обоснование метода, видимо, до настоящего времени еще не получило законченной формы, а первоначально метод имел чисто прагматическое оправдание.
С его помощью были проведены расчеты очень сложных течений. Хотя точность проведенных расчетов признается не всеми и большинство вычислителей согласны с тем, что она не очень высока, качественная картина, полученная в расчетах, выглядит убедительной и правдоподобной. А если учесть, что попытки проведения подобных расчетов другими методами приводят обычно к явно недостоверным результатам, легко понять популярность Р1С-метода. Характерная ситуация, в которой Р1С-метод демонстрирует свои преимущества перед иными, — это течение среды, первоначальна разделенной на простые по форме области, заполненные разными веществами.
В процессе течения эти области сильно деформируются и перемещаются в просгранстве. Для расчета подобных течений неудобны как эйлерова сетка (происходят большие перемещения контактных границ в геометрическом пространстве), так и лагранжева (происходит сильная деформация первоначальной лагранжевой сетки), Р!С-метод — это попытка совместить достоинства эйлерова и лагранжева описаний течений сплошной среды. Перейдем к описанию вычислительной схемы метода. Эйлерова сетка. Область течения покрыта неподвижной в пространстве сеткой. Для простоты будем считать шаг сетки Ь постоянным, одинаковым по х, у. Ячейки сетки занумерованы парами индексов й 1.
Величины, которые в дальнейшем помечиотся индексами г, 1, трактуются либо как относящиеся к ячейке в целом, либо (при разностной аппроксимации уравнений) как относящиеся к центру ячейки. Основные счетные величины В центрах счетных ячеек определены Величины и"... и", — компоненты скорости среды (индекс л показывает их принадлежность,ко времени 1„). Совокупность чисел го, г„..., г„, ... образует временную сетку, вообще говоря неравно- гашение лвхмьтных зхллч гхзовоя линлмнхн 351 меРнУю; г„,, = г„+ г„+ из (шаг т«+ из выбиРаетсЯ на основе сообРа- жений точности и устойчивости в зависимости от состояния среды в момент г„). Кроме того, в ячейках (1, /) определены величины (Е„" .
Р М'„',,) (а = 1, 2, ..., А). Поясним их смысл. Иногда ради удобства мы будем опускать ин- дексы г, /, но не следует забывать, что этн величины — свои в каж- дой ячейке. Напомним, чзо в задаче изучается течение в области, заполненной в разных частях веществами с разными физическими свойствами (а — номер вещества). Вообще говоря, в данной ячейке может быть либо одно вещество, либо несколько разных. В каждой ячейке (й /) в данный момент времени г„индекс а пробегает свой (зависяший от г, /, и) набор значений.
Проще будет считать, что индекс а пробегает все допустимые в данной задаче значения. Од- нако если вещества с номером а в данной ячейке нет, то соответст- вующие значения Е" = М" = О. Физический смысл этих величин таков: Е" .. — удельная внутренняя энергия вещества с номером а, «, ь/ в момент времени г„находящегося в ячейке (й /); М'„' . — масса этого вещества.
Кроме переменных и, щ Е, М, в расчете участвует большое число «частиц». Будем нумеровать эти частицы индексом /г. Чис- ло частиц должно быть много ббльшим числа ячеек: на каждую ячейку в среднем должно приходится, как минимум, пять-десять частиц.
Каждая й-я частица в момент времени г„характеризу- ется следующими величинами: Х", У~ — координаты положения частицы; нзг — масса частицы (не зависящая от г); ૠ— но- мер вещества, из которого «состоит» частица с номером (к=1,2,...,К). Это — основные счетные величины, полностью характеризую- щие (в принятой расчетной модели) состояние среды. Все-остальные величины, которые появятся в дальнейшем, носят вспомогательный характер н выражаются через основные.
Стандартный шаг интегри- рования задачи состоит в переходе от величин (и, и, Е„, М )," (Х, У]~ к величинам (и, и, Е,, М„)«~', (Х, УЦ+' (переход на сле- дующий временной алой), Расщепление уравнений газовой динамики Математическая задача, которую предстоит решать, состоит в интегрировании уравнений газовой динамики, записанных в эйлеровой дивергентной форме (3).
Обозначим плотность полной энергии и = р (е + (и + из)/2), где е — удельная внутренняя энергия, связанная с р и р уравнением состояния р= Р (е, р), своим для каждого типа вещества. 352 пгизлижзнныз мвтоды вычислительной «взнки !ч. и Р1С-метод применяется для расчета быстрых процессов, в которых диффузия не играет заметной роли. Поэтому перемешивання вещества не происходит. Считается, что на протяжении всего времени расчета сохраняются четкие границы, разделяющие разные вещества. Но форма этих границ претерпевает существенные изменения. Может измениться даже нх «топологня»'. если при 1 = О какоето вещество заполняло связную область, в дальнейшем оно может распасться на отдельные части, разделяемые веществом другого типа. Аккуратный расчет таких течений требует знания в кахсдый момент времени положения контактных границ.
Это очень сложная вычислительная задача, и Р!С-метод является попыткой решить ее относительно простыми средствами. Разумеется, эта просюта оплачивается большим обьемом памяти н машинных операций. Мы не будем обсуждать важных вопросов, связанных с граничными условиями, но скажем несколько слов о начальных данных. Физическая постановка задачи обычно связана с заданием начальных данных в виде функций вв(х, у), св(х, у), р~(х, у), ез(х, у) и границ, разделяющих разные вещества. В начальный момент времени эти границы, как правило, имеют простую геометРическУю фоРмУ. Начальные данные длв Расчета из..
во у свЯзаны с ив(х, у), ев(х, у) — это просто значения в центрах ячеек. Если ячейка (1, у) целиком заполнена (допустим, веществом а= 1), то. Ео в ней очевидным образом связано с ез(х, у) в центре. Остальные значения Е« = О. Если через ячейку проходит граница раздела двух веществ (т,е. в ней есть вещества двух типов), то и значения Е~~,. у задаются равными значениям в начальных данных (с точностью до шага сетки Л). Напомним, что Š— это удельная энергия вещества, а не количество энергии данного вещества в ячейке.
Начальные положения частиц (Хы Уо«) задаются так, что, например, в каждую ячейку попадает равное число частиц (впрочем, это ие обязательно, иногда полезно увеличить число частиц в тех ячейках, в которых ожидаются наиболее сложные события). Распределение т должно быть определенным образом согласовано с рз(х, у). Предполагается, что величина М'*, согласована с положениями частиц, а именно: М' =~х" ты Ус: (Х«ы У,",) Е С,, а=а».
(9) В дальнейшем нам часто придется иметь дело с суммамн подобного рода. Условный смысл суммы в (9) — это суммирование величии т„ для частиц с номером к, координаты которых (Х", У") в момент 1„ 353 гвшвние двзмытых злдлч глзовой див»мики 3 231 попали в ячейку (г, г) н которые имеют тип аг = а.'Если известно начальное распределение рв(х, у), то известна масса вещества в каждой ячейке.
И если через ячейку проходит контактная граница, то известно, сколько вещества М",, каждого типа находится в данной ячейке (г, у). Задав положения частиц, нужно приписать им массы т так, чтобм выполнялось соотношение (9) при п = О. Оно будет, как мы увидим, выполняться и в дальнейшем. Уравнения (3) описывают изменения основных физических величин (масса, импульс, полная энергия) за счет процессов двух типов — работы сил давления н перетекания физических величин со скоростью потока ( и, о) . В Р1С-методе переход за малое время т от величин на временном слое и к величинам на слое л + 1 осуществляется в два этапа, на каждом из которых основные физические величины меняются за счет процесса только одного типа. 1.
На первом этапе учитываются изменения основных величин только за счет работы сил давления (процессы перетекания пока исключены). Разностные формулы на этом этапе аппроксимируют (в привычном, наглядном смысле этого алова) следующие уравнения: р,=О, (ри),+ р„=О, (10) (ро),+р =О, в,+(ри) +(ри) =О. 11. На втором этапе используются уравнения, в которых, наоборот, оставлены только процессы перетекания: р, + (ри)„+ (ри) =О, (ри), + (рии)„+ (рии) =О, (ри), + (рио)„+ (рси) =О, ю, + (ыи)„+ (во) = О. Разностная аппроксимация процессов переноса осуществляется с помощью частиц и имеет, как мы увидим, не очень привычный для метода конечных разностей характер.
Скорее, здесь используются методы дискретного моделирования сплошной среды, апеллирующие к основным понятиям механики. Численная реализация лервого этапа. Исходная информация состоит в каждой ячейке (О У) из величин и", о", М'„', Е,". Расчет начинается с вычисления давления рг . Если в ячейке имеются частицы нескольких типов, то прн расчете давления рг . 'используются «физические» соображения о равенстве давлений на границе двух сред.
Введем величины о — части объема лз ячейки, занимаемые веществом типа а. Очевидно, ~ о„= Ь~. Зная массу М"„вещества типа а, находим его плотность р„= М"/о„, а зная его внутреннюю >3 нгивлнженные методы вычислительной визики 1ч. и энергию Е'„' вычисляем давление р = Р (Е",, М"/и). Приравнивая величины Р, друг другу, получаем систему уравнений типа Р,(Еы М1/п~) = Рз(Ез, Мз/пз) = Здесь число уравнений на единицу меньше числа веществ в ячейке. Присоединяя к ним уравнение ,"Р о, = йз, приходим к полной системе уравнений относительно неизвестных о'„.
Она решается итерационным методом. Заметим, что в распространенном случае, когда уравнение состояния имеет вид р /,(е) р, можно выписать явное решение. После того как о, найдены, определяется величина р, которую мы обозначим Р," ., приписав ее центру ячейки. Теперь у нас есть все для того, чтобы рассчитать первый этап по стандартным разностным уравнениям, аппроксимирующим уравнения (10). Результатом будут величины йсл й, Р М,, Е Уравнение р, = 0 аппроксимируется просто: М,,=М:,, Определяем Р",+и =0.5(Р" +Р",+~ ), Р~ .+из=О.5(Р~~ +Р," -+~), Мп ~ч~~ Мв в Мв /ьз а Так как р, = О, то используется аппроксимация (13) Из этих уравнений в явном виде находим йь, йь/ Уравнение изменения энергии аппроксимируегся следующим образом: + — „(Р7~+пз""ь,Дд — Рь/ нее",.+,."знз) =О. (14) Некоторые величины в этой формуле требуют пояснения: и"++Я/ =' (1/4)(й + и", + й„, + и"+ .); аналогично вычисляется н", +/+пни 355 гешение дзхмееных злдлч глзовой динлмики й 331 Сложнее обстоит дело с величиной и",, Напомним, что плотность полной энергии н= р)е+ (из+ вг)/2).