Численное исследование течений около крыла при трансзвуковых скоростях набегающего потока на основе решения уравнений Эйлера тема автореферата и диссертации по механике, 01.02.05 ВАК РФ

ГОСУДАРСТВЕННЫЙ КОМИТЕТ РФ ПО ВЫСШЕМУ ОБРАЗОВАНИЮ МОСКОВСКИЙ ГОСУДАРСТВЕННЫЙ АВИАЦИОННЫЙ ИНСТИТУТ

На правах рукописи УДК 533.Об.011

КРАЙНОВ Виталий Алексеевич

ЧИСЛЕННОЕ ИССЛЕДОВАНИЕ ТЕЧЕНИЙ ОКОЛО КРЫЛА ПРИ ТРАНСЗВУКОВЫХ СКОРОСТЯХ НАБЕГАЮЩЕГО ПОТОКА НА ОСНОВЕ РЕШЕНИЯ УРАВНЕНИЙ ЭЙЛЕРА

Специальность 01.02.05 Кеханика жидкостей,газа и плазмы

Автореферат

диссертации на соискание ученой степени кандидата физико-математических наук

Москва Издательство МАИ 1994

РГ6 од

3.....

Работа выполнена в Московском государственном авиацжо институте.

Научные руководители - доктор технических наук,

академик Ю. А. Рыжов доктор технических наук, профессор С. Б. Свирщевский

Официальные оппоненты :

доктор физико-математических наук, профессор Н. С. Иванов кандидат физико-матемагическях наук Н.Э. Иванов

Ведущая организация - ЦАГИ имени Н.Е.Жуковского

Защита состоится 1994 г. в "/У^^час!

на заседании специализированного совета К 053.18.02 в Московском государственном авиационном институте.

С диссертацией можно ознакомиться в библиотеке маи. Адрес института : 125871, Москва, ГСП, Волоколамское шоссе, 4

Автореферат разослан

аЖу&еи*/' 1994 Г.

Ученый секретарь Совета кандидат физико-математических наук /у Л.Ф. Лобанова

ОБЩАЯ ХАРАКТЕРИСТИКА РАБОТЫ

Диссертация посвящена разработке метода численного решения задач* о стационарной обтекании профиля и крыла конечного размаха идеальной сжинаемой средой в неограниченном потоке.

Актуальность текы. Существенное увеличение лобового сопротивления при трансзвуковых скоростях палата , когда у поверхности летательного аппарата возникают местные сверхзвуковые зоны, замыкающиеся скачками уплотнения , приводит к повышении) экономических затрат на эксплуатацию авиационной техники . Уменьшение лобового сопротивления может быть достигнуто путей изменения обводов поверхности. Быстродействие современных ЭВМ позволяет еще на этапе про-проекткрования использовать численные методы, которые даст возможность рассчитывать обтекание летательных аппаратов потоком сжимаемой среды я выбирать по результатам расчетов более совершенные в аэродинамическом отношении формы. Однако применение таких методов связано с большими затратами ресурсов ЭВМ , а в некоторых случаях наблюдается снижение точности расчетов, поэтому необходимо разрабатывать более экономичные к точные методы.

Целью работы является:

- создание эффективного метода расчета обтекания крыла летательного аппарата потоком сжимаемой среды на основе численного интегрирования системы уравнений Эйлера.

Научная новизна:

- разработан и реализован метод численного интегрирования системы уравнений Эйлера на основе модифицированной неявной Л-схекы;

- предложен метод локализация пространственных скачков уплотнения, позволяющий определять положение фронта скачка с высокой точ-

ностью;

- предложен функционал для оптимизационного метода построения криволинейных сеток, что дало возможность сократить количество вычислений;

- с помощью предложенного метода расчета получены новые результаты, содержащие поля течений около профиля и крыла конечного размаха при различных значениях углов атаки к чисел Маха в трансзвуковом диапазоне.

Практическая ценность.

- на базе предлагаемого метода расчета обтекания тел потокок сжинаемой среды создан программный комплекс, позволяющий генерировать конечно-разностные сетки, рассчитывать аэродинамические характеристики профиля и крыла конечного размаха с высокой точностью на достаточно грубых сетках и проводить обработку результатов расчета с выводом графической информации на экран дисплея, печатающее устройство или графопостроитель. Эти программы могут применяться для проведения исследовательских и проектировочных работ;

- эффективность и экономичность алгоритма позволяют использовать его на персональных компьютерах;

Апробация работы. Результаты работы докладывались и получили положительную оценку:

1. На конференции АН СССР "Юрьевские чтения" (Москва, МАИ,

1989 г.)

2. На Всесоюзном семинаре по проблемам современной аэродинамики под руководством академика Ю. а. Рыжова ( Москва , маи ,

1990 г. )

3. На Всесоюзном семинаре по проблемам современной аэроди-

. намики под руководством академика Ю. А. Рыжова ( Москва , МАИ ,

1991 г.)

4. На конференции ЛИ Р« "Юрьевские чтения" ( Москва , МАИ,

1992 г.)

Публикации. Основные результаты диссертации опубликованы в работах [1-5], список которых лриводон о конце автореферата.

Объем работы. Диссертация состоит из введении, трех глав, заключения, списка литературы и приложения . Изложена на 185 страницах, включая {08 фигур, одну таблицу и список литературы из 98 наименований.

СОДЕРЖАНИЕ РАБОТЫ

Во введении формулируется цель диссертационной работы, отмечается новизна и практическая ценность предлагаемого метода , указывается порядок изложения материала в работе.

В первой главе рассматриваются различные подходы к исследованию течений сжимаемой среды с помощью методов численного моделирования. Отмечается , что точность численных методов существенно зависит от выбора математической модели, описывающей законы движения среды. В § 1. 1 , посвященном обзору литературы по данной тематике, говорится о существующих подходах к решению подобных задач , описываются достоинства и недостатки различных методов, а также рассматриваются вопросы , связанные с построением конечно-раз ноет ной сетки. В § 1, 2 система уравнений Эйлера, имеющая следующий вид-.

5 (а + + оУ'У = О

1 (1).

V + + 5а7а = О

где V - вектор- скорости потока; а - скорость звука; 3=2/( г-1) ,

Гя отношение удельных теплоемкостей, записывается в А-форм

С +(и+а)С +уС +иС =-а(у +и )

I 1 ' X у I у у х'

й + {и-аШ = а(у +« )

I 1 ' х у I у у г'

Е +иЕ + (у-ьа)Е <-уЕ =-а(и +у )

I , 'у х х (2)

Г +аг +(у-а)Е = а(и )

I х у I * х х

С +(и+а)С =-а(и +у )

I х у ' I 1 X у'

II +иН +уН +(и-а)Н = а(и +у )

I х У х * х у

Здесь и, у, V - проекции скорости потока на оси декартовой систе координат; С,0,Е,Г,С,Н - бихарактеристкческие переменные :

С=и+5а; Е=и-5а; Е=у+йа; Г=у-5а; С=и+5а; Н=у-3а (3

Отмечается, что при численном интегрировании системы ( 2 ) кож* возникать неоднозначность решения задачи. Для получения единстве! ного решения система уравнений приводится к следующему виду:

2и +(и+а)С + (и-а)О +у(С +0 )=0

I 1 ' X 1 ' 1 'у у' 1 I I'

2у +и(Е +Г ]+(у+а)Е +(у-а!Р +«(Е )=0

I 4 X X ' * У у 4 X 2 '

(4)

2аа^+(ц+а)С1<-(и-а)0х+(у+а)Е -(у-а)Г +(м+а)01~(*-а)Н1=0

2V +и(С +Н )+у(в +Н )+(у+а)С +(м-а)Н =0

Л 1 X X у у ' I 1 ' 2

В § 1. 3 содержится краткое описание неявного метода численного и> тегрирования системы (4). Рассматривается способ задания граничм условий , параметры потка на верхнем временном слое записывают« через приращения:

ип*1= и% Ди ; уп41= у"+ Ду ;

= Ду ; а"* '= а"+ Да

(5

излагается алгоритм блочного кетсдл. релаксации по точкам (БП-ке-

тод) для определения неизвестных величин приращений.

Во второй главе рассматривается численный алгоритм расчета обтекания профиля крыла сжимаемой средой. Задача решается в физической расчетной области, внешняя граница которой находится на бесконечно большом расстоянии от рассчитываемого объекта. С помощью алгебраического преобразования такая расчетная область отображается на вспомогательную расчетную область, представляющую собой круг единичного радиуса, а которой строится конечно-разностная сетка, обладающая оптимальными свойствами (рис.1). Построение сетки осуществляется с помощью метода сопряженных градиентов, для чего предлагается использовать функционал следующего вида:

♦-Ц с**,.,♦<1(6)-

где - мера ортогональности :

<?г К >г+<?2- гз >2+<гз- г< >г+<*4- >г <7>; - мера гладкости :

к - скалярный весовой коэффициент , задающий соотношение между мерами гладкости и ортогональности; а г - радиус-векторы, соединяющие соседние узлы конечно-разностной сетки и образующие комбинированную ячейку (рис.2) . Отмечается, что применение вариационного метода для построения сеток позволяет получать результаты расчетов обтекания профиля за меньшее число итераций и с большей точностью. В § 2. 2 система уравнений (4), записанная для плоских течений, заменяется конечно-разностнык аналогон, для чего используются кокеч-

рис. 2

но-разностные операторы первого порядка точности по времени и односторонние трехточечные разностные операторы второго порядка точности по пространству :

Г'= { Г"*1- г" )/ Л1 1 ' " (9),

где { - дифференцируемая функция ; п - номер шага по времени ; ш - номер сеточного узла; М - шаг по времени; ДЪ - шаг по криволинейным координатам; б - коэффициент, позволяющий учитывать направление распространения вознущонкй вдоль сеточных линий ( 5=1 о том случае, когда скорость распространения возмущений вдоль бихарактеристик больше нуля, и г=0 - в противном случае). Величина скорости распространения возмущений определяется по формулам :

для С, и? + у£

5 X > V

для Э, цР +уР -а?

5 X ^у X

для Е.

пля Г- и? +у£ -а£

С * у у

для С_ ит) ч-ут) +ат] ; Т) ж у *

для ит) +УТ) -ат) ;

71 хук

для Е ш1 +аг? ;

Ч х у у

для Г ит) +уг) -ат)

Т) х у у

(Ю).

Конечно-разностный аналог системы (4) преобразуется в систему линейных алгебраических уравнений относительно неизвестных величин приращений скоростей, имеющую блочную тринадцатидиагональную матрицу. Отнечается, что приближенное решение такой системы может быть получено с помощью БРТ-метода, для чего на каждом шаге по вренени для каждого сеточного узла необходимо решить систему линейных алгебраических уравнений с матрицей 3x3 . Предлагаются формулы для вычисления элементов этой матрицы.

В § 2.3 рассматривается метод задания граничных условий . На внешней границе области задаются параметры невозмущенного потока, а на

потоком; обтекание профиля NACA0012 дозвуковым потоком ( М-0.72, д-0° х М-0. 63, а-2°), трансзвуковым потоком ( М=0.8, а=0°; М=0.8, а-1.25°: М=0.85, а=0° ; М=0.85, а-Io; М=0.95, а=0°) И сверхзвуковым ПОТОКОМ (М=1.2, сс=0°; М=1.2, а-7°); обтекание профиля Корн (М=0.7, а-0° и М=0.7, а-2°) . На рисунке 3 приведены характерные графики, описывающие изменение рассчитанного коэффициента давления по поверхности профиля . Точками на графики нанесены результаты работ [ 6-9 ]. С целью определения точности и применимости описанного метода были проведены сравнения результатов расчетов с экспериментальными данными. Сформулированы выводы по второй главе.

В третьей главе рассматривается метод расчета обтекания крыла конечного размаха потоком сжимаемой среды. Определяются границы физической расчетной области. Предлагается преобразование этой области на вспомогательную расчетную область (рис. 4) . Во вспомогательной расчетной области строится конечно-разностная сетка с помощью вариационного метода , для чего предлагается использовать функционал вида:

4 «ЕП[* sl)l(+ 11-") v|Jk] (13),

i i *

где s1Jlt~ мера ортогональности •• V(jk- мера гладкости :

v1Jk= ( К +?,>2+< К +?б>2 <15>;

,к - скалярный весовой коэффициент , задающий соотношение между мерами гладкости и ортогональности; а г- радиус-векторы, соединяющие соседние узлы конечно-разностной сетки и образующие комбиняро-

-î.o -o.s л \ o S u-í s \ 0 "8 : 1 4 - ! H я II ce: Cr .0 fv 4 ).S 1.0 1.5 $ i |i|| I|II nj.....» inj mi i ■9 . -1.0 -0.5 ß lu чMl"J11 UU^j^uujffg Z Fi Я NX n : о о : " О pj . : -t- 3 X : » Р О" ТО Я 1 ♦- Or 0 f 0.5 1.0 1.5 " LLj-LlU J.JU1 Jllj-LUj-Li 4 \ t .....fi ... 5 ' ~ О Ï4.' о îO . - ■t- rt X P P: CD: 6 Si —Г- o: .0 / 0.5^ 1.0 l.S лР

* iv>: nj : 5 CÍ 3 -Г-: p x i Ь' o¡ * /1 (U J ; U er! 0 с aKi Í 1 .5 1.0 l.S J ) . -1.0 -0.5 rÖ S N 3 «vi о О -NI С> î H CD Çf Я. V: Л Ф £ $ "l^A Ct 0 / ! 0.5 1.0 1.5 ÍÍ ujHUJ-iHi. fUJ-IJJIU jjUXlUXU > h \ 1 —-—■ * -1.0 -0.5 ¿3 1 pi s О ; M Oj 3 Ц tf >< о сл Oj R. Ci: .0 3.5 1.0 1.5 -P ■XIm "M 1 "I mi4¿»|miíixa.II / /

для

для

для

для

и

ванную ячейку (рис.5) . Эти вопросы рассматриваются в § 3. 1. В § 3. 2 система уравнений (4) записывается в конечных разностях, для чего используются конечно-разностные операторы первого порядка точности по времени и односторонние трехточечные разностные операторы второго порядка точности по пространству (9), Формулы для определения величин скорости распространения возмущений записываются в следующем виде :

и£ +у£ +у£ +а£ ; для С ит) +ут) +иг] +ат) ;

к у I х 7} х у I х

и? +у£ -а£ ; для ич +ут) -ищ -ат) ;

^х ^у ^х чх * Ч X у I X

и£ +у£ +у£ +а£ ; для Е ит) +ут) +ат) ;

^х у х >у ' Т) X у I у

и£ +у£ -а? ; для Г_ ит) +ут) +уп -ат) ;

^Х ^у 2 У 7) X У X У

для С^.О^.Е^.Г^ ; (16).

для Н^ и^+у^+у^+а^ ; для Н^ ит^+ут^+ут^+ат^ ;

для н^ :

для С^ ; для С^ ит^+ут^+ут^-ат^ ;

для С^ }

Записанная в конечно-разкостнок виде система (4) преобразуется в систему линейных алгебраических уравнений относительно неизвестных величин приращений скоростей , имеющую блочную семнадцати-диагональную матрицу. Приближенное решение такой системы получается с помощью БРТ-метода , для чего на каждом шаге по времени для каждого сеточного узла необходимо решить систему линейных алгебраических уравнений с матрицей 4x4 . Предлагаются формулы для вычисления величин элементов этой матрицы.

В § 3.3 рассматривается способ задания граничных и начальных условий . На внешней границе области задаются параметры невозмущенного потока, а на поверхности крыла - условие непротекания:

рис. 5

к

u cosa ♦ v cos/3 + v cos г = О (17),

где cosa , cos0 , cosy - направляющие косинусы нормали к поверх-ноети крыла. Уравнение (17) преобразуется к виду:

Да cosa + Av cos/3 + Av cost = о (18).

В расположенных на поверхности тела сеточных узлах предлагается записать вместо первых трех уравнений системы (4) следующие соотношения: уравнение (18). проекцию на направление касательного к поверхности вектора, лежащего в плоскости сечения. и проекцию на на направление вектора, перпендикулярного к вектору нормали и касательному вектору . Для определения неизвестных величин предлагается алгоритм, позволяющий улучшить характеристики устойчивости и сходимости метода расчета. Для получения однозначного решения задается дополнительное граничное условие - условие Кутта-Жуковского, определяющее на поверхности крыла положение линии схода потока, на разрезе расчетной области . соединяющем заднюю и боковую кромки с внешней границей, задается условие непрерывности параметров потока, а в точке соединения задней х боковой кромок задается точка торможения. В качестве начальных условий используются значения параметров невозмущенного потока.

В § 3. 4 рассматривается нетод локализации скачков уплотнения . Для определения стационарного положения скачка уплотнения вычисляется по предложенному алгоритму величина скорости движения скачка. В стационарном случае скорость движения скачка равна нулю, в противном случае выполняется пересчет поношения скачка , после чего вычисляются значения параметров потока за скачком.

В § 3. 5 приводятся результаты расчетов дозвуковых и трансзвуковых течений. С целью исследования сходимости и устойчивости предлагаемого метода были проведены расчеты обтекания прямоугольного и

иг

стреловидного крыльев , имеющих в поперечном сечении профиль НХСЛ0012, дозвуковым потоком (Н-0. 72, а-0° и К-0.63, а-2°) к трансзвуковым потоком ( Н=0.8,а=0°; М=0. 8, а-1. 25°; М=0.85,а=0°; М=0.85, а-1°) , а также обтекание трансзвуковым потоком крыла ОНЕИА Мб ( М=0.84. а=3.0С°). На рисунке 6 изображено с помощью изолиний поле местных чисел Маха на поверхности крыла ОЫЕЛА Мб , построенное по результатам расчета и имеющее характерную структуру течения [10] . Сформулированы выводы по третьей главе.

В заключении сформулированы основные результаты диссертации.

ОСНОВНЫЕ РЕЗУЛЬТАТЫ

Основные результаты работы сводятся к следующему:

В данной работе разработан и реализован метод численного интегрирования уравнений Эйлера для задач об обтекании профиля и крыла конечного размаха сжимаемым потоком.

Предложен метод локализации скачка уплотнения.

Усовершенствован вариационный метод построения конечно-разностных сеток.

Выполнено численное исследование метода и осуществлены расчеты обтекания профилей и крыльев конечного размаха.

Численные исследования показали достоинства предлагаемого метода : возможность получения результатов при расчете на достаточно грубых конечно-разностных сетках ; простота в задании граничных условий ; обеспечение физически обоснованного учета типа уравнений в дозвуковых и сверхзвуковых областях; проведение расчетов течений со скачками уплотнения методом сквозного счета без введения искусственной диссипации.

Верхняя поверхность

рис. 6

Разработанный алгоритм может быть использован для'Проведения исследовательских и проектировочных работ в авиационной промышленности.

ЛИТЕРАТУРА.

1. Крайнев В. А. Расчет нестационарного обтекания профиля сжинаемым газон//сборник докладов конференции "Юрьевские чтения",-М, МАИ, 1989 г.

2. Крайнов В. А. Расчет потоков сжимаемого газа на оптимальных разностных сетках// трувы МАИ , -Н, МАИ, 1991.

3. Крайнов В.А., Расчет положения скачка уплотнения при трансзвуковом обтекании профиля крыла// труды КАК,-М,КАИ, 19S2.

4. Крайнов В. А. , Метод локализации скачка уплотнения в задаче об обтекании крыла конечного размаха потоком сжимаекой среды// вестник НАЛ,-М,МАИ, 1993.

5. Крайнов В.А., Конечно-разностный метод расчета трансзвукового обтекания крыла конечного размаха с локализацией скачка уплотнения// сборник докладов конференции "Юрьевские чтения", -М, МАИ, 1992 Г.

6. Lock R.C., Test Cases for Humerikal Methods In Tvo-Dimensional Transonic Flows// ACARD-R575-70, 1970.

7. Jameson A., YoonS., Multigrid solution of the Euler Equations Using Implicit Schemes// AIA.A Journal, Vol.24, Nov, 1986.

8. G. Volpe, A. Jameson, Transonic Potential Flow Calculations by Two Artificial Density Methods // AIAA 86 , May 12-14 ,1986 , Atlanta.

9. Г.С. йк , А. Хартен , Неявные схемы TVD для гиперболических систем уравнений , записанных в консервативной форне относительно системы криволинейных координат // Аэрокосническая тех-

и

ника, N 11, -и , Мир, 1987 г. 10. S. Aq Г aval , r.e. Verneland , A. Verhoff , r. в. Lovrie , Euler Transonic Solutions over Finite Wings // AIAA Paper 88-0009, January, 1988.