CC BY
15ВЕСТНИК ПЕРМСКОГО УНИВЕРСИТЕТА
2021
• ФИЗИКА •
Вып. 3
УДК 532.546; 532.5.013.4 PACS 47.56.+Г, 47.55.pd
Численное моделирование процесса высокотемпературного силицирования углерод-углеродного пористого материала
В. А. Демин1'2, Б. С. Марышев1'3, А. И. Меньшиков1*
1 Пермский государственный национальный исследовательский университет, Пермь, Россия
2 Пермский национальный исследовательский политехнический университет, Пермь, Россия
3 Институт механики сплошных сред УрО РАН, Пермь, Россия
* email: aimenshikov@yandex.ru
С целью теоретического описания процесса высокотемпературного силицирования углеродного волокна разработана физико-математическая модель на основе общих уравнений, описывающих закупоривание пористой среды различного рода примесями. В рассматриваемой задаче основным механизмом первичного заполнения пор образца является диффузия газа. Предполагается, что в роли механизма закупорки выступает конденсация газообразного кремния на внутренних стенках пор. Предложены температурные зависимости для коэффициентов конденсации и испарения, характеризующие интенсивность данного процесса. Численное моделирование выполнено методом конечных разностей с использованием явной схемы в случае постоянного значения мобильной концентрации на входе и заданного линейного распределения температуры на вертикальных границах образца, что позволило смоделировать частично силицируемую область и отграничить от нее ту часть изделия, которая остается "сухой" в ходе данного процесса. Проанализирована эволюция основных характеристик композитного материала: пористости и концентрации иммобильного компонента в объеме. Количественно описаны внешние факторы, влияющие на интенсивность процесса конденсации, а именно, проведен расчет концентрации мобильного компонента в объеме в зависимости от времени. Для характерного промежуточного времени расчета вычислено максимальное изменение пористости среды, на основе которого оценено время полной пропитки материала, составившее, как оказалось, порядка 2 ч. Результаты численного моделирования в плане времени силицирования и привеса образца удовлетворительно согласуются с известными количественными экспериментальными данными.
Ключевые слова: пористая среда; силицирование; диффузия; MIM модель
Поступила в редакцию 04.08.2021; после рецензии 11.09.2021; принята к опубликованию 13.09.2021
Numerical modelling of a high temperature silicification process in the carbon-carbon porous medium
V. A. Demin12, B. S. Maryshev13, A. I. Menshikov1*
1 Perm State University, Perm, Russia
2 Perm National Research Polytechnic University, Perm, Russia
3 Institute of Continuous Media Mechanics UB RAS, Perm, Russia * email: aimenshikov@yandex.ru
© Демин В. А., Марышев Б. С., Меньшиков А. И., 2021
распространяется на условиях лицензии
Creative Commons Attribution 4.0 International (CC BY 4.0).
A mathematical model has been developed to describe the process of a high-temperature silicifica-tion process in a carbon-carbon porous material. The statement of the problem is based on general equations describing the clogging of tubules in a porous medium with various kinds of impurities. In the problem under consideration, the cause of blockage is the condensation of gaseous silicon on the inner walls of the pores. Temperature dependences for the coefficients of condensation and evaporation are proposed, which characterize the intensity of this process. It is assumed that the diffusion of gas is the main mechanism of initial filling pores of the material. Numerical modeling was carried out by the finite difference method using an explicit scheme in the case of a constant value of the mobile concentration at the input and for imposed linear temperature distribution on the vertical boundaries of the sample. Such temperature distribution made it possible to describe a partially silicified region and to delimit from it the sample part that remains "dry" during this process. The dynamics of the main physical characteristics, i.e. porosity of the composite material and concentration of the immobile component in the volume, is analyzed. The external factors influencing the intensity of the condensation process are considered quantitatively, namely, the concentration of the mobile component in the volume is calculated as a function of time. For the characteristic time interval, the maximum change in the porosity of the medium is calculated, on the basis of which the time of complete saturation of the material is estimated, which, as it turned out, is equal about 2 hours. The results of numerical simulation for silicification time and weight gain of the product are in qualitative and quantitative agreement with the known experimental data.
Keywords: porous medium; silicification; diffusion; MIM model
Received 04.08.2021; revised 11.09.2021; accepted 13.09.2021
doi: 10.17072/1994-3598-2021-3-56-62
1. Введение
В настоящее время цели и методы силицирова-ния различных материалов весьма разнообразны. В данной работе под силицированием будем понимать специфическую технологию создания композитного материала, заключающуюся в насыщении объема исходной пористой матрицы кремнием при высокой температуре. В качестве пористой среды выступает углеродное волокно с низкой проницаемостью. Кремниевые соединения с углеродом по причине высокой инертности играют роль объемного антиокислительного компонента, наличие которого определяет огнеупорные свойства высокотемпературных камер сгорания. В дополнение, изделия, выполненные из таких композитных материалов, имеют небольшой вес. Производство высококачественных крупногабаритных изделий из композитных материалов на основе углеродного волокна до сих пор является актуальной проблемой современного материаловедения [1-4]. Теоретическое описание в виде прямого численного моделирования подобных процессов отсутствует в литературе.
2. Постановка задачи
Рассмотрим плоский вытянутый вдоль вертикальной оси слой пористого материала толщиной L и высотой H, как показано на рис. 1. Пропитываемая среда представляет собой однородное углеродное волокно с открытой пористостью. Поры имеют регулярную структуру во всех направлениях. Насыщение заготовки атомами газообразного
кремния осуществляется со стороны широких граней, при этом значение концентрации на входе в образец равно Со.
Рис. 1. Геометрия расчетной области
Важной характеристикой в задаче является температура образца, так как считается, что именно она определяет динамику процессов, происходящих внутри пористого материала. Вдоль заготовки поддерживается постоянное линейное распределение температуры, соответствующее
нагреву сверху. Подобные условия нагрева позволяют моделировать естественные условия, которые формируются в реторте в процессе силициро-вания крупногабаритных изделий. Несмотря на то, что в рабочем пространстве создаются условия близкие к среднему вакууму порядка 100 -г- 150 Па, в верхней части реторты формируется тепловая шапка из более горячего остаточного газа. Тем не менее влиянием силы тяжести на сами пропиточные процессы будем пренебрегать, так как предполагается, что основным механизмом проникновения частиц в поры является интенсивная диффузия.
3. Система уравнений
Для моделирования процесса пропитки воспользуемся так называемой MIM-моделью (Mobile - Immobile Model), которая активно используется для описания фильтрации сложной по составу жидкости через пористую среду. Особенностью данной модели является представление примеси в потоке жидкости в виде двух составляющих - мобильной, дрейфующей в потоке, и иммобильной -частиц, осевших на стенках пор материала. Данная модель хорошо подходит для описания фильтрационных процессов наносуспензии [5,6] в пористой среде. Однако эту модель, в силу ее общности, представляется возможным распространить на ситуацию, когда несущей жидкостью является некий инертный газ, а в роли примеси выступает газообразный кремний.
Согласно эксперименту, кремний доставляется в поры материала в виде спутного газа с потоком аргона. Так как концентрация аргона во много раз превышает концентрацию атомов газообразного кремния, то можно сделать допущение, что их смесь представляет собой сплошную среду с задаваемой концентрацией примеси.
Теоретический подход включает в себя закон сохранения массы для атомов кремния, которые сложным образом распределены в несущей среде внутри пористого материала. В то же время давление внутри и снаружи образца предполагается одинаковым, поэтому осредненное напорное фильтрационное движение (типа Дарси) в пористой среде отсутствует. В результате рассматривается процесс исключительно диффузионной пропитки. Таким образом, уравнение, описывающее эволюцию концентрации мобильного компонента, выглядит следующим образом:
ft Фт + Ddiv )• (3.1)
Здесь cm, - концентрации мобильного и иммобильного компонентов, соответственно, ф - значение пористости в текущий момент времени, D -коэффициент диффузии, который для простоты считается постоянным.
Концентрации мобильного и иммобильного компонентов определяются исходя из объемов, занимаемых частицами в каждой из фаз:
Щ = У/У ,
// \ (3.2)
где Ут, Уг - объемы, занимаемые частицами в мобильном и иммобильном состояниях соответственно, Уf - объем несущей газовой фазы, У -общий объем, включающий поровое пространство и твердый скелет пористой среды. Из определения (3.2) видно, что ф не может превышать начальную пористость фф, в то время как Ст лежит в интервале от нуля до единицы.
Кинетическое уравнение, отвечающее за изменение концентрации иммобильной фазы, определяется природой процессов парообразования и конденсации. В качестве основы возьмем обобщенную изотерму Ленгмюра и дополним ее слагаемым, учитывающим обратный процесс испарения:
- = а[ст ( -Щ) -Ка ( -щ)щ] , (3.3)
где аи Ка - соответственно, параметр межфазного обмена и коэффициент распределения примеси, до - концентрация насыщения. Значительную роль в процессах конденсации и испарения играет температура, поэтому для коэффициентов, определяющих интенсивность этих процессов, необходимо адекватно ввести соответствующие от нее зависимости. Определим эти коэффициенты следующим образом:
а = а = А^ТС -Т (Тс > Т),
^ с,__с ; (3.4)
Ь = аКа = ^ТТ-Т; (Т > Тс).
Это феноменологические зависимости, тем не менее они отвечают всем необходимым требованиям, которым должен удовлетворять фазовый переход первого рода. А именно, при Т = Тс имеем нейтральное состояние, когда количество испаряющегося и конденсирующегося кремния одинаково (а = 0, Ь = 0). При Т > Тс в расчетной процедуре в уравнении (3.3) учитывается только второе слагаемое (Ь > 0), и в системе имеет место испарение. Если Т < Тс, то в порах происходит конденсация (а > 0). В этом случае в уравнении (3.3) существенно только первое слагаемое, которое приводит к увеличению концентрации иммобильного компонента. По мере отклонения температуры от критической точки коэффициенты испарения и конденсации возрастают по корневому закону, т.е. интенсивность этих процессов растет при удалении от границы фазового перехода. Амплитуды А и В могут отличаться друг от друга, так
как, из общих соображений, каждый отдельный атом, будучи в газовой среде, т.е. до акта конденсации и собственно в жидкой фазе находятся в разных состояниях. Иными словами, симметрия между прямым и обратным процессами, вероятнее всего, отсутствует.
Еще одной представительной характеристикой при описании процесса пропитки является пористость, которая имеет связь с концентрацией иммобильного компонента в форме:
ф = фо , (3.5)
где ф - начальная пористость среды. Из (3.5) следует, что при полном заполнении пор ф — 0.
Таким образом, полная система уравнений, описывающая диффузионную пропитку, включает в себя уравнения (3.1), (3.3), (3.5). В качестве краевых условий принимается постоянное значение концентрации мобильного компонента на боковых границах и постоянное линейное распределение температуры вдоль заготовки:
T - T T = ^ У +
(3.6)
Здесь T2 = 1800 К и T = 1700 К - температуры верхнего и нижнего торцов, соответственно. Как видно из формулы (3.6), чтобы силицировалась только часть слоя и моделировалась «сухая» область, необходимо подбирать крайние точки T2 и T1 так, чтобы значение Tc оказалось между ними.
4. Методика решения
Система нелинейных дифференциальных уравнений (3.1), (3.3), (3.5) в частных производных решалась численно методом конечных разностей с использованием явной схемы первого порядка точности по времени и второго порядка по пространству [7]. Устойчивость схемы определяется выполнением условия Куранта - Фридрихса - Ле-ви. Задача решалась в двумерной постановке с учетом однородности системы по координате z. Авторский код был реализован на языке программирования FORTRAN - 90.
По условию задачи в отношении начального состояния предполагалось, что заготовка является пустой, что означает равенство нулю значений концентраций мобильной и иммобильной компонент: Cm = ф = 0. Начальное значение пористости принималось равным фф = 0.3. В расчетах использовалась пространственная сетка с числом узлов 25:351 по координатам x и у, соответственно. Расчеты проводились при следующих значениях параметров:
1 "2 L = 0.01м,
A = 4500 с-1К-1/2,
B = 310 с-1К-1/2
c 0 = 4.7 -10-
H = 0.3 м, D = 0.7 м2/с.
(4.1)
Значения параметров А и В выбирались из соображений соответствия общего времени при проведении успешного опыта, в котором выполнялись конденсационные условия и имело место полное заполнение порового пространства. Концентрация кремния в парах аргона и коэффициент диффузии в (4.1) соответствуют действительным значениям для разреженной среды при реализации обсуждаемого технологического процесса. Концентрация парообразного кремния в рассматриваемом диапазоне температур, как и полагается, не превышает значения насыщения [8].
5. Обсуждение результатов
Представим результаты численного моделирования для приведенных выше параметров (4.1). Толщина слоя, от которой зависит время насыщения, тоже максимально приближена к значениям, для которых производятся реальные эксперименты. Сначала проанализируем поведение мобильного компонента в пористой среде. Как видно из рис. 2, концентрация мобильного кремния за очень короткое время t ~ 1 с выходит на установившийся режим и далее со временем меняется незначительно.
2.8Е-09
7
0.0025
0.005
0.0075
Рис. 2. Распределение мобильной фазы в поперечном сечении в следующие моменты времени: 1 - 0.0001 с, 2 - 1 с, 3 - 10 с на срезе 1/4 от общей высоты образца
Столь небольшое время выхода на стационарное состояние является следствием очень интенсивной диффузии, которая определяется двумя факторами: высокой температурой и разреженностью среды. Наиболее сложными для определения являются константы А и В, поэтому их значение было подобранно таким образом, чтобы результаты ре-перных экспериментов согласовались с численны-
ми данными. Критическая температура зависит от концентрации мобильного компонента, т.е. от степени насыщения паров, и подобрана таким образом, что она разделяет расчетную область на две части в соотношении 1:3 (Тс = 1775 К). Как следствие, в верхней части слоя не протекает никаких конденсационных процессов, так как температурные условия отвечают только испарению, а первоначально заготовка предполагается чистой. Конденсационные процессы происходят лишь в нижней половине заготовки.
На рис. 3 приведен профиль распределения концентрации мобильной фазы вдоль оси г в центральном сечении. В верхней части образца, где отсутствуют условия для конденсации, значение концентрации максимально и равно Со. Отметим, что в других сечениях графики концентрации мобильного компонента аналогичны.
время как нижние 2/3 слоя постепенно насыщаются кремнием.
Форма профиля определяется распределением мобильной фазы и температурой. Еще одной необходимой для общего понимания характеристикой является пористость. Как можно увидеть из рис. 6 и 7, максимальное изменение пористости за время равное 10 с составляет порядка 0.2 % от начального значения ф0.
Рис. 3. Зависимость концентрации мобильной фазы от координаты г в момент времени / = 10 с, вдоль оси симметрии слоя
На рис. 4 изображено семейство кривых для концентрации иммобильной фазы вдоль поперечной координаты х в разные моменты времени. Результаты численного моделирования показывают, что для ф, имеет место монотонное увеличение концентрации осевших частиц (кривые 1-5), так что за каждый временной шаг осаждается приблизительно одинаковое количество кремния. При этом нужно отметить, что в поперечном сечении выстраивается приблизительно параболический, симметричный профиль как для иммобильной, так и мобильной фаз (рис. 2, 4). Этот результат необходимо воспринимать как прямое следствие однородной диффузионной пропитки с обеих сторон.
На рис. 5 представлено распределение концентрации иммобильного компонента вдоль оси г для разных моментов времени. Видно, что в верхней части образца имеет место сухая область, в то
Рис. 4. Распределение концентрации иммобильной фазы вдоль оси х в следующие моменты времени: 1 - 2 с, 2 - 4 с, 3 - 6 с, 4 -8 с, 5 - 10 с на срезе 1/4 от общей высоты образца
Рис. 5. Зависимость концентрации иммобильной фазы от координаты г в следующие моменты времени: 1 - 2 с, 2 - 4 с, 3 - 6 с, 4 - 8 с, 5 - 10 с вдоль оси симметрии слоя
С течением времени пористость монотонно уменьшается. Таким образом, на полное закупоривание свободных пор должно потребоваться порядка 2 ч. Масса осевших частиц на единицу площади изделия, вычисленная из полученного значения концентрации иммобильной фазы, со-
с
т
ставляет требуемую технологией величину, что согласуется с массой кремния, полученной во время взвешивания крупногабаритного изделия после удачного эксперимента. По результатам расчетов, выполненных с помощью физико-математической модели (3.1)-(3.5), можно сделать вывод, что интенсивность пропитки кремнием полностью определяется распределением температуры на образце. Если температура некоторой области ниже критической, то процесс силицирования идет и его темпы вполне позволяют реализовать технологический процесс, так как требуемое для заполнения пор время удовлетворяет допустимым нормам (рис. 5 и 7, нижняя часть слоя). В то же время та часть изделия, у которой температура превышает критическую, остается со временим «сухой» (рис. 5 и 7, верхняя часть слоя). Также в модель заложено (см. формулы (3.4) для коэффициентов конденсации и испарения), что скорость прироста имо-бильного компонента в порах зависит от того, как выполнено захолаживание изделия.
6. Заключение
Построена физико-математическая модель, предназначенная для обобщенного теоретического описания пароконденсатной технологии кремниевой пропитки углерод-углеродного композитного материала. В частности, модель объясняет наличие "сухих пятен" на крупногабаритном изделии, отдельные части которого могут быть недостаточно захоложены в ходе технологического процесса. Результаты численного моделирования при правильно подобранных параметрах для коэффициентов конденсации и испарения дают неплохое количественное соответствие экспериментальным данным по времени насыщения, т.е. привесу изделия.
Рис. 6. Распределение пористости вдоль оси х в моменты времени: 1 - 2 с, 2 - 4 с, 3 - 6 с, 4 - 8 с, 5 - 10 с, на высоте 1/4 от общей высоты образца
0.29967
0.1
0.2
Z 0.3
Рис. 7. Продольное распределение пористости в разные моменты времени: 1 - 2 с, 2 -4 с, 3 - 6 с, 4 - 8 с, 5 - 10 с на оси симметрии слоя
Список литературы
1. Гаршин А. П., Карлин В. В., Олейник Г. С. Конструкционные карбидокремниевые материалы. Л.: Машиностроение, 1975. 152 с.
2. Шикунов С. Л., Курлов В. Н. Получение композиционных материалов на основе карбида кремния силицированием углеродных матриц // Журнал технической физики, 2017. Т. 87. Вып. 12. C. 1871-1878. DOI: 10.21883/JTF.2017.12. 45212.2291.
3. Nitin P. Advanced structural ceramics in aerospace propulsion // Nature Materials. 2016. Vol. 15. P.804-809.
4. Narottam P. (Ed). Ceramic matrix composites: materials, modeling and technology. Wiley, 2015. 694 p.
5. Demin V. A., Maryshev B. S., Menshikov A. I. Dynamics of an admixture front during the pumping of a nanofluid through a porous medium // Journal of Porous Media, 2021. Vol. 24, Issue 6, pp. 5367.
6. Demin V. A., Menshikov A. I, Maryshev B. S. On a Threshold of desorption during the pumping of a nanofluid through porous medium // J. of Physics: Conf. Series., 2021, Vol. 1809, 012003. DOI: 10.1088/1742-6596/1809/1/012003.
7. Самарский А. А. Теория разностных схем. М.: Наука, 1989. 616 с.
8. Sevast'yanov V. G., Nosatenko P. Ya. et al. Experimental and theoretical determination of the saturation vapour pressure of silicon in a wide range of temperatures // Russian Journal of Inorganic Chemistry, Vol. 55, No. 13, 2010, pp. 2073-2088.
References
1. Karlin V.V., Garshin A.P., Novikov A.A. Antiabrasion material resistant to erosive and aggressive media. In: Samsonov G.V. (eds) Refractory Carbides. Studies in Soviet Science. Springer, New York, NY: Springer, 1974. DOI: 10.1007/978-1-4684-8598-1_48
2. Shikunov S.L., Kurlov V.N. SiC-based composite materials obtained by siliconizing carbon matrices. Technical Physics, 2017. Vol. 62. No. 12, pp. 1869-1876. DOI: 10.1134/S1063784127120222.
3. Nitin P. Advanced structural ceramics in aerospace propulsion. Nature Materials, 2016, vol. 15, pp. 804-809.
4. Narottam P. (ed.) Ceramic matrix composites: materials, modeling and technology. Wiley, 2015, 694 p.
5. Demin V.A., Maryshev B.S., Menshikov A.I. Dynamics of an admixture front during the pumping of a nanofluid through a porous medium. Journal of Porous Media, 2021, vol. 24, no. 6, pp. 53-67.
6. Demin V.A., Menshikov A.I, Maryshev B.S. On a Threshold of desorption during the pumping of a nanofluid through porous medium. J. of Physics: Conf. Series., 2021, vol. 1809, 012003. DOI: 10.1088/1742-6596/1809/1/012003.
7. Samarskii A.A. The theory of difference schemes. New York: Marcel Dekker, Inc., 2001, 762 p.
8. Sevast'yanov V.G., Nosatenko P.Ya. et al. Experimental and theoretical determination of the saturation vapour pressure of silicon in a wide range of temperatures // Russian Journal of Inorganic Chemistry, Vol. 55, No. 13, 2010, pp. 2073-2088.
Просьба ссылаться на эту статью в русскоязычных источниках следующим образом: Демин В. А., Марышев Б. С., Меньшиков А. И. Численное моделирование процесса высокотемпературного силицирования углерод-углеродного пористого материала // Вестник Пермского университета. Физика. 2021. № 3. С. 56-62. doi: 10.17072/1994-3598-2021-3-56-62
Please cite this article in English as:
Demin V. A., Maryshev B. S., Menshikov A. I. Numerical modelling of a high temperature silicification process in the carbon-carbon porous medium. Bulletin of Perm University. Physics, 2021, No. 3, pp. 56-62. doi: 10.17072/1994-3598-2021-3-56-62
Сведения об авторах
1. Виталий Анатольевич Демин, д-р. физ.-мат. наук, доцент, зав. кафедрой теоретической физики, Пермский государственный национальный исследовательский университет, ул. Букирева, 15, Пермь, 614990.
2. Борис Сергеевич Марышев, канд. физ-мат. наук, н.с., Институт механики сплошных сред УрО РАН, Пермь, ул. Академика Королева, 1, Пермь, 614013; доцент кафедры теоретической физики, Пермский государственный национальный исследовательский университет, ул. Букирева, 15, Пермь, 614990.
3. Александр Игоревич Меньшиков, аспирант, Пермский государственный национальный исследовательский университет, ул. Букирева, 15, Пермь, 614990.
Author information
1. Vitaliy A. Demin, Doctor of Physical and Mathematical Sciences, Docent, Heat of the Theoretical Physics Department, Perm State University, Bukirev str. 15, Perm, Russia, 614990
2. Boris S. Marychev, Candidate of Physical and Mathematical Sciences, Researcher, Institute of Continuous Media Mechanics UB RAS, Akademika Koroleva str. 1, Perm, Russia, 614013; Professor Associate at the Theoretical Physics Department, Perm State University, Bukirev str. 15, Perm, Russia, 614990
3. Alexander I. Menshikov, PhD student, Perm State University, Bukirev str. 15, Perm, Russia, 614990
