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

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

САХАБУТДИНОВ Жавдат Мирсаяпович

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

Специальность: 05.13.16 - применение вычислительной техники, математического моделирования и математических методов в научных исследованиях (механика)

АВТОРЕФЕРАТ диссертации на соискание ученой степени доктора физико-математических наук

Казань - 1997

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

Официальные оппоненты: доктор физико-математических наук,

с.н.с. Ю.Н. Дерюгин, доктор физико-математических наук, профессор В.Н. Игнатьев, доктор физико-математических наук, профессор В.В. Клоков.

Ведущая организация: Институт механики Уфимского научного центра РАН

Защита диссертации состоится " сС^-^^-^Я.—1997 г в

часов на заседании диссертационного совета по защите

докторских диссертаций Д 063.43.03 при Казанском государственном техническом университете имени А.Н. Туполева по адресу: 420111, Казань, ул. Карла Маркса, 10, КГТУ.

С диссертацией можно ознакомиться в библиотеке ЮТУ им. Туполева.

9 С>

Автореферат разослан " ° " ^■л-^-с-'е _ 1997 г.

Ученый секретарь диссертационного совета П.Г. Данилаев

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

Актуальность темы.

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

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

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

1. Создание эффективной методики по интегрированию полной системы трехмерных уравнений газовой динамики. Она подразумевает:

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

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

- реализацию алгоритмов быстрой обработки результатов трехмерных нестационарных расчетов.

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

Новизна работы состоит в

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

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

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

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

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

Диссертационная работа выполнялась на основе Прнкаи-Распоряжсння N 71/32 от 1 марта 1985 года между Президиумом АН СССР и Мин-сельхозмашем; в соответствии с планом научных исследований Инстигута механики и машиностроения по бюджетной теме РАН (номер государственной регистрации 01.9.60. 000596). Проблемно ориентированные варианты программного комплекса были реализованы по заданию ПО ЕлаЗ и АО КамАЗ и нереданы им в виде программ на РС и научно-технических отчетов.

На защиту выносятся следующие результаты:

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

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

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

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

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

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

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

Апробация работы.

Основные положения диссертации докладывались на: V Всесоюзном сьезде по теоретической и прикладной механике (Алма-Ата, 1981); Научно-технической конференции "Механика сплошных сред" (Набережные Челны, 1982); Семинаре ВЦ АН СССР под рук. проф. Шмьвлевского Ю.Д. (Москва, 1984); Семинаре академика Рахмазулина Х.А. (Москва, 1984); Семинаре академика Самарского A.A. (Москва, 1984); Всесоюзной летней школе по теории взаимодействия упрушх оболочек с жидкостью, газом и твердым телом (Казань, 1986); Всесоюзной школе-семинаре "Математическое моделирование в науке и техгппсе", (Пермь, 1986); VIII Всесоюзном симпозиуме по распространению упрушх и упруго-пластических волн, (Новосибирск, 1986); I Всесоюзном рабочем совещашш по методам построения сеток, (Свердловск, 1987); Семинаре им. К.И. Бабенко в ИПМ им. М.В.Келдыша, (Москва, 1987); Всесоюзной научно-технической конференции "Перспективы развития комбинированных двиттелей внутреннего сгорания и двигателей новых схем и на новых тогишвах" (Москва, 1987); VII Всесоюзном семинаре "Теоретические основы и конструирование численных алгоритмов решения задач математической физики" (Кемерово, 1988); II Всесоюзном рабочем совещании по методам построения сеток (Кемерово, 1988); VIII Всесоюзной научно-

технической конференции "Создание компрессорных машин и установок.", (Сумы, 1989); Семинаре академика Сидорова А.Ф. (Свердловск, 1989); III Международном Семинаре "Структура пламени" (Алма-Ата, 1989); захщгге научно-технических отчетов в отделе шавного конструктора ПО Ел A3 в 1986— 1990ir, защите научно-технических отчетов в отделе главного конструктора АО КамАЗ в 1991-1992rr, V Всероссийском совещании "Проблемы построения сеток для решения задач математической физики" (Казань, 1994); Международной конференции "Фундаментальные исследования в аэрокосмической науке" (Жуковский, 1994); Семинаре академика Нигматулина Р.И. (Уфа, 1996); Международной конференции по спектроскопии CSI (Буэнос-Айрос, Аргентина, 1996); Итоговых конференциях Казанского государственного университета (1989 -1996 it.); Итоговых научных конференциях Казанского физико-технического института КФ АН СССР (Казань, 1980- 1996 тт.); Итоговых научных конференциях Института механики и машиностроения КНЦ РАН (Казань, 1993-1996 гг.).

Объем Ii структура диссертации. Публикации. Диссертация изложена на 269 страницах машинописного текста, состоит из введения, пяти глав и выводов. Работа содержит 3 таблицы, 138 рисунков. Список использованной литературы содержит 182 наименования. Основные результаты диссертации опубликованы в 33 работах, список которых содержится в конце автореферата.

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

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

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

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

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

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

Первый этап основан па интерполировании координат внутренних узлов Л/д-гю траничным поверхностям (рис. 1)

Рис. 1.

Рис. 2.

«а = + «Л1 - +И - ка - - с)}п-,л +

+- ка - - «ОК1 - +{к - Кй - ^ - С) --[к-ка-кь-аЫ -К -кМ1-^)!1-^)}^+

+ {к^ - - - С'л) - К - К-й - КлО - -

(1)

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

27 / 27

(2)

где И/ - радиус-векторы узлов, лежащих в 27 вершинах шаблона (рис. 2), щ _ двухзначная функция, равная 1, если /-ьш узел участвует в определении нового положения цсшральной точки (;,_/,£), и 0, если нет. Индекс I дает локальную нумерацию узлов шаблона. В алгоритме используются варианты формулы (2), когда в усреднении участвуют либо все узлы шаблона, либо ближайшие по трем координатным направлениям (/= 5, 11, 13, 15, 17, 23), либо ближайшие по одному из координатных направлений (например, I- 13, 15). Для каждого варианта формулы (2) имеется своя совокупность объемов ячеек, входящих в шаблон. Для них определяются значения

= ппн(ч), Ута( = тах(у,.), (/ = 1,... ,8).

Если для всех вариантов формулы (2) Ут-т<О, то новое положение узла определяется по той из них, для которой \;ПцП имеет наибольшее значение. Иначе из тех вариантов формулы (2), которые дают Ут1П>0, в качестве расчетной берется та, для которой разность (Утах-\'тиг) принимает наименьшее значение. Такой выбор приводит к тому, что в процессе итераций отрицательные объемы устремляются к нулю, а положительные - к относительному их выравниванию. Итерационный процесс длится до тех пор, пока все внутренние узлы не окажутся в "неподвижном" состоянии, которое определяется условием

(3)

где £ - задаваемая малая величина, Ь - характерный размер для данного шаблона, V - номер итерации.

Вопросам оптимизации полученной трехмерной сетки по некоторому набору критериев качества посвящен третий этап. Такими критериями взяты гладкость и ортогональность сеточных линий, а также равномерность сетки по объемам ячеек. Отыскивается экстремум целевой функции /(*•>. г) = Г{х,у,г) + Л,/У(х,у,г) + Ло/°(д:,у,г) + Ль/*(х,у,г) (4)

для шаблона. Здесь Р- мера гладкости криволинейных координатных линий, /V - функция, обеспечивающая примерное равенство объемов соседних ячеек, /> - мера ортогональности линий сетки в узле, А> - мера ортогональности у границ шаблона; Лу, Ад, Х}} - заданные весовые

коэффициенты для соответствующих функций. Поскольку f(x,y,z) является неотрицательной непрерывной функцией, причем

то она достигает своего минимального значения в некоторой точке P(x,y,z). Определив ее координаты, найдем такое положение внутреннего узла шаблона, которое удовлетворяет перечисленным выше требованиям. Обходя все внутренние узлы в процессе итераций, получим такое их взаимное расположение, при котором требование минимума функции (4) приближенно выполняется во всех точках области. Критерием окончания итерационного процесса служит условие вида (3).

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

1) значение Яу, отвечающее за локальную равномерность размеров ячеек, должно изменяться в пределах от 0 до 200 (при нулевых Ло и Л/,). Наиболее быстрые изменения структуры сетки происходят на интервале от 0 до 100. Дальнейшее увеличение Äv мало влияет на движение узлов сетки. При достижении им некоторого значения (в рассмотренных примерах это ÄV-20Q) сеточные линии становятся пилообразными, а итерационный процесс перестает сходиться;

2) контроль за движением узлов в направлении ортогональности сеточных линий осуществляется с помощью параметров Äq и Я/,. Увеличивая /¿) от 0 до 5-7 и не меняя остальные параметры, можно добиться значительного усиления ортогональности. Дальнейшее увеличение Л[, ведет к захлсстам внугренних узлов за границы расчетной области. Изменение параметра Äq в меньшей степени влияет на структуру сетки. Однако при значениях Ä[j>Äq параметр Äq играет роль ограничителя движения внутренних узлов. Увеличивая Äq при постоянном Л¿>, можно добиться "возвращения" внутрь области узлов, попавших за ее пределы;

3) анализ примеров показывает, что наиболее приемлемые разбиения получаются, когда все параметры качества берутся ненулевыми: Лу~ 100-^200, V -10, Л),~ 1+10.

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

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

Рассматривается задача о холодной продувке впускного канала потоком воздуха. Во входном сечении канала задается перепад давлений Лр=р2-р\, где /»[-давление в невозмущенной области, /^-давление во входном сечении. Структура сетки для всей компоновки показана на рис. 3.

В начальный момент времени г=() плоская ударная волна с давлением за фронтом 1.36 атм располагалась вблизи верхнего сечения впускного канала. Затем она проходит цилиндрический участок и начинает взаимодействовать с вогнутой поверхностью канала. Давление в этой области канала увеличивается до 1.6 атм. В дальнейшем у погнутой стенки канала (формируется отраженная ударная волна. Достигая противоположной стенки канала, она создает область повышенного давления -1.7 атм. Ко времени ?=0.0002с основные возмущения проходят нижнее сечение расчетной области. Под верхней крышкой гильзы вблизи стыка с каналом возникает область разрежения (рис. 4).

Изменение ноля скоростей в этом же сечении дано на рис. 5. В местах разрежения в гильзе начинают зарождаться вихри, которые по мере установления течения становятся более ярко выраженными. Характер течения в гильзе камеры сгорания можно изучить по рисункам, представляющих поля скоростей в сечениях г=согШ. Одно из таких сечений на глубине г=-10мм дано на рис. 6. Из анализа полей скоростей в них следует, что в этих сечениях отсутствуют области завихренности. Основными особенностями течения здесь являются движение газа к боковым стенкам гильзы и возвратное течение вблизи верхней крышки гильзы.

На рис. 7 представлена диаграмма зависимости массового расхода в выходном сечении капала от времени для двух перепадов давления 0.36апш и 0.24шш{. Ко времени г=0.0002с возмущение доходит до выходного сечения канала, и расход в нем становится отличным от нуля. Характер роста расхода (его немонотонность) объясняется сложной волновой картиной течения во впускном канале. Ко времени г=0.0006с все основные возмущения проходят выходное сечение канала и массовый расход колеблется около некоторого среднего значения: 0.23кг/с для перепада О.Збатм, 0.15кг/с для перепада 0.24атм.

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

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

клапана. Мгновенное значение площади проходного сечения зависит от величины подъема клапана, геометрий клапанной тарелки, седла и штока.

Рассматривается задача о холодной продувке воздухом полости винтового канала с цилиндрической камерой сгорания при наличии впускного клапана. Начальные условия те же, что и в предыдущей задаче. Кинематические и термодинамические параметры на входной границе полагаются заданными. На выходной границе реализуемся условие продолжения среды как экстраполяция нулевого порядка скорости, плотности и внутренней энергии. На твердых стенках канала, клапана и гильзы камер!,I сгорания компоненты вектора скорости определяются из условия проскальзывания.

При построении трехмерных сеток используется блочный метод. Канат имеет сложную конфшурацию, фаничные поверхности задаются набором координат точек, снятых с чертежа изделия. На рис. В построены расчетная сетка в сечении, проходящем через ось штока и область стыковки канала с крышкой цилиндра дня двигателя Б-216Т. Общее представление о кошуре впускного капала дает рис. 9. Он получен в результате сечения полости канала вертикальной плоскостью, уравнение которой -0.27х+}>=0.

На рис. 10 представлены проекции скоростей на плоскость 2х+ 5^=210 для момента времени /=0.418-10"^с. Возмущения в среде достигли спиралевидной области канала. Здесь основная часть потока слева обходит шток впускного клапана и разворачивается. В то же время менее интенсивная струя газа обтекает шток справа по ходу, ударяется о стенку, разворачивается и сталкивается с основным потоком. В узкой горловине канала формируются локальные вихри, поток в целом устремляется вниз к выходному сечению клапанного отверстия. Максимальное значение скорости в сечении равно 30.1м/с.

Рис. 8. Рис. 9.

Из рис. 11 видно продвижение потока по узкой горловине канала и клапанной щели. Поток растекается по поверхности тюльпана и начинает проходить клапанную щель. Максимальное значение скорости в клапанной щели к моменту времени /=1.98-10"Зс равно 49м/с, а в наиболее узком зазоре между тарелкой клапана и боковой поверхностью скорость достигает 53.47м/с.

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

Особенностью представленного на рис. 12 потока на глубине г=-6мм является радиальная направленность скоростей газа ог тарелки клапана к боковым стенкам гильзы цилиндра. На глубине г=-80мм наблюдается образование двух вихрей, вращающихся в противоположных направлениях (рис. 13).

Рассмотренные выше результаты приведены для максимальной высоты подъема клапана к= 12мм. Вычисления были повторены и для средней высоты подъема клапана /г=8мм. Из сравнения с расчетами для максимальной высоты подъема клапана видно, что картины течений качественно похожи. На основе количественных оценок можно отметить, что для высоты подъема /г-8мм значения скоростей в клапанной щели в разные моменты времени на 20-23% больше, чем для 1г= 12мм.

Проведение расчетов течения газа для малой высоты подъема клапана (/¡=4мм и меньше) из-за сильной нерегулярности ячеек наталкивается на значительные математические трудности. При используемой логической системе координат получаются слишком вытянутые ячейки. Это обстоятельство при работе газодинамического комплекса приводит к достижению максимального числа итераций в неявной фазе расчетов. Шаг интегрирования <4/ становится очень малым для проведения реальных расчетов и программа прерывает вычисления. Используемый подход к решению поставленной задачи для малых высот подъема клапана результатов не дает. Проведение расчетов для высоты клапана /г=4мм и менее требует иною подхода к разбиению на подобласти и установлению логических связей между ними.

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

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

MARS FLOW «а/SOC

Рис. 14.

В следующем параграфе главы рассматривается задача о движении газа в цилиндре двигателя с камерой сгорания в днище поршня на тактах сжатия и расширения. Угол поворота и геометрические характеристики крнвошипно-шатунного механизма связаны с мгновенной скоростью поршня Wp соотношением:

cos в

W = - IV

п 2

1

(/2 / а2 - sin2 6>)"2

sin О,

где И'о-средняя скорость поршня, На - отношение длины шатуна к длине кривошипа. Схема поршневого двигателя представлена на рис. 15, на котором зависящая от времени расчетная область окаймлена штриховкой.

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

головкой цилиндра. Для снятия его применяется полностью консервативная процедура перестройки сетки, которая при достижении некоторою значения г1 начинает удалять слой ячеек, увеличивая тем самым размеры ячеек в сжимаемом направлении. В терминах угла поворота кривошипа началу перестройки сетки соответсвует значение 9= -45°. Для углов, значения которых больше #=-27°, число ячеек вдоль оси Ог принимает минимальное значение, равное 3. При обратном движении поршня (ход расширения) происходит постепенное добавление ячеек вдоль 01 до полного восстановления прежнего значения.

Начальные условия задачи задаются однородными но всей области: плотность Д)=1.23кг/м3, давление /?0=1атм, внутренняя энергия вычисляется на основе уравнения состояния. Размеры и характеристики кривошипно-шатун-нош механизма выбирались одинаковыми для всех рассмотренных ниже вариантов задачи: диаметр цилиндра Д^-=9.843см, ход поршня ¿=9.95см, дайна шатуна /= 16.269см, диаметр и глубина камеры сгорания Э^Зсм и Я/Ь=3см соответственно. Минимальное расстояшю между головкой цилиндра и поршнем в верхней мертвой точке (ВМТ) равнялось 0.18см. На поверхности подвижного поршня выполнялись кинематические условия - равенство нормальной составляющей газа скорости поршня. На остальных неподвижных поверхностях области выполнялись условия скольжения. Ниже представлены результаты расчетов, выполненных для разных компоновок "цилиндр плюс камера сгорания" с анализом влияния того или иного фактора на структуру потока.

Рассматривается цилиндрическая симметрично размещенная камера сгорания при отсутствии начального вращения потока. Расчетная область состоит из двух соосных цилиндров - цилиндра двигателя и симметрично расположенной камеры сгорания в днище поршня. Вычисления начинаются дтя угла поворота кривошипа в 90° до верхней мертвой точки (ВМТ), что соответствует среднему положению поршня и максимальной его скорости. На рис. 16 для 6=-56° до ВМТ представлено поле скоростей в вертикальном сечении, проходящем через ось цилиндра. Максимальное значение скорости газа в сечении практически совпадает с мгновенной скоростью поршня.

Развившееся движение затекания газа из сдавленной области над поршнем в камеру сгорания видно в сечении 1=сош при #=-23° до ВМТ (рис. 17). Радиальная скорость над камерой сгорания достигает 2016см/с. что примерно в три раза больше максимальной скорости движения поршня. Движение затекания газа в вертикальной плоскости для этого же момента времени представлено на рис. 18.

II /«

II г. ■ /,

IX

1

I 1

/

Г I

II I 1 м

11

Рис. 15.

Рис. 16

При &=0° в камере сгорания формируется кольцевой тороидальный вихрь, поперечное сечение которого представлено на рис. 19. Поток в ВМТ характеризуется значительной радиальной скоростью затекания в область над камерой сгорания. При движении поршня вниз от головки цилиндра (#=-24° после ВМТ, рис. 20) можно наблюдать механизм растекания, коща воздух из внешней окрестности камеры сгорания "выплескивается", падает на поверхность головки цилиндра, проникает в зазор между плоскостью днища поршня и головки, разворачивается и увлекается нисходящим движением поршня. Вихрь какое-то время остается в камере сгорания, затем вытягивается и постепенно исчезает.

Рис. 17.

Рис. 18.

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

< « « II

т—■-1-1-г

Рис. 19. Рис. 20.

Исследовалось влияние на структуру потока нессимметричного раз мещения камеры сгорания. Подобные конструкции со смещенной камеро! сгорания имеют место в двигателях с наклонно расположенной форсункой (двигатели "косого впрыска"). В рассматриваемом варианте камера сгораши сдвинута по оси (к на величину л£=0.75см. Все характерные фазы, выявленные для симметрично расположенной камеры сгорания, присутствуют и в это? задаче, но в несколько видоизмененном виде. Структура движения затекания существенно отличается от симметричного варианта расположения камеры Поток оказывается несимметричным, основное движение газа смещается I сторону смещения камеры. Преобладающее влияние объема сдавленно? области в одной части цилиндра в дальнейшем приводит к организации вихря ближе к противоположной части камеры сгорания. На рис. 21 показан момет времени, соответствующий верхней мертвой точке. Образовался отчетливый вихрь в правой половине камеры сгорания. При движении поршня вниз от ВМТ структура течения фактически не отличается от полученной для задачи с несмещенной камерой сгорания.

Влияние геометрии камеры сгорания на внугрицилищровые течения рассматривается на примере квазиквадратной смещенной камеры. На рис. 22 представлена вычислительная сетка для такой задачи. Видно, что качество расчетных сеток является вполне удовлетворительным. Размеры камеры подобраны так, чтобы степень сжатия (отношение максимального объема цилиндра к минимальному) совпадала с рассмотренными выше случаями.

Рис. 21. Рис. 22.

Все основные этапы в развитии течения на такте сжатия для квадратной в плане камеры сгорания качественно совпадают с процессами в

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

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

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

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

Внутри трубки при помощи специального крепления "подвешена" платформа, являющаяся элементом цилиндра. Боковые кромки платформы располагаются несимметрично относительно плоскости х=0, торцы ее снабжены специальными бортиками. На рис. 24 дано сечение расчетной области плоскостью у=0.

у,

Рис. 23.

о '

------------- ,:21-----------

Рис. 24.

Стенки графитовой трубки нагреваются электрическим током от 300А' до ЗОООК за 2.5-3 секунды. Для описания процессов в газе в уравнении энергии учитываются процессы переноса теша. Вектор потока тепла моделируется на основе гипотезы Фурье. Интенсивный подвод теши предполагает сильную зависимость коэффициента удельной теплоемкости графита, коэффициентов вязкости и теплопроводности газа от температуры. Описывается процедура их расчета в случае одно- и многоатомного ¡азов. Обсуждаются детали формулировки граничных условий, задания начального состояния.

Нагрев платформы происходит за счет теплообмена излучением и взаимодействия с газом. Записываются уравнения теплового излучения. Температура платформы определяется из решения уравнения

2н'Цге + г){ге - #;.) р с,* с1Тр / Л = Рх + Р2 + • м), (5)

где скудельная теплоемкость материала платформы, Бр - площадь поверхности платформы,

Р> = фаР,Х - Р,Х)> Р2 = фа(р=Х + + - Р^Тг\

Последнее слагаемое в правой части (5) определяет передачу тепла от платформы к газу.

Есгш температуры стенок и платформы не зависят от пространственных координат и не учитывается теплообмен с газом, то уравнение (5) при заданной 7И(0 даегг зависимость Тр((). На рис. 25 предполагаемый закон на!репа стенок трубки представлен сплошной линией, а рассчитанные кривые нагрева платформы для двух значений коэффициеггга поглощения Са отмечены индексами 1 и 2. В интервале 0<?<2с температура платформы значительно ниже температуры стенок трубки, что полностью согласуется с требованиями к технологическому процессу атомизации исследуемого вещества

ч 3500 К

К 3000

2500 2000 1500 1000 500 0

0 1 2 3 4 5 6

t tima

Рис. 25.

Изучаются процессы, протекающие в атомизаторе THGA фирмы Perkin-Elmer1. Рассматривается случай, когда течение симметрично относительно плоскости г=0. В таком приближении пренебрегается влиянием узла крепления платформы на поток газа. Полная задача рассматривается как сумма идентичных потоков по обе стороны от плоскости В качестве несущего газа используется азот.

Закон изменения температуры стенок трубки в интервале от 0 до 6с задается дискретно через каждые 0.1с в соответствии с экспериментальными данными1. Под действием интенсивного притока тепла азот расширяется, возникает истечение газа из внутренней полости атомизатора через круговой вырез и кольцевой зазор. Горячий газ из активной зоны нагрева обменивается энергией с холодным газом в буферной зоне.

На рис. 26 к моменту времени ?=0.1с представлены поля скоростей и температур в плоскости х=0. Видно, что скорость в буферной зоне мала (рис. 26-а). Градиенты скоростей значительны в области кольцевого зазора, в окрестности бортиков платформы и на выходе из кругового отверстия в трубке. Максимальные значения скорости истечения азота достигаются в круговом вырезе. Отметим, что примыкающая к платформе масса газа малоподвижна, поток тормозится из-за наличия сравнительно высоких бортиков. Ко времени í=0.1c температура газа в буферной зоне практически не изменилась (рис. 26-Ь). Максимально разогретые ячейки

1 M. Sperling, B. Welz, J. Hertzberg, C. Ricck and G. Marowsky. Temporal and spatial temperature distributions in transversely heated graphite tube atomizers and their analytical characteristics for atomic absorption spectrometry. Spectrochimica Acta, Part B, CSI XXIX Post Symposium, University of Ulm, Germany, September 1995.

(Т=5ПК) располагаются па внутренней поверхности трубки 4.9<г<9мм, температура платформы поднялась лишь до 324К.

Пространственное распределение температуры при /=0.1с дают изолинии в трех сечениях трубки: г=0.2, 5.1 и 7мм, представленные на рис. 27 -а, -Ь и -с. Поле температур в сечении г=0.2мм, перпендикулярном образующей платформы, неоднородно (рис. 27-а, 0-3072К, 9-519.ПА', ЛТ=2А.03К). Температура азота максимальна у стенок трубки. Небольшое остывание потока имеет место на выходе из трубки (отверстие расположено в шжней части крута), менее всего прогрет газ над поверхностью платформы.

Влияние платформы на распределение температуры в газе видно из рис. 27-Ь, (0-343.12АГ, 9-520.33К, АТ=20.09К)- В этом случае сечение г=5.1мм располагается полностью в газе и проходит в непосредственной близости от бортиков. При смещении секущей плоскости вдоль оси 0; ближе к концу трубы г=7мм температурное поле в газе становится однородным и определяется только нагревом цилиндрических стенок трубки (рис. 27-с, 0-374.43Л", 9-520.93К, АТ=16.6\К).

Наиболее интенсивно процесс нагрева стенок трубки протекает до времени г=1.2с. На рис. 28-а представлена картина распределения температуры в сечении у=0 при г=1с. Из нее видно, что газ хорошо профелся в области 5.1<г<9мм (Т~2500К). Над платформой температура газа на 1000—1200АГ ниже, так как значительная часть энергии уходит на прогрев платформы. Тепло медленно распространяется вглубь буферной зоны.

0= 302.1380 1= 324.2610 2= 348.5730 3= 372.8840 4= 397.1950 5= 421.5060 6= 445.8170 7= 470.1280 8= 494.4390 9= 516.5630

С

Рис. 26.

Рис. 27.

Сложный существенно трехмерный характер потока виден из рис. 28-Ь, на котором для /=1с построено поле скоростей в сечении 2=5.1 мм. В окрестности взаимодействия потока с кромками бортиков платформы образуются вихревые движения газа. Несимметричное расположение центров вихрей относительно плоскости х=0 согласуется с тем, что платформа развернута вокруг оси Ог на угол в 45°. Заметам, что величины составляющих скоростей в представленной плоскости значительно меньше составляющей вдооль оси Ог.

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

Из рис. 29-Ь видно, что вся внутренняя область трубки 0^г^9мм прогрета до 2500К. Нагретый газ более чем на половину длины проник в буферную зону. Поле скоростей при г=3с в сечении г=5.1мм уже не содержит вихрей, распределение температуры в сечении г=5.1мм становится однородным.

На рис. 30 представлены зависимости массового расхода от времени в характерных сечениях. Сплошной линией отмечен расход на входе в трубку в сечении г=9мм, штриховой линией - в круговом отверстии, штрихпунктирной линией - в кольцевом зазоре. Из рисунка видно, что наиболее интенсивно процессы развиваются в интервале 0йГ^0.5с. Начиная с момента />1.5с кривые слабо осциллируют около некоторого среднего значения и постепенно выправляются в горизонтальную прямую.

Оа 322.0990

1= 544.3650

2= 788.5120

3=1032.8610

4=1277.1080

5=1521.3560

6=1765.6040

7=2009.8510

8=2254.0990

9=2476.3650

а)

Рис. 28.

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

0= 322.9270 553.3030 2» 806.4640 3=1059.6250 4=1312.7860 5=1565.9460 6=1819.1070 7=2072.2680 8=2325.4280 9=2555.8050

а) Ь)

Рис. 29.

Для рассчитанного выше атомизатора ТСЮА большая серия эксперимент по прсхпранственному и временному распределению температуры в газе, I стенках печи и платформы была опубликована в упомянутой в сноске 1 рабоп На рис. 31 приведена схема пространственного расположения тех точек, в кот рых проходились экспериментальные замеры и результаты сравнивались с рг четами. В прилагаемой таблице перечислены точки замеров и их координаты.

00 05 1.0 15 20 25

TIME (s)

Рис. 30.

Реальный закон нагрева стенок графитовой трубки показан на рис. 32 закрашенными треугольниками, соединенными штриховой линией. Незакрашенными перевернутыми треугольниками отмечена измеренная температура платформы в зависимости от времени. Максимальная разница температур трубки и платформы п 860К достигается при t- 1.2с. С течением времени (f>3c) обе зависимости сливаются в одну кривую.

z=12

точка хоордЕяаты

к 0;0;12

В 0;0;7

С 0; В;5

D 0;0;2.5

Е 0; 0; 0

G 0,-1.9;0

I 0;1.7;0

J 1.7;0;0

К 0.85;0;0

L -1.4;0;0

И 0;1.7;5

Я 0;-1.7;5

р 0;1.7;7

Рис. 31.

Расчетная кривая нагрева платформы от времени построена для двух вариантов: 1) нагрев происходит только за счет лучистого теплообмена и при моделировании не учитывается передача тепла между газом и платформой; 2) при нагреве принимаются во внимание оба механизма теплопередачи. На рис. 32 эти кривые отмечены индексами 1 и 2 соответственно. Видно, что кривая 1 на всем расчетном интервале времени дает сильно заниженные значения температур по сравнению с экспериментальной кривой (незакрашенные треугольники). Это говорит о том, что при моделировании нагрева платформы недостаточно учитывать

только лучистый теплообмен. Кривая 2 более реалистично отражает нагрев платформы. Она совпадает с экспериментальной кривой в момент времени г=1.6с при температуре 2250А-. При /< 1,6с расчетная зависимость проходит ниже экспериментальной, разница в значениях при ;=1.2с достигает 380АГ. При Г> 1,6с численные расчеты хорошо согласуются с опытными данными, максимальное расхождение не превышает 90К.

Штрихпунктирной линией на рис. 32 представлена температура в газе в пространстве между трубкой и платформой в сечении г=0. Как и следовало ожидать, эта зависимость наиболее близка к кривой нагрева трубки. Температура азота от времени на оси цилиндра в этом же сечении на рисунке представлена пунктирной линией. Характер изменения температуры в этих точках полностью согласуется с общей картиной нагрева трубки и платформы.

Рис. 32.

На рис. 33 приведены данные о температуре азота в трех точках С, Е и I, расположенных в плоскости г=0 на оси 0у на высотах >=—1.9, 0 и \.1мм соответственно. В точке С экспериментальные значения представлены закрашенными кружочками, а расчетная кривая показана сплошной линией. До времени г=0.9с выполняется количественное согласие расчетных и опытных данных. В интервале 1</<2.7с численные расчеты лежат выше экспериментальных значений, максимальное различие в 14% наблюдается при /=1.7с. В точке Е замеры представлены ромбиками, а расчетные кривые - штриховой линией. Для точки I, расположенной вблизи кругового выреза в трубке, закрашенные треугольники и штрих-пунктирная линия есть данные эксперимента и расчетов соответственно.

Для этих двух точек газа в интервале 1<г<2с расчетные значения температуры выше экспериментальных на 10-14%.

4

!ипе в

Рис. 33.

На основе сравнения данных эксперимента с результатами численных расчетов можно сделать следующие выводы: 1) зависимости температуры от времени в замеряемых точках для численного и натурного экспериментов качественно полностью согласуются; 2) результаты численных расчетов в интервале времени, соответствующем интенсивному тештоподводу, в среднем на 340Л- выше экспериментальных данных. С течением времени (/>2с) различие в данных уменьшается до 90К, что составляет менее 4% и можно говорить о количественном согласии сравниваемых величин.

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

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

N N N

р = - ^ди,, >у,=и,-и, ( /=1, 2, ..., /V). (6)

1-1 1=1 /-1

Для скорости диффузии используется гипотеза Фика =-£) У1п(д /р). (7)

Уравнение неразрывности для всей смеси в целом берется в виде др!д1 + 4-{ри) = р'. (8)

Здесь р' -скорость образования газовой компоненты в выделенном объеме в результате испарения с поверхности жидких частиц. Изменение массы ¿-й компоненты в объеме происходит либо за счет рождения ее в результате химической реакции, либо за счет внедрения извне.

С учетом (7) уравнение неразрывности дня /-и компонент берегся в виде др., / л + V • (р, и) = V ■ [р й Цр, / р)] + р\ + р;8ч, (9)

где р- - характеризует скорость рождения / й компоненты за счет химических реакций. Дельта функция в (9) позволяет учесть изменение плотности только той компоненты, для которой приращение отлично от нуля.

В векторном виде уравнение изменения импульса для смеси по форме совпадает с соответствующим уравнением для простого газа д(р и) / Л + V • (р ии) = р ё - Ур + V • В + ^. (10)

В диффузионном приближении для гомогенных смесей в (10) пренебрегают слагаемыми, содержащими величины второго порядка малости относительно скоростей \у/. Тензор вязких напряжений записывается в форме Навье-Стокса

а = //(Уи) + (Уи)Г

(И)

где ¡1 коэффициент вязкости, верхний индекс Т - знак транспонирования.

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

д{р I) I д г + V • (р I и) = -рУ ■ и + <5 : Уи - V ■ Л + £>. + (12)

где 3 - вектор теплового потока, <2, - тепло, связанное с взаимодействием жидких капель и газа, ()с — скорость выделения тепла в химических реакциях. Тепловой поток в случае смеси газов определяется уравнением

1 = р),

/ р). ■ (13)

где Л,- - парциальная удельная энталышя компоненты.

Для модели идеальной газовой смеси каждая компонента занимает весь объем V и является совершенным газом. Термодинамические параметры смеси определяются через суммы

n n n

р = «т£(р, / Л/,.),/(г) = £(р, / р) 1,{Т)МТ) = / Р) си(Г), 1=1 ¡=1 ¿=1

где // (7) - удельная внутренняя энергия г-й компоненты, сУ](7), су(Т) - удельные теплоемкости при постоянном объеме для отдельной компоненты и смеси в целом, А/^-молярная масса. Из определения энтальпии следует Л1(Г) = /((Г)+ЛГ/М1.. (15)

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

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

Вводится в рассмотрение специальная функция распределения /. По физическому смыслу она является вероятным числом капель, находящихся в единице объема в точке х, у и г к моменту времени Г со скоростями в интервале v, v+с1х, радиусами г, г+с1г, температурой Т, Т+(1Т, а также параметрами отклонения от сферичности капли 7/, г/+сЬ] и скоростью Г], т]+с1т].

Уравнение для функции распределения / аналогично уравнению Больцмана в кинетической теории газов. Предполатется, что фуппа капель с идентичными свойствами объединены в частицы. Капли моделируются упру! ими сферами радиуса г, с температурой Т([ и скоростью v. Два параметра ц \\ Т] учитывают степень отклонения капли от сферической формы в зависимости от скоростного напора газового потока.

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

Т1+ + ^ М ■-■♦ 4 М + + •>> + >>■

В (16) введены обозначения /' = V. К =г. Функции /г и /ь строятся аналогично интегралам столкновений в молекулярной теории (с учетом отсутствия сил притяжения между каплями).

Для каждой капли с индексом / решается дифференциальное уравнение взаимодействия с газом

ТП:

т,

= ЩЪ-~±Чр + П1[ и-у,.). (17)

А " р

р

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

£>р = 6л^г + ш2р Сд|и- \к\!2. (18)

Первое слагаемое в (18) представляет собой сопротивление Стокса, и оно учитывается при малых значениях числа Яе, второе слагаемое справедливо при больших значениях числа Не. Через С([ обозначен коэффициенг аэродинамическою сопротивлеиия, который задается эмпирически для твердой сферы. Выражение для него в зависимости от границ изменения числа Рейнолщса имеет вцц

[24(1 + Яе2/3/ б) / Ле для Яе < 1000, СЙ=М > (19)

[0.424 для Ле > 1000.

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

Скорость испарения жидких капель находится в результате интегрирования уравнения стефановского потока по толщине диффузионного слоя. Изменение радиуса капли от времени в зависимости от давления испарения находится из уравнения

¿г _ (/*>)„■>% у,' - к, 7 =_____>

Л 2грр 1 - У,* ' '' Ц + щ(р / Ру(т,) - 1))'

ще /ид-средний молекулярный вес газовой смеси, т [-молекулярный вес паров октанового топлива, Ру(7^>-равновесное давление испарения капель.

Для определения скорости изменения температуры капли Т^ рассматривается изменение баланса энергии

Рр\* ~ г'ЯИ^Т,) = 4пггйл. (21)

V' = т-^___, , (20)

Первое слагаемое в левой части (21) дает приток теплоты к капле за счет изменения температуры, второе слагаемое учитывает изменение теплоты за счет испарения. Через ¡¡¿(Trf) обозначена величина скрытой теплоты парообразования. В правой части (21) стоит теплота, передаваемая капле ог газа. Выражение для нее берется в виде

2 г

Здесь Kah. - коэффициент теплопроводности воздушной смеси

Kair{f} = К,ТУ2/(т + К2у i = {Т+ 2Ти)/Ъ. (23)

а АГ,

и Kj — заданные константы, Niij — число Нуссельта, определяющее отношение характерного размера частицы к толщине температурного пограничного слоя с учетом испарения.

Скрытая теплота парообразования находится из допущения, что плотность капель жидкости постоянна. В этом случае внутренняя энергия жидкой фазы зависит только от температуры. Энергия для преобразования единицы массы жидкости в пар при постоянном давлении определяется через разность энтальпии пара и жидкости л,. =lh{Td)-h(Td^(T)). (24)

Для газа имеем = /,(■/;) + а для жидкости энтальпия равна

h, =l,{Tj) + p/p„.

Подставим эти представления в (24) и получим формулу для /г^

Л, =/.fa) +/?g V"' -(l, + Pv/pp). (25)

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

I

д^И^Г^Д./А), (26)

7 = 1

причем А',- - число молей 1-й компоненты в единице объема смеси, ¿уЛ',-- изменение числа молей г'-й компоненты в /-ой реакции.

Если ввести в рассмотрение скорость химической реакции а>1, то

¿ДЫг =(ь0-ац)а>г (27)

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

N N

ч=1 '"<)"" - ' ' (28)

/=1 ы

где Кщ -коэффициенты прямой и обратной реакций. Выражения для них берутся в обобщенной форме Аррениуса

К(1 = А{} Г*' ехр(- £г; / 77?), К^ = А„Т'« ехр(- Е^ / ТЯ). (29)

Здесь Е^, Е}у - энергии активации прямой и обратной реакций. Предэксионенциальные множители А^у Л ¿у и показатели степени ^¿у, как и энергия активации, задаются для каждой реакции.

Для реакций, идущих с поглощеггием тепла, константа скорости обратной реакции (29) часто гге завиеггт от экспонеггциальньгх множителей и является медлегпю меняющейся функцией температуры Кц-ТМ. Перепишем соотношение (28), выделив слагаемое, не зависящее от концентрации

(30)

^ = 'т)Ь" К'{Т) ~ П<*; ,П')(Ь"

1=1 I. ¡-1

В (30) величгггга К0^Т)-Кд/К^ - константа равновесия у'-ой реакции.

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

N

= (31)

ы

Константа равновесия зависит только от температурьг и представляется в виде разложения

К°(Т) = ехр(д, 1п ТА + В] / ТА + С1 + И ¡Г, + Е^). (32)

В (32) введено обозначение Тд -Г/1000. Постоянные в (32) задаются экспериментально. Уравнение (31) является одной из форм записи законов действующих масс, определяющих равновесный состав смеси.

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

При описании химических реакций исследователи сталкиваются с системой обыкновенных дифференциальных уравнений (26Н30), которую относят к "жесткой" системе. Основная особенность этих уравнений состоит в том, что вектор действующих am в правой части (26) существенно нелинеен, и коэффициенты перед слагаемыми различаются между собой на несколько порядков величины. Если решать систему обычными численными методами, например широко распространенным методом Рун-I е-Купа, то, как правило, получаются неверные результаты. В монографии 132] подробно исследуются особенности такого вида уравнении и обсуждаются методы их решения. Основная идея метода заключается в использовании неявной центрированной по времени разностной схемы для решения системы. В диссертации приводятся примеры, демонстрирующие трудности при решении уравнений химических реакций.

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

Рассматривается цикл сжатия воздуха в цилиндре двигателя внутреннего сгорания. Камера сгорания выполнена в вице чашеобразной выемки в днище иоршш, и расположена несимметрично относительно оси цилиндра. Степки камеры ir щиппщра теплоизолированы, перераспределение тепла происходит в замкнутой области. Все термодинамические параметры задачи, математическая модель пщрокарбошюго горения, необходимые константы взяты в соответствии с тестовой задачей2. Существеннььм отличием рассмотренного ниже случая

? Amsden A.A., Ramshaw J.D., O'Rourke P.J., and Dukowicz J.K. KIVA: A computer programn for two- and three-dimensional fluid flows with chemical reactions and fuel sprays //Los Alamos National Laboratory report LA-10245-MS, 1985

,v

(33)

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

При достижении определенного угла кривошипа в цилиндр двигателя впрыскивается жидкий октан. Впрыск производится через форсунку, наклоненную иод углом фг-Иу=15° к оси цилиндра. Начинается он при угле кривошипа в=-52° и заканчивается при 9=-39.3° до ВМТ. За это время в камеру сгорания в воде полого конуса вводится 0.0116г топлива, составляющих Л^=2000 дискретных частиц. Начальные значения: плотности -рр=0.7436г/см3, скорости - 4000см/с, среднего радиуса Саугера -5-10~4см, температуры - Т^-350К, поверхностного натяжения капельки при Т=350К - ор= 16.23г/с2. Критическая температура капель - 569.36А'.

В результате взаимодействия с потоком газа частицы топлива смешиваются, испаряются и при достижении некоторого значения температуры начинают вступать в реакции. Зажигание смеси газов наров топлива производится при 6=-21° в двух центральных ячейках.

С точки зрения химии рассматриваемый процесс включает 12 компонентов. Из десяти реакций выделяются четыре кинетические и шесть равновесных.

Проекции координат жидких частиц на плоскость у=0 к моменту 6>-37.41°, близкому к завершеншо впрыска, представлены на рис. 34. Из впрыснутых в газ 2000 жидких частиц за счет взаимодействия между собой и с движущейся смесью число их сократилось до 1110. Для этого же угла на рис. 35 построено поле скоростей. Максимальная скорость I аза в сечении и,)ШХ= 1225.4см/с, и>?МЛ-=1522.4см/с достигается в ячи'ках, примыкающих к форсунке. Такую скорость приобретает газ за счет дополнительного импульса со стороны жидких капель. С течением времени начинает преобладать движение "затекания". Максимальные значения скорости достигаются в сдавленной части объема между поршнем и крышкой цилиндра.

Распределение температуры в сечении у=0 при 6Ъ-25.37° дано па рис. 36. Изолинии с индексом Ь соответствует температура 7",Ш';;=558.УА', с индексом Н - Ттах=6%2.\К. Из рисунка видно, что наиболее холодная часть в рассматриваемый момент времени расположена в центре камеры сгорания. Физически это явление объясняется введением в область холодного топлива и интенсивным испарением капель. Со временем холодное ядро газа прогревается и дробится. Это видно из рис. 37, на котором для 6^-12.22° построены изотермы в том же сечении. Предельные значения температуры в области равны Тт[п-121К и Ттах=1%5К соответственно.

> I.

Рис. 34.

Рис. 35.

Рис. 36.

Рис. 37.

Распределение массовой фракции октана в сечении у=0 при 6=-2.22° показано на рис. 38. Наибольшего значения она достигает в центре камеры сгорания. Со временем область максимальных значений смещается, изолинии "разрываются" на боковой поверхности камеры сгорания. На рис. 39 представлено распределение массовой фракции октана в сечении у=0 для меньшего угла кривошипа 0= -11.84° Смещение вправо "центра" контурных ли-шш связано с углом наклона форсунки к оси и несимметричным положением камеры сгорания в цилиндре. С точки зретптя процессов горения полученная картина распределения массовой фракции топлива является вполне

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

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

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

Рис. 38.

Рис. 39.

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

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

Основные результаты и выводы.

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

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

-алгоритмы построенггя гг оптимизации криволинейных сеток в пространственных областях;

-возможность компоновки сложных областей с разветвлениями из нескольких простых подобластей;

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

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

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

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

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

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

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

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

6. Рассмотрены основные соотношения кинетики гомогенных химических реакций. На модельных примерах исследованы особенности решения этих уравнений.

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

Результаты диссертации изложены н следующих публикациях:

1. Сахабушшов Ж.М. Нелинейные задачи аэропщроупругости п ла1ранжсвых координатах // Труды семинара по теории оболочек. / Казанск. физ.-техн. ин-т. 1971, N 2. С. 165-187.

2. Сахабугдинов Ж.М. Уравнения вязкой жидкости в лагранжевых переменных как. частый случай соотношений нелинейной теории упругости // Труды семинара по теории оболочек. / Казанск. физ.-техн. ин-т. 1975, N 6. С. 286-295.

3. Сахабутдннов Ж.М. Приложение вариационного принципа к задаче взаимодействия упругого тела со сжимаемой жидкостью. // Труды семинара но теории оболочек./Казанск. физ.-техн. ин-т. 1976, N 7. С. 52-62.

4. Сахабутдннов Ж.М. Динамика газа в трубе с эластичным участком // В сб.: Нелинейные проблемы аэрогидроупругости. - Казань, 1979, вып. XI.- С. 82-97.

5. Гильманов А.Н., Сахабутдннов Ж.М. Динамика упругих мембран из несжимаемого материала. Сб. "Нелинейная теория оболочек и пластин". - Тез. докл. Казань, 1980.

6. Гильманов А.Н., Сахабугдинов Ж.М. Численное решение задачи динамики мягкой оболочки в потоке пш// Препргагг, Новосибирск: ИТ11М СО АН СССР, 1980, N 47. С. 4-6.

7. Гильманон А.Н., Сахабутдннов Ж.М. Произвольный лагранжево-эйлеров метод в нелинейных задачах взаимодействия упругого тела потоком газа // Взаимодействие оболочек с жидкостью. Труды семинара/ Казанск. физ.-тех. ин-т. 1981, N 14. С. 127-145.

8. Гильманов А.Н., Нльгамов М.А., Сахабутдннов Ж.М. Взаимодействие мягких оболочек со сжимаемой жидкостью. //Пятый Всесоюзный съезд по теоретической и прикладной механике. Аннотация докладов. Алма-Ата: Наука, 1981, С. 111.

9. Гильманов А.Н., Сахабутдинов Ж.М. Задачи динамики упругих мембран из несжимаемого материала. // Взаимодействие оболочек с жидкостью. Труды семинара/ Казанск. физ.-тех. ин-т. 1981, N 14.

10. Гильманов A.II., Сахабутдинов Ж.М. Произвольный лшранжево-эйлеров метод в нелинейных задачах азропщроунругостн.-Сб. ЧММСС. - Новосибирск, т. 12, N 6. 1981.

11. Гильмапои А.Н., Сахабугдпнов Ж.М. Нестационарное взаимодействие мягкой оболочки с потоком газа. - Сб. "Механика сплошных сред" .- Тез. докл. -Наб. Челны, 1982 г. •

12. Сахабутдинов Ж.М. Прямое численное моделирование в газовой динамике. // В сб.: Гидроупругость оболочек. - Казань, 1983, вып. XVI,- С. 89-106.

13. Сахабутдинов Ж.М. Анализ разностных схем для уравнений движения материальной точки. Дсп. ВИНИТИ N 3230, 1984, 25 с.

14. Сахабугдпнов Ж.М. Прямое численное моделирование в механике сплошной среды. // Числен, методы механ. сплошной среды /-Новосибирск: Изд. ВЦ и Ин-та теорет. и приют, механ. СО АН СССР. 1984, Т. 15, N 1. С. 151-152.

15. Кулачкова H.A., Сахабугдпнов Ж.М. Построение расчетных сеток дня областей сложной конфигурации. // Числен, методы механ. сплошной среда / Новосибирск: Изд. ВЦ и Ин-та теорег. и ирикл. механ. СО АН СССР. 1985, Т. 16, N 3. С. 68-76.

16. Петров Г.А., Сахабутдинов Ж.М. Произвольный лафаижево-энлеров численный метод интегрирования полной системы уравнений газовой динамики в пространственном случае // Тезнсы докладов Всесоюзной школы-семинара "Магматическое моделирование в науке и технике". - Пермь, 1986.- С.235-236.

17. Петров Г.А., Сахабутдинов Ж.М. Построение трехмерных расчетных сеток для областей произвольной конфигурации. // Тезисы докладов конференции молодых ученых. КФТИ КФ АН СССР, Казань. 1986.- С.71-72.

18. Сахабутдинов Ж.М., Петров Г. А. Применение произвольного лагранжево-эйлерова численного методы к трехмерным задачам газовой динамики / Казан, фи з. -техн. ин-т. КФ АН СССР.-1986, 69с. - Дсп. в ВИНИТИ 23.07.86, 3692-В86.

19. Зарипов Р.Р., Петров Г.А., Сахабутдинов Ж.М. Графическая обработка результатов расчета трехмерных задач механики. // Тезисы докладов И Республиканской научно-технической конференции. Секция механики жидкости, газа и плазмы. - Наб. Челны, 1987.- С.27.

20. Петров Г.А., Сахабутдинов Ж.М. Нестационарное взаимодействие сильной ударной волны с пространственной конструкцией // Сб.: Динамика неоднород. сред и взаимодействие волн с элементами конструкции.-Новосибирск: ИГД СО АН СССР, 1987. - С. 137-138.

21. Сахабутдинов Ж.М., Негров Г.А., Мангурова С. В. Геометрический подход к построению и оптимизации крива»шейных ссток в пространственных областях. // Казань, физ - техн. ин-т. КФАН СССР.-1987.-90с.-Деп. в ВИНИТИ 31.07.87, 6976-В87.

22. Сахабутдинов Ж.М., Петров Г.А., Машурова С. В. Геометрический метод построения сеток для пространственных областей произвольной формы // Взаимодействие оболочек со средой. Труды семинара.-Казань, КФТИ КФ АН СССР, 1987, вып. ХХ.-С. 209-222.

23. Петров Г.А., Сахабутдинов Ж.М. Применение произвольного лафанжево-эйлерова метода к расчету пространственных течений газа в каналах сложной формы // Тезисы докладов VII Всесоюзного семинара "Теор. основы и конструирование численных алгоритмов решения задач математической физики".-Кемерово, 1988.-С. 92.

24. Сахабутдинов Ж.М., Петров Г.А., Майгурова С. В. Построение и оптимизация трехмерных криволииейных сеток // Вопросы атомной науки и техники, серия: "Мат. модели физ. процессов", 1988, N 3. С. 76-85.

25. Сахабутдинов Ж.М. Алгоритм перестройки произвольной трехмерной сетки // Тезисы докладов VII Всесоюзного семинара "Теор. основы и конструирование числа пил х алгоритмов решения задач математической физики". - Кемерово, 1988.-С. 109.

26. Бекбулагов И.Г., Сахабутдинов Ж.М. Течение газа в межлопаточном канале рабочего колеса центробежного компрессора // Моделирование нелинейных процессов в механике и теплотехнике. Труды семинара / Казанск. физ.-техн. ш-т. 1989, N 24. С. 6-16.

27. Бекбулатов И.Г., Сахабутдинов Ж.М. Динамика газа в межлопаточном канале рабочего колеса компрессора // Аннотации докладов VIII Всесоюзной научно-технической конференции "Создание компрессорных машин и установок.", Сумы, 1989.

28. Зарипов Р. Р., Сахабутдинов Ж. М. Моделирование трехмерных нестационарных течений газа в цилиндре двигателя внутреннего сгорания // Моделирование нелинейных процессов в механике и теплотехнике. Труды семинара / Казанск. физ.-техн. ин-т. 1989, N 24. С. 17-26.

29. Kaminsky V. N.. Olisov P. A., Warshavsky Y. M., Sakhabutdinov Zh. M., Petrov G. A. Simulation of Mixing and Combustion Processes in Internal Combustion Engines // Abstracts of the III International Seminar on Flame Structure. Alma-Ata, USSR. September 18-22. 1989, P. 57.

30. Sakhabutdinov Zh.M., Baranova E.N., Belaeva U.A. Interaction of Fluid Cloud of Droplet with Multicomponent Gas Mixture // Abstracts of the International Conference: Fundamental Research in Aerospace Science. Zhukovsky, Russia. 22-24 September 1994, Section 4. P. 63-64.

31. Гильманов A.H., Сахабутдинов Ж.М. Численное решение задач аэрогидроупругости // Обзоры исследований по механике сплошном среды. Труды семинара / Институт механики и машиностроения КНЦ РАН. 1995г. С. 131-145.

32. Сахабутдинов Ж.М. Анализ дискретных моделей движения точки. ИММ КНЦ РАН, Казань, Полиграфический комбинат им. К. Якуба, 1995. С. 196.

33. Gilmutdinov A.Kh., Araslanov Sh. F., Sakhabutdinov Zh.M. Radziuk В, M. Sperling, B. Welz. Computer modeling of a transversely heated tube electrotermal atomizer. Fourth Symposium on Atomic Spectrometry, 24-30 November, Buenos-Aires, Argentina, 1996.

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