автореферат диссертации по информатике, вычислительной технике и управлению, 05.13.18, диссертация на тему:Метод псевдохарактеристик для расчета сверхзвуковых пространственных течений газа

На правах рукописи

НУРБАЕВ УЛУКБЕК ДЖАНЫБЕКОВИЧ

МЕТОД ПСЕВДОХАРАКТЕРИСТИК ДЛЯ РАСЧЕТА СВЕРХЗВУКОВЫХ ПРОСТРАНСТВЕННЫХ ТЕЧЕНИЙ ГАЗА

05.13.18 - Математическое моделирование, численные методы и комплексы программ

Автореферат диссертации на соискание ученой степени кандидата технических наук

г "т ?015

005570758

Москва-2015

005570758

Работа выполнена в ФГОУ ВПО «Московский физико-технический институт (государственный университет)»

Научный руководитель: Михайлов Игорь Ефимович

доктор физико-математических наук ФИЦ «Информатика и управление» РАН профессор, ведущий научный сотрудник

Официальные оппоненты: Сидняев Николай Иванович

доктор технических наук ФГБОУ ВПО «Московский государственный технический университет имени Н.Э. Баумана» профессор, заведующий кафедрой высшей математики ФН-1

Грудннцкий Виктор Георгиевич

доктор физико-математических наук ФГБОУ ВПО «Государственный университет управления»

профессор кафедры математики

Ведущая организация: ФГУП «Центральный научно-

исследовательский институт машиностроения»

Защита состоится « 17 » сентября 2015 г. в 14 ч. 00 мин. на заседании диссертационного совета Д 212.110.08 при ФГБОУ ВПО «МАТИ - Российский государственный технологический университет имени К.Э. Циолковского» по адресу: 121552, г. Москва, ул. Оршанская, д. 3, ауд. 612А.

С диссертацией можно ознакомиться в библиотеке и на сайте ФГБОУ ВПО «МАТИ - Российский государственный технологический университет имени КЗ. Циолковского» http://www.mati.ru

Автореферат разослан члО Л. л. го (Гг.

(дата)

Ученый секретарь диссертационного совета Д 212.110Л1^ /

кандидат физико-математических наук л/1/ Мокряков A.B.

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

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

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

Таким образом, имеются различные по своему происхождению запросы на разработку эффективных численных методов при изучении движения летательных аппаратов. При этом движение летательных аппаратов со сверхзвуковой скоростью является одним из самых востребованных. Поэтому разработка численных методов для расчета движения летательных аппаратов со сверхзвуковыми скоростями имеет важное практическое значение.

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

Неявные сеточные конечно-разностные методы (см., например: Бабенко К.И., Воскресенский Г.П., Любимов А.Н., Русанов В.В., 1964 г.) имеют простую реализацию, но требуют для получения нужной точности много машинного времени. В методах такого типа при решении задачи о сверхзвуковом трехмерном течении газа сложно удовлетворять краевым условиям на поверхности тела и на ударной волне. В 1960-е годы Годунов С.К. предложил явный сеточный метод решения нестационарных задач газовой динамики, записывая уравнения движения в интегральной форме и используя автомодельные задачи о распаде разрыва. В 1972 г. Иванов М.Я. и Крайко А.Н. распространили метод на решение пространственных стационарных задач газовой динамики. (Годунов С.К., Забродин A.B., Иванов М.Я., Крайко А.Н., Прокопов Г.П., 1976 г.).

В обратных характеристических методах (Кацкова О.Н., Чушкин П.И., 1968 г.) необходима двумерная интерполяция искомых величин на предыдущем расчетном слое. Прямые характеристические методы (Борисов В.М., Михайлов

И.Е., 1978 г.) обладают высокой точностью и быстродействием, но имеют плавающую сетку, что требует дополнительных вычислений для определения газодинамических функций в фиксированных точках.

Сеточно-характеристические методы (Магомедов K.M., Холодов A.C., 2988 г.) опираются на достоинства сеточных и характеристических методов, однако для его использования необходимы одномерные интерполяции газодинамических величин на предыдущем расчетном слое.

Актуальность диссертационной работы объясняется тем, что необходим численный метод расчета сверхзвуковых течений, объединяющий в себе достоинства существующих численных методов, таких как высокая точность, быстродействие и фиксированная сетка.

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

1. Разработка численного метода, обладающего фиксированной сеткой, быстродействием и точностью

2. Исследование устойчивости разработанного метода

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

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

Научная новизна. В диссертационной работе предлагается маршевый численный метод для расчета нелинейных уравнений в частных производных гиперболического типа с тремя независимыми переменными. В рассчитываемой области течения газа рассмотрим два семейства заданных поверхностей противоположного типа. На заданных поверхностях запишем исходные дифференциальные уравнения газовой динамики определенным образом - так, чтобы дифференциальные уравнения содержали две касательные производные от газодинамических функций к заданным поверхностям и одну выводящую производную, в направлении, общем к обеим поверхностям. Далее выводящая производная берется в качестве дополнительной искомой функции, наряду с газодинамическими функциями (следует отметить, что идея использовать производную искомой функции в качестве дополнительной искомой функции для конечно-разностных методов рассматривалась в работе Грудницкого В.Г. и Прохорчука Ю.А. в 1977 г.). Таким образом, мы имеем

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

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

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

Апробация работы. Научные результаты были доложены, обсуждены и получили одобрение специалистов на следующих научных конференциях и семинарах:

1. Научно-практический семинар "Теория, численные методы и математический эксперимент в газовой динамике", ЦИАМ, Москва, 2009 г.

2. 53 научная конференция МФТИ, Москва-Долгопрудный, 2010 г.

3. Седьмой международный аэрокосмический конгресс, Москва, 2012 г.

4. Научно-практический семинар в Институте автоматизации проектирования, Москва, 2013 г.

В диссертационной работе использованы результаты, полученные автором в ходе исследований, проводимых в рамках научно-исследовательской работы по проекту РФФИ № 11-01-0083 5-а.

Публикации и личный вклад. По теме диссертации опубликованы 5 работ [1-5], в том числе 2 работы в издании, входящем в перечень ведущих журналов и изданий, рекомендованных ВАК для публикации основных результатов диссертаций на соискание ученой степени кандидата наук. Получено свидетельство о государственной регистрации программы для ЭВМ [6]. Личный вклад автора заключается в разработке и реализации численного метода, разработке комплекса программ, адаптированного для расчета течения газа при пространственном обтекании тел, анализе результатов расчетов, построении математической модели для исследования устойчивости метода, а также подготовке публикаций.

Структура и объем диссертации. Диссертационная работа состоит из введения, обзора литературы, пяти глав, заключения и списка использованной литературы. Диссертация изложена на 93 страницах, включает 10 таблиц, 37 рисунков и 77 наименований использованной литературы.

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

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

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

В первой главе описывается идея метода псевдохарактеристик на примере уравнения переноса

и, + аи,= 0 (1)

где а = сош[, их, и, - частные производные по х и г, и(х,г)- искомая функция. Добавив в левую и правую части уравнения выражение Ьиг, выбрав Ь так,

чтобы выполнялось выражение (а+Ь) = —, где — - угол наклона некоторой

сЬс с!х

заданной кривой г = г(х), и записав (1) вдоль этой кривой, получим — = (— -а)иг. Выбирая Ь другим образом, например а+Ь = -—, получим запись

с1х дх а!г

сЫ , Лг .

уравнения переноса — = (---а)иг на другой заданной кривои

<3х сЫ

(противоположного типа).

Далее введем сетку в рассматриваемой области течения газа и зададим координаты узлов 1, 2 и 3. По известным данным в узлах 1 и 2 необходимо определить искомые функции в узле 3. Расчетный шаблон представляет собой треугольник 1-3-2 (см. рис. 1).

<>г

а + Ь =--

с!х

Рис. 1. Расчетный шаблон метода псевдохарактеристик Рассматривая линии 1-3 и 2-3 как семейство заданных кривых противоположного типа, получим систему уравнений для определения искомой функции в точке 3

сЬс с!х

), ^ ск '

где индексы "1" и "2" означают, что полная производная функции и(х,г) записана вдоль линий 1-3 и 2-3 соответственно. Тогда левая часть уравнений не содержит выводящих производных с линий 1-3 и 2-3, при этом правая часть -содержит (в данном случае иг). Поэтому введем дополнительную искомую переменную Я = иг. Тогда получим следующую систему, где искомыми функциями являются и(х,г) и Я

т)

Ох ^ ах й.х )1 ах

Полученная система, как и условия совместности — = 0 в методе характеристик

с!х

для уравнения переноса, уже не содержит выводящих производных с линий 1-3 и 2-3. Но поскольку линии 1-3 и 2-3 могут, вообще говоря, не совпадать с Л-

характеристиками — = а, то мы их назвали псевдохарактеристиками, а с1х

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

Разностные уравнения, аппроксимирующие систему на расчетном шаблоне 1-3-2, примут вид

Л, \ 2

_ / К К К 2

«.-и.

Значения искомой функции и(х,г) в фиксированной точке 3 находится с помощью известных методов решения систем алгебраических уравнений.

Заметим, что идея использовать производную искомой функции в качестве дополнительной искомой функции рассматривалась Грудницким В.Г. и Прохорчуком Ю.А. в 1977 г. Однако, авторы записывали исходные уравнения в дивергентном виде и выписывали отдельное дифференциальное уравнение для дополнительной искомой функции. Авторы использовали данную идею для расчета двумерного нестационарного сверхзвукового обтекания источника энергии.

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

Рассмотрим задачу о внешнем обтекании пространственных тел сверхзвуковым стационарным потоком идеального газа. При таком обтекании течение газа будет смешанным - между поверхностями ударной волны (УВ) и тела имеются дозвуковые (I), трансзвуковые (ЗЛ) и сверхзвуковые (II) области (рис. 2). В дозвуковой области течение газа описывается уравнениями в частных производных эллиптического типа, а в сверхзвуковой -гиперболического типа.

Рис. 2. Внешнее обтекание пространственного тела сверхзвуковым потоком газа.

УВ

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

Введем цилиндрическую систему координат (х,г,<р). Тело будем считать неподвижным, а поток на бесконечности равномерным. Система уравнений, описывающая течение стационарного идеального газа, представляет собой уравнение неразрывности, уравнение импульса и уравнение Бернулли

сНу(р6Р6) = 0

= 0 (2)

Рв

2

2 (к-\)ре 2 " к-1

где р5 - плотность газа, Ус =(«,у,и<)'- скорость (верхний индекс штрих - знак транспонирования), р6 - давление газа, к - показатель адиабаты, м„ - число Маха набегающего потока и V— оператор Гамильтона. Все газодинамические функции взяты в безразмерном виде следующим образом

Р р V

Р*=—. Ре= — Р.

где р„ и а

Ус =-

Р. " Р. у1Р„ / давление и плотность набегающего потока газа. Опуская в дальнейшем индекс "б", указывающий на безразмерные величины, запишем систему (2), которая содержит четыре уравнения для четырех искомых

газодинамических функций V, р в следующем матричном виде

ли,+ ви,+сиг = /

где

(3)

А =

и' а' иу ми> 'V 0 иу с? У! а' 0 ИМ* а1 1Л> а' 0

ри 0 0 1 , в = Р» 0 0 0 . с=1 рю 0 0 0

0 ри 0 0 0 ру 0 1 г 0 р\у 0 0

0 0 ри . о 0 рV , 0 0 рм К

/ = -(-у,0,р1У:,-1>ри>У, а1 = к— - местная скорость звука, и=(и,\,ю,р)' - вектор г4 ' р

искомых переменных, и,, V,, - частные производные по соответствующей координате.

Начальные условия задаются на поверхности П. Граничными условиями являются условие непротекания на поверхности тела

(?,«») = О (4)

и соотношения Ренкина-Гюгонио на ударной волне

К=К, РК

Р + РК'=Р.+Р,С (5)

к-\р 2" к-\р„ 2 ™

где V - скорость в точке за ударной волной, - скорость в точке перед ударной волной, - проекции на нормаль п к ударной волне и

—» ч

касательный вектор / к ударной волне, л» - нормаль к поверхности тела.

Система (3), начальные данные на поверхности П и граничные условия на теле (4) и ударной волне (5) позволяют определить искомые величины во всей рассматриваемой сверхзвуковой области движения газа.

В третьей главе описывается численный метод псевдохарактеристик для решения системы (3). Для простоты изложения в качестве поверхности П будем рассматривать плоскость х = хи. Для численного решения задачи введем сетку в рассматриваемой области течения. На плоскостях х, =*„+.*-Л,, где 5 = 1,2,3,... и Иж — шаг сетки по координате х, введем семейства непересекающихся замкнутых линий г1=г/(х„<р), у = 0,7,2,...,У. В дальнейшем эти линии будем называть венками. При этом венок (*,,$>) лежит на ударной волне, а венок г0(*,,?>) — на обтекаемом теле. Точки пересечения венков г, =/■(*,,р) и меридиональных плоскостей <р = (р,, / = 0,1,2,...,N. образуют фиксированную сетку в плоскости х = х, (см. рис. 3). Положение узлов сетки в плоскостях х = х, при ^ > 1 определяются в процессе расчета после нахождения точек о=0(*..<г'|) на ударной волне.

Меридиональные плоскости

Рис. 3. Расчетная сетка в плоскости х =

Расчет течения в рассматриваемой области разобьем на три элементарные задачи:

1. расчет искомых функций во внутренней точке области течения

2. расчет искомых функций на поверхности обтекаемого тела

3. расчет искомых функций на поверхности ударной волны.

\

характеристика

X, X,., X X, У,„ X X, Х3.1 X

Рис. 4. Шаблоны для расчета искомых величин в некоторой меридиональной плоскости: а) в поле течения, б) на поверхности тела, в) на поверхности ударной волны.

Расчет искомых функций во внутренней точке области течения

На линии пересечения плоскостей х,и р, зададим точки 1 и 2, а на линии пересечения и - точку 3 (см. рис. 4а). Положение точек в поле течения задается по следующим формулам (предполагается, что координаты венка О (■*•+!> й-) известны, то есть известно положение ударной волны в плоскостях

■*=■*,♦! и <Р = <Р,У-

где 0<7|,уг,у3 <J. Кроме того, j2-j, =2 и j2-j, =j\-j,=\.

Производная произвольной функции g(x,r,<p) на некоторой линии г = г(х) в

меридиональной плоскости (р = <р, имеет вид — = gx+gr — ■ Производная

¿х 8х

функции g(x,r,<p) на некотором венке г = г{х,,<р) в плоскости х = х, имеет вид ^- = SrС Учетом этого перепишем исходную систему (3) в точке 3 с использованием линий 1-3 и 2-3, получим

'Sr

(

SU в

А +

5х

Ч /13

Г \

SU в

А +

Sx

К Л 3

1 -с-

Лз 5<Р.

Ur = f-C

SU S(p

Scp

и,-f-C™

Sq>

где индекс 1 указывает, что равенство записано на линии 1-3, а индекс 2 - на линии 2-3.

Поверхности, проходящие через точки 1 и 3 расчетного шаблона, будем называть псевдохарактеристическими поверхностями первого семейства, и поверхности, проходящие через точки 2 и 3 — псевдохарактеристическими поверхностями второго семейства. Вводя новую искомую вектор-функцию

Я = иг, получим в точке 3 следующую систему уравнений

( ^

SU

А Sx +

V /13

( ^

5U

А Sx +

У /23

В-Л\£] -с — SxJn Scp

В-А\^) -С — öx )и S<p

Uj-C™

S<p

(7)

Система (7) является замкнутой системой из восьми уравнений относительно восьмимерного вектор-столбца (и,Я)'. Показано, что функция и,

являющаяся решением системы (7) в точке 3, также является решением исходной системы (3).

Производные по х аппроксимируем центральными разностями, а производные по <э дифференцированием кубического сплайна, построенного по всем точкам соответствующего венка (Борисов В.М., Михайлов И.Е. 1978 г.). Разностная схема, аппроксимирующая систему (7) на шаблоне (см. рис. 4а) со вторым порядком точности, примет вид:

X»

Х31 1

Л

я =

(8)

где цк и //<'> - числа при нахождении производной кубического сплайна, 41' = ^ , 4* = > ' =2 _ точки шаблона.

Отметим два свойства системы (8). Во-первых, система имеет 8Ы уравнений и 8Ы переменных, поэтому решение ищется сразу для всех меридиональных плоскостей. Во-вторых, система (8) является нелинейной, так как матрицы содержат искомые функции, и решение ищется с помощью метода Рунге-Кутта с двойным пересчетом, при этом число итераций бралось равным трем (Кацкова О.Н., Наумова Н.И., Шмыглевский Ю.Д., Шулишнина Н.П., 1961 г.).

Расчет искомых функций в точке на поверхности обтекаемого тела

Для расчета искомых функций в точке на твердой стенке используем шаблон, изображенный на рис. 4.6. Координаты точек задаются по формуле (6) при условии, что у, = 0, у, =1 и Уз =0. Разностная система уравнений совпадает с (8).

На поверхности тела выполняется граничное условие непротекания (4), которое является дополнительным уравнением к системе (8). Записывая вместо одного из уравнений системы (8) условие непротекания, получим модифицированную систему, которую также решаем с помощью трех итераций. Опыт расчетов показал, что наилучшая точность достигается при замене условия непротекания на первое уравнение системы (8) вдоль псевдохарактеристики 1-3.

Расчет искомых функций в точке на поверхности ударной волны

Искомые величины за ударной волной определяются из соотношений Ренкина-Гюгонио (5). Однако сложность расчета точки на ударной волне

заключается в том, что угловой коэффициент наклона и положение

ударной волны г,в искомой точке заранее неизвестны и также подлежат определению. Расчет углового коэффициента наклона и искомых газодинамических функций на ударной волне проводился с помощью итераций (как и в методе характеристик Кацкова О.Н., Наумова Н.И., Шмыглевский Ю.Д., Шулишнина Н.П., 1961 г.). Расчетный шаблон показан на рис. 4.в.

Расчет газодинамических величин в рассматриваемой области течения газа осуществляется по следующему алгоритму. Пусть известны газодинамические функции в плоскости дг,. Сначала происходит расчет точки на ударной волне, определяются газодинамические величины и положение ударной волны оС*.»1 >й) Для всех <р: в плоскости . По найденному положению ударной волны и заданным координатам обтекаемого тела ^(х,^,^) находятся координаты искомых точек г,{хм,<р,) для всех р, в плоскости по формулам (6). Далее рассчитываются газодинамические величины в поле течения при \<.J<.J-\ и на обтекаемом теле при у = 0 на плоскости . Таким образом определены искомые величины на всех венках в плоскости Аналогично проводится расчет плоскости дг,,; по известным величинам в плоскости дг„, и так далее.

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

лих+ви, + сит= 0 (9)

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

Р = Ро+Р„, где У0 »

У0 = сот1, ра = сот1, V„ и рм - новые

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

А„ =

V}

--у 0

<

лЛ 0

0 АЛ

0 0

О 1

О О

'0 1 0 0~ '0 0 1 о4

0 0 0 0 С Л > ^ л Г 0 0 0 0

0 0 0 1 0 0 0 0

0 0 о; ,0 0 0

где ра=кЕ± И +

а„ 2 2

Рассмотрим шаблон, изображенный на рис. 5 и введем новый базис единичных векторов (е|,г2,е3), где е, = (е, + гх е,)т, е2 =(ех-гх ег)т,

и ег, е<р - орты вдоль осей г и р.

3" 3 * е. \3'

\ «и/ \ в!

\ еч

V \1 г

/ / 3 у^Н /

2" 2 Г

Рис. 5. Шаблон для расчета точки методом псевдохарактеристик

Записывая систему (9) в новом базисе, приводя подобные слагаемые и вводя новую переменную -(е2У)С/ = Д, получим систему уравнений, не содержащую выводящих производных

2г1 ' г,) 5?, 2[ гх) ' дд,

1 дЦ. г г?2

-л = о

где с/?, - длина дуги линии, к которой вектор а является касательным (/ = 1,2,3). Боковые грани шаблона будем называть псевдохарактеристическими поверхностями

Производные по х и у аппроксимируем центральными разностями, тогда разностная схема примет вид:

Г, J h.

2 2

U Mj - Ui-i) UMJ-I - Ui-ij-i

2h

2h

= 0

Uij -Uij*i Rij + Rij+\

= 0

Для исследования устойчивости воспользуемся спектральным признаком

■ У -»

устойчивости Неймана. Представим искомые вектор-столбцы (U„,Ry в виде

рядов Фурье: U,j = YjY,U е'и{е*( и R,j = e'^e1'1', где U , R - вектор

коэффициентов Фурье, i^kji^, j = krhr, kr и kr - целые волновые числа и /, -мнимая единица.

Для устойчивости метода псевдохарактеристик необходимо выполнение условия (Годунов С.К., Рябенький B.C., 1977 г.)

|Л|<;1 + Оф„) (10)

где X - собственное число матрицы перехода со слоя s на слой s + 1 для

восьмимерного вектор-столбца коэффициентов Фурье (U,R). Собственные числа определяются из уравнения

{W,-D, W4-DJ

/Ш,-05=Л,+——(Л-е + —--—(Л + е ''') и £4 - единичная матрица

I К 2

размера 4x4.

Определитель матрицы перехода представляет собой уравнение восьмой степени относительно Л. Подставляя выражения для определителя и приводя подобные, получим, что оператор перехода имеет два собственных числа Л = ±1 кратности два. Таким образом, модуль первых четырех собственных чисел равен единице

Ш=КИ> ИзЬКИ

Остальные четыре собственных числа определяются из уравнения

q4A* +д2Л2 +q¡A + q0 =0

(П)

К l6hÍ U гг) Л, h) 16 h\ l^r, г,

г) sin2ífl lY ^Г, rj

Решая данное уравнение с помощью замены у = л + —, получим

Л

где у, =

( V М -4 / \ —-21

J

^ И

;

4i

Г \ h. 2 -4 / \ 1--А

U<J J

Определим соотношение между параметрами течения, при которых выполняется (10). Для удобства введем следующие коэффициенты А,>0 и к2> 0, связывающие шаги по пространственным переменным между собой:

К = КК

к=кк

(12)

Рассмотрим течение газа в сверхзвуковой области, то есть такой, что в любой

V1

точке области выполняется неравенство М]ох > 1, где М\ак = -у. Значения

Оо

собственных чисел, являющихся решениями уравнения (11) при различных значениях входных параметров, вычислялись на компьютере. Анализ полученных результатов позволил сделать следующие выводы относительно устойчивости метода псевдхарактеристик.

1. Метод псевдохарактеристик удовлетворяет условию устойчивости при выполнении определенных ограничений на шаги Их и А,, При выполнении этих ограничений, дополнительные ограничения на шаги по И9 не требуются.

2. Метод псевдохарактеристик остается устойчивым при г->0. Поэтому метод может быть также применении при расчете окрестности оси симметрии при внутренних течениях газа.

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

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

(И)

Графически это означает, что область устойчивости метода псевдохарактеристик находится под графиком к2 тах. Ниже на рис. 6 показан

график А2тах{М2тк) при 0(/гг) = 0,05. «Квадратами» указаны значения кг П1ах, вычисленные автором с помощью компьютера.

м2„.

Рис. 6. Зависимость коэффициента к2та от М]ОК.

Отметим, что спектральный признак оценивает устойчивость только по начальным данным. Хотя данный аппарат не дает точного ответа на вопрос об устойчивости, тем не менее схемы, признанные на его основе устойчивыми, как правило, на самом деле являются таковыми, и в более или менее общем случае из такой устойчивости следует устойчивость по правой части. Анализ устойчивости по краевым условиям выходит за рамки данной диссертации.

В пятой главе приводятся результаты расчетов и сравниваются результаты с существующими табличными. В качестве первой задачи рассматривался расчет осессимметричного течения идеального газа (¿ = 1,4)

около кругового параболоида г1 под нулевым углом атаки.

Осесимметричное течение характеризуется тем, что газодинамические величины равны во всех меридиональных плоскостях, таким образом Uv = 0. Число Маха набегающего потока = 2. Начальные данные брались из таблиц1 на плоскости х0=1,5. Число венков J-9, шаг hx= 0,1. Ниже на рисунках 7-9

1 Любимов А.Н., Русанов В.В. Течение газа около тупых тел. М.: Наука, 1970. Т.2. Таблицы газодинамических функций. 379 с.

представлены расчеты методом псевдохарактеристик искомых газодинамических величин на теле в зависимости от координаты х.

1,450 1,400 1,350 1,300 1,250 1,200 1,160 1.ЮО 1,050 1,000

Рис. 7. Табличные (квадраты) и расчетные (ромбы) значения давления газа в точке на обтекаемом теле

1,950

1,900 ...................-..................:...........................................................................:..................-...................

23456789 10

Рис. 8. Табличные (квадраты) и расчетные (ромбы) значения скорости и газа в точке на обтекаемом теле

0,550 0,500 0,450 0,400 0,350 0,300 0,250 0,200

Рис. 9. Табличные (квадраты) и расчетные (ромбы) значения скорости V газа в точке на обтекаемом теле

В качестве второй задачи рассматривался расчет внешнего обтекания эллиптического параболоида 2х-(5со$2р+Зь'т2р)-г2 =0 идеальным газом {к = 1,4)

под различными углами атаки 5° и 0°. Число Маха набегающего потока = 4. Начальные данные брались из вышеуказанных таблиц для х0 = 2. Число меридиональных плоскостей N = 10, число венков 7 = 9, шаг Их=0,\. При расчете внешнего обтекания сверхзвуковым потоком эллиптического параболоида методом псевдохарактеристик были получены результаты, которые представлены на рис. 10-13.

3,20 <5, 2,70 2,20 1,70 1,20 2

N.

Рис. 10. Табличные (квадраты) и расчетные (ромбы) значения давления на теле в меридиональной плоскости ф=20° при угле атаки а) 5е и б) 0°.

4,30

4,10

3,90

3,70

3,50 ;

3,30

4,10 3,90 3,70 3,50 3,30

6 6)

Рис. 11. Табличные (квадраты) и расчетные (ромбы) значения скорости и на теле в меридиональной плоскости ф=20° при угле атаки а) 5° и б) 0°.

0,90 0,70 0,50 0,30

Рис. 12. Табличные (квадраты) и расчетные (ромбы) значения скорости V на теле в меридиональной плоскости ф=20° при угле атаки а) 5° и б) 0°.

-0,08 -одо -0,12 -0,14 -0,16

Рис. 13. Табличные (квадраты) и расчетные (ромбы) значения скорости м> на теле в меридиональной плоскости ф=20° при угле атаки а) 5° и б) 0°.

Как видно из вышеприведенных рисунков, результаты расчета методом псевдохарактеристик хорошо согласуются с табличными данными.

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

1. Разработан эффективный вычислительный маршевый метод псевдохарактеристик для расчета нелинейных дифференциальных уравнений в частных производных гиперболического типа с тремя независимыми переменными. Метод использует фиксированные шаги по маршевой переменной и имеет второй порядок аппроксимации. Метод реализован для расчетов стационарного пространственного обтекания тел сверхзвуковым потоком газа под различными углами атаки. Результаты расчетов с высокой точностью совпадают с данными таблиц Любимова А.Н. и Русанова В.В., которые, в свою очередь, хорошо согласуются с экспериментальными данными.

2. Разработана математическая модель стационарного течения сверхзвукового газа в цилиндрических координатах (х,г,<р), линеаризованная вдоль маршевой переменной х. Модель используется для исследования устойчивости разработанного вычислительного алгоритма и его области применимости, в том числе при г-*0. Получены ограничения на шаги по пространственным переменным, обеспечивающие устойчивость алгоритма. Показано, что при г-*О алгоритм сохраняет устойчивость, что позволяет его использовать и для расчетов течений газа в окрестности оси симметрии. Отметим, что ранее для исследования устойчивости разностных схем для стационарных пространственных течений на основе характеристических соотношений использовались линеаризованные уравнения только в декартовых переменных (х,у,:). Такой подход не позволял исследовать устойчивость в окрестности оси симметрии.

3. Алгоритм метода псевдохарактеристик реализован в виде комплекса программ, написанного с использованием пакета прикладных программ МайаЬ. Комплекс программ адаптирован к расчету обтекания стационарным сверхзвуковым потоком газа пространственных тел под различными углами атаки.

СПИСОК ПУБЛИКАЦИЙ ПО ТЕМЕ ДИССЕРТАЦИИ

1. Михайлов И.Е., Нурбаев УД. Метод псевдохарактеристик для расчета стационарных пространственных сверхзвуковых течений газа // Теория, численные методы и математический эксперимент в газовой динамике. Материалы Научно-практического семинара. Москва. 25-26 августа 2009 г.-М.: ЦИАМ, 2009. С. 91-92.

2. Нурбаев УД. Об одном методе расчета стационарных пространственных сверхзвуковых течений газа // Современные проблемы фундаментальных и прикладных наук. Часть VII. Управление и прикладная математика. Труды LIII Научной конференции. / МФТИ. - М. - Долгопрудный, 2010. С. 104-106

3. Mikhailov /., Nurbaev U. On method of pseudocharacteristics for solving 3D stationary supersonic flow of ideal gas // VII International Aerospace Congress IAC'12, August 25-31, 2012, Moscow, Russia, Proceedings. Электронный вид. Зарегистрировано в ВГУП НГЦ в ИНФОРМ-РЕГИСТР. Гос. per. № 0321303652. 2013. Р. 224-229.

4. Михайлов И.Е., Нурбаев УД. Метод псевдохарактеристик для расчета стационарных пространственных сверхзвуковых течений газа // Вестник Тверского Государственного Университета. Серия: Прикладная математика. Тверь, 2013. № 15. С. 19 — 27.

5. Нурбаев УД. Об устойчивости метода псевдохарактеристик // Вестник Тверского Государственного Университета. Серия: Прикладная математика. Тверь, 2014. №2. С. 21-28.

6. Нурбаев УД. Метод псевдохарактеристик для расчета пространственных сверхзвуковых течений газа // Свидетельство о государственной регистрации программы для ЭВМ № 2015610150 от 12 января 2015 г.

Подписано в печать:

02.07.2015

Заказ № 10812 Тираж -120 экз. Печать трафаретная. Типография «11-й ФОРМАТ» ИНН 7726330900 115230, Москва, Варшавское ш., 36 (499) 788-78-56 www.autoreferat.ru