МОДЕЛИРОВАНИЕ ОДНОМЕРНОЙ ДИФФУЗИИ ВОДОРОДА В МЕТАЛЛАХ. III. ПРОНИЦАЕМОСТЬ ЧЕРЕЗ ПЛОСКОПАРАЛЛЕЛЬНУЮ ПЛАСТИНУ, КОНТАКТИРУЮЩУЮ С ЗАМКНУТЫМИ ОБЪЕМАМИ Текст научной статьи по специальности «Математика»

ВОДОРОДНАЯ ЭКОНОМИКА

ti

HYDROGEN ECONOMY

Статья поступила в редакцию 17.04.11. Ред. рег. № 975 The article has entered in publishing office 17.04.11. Ed. reg. No. 975

Части I и II опубликованы в № 10, 11 за 2010 г.

УДК 539.219.3:669.788

МОДЕЛИРОВАНИЕ ОДНОМЕРНОЙ ДИФФУЗИИ ВОДОРОДА В МЕТАЛЛАХ. III. ПРОНИЦАЕМОСТЬ ЧЕРЕЗ ПЛОСКОПАРАЛЛЕЛЬНУЮ ПЛАСТИНУ, КОНТАКТИРУЮЩУЮ С ЗАМКНУТЫМИ ОБЪЕМАМИ

1 2 В.Н. Лобко , И.Н. Бекман

владимирский государственный университет 600000 Владимир, ул. Горького, д. 87 Тел. (4922)47-98-67, e-mail: lobko_vn@laser-2.vpti.vladimir.ru 2МГУ им. М.В. Ломоносова 119991 Москва, ГСП-1, Ленинские горы, д. 1, стр. 3 Тел. (495)939-32-12

Заключение совета рецензентов: 27.04.11 Заключение совета экспертов: 28.04.11 Принято к публикации: 30.04.11

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

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

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

проницаемость

NUMERICAL MODELING OF ONE-DIMENSIONAL DIFFUSION OF HYDROGEN IN METALS. III. PERMEABILITY THROUGH THE PLANE-PARALLEL PLATE CONTACTING

TO CLOSED VOLUMES

V.N. Lobko1, I.N. Beckman2

Vladimir State University 87 Gorky st., Vladimir, 600000, Russia Tel. (4922)47-98-67, e-mail: lobko_vn@laser-2.vpti.vladimir.ru 2Moscow State University 1 Leninskiye Gori, Moscow, GSP-1, 119991, Russia Tel. (495)939-32-12

Referred: 27.04.11 Expertise: 28.04.11 Accepted: 30.04.11

In the first article the mathematical foundations of the method for cases of constant diffusion coefficient and its dependence on the arbitrary concentration of hydrogen were presented. In this work we describe the method for the case of permeability through a plane-parallel plate. The influence of geometric dimensions of the diffusion system was considered. The modeling of integral and differential methods of permeability was held. The correctness of the results was proved.

Keywords: diffusion, hydrogen, metals, modeling, coefficient of diffusion, permeability, finite difference method, concentration dependent of diffusion coefficient, profile of concentration of hydrogen.

Игорь Николаевич Бекман

Сведения об авторе: МГУ им. М.В.Ломоносова, Химический факультет, кафедра радиохимии, д-р хим. наук, профессор.

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

Публикации: более 200.

Владимир Николаевич Лобко

Сведения об авторе: доцент Владимирского гос. университета, Факультет химии и экологии, кафедра химии, канд. хим. наук.

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

Введение

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

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

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

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

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

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

В данной работе общими параметрами были: толщина пластины, 1,5-10-3 м; температура, 680,7 К; площадь пластины, 3,14-10-4 м2; коэффициент диффузии, 5,66-10-10 м2/с; константа Генри, 2,45-10-2 моль/м3/Па; константа Сивертса, 2,45-10-2 моль/м3/Па0,5; точность метода итераций Ньютона, 10-6.

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

Влияние ограниченной емкости резервуара и приемника на характер изменения во времени концентрационного профиля диффузанта продемонстрировано на примере, представленном следующими параметрами: объем резервуара, 5-10-7 м3; объем приемника, 5-10-6 м3; давление в резервуаре, 8,95104 Па; давление в приемнике, = 0 Па; начальная концентрация в пластине, 1-10-12 моль/м3; число разбиений по оси х, 500; конечное время 1000 с; число разбиений по оси /, 100. (Отношение объема резервуара к объему пластины приемника = 10, то же для приемника = 1.)

/, м, х10 "

Рис. 1. Концентрационные профили водорода в пластине металла при различных временах диффузии (дифференциальный вариант метода проницаемости, реализованный при малых величинах объемов резервуара и приемника).

Стрелкой показан ход времени t. Шаг по времени 100 с Fig. 1. Concentration profiles of hydrogen in a metal plate at various times of diffusion (a differential version of method of permeability, realized at small sizes of volumes of the reservoir and the receiver). The course of time t is shown by an arrow.

A step of time 100 s

International Scientific Journal for Alternative Energy and Ecology № 4 (96) 2011

© Scientific Technical Centre «TATA», 2011

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

IIIII

I Ifltt HI ttt

- . X' • ! :! Ülílj № jj i[ :¡ fjjfj If 1 1

Spf fit if U

Г'«.*'j././

J_li-_1 . **« 8рШ

—4— i-

2-

- ц ......[......- .......

I, м, *1 о4

Рис. 2. Концентрационные профили диффузанта в пластине при различных временах диффузии (интегральный вариант метода проницаемости, реализованный при малом объеме приемника и большом - резервуара). Стрелкой показан ход времени t

Fig. 2. Concentration profiles of diffusant in a plate at various times of diffusion (an integral version of method of permeability, realized at small volume of the receiver and big volume of the reservoir). The course of time t is shown by an arrow

Рис. 3. Выходной поток J в зависимости от времени t для одноатомного газа (интегральный вариант метода проницаемости, реализованный при малом объеме приемника и большом - резервуара) Fig. 3. Output flux J depending on time t for one-atom gas (an integral version of method of permeability, realized at small volume of the receiver and big volume of the reservoir)

Моделирование диффузии в замкнутый малый приемный объем в интегральном варианте метода проницаемости для одноатомного газа проводили при параметрах: объем резервуара, 5-10-2 м3; объем приемника, 5-10-5 м3; давление в резервуаре, 8,95-104 Па; давление в приемнике, = 0 Па; начальная концентрация в пластине, 1-10-12 моль/м2; число разбиений по оси х, 100; конечное время 10000 с; число разбиений по оси t, 500. Эволюция во времени концентрационных профилей диффузанта по толщине пласти-

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

Моделирование интегрального и дифференциального метода проницаемости

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

В этом случае происходит замедление роста давления газа в приемнике, хотя на рис. 4, а это мало заметно. Соответствующая этому давлению приповерхностная концентрация газа (рис. 4, Ь) определяется законом Сивертса.

Рис. 4. Изменение во времени давления водорода в приемнике (а) и его поверхностной концентрации со стороны выхода (b) в интегральном варианте метода проницаемости Fig. 4. Change of pressure of hydrogen in the receiver depending on time (a) and its superficial concentration from an output side (b) in an integral version of method of permeability

Рис. 5. Эволюция концентрационного профиля в интегральном варианте метода проницаемости . Стрелкой показан ход времени t Fig. 5. Evolution of a concentration profile in an integral version of method of permeability. Time course t is shown by an arrow

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

При помощи моделирования можно решить практически важную задачу реализации эксперимента по интегральному методу проницаемости - подбор или проверку величины приемного объема и времени окончания опыта по достижении квазистационарного потока.

Рис. 6. Выходной поток J в зависимости от времени t для водорода. Интегральный (1) и дифференциальный (2)

варианты метода проницаемости Fig. 6. Output flux J depending on time t for hydrogen. Integral version (1) and differential version (2) of method of permeability

На рис. 5, представляющем концентрационные профили по толщине пластины, заметен рост концентрации диффузанта на выходе, со стороны приемника, объем которого достаточно мал. В этом состоит отличие интегрального метода от дифференциального (см. рис. 7), где выходная концентрация близка нулю. Хотя, как уже отмечалось, кривизна на линии выходного давления (загиб вниз) очень мала (рис. 4, а), снижение выходного потока на заключительной стадии диффузии видно отчетливо (рис. 6, кривая 1).

На рис. 6, 7 показаны результаты моделирования метода установления стационарного потока (дифференциальный метод проницаемости) для водорода при следующих условиях: объем резервуара, 5-10-2 м3; объем приемника, 5 м3; давление в резервуаре, 8,95-104 Па; давление в приемнике, = 0 Па; начальная концентрация в пластине, 1-10-12 моль/м2; число разбиений по оси х, 100; конечное время 10000 с; число разбиений по оси /, 500.

Рис. 7. Эволюция концентрационного профиля в дифференциальном варианте метода проницаемости. Стрелкой показан ход времени t Fig. 7. Evolution of a concentration profile in a differential version of method of permeability. Time course t is shown by an arrow

«Особые точки» на интегральных кривых проницаемости

Как известно, на кривых выходного потока дифференциального варианта метода проницаемости существуют «особые точки»: а) время достижения половинного значения стационарного потока Ту2; б) время точки перегиба тпер; в) время запаздывания тзап; г) время достижения максимального значения потока ттах (если выходной поток проходит через максимум). Наиболее показательна в этом плане кривая выходного потока (рис. 3), в меньшей степени эти точки заметны на графике рис. 6, кривая 1. По результатам моделирования указанные особые точки могут быть довольно точно определены обработкой полученных зависимостей выходного потока от времени ДО (например, рис. 3) обычными математическими методами (рис. 8). Время точки перегиба (рис. 3) соответствует максимуму на кривой рис. 8, а время максимального потока (рис. 3) - пересечению с осью абсцисс (рис. 8). По особым точкам можно рассчитать значение коэффициента диффузии Б.

Из рис. 3 и 8 видно, что на кривой потока существует еще одна «особая точка». После прохождения

International Scientific Journal for Alternative Energy and Ecology № 4 (96) 2011

© Scientific Technical Centre «TATA», 2011

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

Рис. 8. Первая производная выходного потока J (рис. 3) по времени t

Fig. 8. The first derivative of output flux J (Fig. 3) on time t

Рис. 9. Концентрационные профили по толщине l. Средняя линия - в момент прохождения максимума на кривой выходного потока. Стрелкой показан ход времени t Fig. 9. Concentration profiles on a thickness l. An average line -at the moment of maximum passage on a curve of output flux. Time course t is shown by an arrow

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

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

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

Заключение

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

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

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

Анализ «особых точек» на кривой проницаемости (например, рис. 3 и рис. 6, кривая 1) показывает, что экспериментальный метод нестационарной проницаемости с приемником ограниченного объема объективно более информативен, чем метод установления стационарного потока. В первом случае снимается временная зависимость выходного давления р(/), из которой затем рассчитывается выходной поток ,/(/) (нестационарный). Во втором случае снимается зависимость •/(/), где поток достигает стационарного состояния (см. например, рис. 6, кривая 2). При этом полученная кривая имеет одну точку перегиба и точку максимального (стационарного) потока. В первом же случае кривая .(О имеет две точки перегиба и точку максимального потока.

Список литературы

1. Взаимодействие водорода с металлами / Под ред. А.П. Захарова. М.: Наука, 1987.

2. Лобко В.Н., Бекман И.Н. Моделирование одномерной диффузии водорода в металлах. I. Метод расчета изотермической диффузии при постоянном и переменном коэффициенте диффузии // Альтернативная энергетика и экология - ШАБЕ. 2010. № 10. С. 36-42.

3. Лобко В.Н., Бекман И.Н. Моделирование одномерной диффузии водорода в металлах. II. Сорбция и десорбция в плоскопараллельной пластине, находящейся в замкнутом объеме // Там же. № 11. С. 38-42.

Г'-": - TATA — LXJ