Fedorenko-RP-Vvedenie-v-vychislitelnuyu-fiziku (810773), страница 77
Текст из файла (страница 77)
Слабым местом такого подхода является то обстоятельство, что часто рассчитываемое течение имеет разные характеристики гладкости в разных частях области. Адаптирующиеся сетки могут это учитывать (хотя и не без определенных трудностей). Расчет же на равноь ерной сетке, размер которой диктуется наиболее «узким» местом, требует слишком малого шага. Поэтому, хотя это направление начало развиваться в нашей стране с первых лет работы на ЭВМ, на некоторое время оно 9 гз( гашение двгмкгных злдлч глзовой динлмики 369 был6 оставлено. Это было связано, видимо, с малым объемом оперативной памяти ЭВМ того времени.
ДО недавнего времени наиболее распространенной ЭВМ была БЭСМ-б, имевшая оперативную память в 32 000 слов. После введения в эксплуатацию новых ЭВМ с большими ресурсами оперативной памяти (порядка 10ь) интерес к этому направлению стал возрождаться. Конечно, число узлов равномерной сетки должно быть существенно больше числа узлов адаптирующейся, но расчетные формулы оказываются намного проще. По оценкам специалистов, программа трехмерной газовой динамики с адаптирующейся сеткой тратит порядка 104 операций на один узел (прн расчете одного временного шага); программа же, основанная на прямоугольной декартовой сетке, около 10з операций.
Смогут ли такие простые программы конкурировать с расчетами на адаптирующихся сетках — покажет будущее. Работа здесь в сущности только начинается. Ниже мы опишем некоторый способ учета криволинейной границы на прямоугольной декартовом сетке, который продемонстрирует, какого сорта проблемы здесь возникают. Аппроксимация около криволинейной границы. Рассмотрим течение газа, описываемое уравнениями в форме Эйлера в декартовой системе координат. Граница, на которой поставлено условие не- протекания (контур обтекаемого тела), является относительно гладкой н проходит более илн менее произвольно относительно системы координат.
Введем 'прямоугольную равномерную сетку с шагом Ь по обоим направлениям (для простоты). Вертикальные и горизонтальные линии сетки занумеруем целыми индексами 1 и 1, а различные счетные величины будем вводить в разных точках сетки. Определим в центрах ячеек сетки термодинамические неизвестные р",+„„,.+„, р", х,. „, Определим компоненты скорости иа серединах сторон ячейки и",... е",+,, Такая «шахматная» сетка удобна для аппроксимации оператора дивергенция в точках (1+ 1/2, /+ 1/2) и йтад р, причем х-компоненту удобно аппроксимировать в и-точке (где она только и нужна), а у-компоненту— аютвегственно в в-точке.
Среди введенных формально узлов квадратной сетки выделим счетные узлы, т.е. те, в которых соответствующая величина считается определенной. Узел сетки считаем счетным, если квадрат Ь х й с центром в этом узле целиком помещается в области течения (не пересекается с обтекаемым телом). Если такой квадрат частично находится в области потока, частично — внутри тела, он считается фиктнвныло ему не соответствует никакая величина, входящая в основные массивы счетных величин. Однако такая величина в фиктивном узле может появиться как вспомогательная, вычисляемая через основные, зто пгнвднженные методы вычислительной «нвнкн ч.п Рисунок 41 поясняет сказанное.
На нем в части области показаны сетка, счетные (р, р)-узлы (темные кружки), счетные и-узлы (Знаки «минус») и о-узлы (крестики). Указана граница и расположенрые на ней нестандартные счетные узлы. Будем считать, что рассматриваемый участок А границы задается кривой У(х), причем :«. ° 1У,(х)! и1. (При ~У,(х)~ >1 контур * ° следует задавать функцией Х(у) и в последующем поменять роли х и у.) Итак, на контуре вводятся (р, р)-узлы с координатами хв, з = (1+ 1/2)Л, У(х1+ 1п), в которых определены счетные Гнс. 41 величины р"; 1лн рв+нз (отсутствие второго индекса — признак принадле;кносги границе). Кроме тою, на контуре вводятся н узлы с координатами х1, У(хв), в которых определена счетная величина нл1, имеющая смысл касательной к контуру компоненты скорости.
Опишем структуру одною шага интегрирования по времени в явной схеме. Шаг состоит из трех основных операций. 1. По значениям величин в счетных узлах (внутренних и граничных) производится интерполяция соответствующих величин в фиктивные узлы. Тем самым создается ситуация, позволяющая во внутренних счетных узлах рассчитать величины на следующем шаге по стандартным формулам. 2. Во внутренних счетных узлах вычисляются значения на верхнем, (л + 1)-м, слое по стандартным формулам. 3. Расчет граничных значений р"; 111, рв+11ди 1«", ' производится по специальным формулам.
Конкретизируем в1ш расчетных формул. Вычисление какой-либо величины в фиктивном узле производится с помощью линейной интерполяции по значениям соответствующей величины на той же вертикальной линии сетки. Используются последняя внутренняя точка и граничная. Если в граничной точке нет нужной величины, например и",+,, она вычисляется как 0.5(нл + нв,) соз ив+и, где а,»пз — угол наклона границы. Уравнение для плотности р аппроксимируется так (за скобки выносятся общие индексы): .г '(р"+1 — р")1, +, + (и411р),". 1 ., + + (нАР)1»11хг«11з+ (Р б1» )1+1р,1+11з = 11.
гашения днгмьтных злдлч глзовой динлмики й гз1 377 Поясним смысл некоторых величин. Если в указанной точке отсутствует, например, величина и";~цкгчнг, она вычисляется линейной интерполяцией по ближайшим точкам: и";.„цх/,,цг -— (и," „, + и".. н,.~нг)/2. Символом Н, обозначена односторонняя разностная производная по х, ориентированная против потока. Так, Н (~/~Р)гчнхгг цг = '" (Р;+цкг+нг Р— цгт.
цг) если и",,н чц > О. Так же строится Иг — аппроксимация производной по у. Оператор дивергенции аппроксимируется естественно: (О1т)Нц2,/чцг=" (("гн,,чцг ",,+цг) +,("ННХ7+г "Нику)). Аналогично строится аппроксимация уравнения для давления р;: р,+яр +яр +[ур+(у — 1)я1(и„+о )=О (в случае уравнения состояния идеального газа). Здесь величина д = С/гг р(и„+ н ) г — искусственная вязкость.
ФоРмУлы Расчета гРаничных значений Р,".,+717г, Рг++,'гг бУдУт 'описаны отдельно. После вычисления всех р"+', р"+' величины и"~', о" ч' находятся по формулам типа (" " ) Ь та Нг + (Н'(гц) С 7~ Нг + (™гя) Ь гз-Цг + + " '((Рл '+ Я")овца,/ Нг (Р"+'+ Я")7-НХ/+Нг)/Р",/+Нг.
Разностное уравнение для о строится по такому же принципу. Подчеркнем, что р берется уже с верхнею слоя. Перейдем к аппроксимации уравнений на границе. С учетом условия непротекания они имеют вид р, + ггр, + р(п1т) = О, где д1т= и„+ о„, х — длина дуги на контуре. Члены р, + чр, аппроксимируются по тем же принципам, но с учетом знака и.
Пояснения требует только вычисление (п1у )гчнг. Каждой (р, р)-точке ставится в соответствие свой элементарный объем (см. рис. 41). Если точка внутренняя, элементарный объем есть ячейка И х Ь с центром в точке (г + 1/2„ / + 1/2). Если точка граничная, это — четырехугольник (хг а х а хг„( х (у, а у а г (х) ), где уг — уровенр верхней границы элементарной ячейки, соответствующей последней на зтг ПРИБЛИЖЕННЫЕ МЕТОДЫ ВЫЧИСЛИТЕЛЬНОЙ ФИЗИКИ 1ч, и линии х = хг»пг внутренней счетной (р, р)-точке, таким образом, элементарные объемы, соответствующие (р, р)-точкам (как внутренним, так и граничным), покрывают всю область течения без пустот и перекрытий. Теперь для построения аппроксимации д1т используем известную формулу о!тв»»йш (1/о) ~> (в, п) Б(л, в»» (и, н), «В где о' — малая, стягивающаяся к точке (х, у) область, п — нормаль к ее границе (внешняя), Мы используем допредельный аналог этой формулы, беря в качестве О элементарный объем, соответствуюший данной граничной (р, р)-точке.
Интеграл аппроксимируется почти очевидным образом: на сторонах АВ и С.0 (см, рис. 41) по значениям м на границе и ид П»нг интерполируется и, и эта величина интегрируется. На стороне ВС известно и, пг г, и с ~ (в, и) ~/зж — АР»пх . в Остальные детали аппроксимации не уточняем; они аналогичны тем, которые используются в стандартных счетных точках.
Опыт показал, что при вышеописанном способе формирования счетной сетки могут образоваться ячейки с маленькими линейными размерами. Наличие подобных ячеек вынуждает уменьшить пгаг по времени т по сравнению с тем, который обусловлен требованиями вычислительной устойчивости Во внутренних счетных ячейках. В противном случае в Области «уменьшенных» ячеек проявляется вычислительная неустойчивость. Чтобы избежать этого, следует вводить более крупные приграничные ячейки, присоединяя к ним примыкающие внутренние ячейки /Б х /з и исключая оютветсгвующие внутренние счетные точки. Изложенная выше схема была испытана расчетом следующей задачи.
Пусть начальные данные имеют вид (Х„берется несколько левее острия клина; см. рис. 42) ~ О, 1, О, 0 прих>ХВ, ~!О, 4, 25, 0 при х< ХВ. Они удовлетворяют (при у = 5/3) условиям Гюгонио. Ударная волна, движущаяся вправо, падает на острый клин, возникает сложное течение, известное в газовой динамике под названием «гришокм Образуется точка, в которой «сходятся» три ударных волны: исходная ударная волна (ее фронт — вертикаль), «преломленная прямая волна» (ее фронт ортогонален линии клина) и отраженная ударная зтз 1 гз1 гвшвнив двгмвгных злдлч газовой хин»михи волна.
Эти три линии сильного разрыва сходятся в «тройной точке», движущейся по прямой линии под некоторым углом к линии клина. На рис. 42 представлены линии уровня давления, которые соответствуют (справа — налево) значениям р, равным 2, 4, б, 8, 10, 12, 13, 14, 15, 1б, 17. Перед фронтом волны р = 1. Хорошо видна характерная конфигурация тришока, хотя число узлов на рис. 42 не очень велико для столь сложного течения.