Введение в вычислительную физику - Федоренко Р.П.
ISBN 5-7417-0002-0
Скачать (прямая ссылка):
Аппроксимация около криволинейной границы. Рассмотрим течение газа, описываемое уравнениями в форме Эйлера в декартовой системе координат. Граница, на которой поставлено условие не-протекания (контур обтекаемого тела), является относительно гладкой и проходит более или менее произвольно относительно системы координат. Введем прямоугольную равномерную сетку с шагом h по обоим направлениям (для простоты). Вертикальные и горизонтальные линии сетки занумеруем целыми индексами і и /, а различные счетные величины будем вводить в разных точках сетки.
Определим в центрах ячеек сетки термодинамические неизвестные p?+i/2,y+i/2> ^+1/2,у+1/2- Определим компоненты скорости на серединах сторон ячейки Uni j+ll2, f"+1/2y Такая «шахматная» сетка удобна для аппроксимации оператора дивергенции в точках (г + 1/2, j + 1/2) и grad р, причем х-компоненту удобно аппроксимировать в и-точке (где она только и нужна), а _у-компоненту — соответственно В D-ТОЧКЄ.
Среди введенных формально узлов квадратной сетки выделим счетные узлы, т.е. те, в которых соответствующая величина считается определенной. Узел сетки считаем счетным, если квадрат А х Л с центром в этом узле целиком помещается в области течения (не пересекается с обтекаемым телом). Если такой квадрат частично находится в области потока, частично — внутри тела, он считается фиктивным: ему не соответствует никакая величина, входящая в основные массивы счетных величин. Однако такая величина в фиктивном узле может появиться как вспомогательная, вычисляемая через основные.
370 ПРИБЛИЖЕННЫЕ МЕТОДЫ ВЫЧИСЛИТЕЛЬНОЙ ФИЗИКИ /Ч. II
Рисунок 41 поясняет сказанное. На нем в части области показаны сетка, счетные (р, р)-узлы (темные кружки), счетные м-узлы (знаки «минус») и и-узлы (крестики). Указана граница и расположенные на
ней нестандартные счетные узлы. Будем считать, что рассматриваемый участок границы задается кривой У(х), причем IXx(^)I ^ 1. (При I Ух(х) I > 1 контур следует задавать функцией X (у) ив последующем поменять роли хну.)
Итак, на контуре вводятся (р, р)-уз-лы с координатами хі+ІІ2 — (і + 1/2)h, Y(xi +1/2), в которых определены счетные
величины р"+1/2> P?+1/2 (отсутствие второго индекса — признак принадлежности границе). Кроме того, на контуре вводятся w-узлы с координатами xt, Y(xt), в которых определена
счетная величина tv", имеющая смысл касательной к контуру компоненты скорости.
Опишем структуру одного шага интегрирования по времени в явной схеме. Шаг состоит из трех основных операций.
1. По значениям величин в счетных узлах (внутренних и граничных) производится интерполяция соответствующих величин в фиктивные узлы. Тем самым создается ситуация, позволяющая во внутренних счетных узлах рассчитать величины на следующем шаге по стандартным формулам.
2. Во внутренних счетных узлах вычисляются значения на верхнем, (л + 1)-м, слое по стандартным формулам.
3. Расчет граничных значений р"+//2> РІ+и2> wTl производится по специальным формулам.
Конкретизируем вид расчетных формул. Вычисление какой-либо величины в фиктивном узле производится с помощью линейной интерполяции по значениям соответствующей величины на той же вертикальной линии сетки. Используются последняя внутренняя точка и граничная. Если в граничной точке нет нужной величины, например м"+1/2, она вычисляется как 0.5(w” + wni+l) cos aj+1/2, где Cti+ід — угол наклона границы.
Уравнение для плотности р аппроксимируется так (за скобки выносятся общие индексы):
T Чр" + 1 — рП)Ї + 1/2, / + 1/2 4" (й^іР)н-1/2,у + 1/2 +
+ (^^2р)?+1/2,/ + 1/2 + (Р )7+1/2, у+ 1/2 = 0*
Рис. 41
§ 23] РЕШЕНИЕ ДВУМЕРНЫХ ЗАДАЧ ГАЗОВОЙ ДИНАМИКИ 371
Поясним смысл некоторых величин. Если в указанной точке отсутствует, например, величина м"+1/2 у+і/2> она вычисляется линейной интерполяцией по ближайшим точкам:
U/ + l/2,y + l/2 = (uIj +1/2 + Ul+lJ + l/2
Символом d{ обозначена односторонняя разностная производная по х, ориентированная против потока. Так,
(^lP)?+l/2,y + l/2 = h l(pl+ll2,j + ll2~ P?-1/2, у + І/г)’
если м"+і/2 у+1/2 > 0- Так же строится d2 — аппроксимация производной по у.
Оператор дивергенции аппроксимируется естественно:
(div)"+ll2,j + H2 = h 1i(Ul+l,j + U2~ uIj+ 112) +¦ (^+1/2, у + 1 ~ V?+U2, j)]-Аналогично строится аппроксимация уравнения для давления р;: Pt + uPx + vPy + Ьр + (У - 1)<?]Ох + vy) = О
(в случае уравнения состояния идеального газа). Здесь величина q = Ch2 р(их + Vy)2 — искусственная вязкость.
Формулы расчета граничных значений p?+i/2> р?+//2 будут описаны отдельно. После вычисления всех pn+i, p"+1 величины un+l, vn+l находятся по формулам типа
т~l(un+l - un)iJ+ll2 + (UdlU)I J+l/2 + (vd2u)nt J+in +
+ h '[(p"*1 + ?")( + 1/2,7 + 1/2- (Р" + 1 + ?")/-1/2, у+ 1/2І/P?, у +1/2-
Разностное уравнение для v строится по такому же принципу. Подчеркнем, что р берется уже с верхнего слоя.
Перейдем к аппроксимации уравнений на границе. С учетом условия непротекания они имеют вид
р, +wps + p(div) =0,
где div = их + vy, s — длина дуги на контуре. Члены р( + Wpj аппроксимируются по тем же принципам, но с учетом знака w. Пояснения требует только вычисление (div )/+1/2- Каждой (р, р)-точке ставится в соответствие свой элементарный объем (см. рис. 41). Если точка внутренняя, элементарный объем есть ячейка h х h с центром в точке (і + 1/2, j + 1/2). Если точка граничная, это — четырехугольник [х, ^ х ^ х/+1] х [yt ^ у < У(х)], ще Уі — уровень верхней границы элементарной ячейки, соответствующей последней на