Содержание к диссертации
Введение
1. Численное моделирование тепломассопереноса в композиционных материалах 14
1.1. Идентификация закона разложения связующих композиционных материалов 15
1.2. Идентификация нелинейного закона фильтрации продуктов разложения связующих через пористый остаток 18
1.3. Физико-математическая модель теплового состояния композиционных материалов при высокоинтенсивном тепловом нагреве 23
1.4. Методология численного решения 31
1.4.1. Базовая конечно-разностная схема для незатронутой области 32
1.4.2. Момент появления подвижной гранж\ы хн (/) npuTwl >ТН и (или) границы x*K\t) при Tw] > Т*к* (Т** > Т**) 37
1.4.3. Определение скоростей движения границ xH(t) и хк (/) 43
1.4.4. Решение задачи в области разложения связующего 44
1.4.5. Определение температурного поля в области фильтрации пиролизных газов 46
1.5. Алгоритм и программный комплекс 49
2. Моделирование тепломассопереноса в анизотропных композиционных материалах 55
2.1. Математическая модель двумерного нестационарного тепломассопереноса в анизотропных композиционных материалах 55
2.2. Метод переменных направлений с экстраполяцией численного решения многомерных задач анизотропного тепломассопереноса 67
2.2.1. Базовая конечно-разностная схема метода МПНЭ с учетом фильтрации 68
2.2.2. Модификация метода МПНЭ с использованием интегроинтерполяционного метода Самарского А. А. в многомерных задачах анизотропного тепломассопереноса 71
2.3. Методология численного решения многомерных задач тепломассопереноса в анизотропных композиционных материалах 85
2.4. Программный комплекс и тестовые результаты 90
3. Аналитическое и численное исследование тепломассопереноса в композиционных материалах 94
3.1. Аналитическое решение задачи типа Стефана с двумя нестационарно подвижными границами 94
3.2. Численный анализ тепломассопереноса в композиционных материалах 106
3.2.1. Выбор шага численного интегрирования по времени 106
3.2.2. Результаты численных экспериментов' 107
3.3. Анализ результатов численных экспериментов по определению теплового состояния анизотропных материалов 113
Заключение 126
Список использованной литературы
- Идентификация нелинейного закона фильтрации продуктов разложения связующих через пористый остаток
- Базовая конечно-разностная схема для незатронутой области
- Метод переменных направлений с экстраполяцией численного решения многомерных задач анизотропного тепломассопереноса
- Численный анализ тепломассопереноса в композиционных материалах
Введение к работе
Композиционные материалы широко применяются в авиации и космонавтике в качестве конструкционных и теплозащитных материалов (ТЗМ) благодаря своим уникальным свойствам, вытекающим из технологии их изготовления. Матрица из тонковолокнистого наполнителя (стекло-, асбо-, угле- и т.д. волокна) пропитывается связующими (легко разлагающимися при умеренных температурах смолами). В результате получаются различные пластики: стеклопластики, углепластики, асбопластики и т.п. При использовании таких материалов в качестве теплозащитных при гиперзвуковых скоростях полета (числа Маха М>5-6) летательных аппаратов (ЛА) тепловые потоки от высокотемпературных пограничных слоев поглощаются за счет следующих факторов: при температурах до ~ 600К теплозащитный материал работает за счет своих теплофизических характеристик, а именно, объемной теплоемкости (произведение плотности на теплоемкость) и отвода тепла вглубь материала за счет теплопроводности; при достижении температуры ~ 600К начинается разложение (пиролиз) связующего под действием эндотермических химических реакций с образованием газообразных продуктов разложения и пористого остатка, состоящего из волокон наполнителя и коксового остатка от разложения связующего; процесс разложения связующего заканчивается при температурах ~ 900 - 1100К, в результате чего область разложения, ограниченная координатами начала и окончания разложения, является очень тонкой и разделяет незатронутый разложением материал и пористый остаток, в котором разложение связующего отсутствует; через пористый остаток под действием градиента давления между областью разложения, где давление из-за низкой скорости фильтрации пиролизных газов считается давлением торможения, и наружной границей, где давление равно давлению окружающей среды, пиролизные газы фильтруются, поглощая тепло за счет конвекции; при этом для реализации процесса фильтрации давление
5 торможения должно превышать давление окружающей среды на величину гидравлического сопротивления пористого остатка; пиролизные газы вдуваются в высокотемпературный пограничный слой, оттесняя его от наружной границы и уменьшая его температуру, что влечет за собой уменьшение конвективных тепловых потоков к наружной границе; при достижении пористым остатком температуры уноса массы, на наружной границе уносится масса за счет физико-химических превращений (оплавления, испарения, возгонки , гомо- и гетерогенных химических реакций), поглощающих значительное количество тепловой энергии с суммарным тепловым эффектом, называемым теплотой уноса массы.
Как видно из этого перечня композиционные материлы, используемые в качестве теплозащитных для гиперзвуковых ЛА, поглощают значительное количество тепловой энергии при аэродинамическом нагреве за счет различных физико-химических превращений, происходящих в композиционных материалах.
Отсюда видно, что математическое моделирование процессов тепло- и массопереноса в композиционных материалах представляет собой сложную комплексную проблему, которая в полной мере не решена до сих пор и которой уделяется значительное внимание не только теплофизиками, аэродинамиками, материаловедами, но и математиками.
При математическом моделировании . тепломассопереноса в композиционных материалах необходимо в комплексе моделировать следующие процессы: разложение связующих с образованием пиролизных газов и пористого остатка; фильтрацию пиролизных газов через пористый остаток и учет этой фильтрации в теплопереносе; тепломассоперенос в области разложения связующего; унос массы и его влияние на нестационарное температурное поле; вдув пиролизных газов в газодинамический пограничный слой. При этом необходимо учитывать различные явления, приводящие к существенной нелинейности и нестационарности математических моделей при высоких температурах. К этим явлениям можно отнести учет излучения, зависимость теплофизических характеристик (ТФХ) материалов от температуры, разрывы ТФХ, анизотропию и многомерность распространения тепла и др.
Таким образом, математическое моделирование тепломассопереноса в композиционных материалах при высокоинтенсивном нагреве является актуальной проблемой.
В диссертационной работе исследуются вопросы численного моделирования процессов совместной нелинейной теплопроводности, нелинейной фильтрации, разложения связующего в нестационарно подвижной области, в том числе в условиях уноса массы, а также с учетом анизотропии свойств переноса тепла и массы в композиционных материалах при сложном теплообмене (конвективно-кондуктивном и лучистом видах теплообмена).
В виду комплексности проблемы тепломассопереноса в композиционных материалах при обзоре литературы следует проанализировать публикации по композиционным материалам, по совместной теплопроводности и фильтрации, по кинетике разложения связующих и законах ее описывающих, по уносу массы, по задачам типа Стефана (задач теплопереноса при фазовых превращениях) в том числе при наличии многих (а не одной) подвижных границ фазовых превращений.
Насчитываются ~ десятки работ по совместной теплопроводности и фильтрации. Например, работы отечественных ученых таких как Лыков А.В. со своей школой [3], [65-68], Полежаев Ю.В. [5,79-81,27], Боровой В.Я. [16], Зарубин B.C. [38, 39], Формалев В.Ф. [96-98, 100, 104, 105], Баренблатт Г.И. [10] и др., в которых фильтрация моделировалась на основе линейного закона Дарси.
По теории фильтрации через пористые среды можно назвать работы Бана А. и др. [9] для сжимаемых пористых сред , Бернадинера М.Г. [14] в задачах нелинейной фильтрации, Глущенко А.А. для многомерных задач фильтрации [29], монографии Коллинза Р. [50] и Шейдеггера А.Э. [112]. В этих работах исследовалась изотропная фильтрация. Публикации по анизотропной фильтрации с учетом неизотермичности автору не известны.
Значительный вклад в развитие теории теплопроводности внесли советские и российские теплофизики Лыков А.В. [67, 68], Леонтьев А.И. [64], Зарубин B.C. [38, 39], Дульнев Г.Н. [35], Карташов Э.М. [43,44], Формалев В.Ф. [97, 100, 101] и многие другие, а также зарубежные теплофизики Карслоу Г. и Егер Д. [41], Берман Р. [12], Гребер Г, Эрк С, Григулль У. [28], Шнейдер П. [114].
7 Различные методы решения задач, моделирующих теорию теплопроводности в условиях уноса массы, когда возникают подвижные границы, разрабатывались в работах Самарского А.А. [90, 91], Будака Б.М. и Гольдмана Н.Л. [15], Карслоу Г. и Егера Д. [41],Лыкова А.В. [65], Адамса М.К. [2], Полежаева Ю.В. [79], Скала СМ. и Гильберта Л.М. [93], Розенсвеига Р.Е. и Бичера Н. [85], Формалева В.Ф. [97, 100, 104], а также в работах [20, 73, ПО, 121, 128-130]. В этих работах рассматривались в основном одномерные по пространственным переменным задачи.
Большой вклад в разработку методологии моделирования теплового состояния композиционных материалов при высоких температурах с учетом их механического состояния внес Димитриенко Ю.И. [30-34], а также Головин Н.Н. и Кувыркин Г.Н. [25, 26]. Исследованию тепломассопереноса в композиционных материалах при высоких температурах и давлениях посвящены работы Кузнецова Г.В. [62] и Кузнецова Г.В. и Рудзинского В.П. [63]. Развернутый обзор по - материалам и теплозащитным покрытиям в экстремальных условиях сделан в работе Резника СВ. [73].
Следует отметить, что практически все композиционные материалы в той или иной степени являются анизотропными. При моделировании теплового состояния таких материалов возникают значительные трудности, связанные, прежде всего, с тем, что теплопроводность и проницаемость в них. являются не скалярными, а тензорными характеристиками, следствием чего является многомерность задач и усложнение краевых условий при сложном теплообмене. Этими сложностями обусловлено незначительное число работ по анизотропной теплопроводности и анизотропной фильтрации. Здесь, прежде всего, следует отметить работы Пэдовена Д. [82], Пуня К.С., Цзоу Р.С. и Чжана Ю.П. [83, 107, 108], Формалева В.Ф. [100], в которых получены аналитические решения некоторых простейших задач анизотропной многомерной теплопроводности операционными методами.
В работе Кима Л.В. и Микова В.Л. [49] рассмотрен полуаналитический метод решения линейных задач анизотропной теплопроводности в полупространстве, основанный на численном обращении преобразования
8 Лапласа, а в работе Чепрасова А.И. [106] численно моделируется трехмерная анизотропная теплопроводность в аномальной подобласти, принадлежащей бесконечной расчетной области с использованием метода дробных шагов Яненко Н.Н. [117]. Таким образом, работ по анизотропной теплопроводности очень мало, а работ по совместной анизотропной теплопроводности и анизотропной фильтрации в многомерных телах автору не известны.
Из проведенного обзора по тепломассопереносу в композиционных материалах ясно что физико-математическая модель разрабатывалась для каждого композита с учетом химического состава его связующего и наполнителя, химических реакций и физико-химических явлений, специфических для каждого материала. Ясно, что для математического моделирования тепломассопереноса в композитах такой подход не приемлем. Поэтому в диссертации предложен подход, являющийся универсальным для любого композиционного материала и позволяющий обойти химическую кинетику разложения связующего. Существо подхода заключается в том, что для большинства композиционных материалов экспериментально определены температуры и плотности начала и окончания разложения связующих, что позволяет на основе этих данных идентифицировать закон разложения связующих и включить его в общую физико-математическую модель тепломассопереноса.
На основе изложенного формулируются следующие цели диссертационной работы:
Идентификация закона разложения связующего (типа закона Аррениуса), позволяющего обойти химическую кинетику и пригодного для любого композиционного материала.
Идентификация нелинейного закона фильтрации пиролизных газов через пористый остаток, позволяющий, в отличие от линейного закона Дарси, учитывать не только вязкостные, но и инерционные силы.
Разработка на основе этих законов физико-математической модели тепломассопереноса в композиционных материалах в виде задачи типа Стефана с учетом фильтрации с тремя нестационарно подвижными границами фазовых
9 превращений - границ начала и окончания разложения связующего и границы уноса массы.
Применение полученных законов для моделирования многомерного нестационарного тепломассопереноса в анизотропных композиционных материалах при высокоинтенсивном нагреве.
Модификация и _ использование численного метода переменных направлений с экспраполяцией с позиций увеличения точности при аппроксимации краевых условий при сложном теплообмене.
Разработка программных комплексов, исследование температурных полей и полей давления фильтрации, получения аналитического решения задачи типа Стефана с двумя подвижными границами и тестирование с его помощью численных результатов.
Для численной реализации комплексов задач, перечисленных в этих целях, необходимо выбрать численный метод, удовлетворяющий следующим требованиям: экономичности, абсолютной устойчивости, простоте алгоритмизации. К настоящему времени существует огромное число работ, посвященных численному решению параболических задач. Среди экономичных численных схем необходимо отметить следующие: метод дробных шагов Яненко Н.Н. [117], метод переменных направлений (МПН) Писмена-Рэчфорда [126], центрально-симметричный метод Самарского А.А. [86], метод переменных направлений с экстраполяцией (МПНЭ) Формалева В.Ф. [99, 101], метод полного расщепления (МПР) Формалева-Тюкина [101]. Метод стабилизирующей поправки, являясь разновидностью методов расщепления, был предложен Дугласом Д. И Ганном Д. [123]. Метод предиктор-корректор для уравнений параболического типа был развит в работе Дугласа Д. и Джонса Б. [124].
Все эти методы являются экономичными, поскольку сводятся к скалярным прогонкам по координатным направлениям и абсолютно устойчивыми, если дифференциальные уравнения не содержат смешанных дифференциальных операторов. При наличии смешанных производных все перечисленные методы за исключением методов МПНЭ и МПР [101] являются условно устойчивыми, даже такой метод как метод дробных шагов Яненко Н.Н. Поэтому для программной
10 реализации многомерных задач анизотропной теплопроводности и фильтрации в качестве базового выбран метод МПНЭ, разработанный научным руководителем.
По точности и порядку аппроксимации он уступает таким методам, как МПН, метод стабилизирующей поправки, но по запасу устойчивости он не имеет себе аналогов. Кроме этого, МПНЭ обладает полной аппроксимацией на каждом дробном временном подслое.
Ниже излагается краткое содержание диссертации. Она состоит из введения с обзором литературы, трех глав с выводами, заключения и списка литературы.
В первой главе в одномерной по пространственным переменным постановке моделируется нестационарный тепломассоперенос в композиционных теплозащитных материалах в виде задачи типа Стефана с тремя нестационарно подвижными границами, определяемыми из теплового состояния тела.
В параграфе 1.1 на основе известных температур и плотностей начала и окончания разложения связующих идентифицирован закон разложения связующих (аналог закона Аррениуса), позволяющий в пространстве и времени определять плотность связующего в зоне разложения, генерацию пиролизных газов и давление торможения газов.
В параграфе 1.2 на основе анализа течения вязкого газа в капилляре, введения функции от фильтрационного числа Рейнольдса, учитывающей нелинейность фильтрации, и ее дальнейшего определения идентифицирован нелинейный закон фильтрации пиролизных газов, существенно уменьшающий скорость фильтрации при больших градиентах давления по сравнению с законом фильтрации Дарси.
В параграфе 1.3 на основе двух идентифицированных законов предложена физико-математическая модель определения теплового состояния композиционных материалов при выскоинтенсивном тепловом нагреве. Модель учитывает возникновение и нестационарное движение вглубь композиционного материала трех нестационарных подвижных границ - двух границ начала и окончания разложения связующего и границы уноса массы, квазистационарную фильтрацию, определение массы пиролизных газов и давление торможения их, а также различные виды нелинейностей.
В параграфе 1.4 изложена общая методология численного решения всей комплексной проблемы, определяются моменты начала движения подвижных границ и итерационного уточнения их положения, наложены ограничения на шаг численного решения по времени из-за наличия подвижных границ, коротко описан программный комплекс и тестовые результаты, полученные с помощью него.
В главе 2 численно моделируется многомерный тепломассоперенос в анизотропных композиционных материалах на основе тех же двух идентифицированных законов.
В параграфе 2.1. предложена математическая модель двумерного нестационарного тепломассопереноса в анизотропных композиционных материалах в условиях конвективно-кондуктивного и лучистого видов теплообмена на свободных границах. Краевые условия типа Стефана на подвижных границах представлены в виде, пригодном для применения численных методов расщепления по координатным направлениям. Модель предполагает определение двумерных нестационарных температурных полей в различных фазах анизотропного композита, полей давления фильтрации в анизотропном пористом остатке, полей компонентов вектора скорости фильтрации и скорости движения многомерной границы фазовых превращений.
В параграфе 2.2 предложена модификация метода переменных направлений с экстраполяцией (МПНЭ) численного решения многомерных задач анизотропного тепломассопереноса. Модификация основана на использовании в МПНЭ интегроинтерполяционного метода Самарского А.А., позволяющая сохранять второй порядок аппроксимации в краевых условиях, содержащх производные искомых функций. Рассмотрена конечно-разностная аппроксимация в различных нерегулярных узлах.
В параграфе 2.3 изложена методология численного решения многомерных задач тепломассопереноса в анизотропных композиционных материалах. Особое внимание уделяется моделированию задачи многомерной фильтрации пиролизных газов в анизотропном пористом остатке.
В параграфе 2.4 коротко описан программный комплекс и тестовые результаты, полученные с его помощью для нестационарного температурного поля и координат двумерной нестационарно подвижной границы, разделяющей исходную фазу и пористый остаток. —
Идентификация нелинейного закона фильтрации продуктов разложения связующих через пористый остаток
Таким образом, линейный закон фильтрации Дарси (1.5) характеризуется полным отсутствием инерционных членов в уравнении сохранения импульса (1.10). Сравнивая (1.13) с (1.5), можно записать нелинейный закон фильтрации в следующем виде: поскольку с падением давления (dp/dx 0) наблюдается рост продольного компонента скорости фильтрации (du/dx 0). Результат, аналогичный (1.15), можно получить и путем введения фильтрационного числа Рейнольдса ReL = pniL/p. М.Д. Миллионщиков [68] предложил в качестве характерного размера в ReL принять величину L = лД/П (П - пористость), а в качестве характерной скорости - среднюю скорость и = ит/П, um=G/p (G - секундный расход охладителя в капилляре). Тогда ReL=pumJk/{n3/2p). Рассмотрим закон нелинейной фильтрации в виде, когда антиградиент давления пропорционален некоторой функции фильтрационного числа Рейнольдса ReL и в соответствии с (1.15) - величинам р и и следовательно, уменьшает скорость фильтрации. С другой стороны, увеличение градиента давления ведет к возрастанию скорости фильтрации и числа Re , что уменьшает скорость фильтрации. Проанализируем, начиная с каких значений динамической вязкости и градиента давления необходимо учитывать нелинейность закона фильтрации. На рис. 1.2 представлены результаты сравнения скоростей линейной и нелинейной фильтрации в зависимости от градиента давления и динамической вязкости. Из рисунка видно, что нелинейность закона фильтрации начинает оказывать влияние на скорость фильтрации при значениях динамической вязкости Г О охладителя от // = 8x10 кг/м-с и ниже при градиенте давления равном -1.7x10
Па/м. При уменьшении вязкости до значения 1х10 5 кг/м-с величина градиента давления, при котором необходимо учитывать нелинейную фильтрацию, составляет -0.3x108 Па/м.
Для многослойной области (рис. 1.3), в которой наружный материал является композиционным, рассматривается следующая задача тепломассопереноса при сложном теплообмене на свободных границах w\, w2.
Баланс конвективно - кондуктивных и лучистых тепловых потоков на наружной границе wl, соприкасающейся с высокотемпературным газодинамическим пограничным слоем, с учетом уноса массы
Индексы: wl, w2-границы x = xQ и x x0+2_JSs соответственно; наружная подвижная граница, - подвижные границы разложения связующего; н, к-начальные и конечные значения; s -граница разрыва теплофизических характеристик; г-число слоев, номер наружного слоя; Г-газ; eff-эффективное значение; е\, е2 - высокотемпературная среда, соприкасающаяся с границей w\, и окружающая среда, соприкасающаяся с границей w2 соответственно; кс - коксовый остаток.
Математическая модель (1.24)-(1.42) является комплексной. Она состоит из следующих частных моделей. Соотношения (1.24)-(1.27), (1.41) представляют собой задачу теплопроводности в многослойных областях с возможным уносом массы на наружной границе xwi. Задача (1.24), (1.25), (1.29), (1.30), (1.31), (1.36), (1.41) является задачей теплопроводности с конвективным членом в пористой области х к (/) х x wl (?), ограниченной двумя нестационарно подвижными границами x K (t) и x w](t). Массовые и линейные скорости движения этих границ описываются соотношениями (1.36) и (1.25) соответственно. Массовая скорость фильтрации (ри)г, входящая в уравнение (1.29), определяется математической моделью (1.37)-(1.40). Распределение температуры и плотности композита в подвижной области х (t) х х к (t) разложения связующего моделируется соотношениями (1.31) (1.36). При этом координаты подвижных границ х (?) и x"(t) определяются соотношениями (1.35), (1.36).
Следует отметить, что коэффициент теплоотдачи задается с учетом вдува в окружающую газообразную среду продуктов разложения связующего и ={a"\)o fl(PrUr)w] где коэффициент вдува /?-определяется экспериментально [1J.
Решение подобных комплексных, существенно нелинейных задач тепломассопереноса в композиционных материалах при наличии трех нестационарно подвижных границ представляет собой значительные трудности.
Математическая модель (1.24-1.42) является комплексной, включающей в себя следующие частные задачи. 1) Обычного прогрева с граничными условиями различных родов (1-го, 2-го, 3-го, 4-ого родов) при T(x,t) Т . 2) Первоначального определения координат подвижных границ. 3) Теплопроводности в пористом остатке с учетом фильтрации и подвижных границ - наружной границы уноса массы и границы окончания разложения связующего. 4) Генерации пиролизных газов в области разложения связующего и определения плотности давления торможения рат,р(х"). 5) Распределения давления газов в пористом остатке, поскольку фильтрация возможна в случае превышения давления в пористом остатке над давлением /?„., окружающей среды. 6) Возможного уноса массы с температурой Ґ и эффектом Q . 7) Определения скоростей движения границ (/) и x {t)
При решении этих задач делаются следующие Допущения: 1) Наружный композиционный материал состоит из наполнителя (стекло-, угле-, и т.д. волокон) и связующего (например, фенольно-формальдегидная смола), начинающего разлагаться при температуре Т" с образованием твердого (коксового) остатка и газового компонента и заканчивающего разлагаться при температуре Т . При этом через образовавшуюся пористую структуру газовый компонент фильтруется от границы А- к наружной границе w\ под действием перепада давления.
Базовая конечно-разностная схема для незатронутой области
Зафиксируем начало_ координат по пространственной переменной х и введем фиксированную пространственно-временную сетку с послойными пространственными шагами /2,,/2 ,...,/2,. и числом шагов в каждом слое ml,m2,...,mr и шагом т по времени
Здесь г-число слоев с различными теплофизическими характеристиками (ТФХ) материалов. При этом узлы сетки попадают на границы разрыва ТФХ.
Рассмотрим вначале конечно-разностную схему решения для случая, когда наружный материал не разлагается 7] Гн , т.е. осуществляется обычный прогрев. Тогда конечно-разностная аппроксимация соотношений (1.24), (1.26) - (1.28) с использованием начального условия для температуры (1.41) принимает вид {Т ,Т, +,і = \,N,к = 0,1,2,... - сеточные функции на нижнем /К=К2" и верхнем tK+ = (к + 1)г временных слоях соответственно),
Для сохранения устойчивости по граничным условиям при высоких температурах Тт , аппроксимацию лучистого теплового потока на границе w\ будем осуществлять неявно (на верхнем временном слое), а на внутренней (низкотемпературной) границе w2 - явно (на нижнем временном слое).
В противном случае, т.е. если Х„ »к+1 - X. к+1 hr, в промежутке х є хн ,хк решается первая краевая задача для уравнения (1.31) со стоками энергии вследствие разложения связующего (и -возможно фильтрации с переменной пористостью).Решение этой задачи см. в п. 1.4.4.
После того, как подвижная граница x = x"(t) сформирована, в области л єГх0,дг (/)1 решается следующая задача теплопроводности с граничными условиями различных родов (см. рис. 1.4).
Аппроксимация в остальных узлах q-1,...,0 аналогична аппроксимации базовой модели (1.44) - (1.48). Прямой ход прогонки осуществляется в направлении от узла 0 к узлу q с определением Т + , а обратный ход - в направлении от узла q к узлу 0 с определением распределения температур 7)к+ , і = q,q-1,...,0. После этого по распределению температур в области х є х0,х и можно определить плотность теплового потока на границе хн (t) со стороны незатронутого материала 2( „ - ,)
причем, если подвижная граница хн {t) войдет в є - окрестность узла q, то плотность теплового потока (1.63) определяется по значениям температур Гк+ и 7?t. За є - окрестность узла хп можно принять окрестность в объеме, например, Ч Ч hr /100. Если этого не сделать, то возможно деление на нуль!
Далее необходимо определить плотность теплового потока на границе хн {t) со стороны прореагировавшего материала. Для этого осуществляется следующий анализ. 2) xmi-\ хи xmZ-o 5 этом случае в конечно-разностной аппроксимации краевого условия (1.24) аккумуляция энергии при температуре Т = Т будет в объеме hr/2, т.е. в выражении (1.65) в правой части вместо сомножителя hr-% будет сомножитель hr/2. Решение его также дается выражением (1.53), а тепловой поток qn - выражением (1.66).
В случае, если х0 + 2_jSs хн xmj-_{, то в области хє[х ,х1(., 1 решается" 5=1 задача с разложением связующего (см.п. 1.4.4), после чего определяется тепловой поток на границе хн (t) со стороны разлагающегося материала.
После определения тепловых потоков на границе \хн 1 со стороны .+ разлагающегося материала qH и со стороны незатронутого материала qn массовую тн , а затем и линейную хн (t) скорости движения границы хн (/) можно определить из граничного условия (1.35) ,к+1\" 1 PrQ х« {ti-чі), (L67) откуда получаем новое значение границы начала разложения связующего на п -ои итерации к+1 \л к Ы«"- ) -б8) PrQ где плотности тепловых потоков qH ,qH определяются выражениями (1.63)-(1.66), или подобными им выражениями. Аналогично с использованием краевых условий (1.36) определяется скорость движения границы 1хк I и значение координаты этой границы (С)"—(?;-?;) " (1.69) у prQ { »«к+і\ к Т PrQ (?;-«;) (1-70) 1.4.4. Решение задачи в области разложения связующего Эта область ограничена слева координатой хн (/), а справа или координатой хк (t), если граница хк (t) появилась, или границей xwl, если граница хк {t) еще не возникла.
В первом случае в этой области (если она больше hr) решается следующая задача:
Как уже отмечалось выше, в области разложения линейная скорость фильтрации иг мала и фильтрационным членом в" уравнении (1.71) можно пренебречь.
Распределение плотности материала ркг{х), х н х х" можно определить из идентифицированного в п. 1.1. закона разложения связующего (1.74).
После определения рг(хрҐ+]) в области разложения связующего конечно-разностная аппроксимация и решение первой начально-краевой задачи не вызывают затруднений.
Метод переменных направлений с экстраполяцией численного решения многомерных задач анизотропного тепломассопереноса
Все известные конечно-разностные методы с расщеплением разностных операторов по координатным направлениям являются экономичными и используют либо неявно-явную аппроксимацию дифференциальных операторов на каждом дробном временном слое, обладая при этом свойством полной аппроксимации (например ,схема Писмена -Рэчфорда [126] ), либо неявную аппроксимацию одного из дифференциальных операторов по пространственным переменным на каждом дробном шаге по времени, обладая частичной аппроксимацией (например, схема Яненко Н.Щ117]).
Однако при наличии в задаче смешанных дифференциальных операторов во всех схемах расщепления , последние аппроксимируются явно, что влечет за собой условную устойчивость этих схем. Например, схема Писмена -Рэчфорда условно устойчива уже в двумерных по пространственным переменным задачах , содержащих смешанные производные , а схема Лненко Н.Н. условно устойчива в случае недиагонального преобладания тензора теплопроводности.
Поэтому для численного решения задач, сформулиированных в математической модели (2.1)-(2.24), выбран экономичный абсолютно устойчивый метод переменных направлений с экстраполяцией (МПНЭ), разработанный научным руководителем [і01]. Экономичность в нем достигается использованием только скалярных прогонок по координатным направлениям , а устойчивость -использованием полностью неявных схем по координатным направлениям, включая неявную аппроксимацию смешанных дифференциальных операторов. При этом на каждом дробном шаге достигается полная аппроксимация.
Для достижения второго порядка при аппроксимации краевых условий , содержащих производные, метод МПНЭ в данном параграфе модифицирован путем использования интегро-интерполяционного метода Самарского А.А.
Желание получить схемы экономичные и одновременно абсолютно устойчивые с большим запасом устойчивости приводит к необходимости использовать апостериорную информацию о сеточной функции, которая получается в процессе счета. К таким схемам относится схема метода переменных направлений с экстраполяцией (МПНЭ). Эта схема использует распределение сеточной функции в левом (от расчетного) пространственном сечении, уже полученное на верхнем временном слое, т.е. неявно. В правом пространственном сечении значения сеточной функции в первом приближении и с контролируемой точностью определяются линейной экстраполяцией по времени с нижних временных слоев на верхний. Затем эти значения уточняются с помощью скалярных прогонок.
В подсхеме (2.34) значения T f, Ту+]/2, Т$2 являются искомыми, определяемыми из скалярных прогонок в направлении переменной х, значения -\%, Т 2, Т [Х!2Х уже известны на верхнем временном полуслое из прогонки вдоль координатной линии ун= (; -1)/ , а значения f$%, 7 I/2, 72 с порядком 0\т2) определяются экстраполяцией (2.36) по распределениям функции на двух предыдущих временных полуслоях. При этом все конечно В подсхеме (2.35) значения 7У, , Т! , TIJ+l являются искомыми, определяемыми из скалярных прогонок вдоль переменной у, значения 7] ,+!_,, разностные операторы по пространственным переменным, за исключением оператора по переменной JC, переводятся в правые части (хотя они и являются практически полностью неявными). ik+\ rpk+\ T i + l T j, 7] !+, известны как значения сеточной функции в левом пространственном сечении, а значения 7J J"J_,, 7 ], 7J +1 с порядком 0\т2) определяются экстраполяцией (2.37) по двум предыдущим временным полуслоям. При этом все пространственные операторы, за исключением оператора по переменной у, переводятся в правые части, хотя они практически являются полностью неявными.
Для получения второго порядка при аппроксимации краевых условий, содержащих производные , в задачах анизотропного тепломассопереноса модифицируем метод МПНЭ с использованием интегроинтерполяционного метода Самарского А.А. и примененного к задачам для обыкновенных дифференциальных уравнений в книге [89].
Рассмотрим для этого некоторые типовые граничные узлы, а именно узлы (/,0),(0,0), (0, MI,J). Проинтегрируем уравнение (2.38) по переменной х в пределах от х = _2 = xt.-— до х = х = х;. + - и по _у в пределах _у = 0 до у = — , получим
С помощью квадратурной формулы прямоугольников интеграл в левой части этого выражения равен значению подынтегральной функции в узле (/,0), h -к домноженной на площадь прямоугольника ——-. Первый интеграл в правой части проинтегрируем по переменной х, а второй- по переменной у, получим Первый определенный интеграл в правой части этого выражения со значением подинтегральной функции в нижнем пределе вычисляется по квадратурной формуле прямоугольников, а второй - по переменной х по той же формуле со значением подынтегральной функции при х = хп получим :
Численный анализ тепломассопереноса в композиционных материалах
В задачах тепломассопереноса, вследствие существенной нелинейности комплексной проблемы, используются только неявные численные схемы, вследствие чего можно ожидать выполнения устойчивости численного решения и, при условии аппроксимации, - сходимости численного решения к точному. При этом известно, что при использовании неявных схем допускается произвол в выборе сеточных характеристик вообще и шага численного интегрирования по времени в частности.
Однако это не так в задачах тепломассопереноса с фазовыми превращениями, поскольку распределение температур в существенной степени влияет на скорость фазовых превращений, а фазовые превращения, в свою очередь, значительно изменяют температурные поля. Вследствие этого в таких задачах шаг численного интегрирования по времени должен быть связан со скоростью фазовых превращений или со скоростями движения границ фазовых превращений даже при использовании полностью неявных численных схем. Это подтверждают численные эксперименты. Если шаг численного интегрирования по времени не согласован со скоростями движения границ фазовых превращений , то решение разваливается вследствие того, что при произвольном шаге по времени границы фазовых превращений могут продвинуться вглубь композита так глубоко, что температурное поле , зависящее от координат этих границ, не будет связано с температурным полем на предыдущих временных слоях.
Вследствие изложенного необходимо a priori до начала расчета оценить шаг численного интегрирования по времени, используя предварительную оценку скоростей движения одной из границ фазовых превращений. Если при решении конкретной задачи имеется подвижная граница уноса массы w\ ( TM(t) может быть выше Ґ ) с известной температурой Ґ, теплотой Q (Te1) фазовых превращений и параметрами теплообмена aw](t), Te](t), то из баланса (1.25) шаг г численного интегрирования по времени можно оценить следующим образом:
С помощью разработанного программного комплекса был проведен параметрический анализ результатов численного решения по определению теплового состояния двухслойного пакета состоящего из композиционного материала толщиной 82 и конструкционного материала толщиной 5l. В расчетах использовались следующие входные данные: Q =800-1200 кДж/кг; Q" = 600 кДж/кг; Ґ = \100К; Ґн = 600К; Г" = 1000АГ; рн = 2000кг/м3; рн = 2000кг/м3; рк=\850кг/м3; tK=A0c; т = 0,2с; 52=0,0\м; ,=0,0Ъи; .irI=w2 = 0,8; ТН=300К; рм = 105Па; с \,2кДж/кг-К; с2 = \,6кДж/кг-К; рх = 6000кг/V; р2 = 2000кг/м3; awl =\кВт/м2 -К; aw2 = 0,2кВт/У -К; Те1 = 3000 - 4000Я; Те2 = 500К; Лкс = 0,02кВт/м -К; сГ= \,5кДж/кг К; рГ = 10 кг/м3; М = 30 кг/кмолъ.
На рис. 3.5. представлены зависимости координат трех подвижных границ л- (О» Г(0 Х А0 от времени , построенные для трех значений теплоты уноса массы Q композиционного материала (Q = 800,1000,1200кДж/кг). Результаты показывают линейные зависимости координат подвижных границ от времени за исключением начала и окончания фазовых превращений. Некоторые колебания подвижных границ около линейных зависимостей связаны с двусторонней сходимостью в разностно-итерационном процессе.
Кроме этого, из рисунка видна эквидистантность подвижных границ хн и х , когда их движение вышло на маршевый режим , ( т.е. исключая небольшие промежутки времени в окрестности начала и окончания движения границ ). Это означает, что если теплота Q разложения связующих на границах x (t) и x"(t), одинакова, то скорость движения этих границ можно рассчитывать не по двум условиям Стефана, а по одному, например, на границе x {t), а граница x"(t) будет формироваться на некотором примерно постоянном по времени расстоянии. На рис. 3.6. приведены результаты расчета зависимости температуры от времени в точке x = xwl Sr/2 (срединная точка в слое композиционного материала) для различных значений температуры Те] окружающей среды на границе wl (ТеХ =3000К,3500К,4000К) .
