Д.В. Маклаков - Нелинейные задачи гидродинамики потенциальных течений с неизвестными границами (1163230), страница 30
Текст из файла (страница 30)
Тем самым мы негласно предположили, что функция Ф(1), определяемая по формуле 11.23, ведет себя при больших 1 как 1/1, то есть 1пп Ф(1)1 = сопв1. Введенная аппроксимация функции Л(и) позволяет в аналитическом виде выразить функцию Ф(1), а затем получить аналитические 192 Глава 4. Препхтствнн вблизи границы раздела сред выражения для значений интегральных операторов 8[Л, 0], Бг [Л, О] в точках з. и подсчитать функционал Б'[Л, О]: Ф(1) = — / 1 ОЛ(и)ди 1 2те',/ и — 1 2н1 Е, п (Аг 1 — А1)(1 — иг) 1п(1 — и.)+ т=г +(Авг + „— — ) 1п(М вЂ” и„~.г) — — (Ад С + Вг — — ) 1п(8 — иг)+ А„ег Ао А„ег Ао А„+г Ао] + — 1п(-и„е1) — — 1п(-иг ) + — — — ~, не+г нг (12.2) В[Л,О](е.) = Ф (и ) т — ~ 2 Я,[Л,О](е,) = ф(1/иг) — ф(0), (12.3) 8'[Л,О] = 2 1тп Ф(0). Алгоритм последовательных приближений выглядит следующим образом.
Зададим значения параметра ~0, а в узлах з — впадения функций Л(е;), 0(е.). С помощью формулы (12.1) определим точки и и значения операторов Бг[0,11], Из[0] в зтих точках. С помощью формул (12.2), (12.2) найдем значения операторов Б[Л, О], Вг [Л, 0] и вычислим функционал Б'[Л,О].
Полученные числа подставим в правые части уравнений (11.29), (11.32), (11.33) и определим новые значения функций Л(е1), 0(е1 ) и параметра 0 и т. д. За нулевое приближение при больших йе можно взять Л(з ) = 0(е ) = О, г = 1, и, О = и — Д. Численное решение для конкретного профиля сначала найдем для больших пе. Затем величину ле постепенно уменьшим с некоторым шагом и за нулевое приближение для следующего т1е выберем решение, найденное для предыдущего значения пе.
Таким способом практически всегда удается добиться сходимости метода вплоть до отстояний линии раздела от профиля порядка сотых долей хорды. Однако для малых отстояний число итераций для достижения хорошей точности может быть довольно большим: иногда более пятидесяти. Немаловажное значение для сходимости и точности вычислений имеет также выбор сетки узлов е . Для достижения хорошей точности сетку выгодно брать неравномерной так, чтобы в результате получилось, большее количество точек и; на участках большой кри- 112.
Численный метод визны контура Х„. Последние образуются, если отстояние профиля от линии раздела мало. В этом случае контур Ь, вблизи параметрического круга практически повторяет форму последнего, а затем его кривизна резко падает. Для получения неравномерной сетки в этом и в следующем параграфах настоящей главы был использован следующий прием. Строилась функция а(а), которая взаимно однозначно переводит интервал ( — х,и) на числовую прямую Я, Равномерная сетка а на отрезке (-х, и) преобразовывалась по формуле и = а(ву) в неравномерную сетку и . В качестве а(а) выбиралась функция вида: / )а — Ь соа а~ савау а где а и 6 — положительные числа, подбираемые в результате численных экспериментов.
Возможны два случая: 6 < а и 6 > а. В первом случае а(а) = (а + 6) ай — — 6 2 (12.4) Производная функции а(а) нигде не обращается в нуль, и расстояние между узлами сетки а; будет монотонно увеличиваться в обе стороны от нуля. При увеличении величины 6/а сгущение сетки вблизи нуля будет усиливаться, Функцию а(гт) в виде (12.4) выгодно применять при больших отстояниях профиля от линии раздела 1 г. Во втором случае Ь > а и (а+Ь) ак — — Ьа — с, — тг < а < — ха, а агт — (а + Ь) аб —, -ха «т < хо, а (а + Ь) 1к — — Ьа + с, ха < а < х, а(а) = (12. 5) где хо — — агссоа(а/Ь), с = 2(Ьхо — ~/Ьа — аа). Производная функции а(а) обращается в нуль в точках а = ххо. Поэтому сетка а будет сгущаться вблизи точек а = х(Ьхо — ~/Ьа — аа).
Такую сетку выгодно использовать на малых отстояниях профиля от линии раздела, добиваясь за счет выбора а,6 максимально информативного описания участков резкого изменения кривизны на контуре Ьо Глава 4. Прелятствяя вблизи границы раздела сред 194 12.2. Вычисление гидродинамических характеристик. После определения функции Х(() силу В,+(Вю действующую на контур профиля, можно вычислить по формуле Блазиуса-Чаплыгина [67): 'Р ((,— г 1' [' '~Ы~)),-гх(0,~( 2 ,/ ~'(() с где Сн — окружность радиуса В > 1 с центром в начале координат плоскости г.
Поскольку для рассматриваемого течения справедлив парадокс Даламбера, то Вг — — 0. Если передняя кромка профиля острая или кривизна ее велика, непосредственное применение формулы (12.6), даже в случае интегрирования по окружности радиуса В > 1, приводит к большим погрешностям при расчетах. Пусть е'Р" — образ точки г„(носика профиля) с максимальной кривизной.
Из формулы (12.6) нетрудно получить, что 8т соя е гх(е )Ке-(ао ( (7 — (уг .з 2 + [е ' 'КХО(()3 [ -гх(0 -гх(е'з")[6( У'(() В .В Р (р-)г (12.7) В формуле (12.7) поцынтегральное выражение не обращается в бесконечность при ( > 1 даже цля профилей, имеющих в носике точку возврата. Для главного момента сил давления имеем: где момент вычислен относительно точки г„с максимальной кривизной. Подынтегральное выражение в (12.8) также не имеет особенностей. Формулы (12.7), (12.8) использовались при расчетах подъемной силы и момента.
При этом помимо В вычислялось также значение В,. Теоретический результат — В, = 0 — служит хорошей проверкой точности вычислений. Можно ожидать, что порядок полученного при расчетах значения В дает порядок точности вычисленного значения В„. Результат считался удовлетворительным, если [В,[ < 0.001[Вт[. 190 112. Численный метод Другой внутренней проверкой точности вычислений является способность алгоритма генерировать прямолинейный экран при 7 = оо. Здесь результат считался удовлетворительным, если ордината линии раздела отличалась от некоторой постоянной величины на 1% характерного размера. В рассматриваемой задаче обтекания профиля вблизи границы раздела двух невесомых жидкостей линия раздела будет либо бесконечно подниматься вверх при удалении в обе стороны от профиля, либо бесконечно опускаться вниз. Исключение составляет особый случай 7 = оо — движение профиля вблизи твердой стенки.
Ввиду отсутствия характерной горизонтальной линии для отсчета глубины погружения в качестве параметра, характеризующего эту глубину, возьмем величину 6, равную расстоянию на бесконечности между разделяющейся линией тока и линией раздела (см. рис. 58). Этот параметр выражается через 1э по формуле: 6 Ке-$ао~ (1)е-хр) 01 где г1 — любая точка на параметрической окружности. 12.3. Результаты расчетов модельных и тестовых задач Обтекание вихря.
Рассмотрим задачу об обтекании вихря двухслойным потоком невесомых жидкостей, являющуюся простейшей из рассмотренных в настоящей главе. Она была исторически первой из решенных автором задач обтекания препятствий вблизи границы раздела сред (см, [79) ) и на ней была апробирована изложенная в разделе 12,1 техника решения систем нелинейных интегральных уравнений вида (11.29), (11.32). В частных случаях обтекания вихря вблизи твердой стенки (7 = со) и свободной поверхности (7 = О) легко построить точные аналитические решения, что делает возможным проведение тестовых расчетов. Схема течения показана на рис.
02. Здесь мы считаем, что свободная поверхность проходит через начало координат, вихрь расположен на мнимой оси в точке — гй(6 > О), а — Ц(0 > 0) — критическая точка вихря. При ф < 6 будем иметь обтекание вихря с положительной циркуляцией (отрицательной подъемной силой), при 8 > 6 — с отрицательной циркуляцией (положительной подъемной силой). В данной задаче нет необходимости использовать для решения параметрическую плоскость, и систему уравнений можно вывести 196 Глава 4, Препятствия вблизи травины раздела сред Рис. 62. Схема обтехавии вихря непосредственно в области течения. Введем функцию у(4), которая связана с комплексным потенциалом течения Ит соотношением: 'т' 'Ь(т)е ~'~, б а 1, ~' уо(т)е х1'1, г Е С где Уо(г) = (х + 1Ф)/(в + й) — комплексно-сопряженная скорость обтекания вихря безграничным потоком.
Систему уравнений выведем относительно функций Л(в) = т (и) — т (и) и в(в), где в — длина дуги контура Ь, отсчитываемая от начала координат, д — угол наклона к оси абсцисс касательной к контуру Е. В результате мы придем к уравнениям (11.29), (11.32), в которых Яз[д) = О. После решения этой системы функцию х(в) выразим через А по формуле (11.23), где 1 = ю Подъемную силу Яю действующую на вихрь, вычислим по формуле Влаэиуса-Чаплыгина.
В табл. 10, 11 приведены значения коэффициента подъемной силы с р (р~)зб' ПОЛУЧЕННЫЕ ПО ПРЕДЛаГаЕМОМУ МЕТОДУ, И ЗНаЧЕНИЯ 4 'вт, НайДЕННЫЕ ПО точным аналитическим формулам из решения задач о вихре вблизи твердой стенки (1 = оо) и свободной поверхности (1 = О).
Из 197 512. Численный метод Таблица 10 Звачеиия Си и Ст для вихря при т = со (твердая стелла) Д/Ь 0.1 0.3 0.5 0.7 0.9 1.1 1.3 1.5 1.7 1.9 С„ -4.692 -4.420 -3.832 -2.799 -1.138 1.386 5.078 10.30 17.48 27.11 С -4.681 -4.417 -3.829 -2.796 -1.137 1.389 5.083 10.31 17.49 27.10 Таблица Ы Зиачеиия Сх и СГ для вихря при ч = О (свободиая поверхяоств) 13/Ь 0.1 0.3 0.5 0.7 0.9 1.1 ' 1.3 1.5 1.7 1.9 Ст -42.17 -20.63 -10.95 -5.161 -1.388 1.140 2.809 3.840 С -42.37 -20.67 -10.96 -5.168 -1.391 1.138 2.810 3.846 4.370 4.464 4.377 4.470 Пластина под свободной поверхностью. Обтекание пластины под свободной поверхностью было рассчитано в качестве тестового примера.
Аналитическое решение этой задачи получено рядом авторов (см. [124], [56], [57], [20]). Наиболее полные числовые расчеты выполнены в статье [58]. В табл. 11 при различных Ь/( и угле атаки од = 20' даны значения С„, полученные настояшим методом, и значения С„из [58). Здесь ф— коэффициент подъемной силы пластины, отнесенный к ее длине (, С„" = 2лейп ел — коэффициент подъемной силы в безграничной жидкости, Ь вЂ” расстояние на бесконечности между разделяюшейся линией тока и линией раздела. Из сравнения сопоставления результатов видно, что погрешность вычислений не превосходит 0.5%. Данный расчет был проведен при 31 точке коллокаций на линии раздела.