Теоретические основы химической технологии, 2020, T. 54, № 5, стр. 565-571
Расчет выхода вещества из пористой системы при переменном коэффициенте диффузии
П. Г. Ганин a, А. В. Маркова a, А. И. Мошинский a, Л. Н. Рубцова a, В. В. Сорокин a, *
a Санкт-Петербургский химико-фармацевтический университет
Санкт-Петербург, Россия
* E-mail: spcpa@outiook.com
Поступила в редакцию 17.12.2019
После доработки 19.02.2020
Принята к публикации 15.05.2020
Аннотация
Рассматривается задача экстрагирования целевого вещества из пористой системы, когда коэффициент диффузии зависит от координаты по экспоненциальному закону. Разыскивается выход целевого вещества через границу области как функция времени. Получены асимптотические при малых и больших значениях времени формулы. Построена кинетическая кривая извлечения целевого компонента в широком диапазоне времени. Предложена аппроксимационная формула достаточно простого вида.
ВВЕДЕНИЕ
При описании процессов извлечения целевого компонента (ЦК) из пористой системы часто используют уравнение диффузии. Специфику конкретной задачи при этом учитывают, выдвигая необходимые дополнительные условия. Обсуждение различных вариантов задачи экстрагирования и родственной ей задачи пропитки пористого тела представлено, например, в работах [1–12]. Как правило, для данного типа задач основной интерес представляет определение выхода ЦК из системы как функции времени. В ряде задач для этой величины удается получить представление в виде отрезков рядов, описывающих решение при малых и больших значениях времени. Возможно, в некоторых случаях, получить решение в виде интегралов, зависящих от параметров.
Развитие и распространение программных средств типа Mathcad, MathLab и др. и рост возможностей ПЭВМ позволяют легко составлять программы расчетов и получать надежные результаты. Цель настоящей работы – используя программную среду Mathcad, построить графики функций, описывающих извлечение ЦК из пористой системы при экспоненциальной зависимости коэффициента диффузии. В работах [13, 14] представлены определенные аргументы необходимости рассмотрения процесса экстрагирования в двухконтинуальной пористой системе с переменным по координате коэффициентом диффузии. Задача при экспоненциальной зависимости коэффициента диффузии была исследована в работах [13, 14] и получены несколько слагаемых разложения решения в ряды определенного вида. Эта задача нами выбрана потому, что найденных в [13, 14] слагаемых явно недостаточно для исчерпывающего представления интересных для практики интегральных характеристик решения во всем диапазоне времени. Построены необходимые кинетические кривые, которые позволяют оценить диапазоны пригодности для описания процесса результатами [13, 14].
ПОСТАНОВКА ЗАДАЧИ
Экспоненциальная зависимость коэффициента переноса (диффузии) от координаты при описании процессов тепломассопереноса использовалась ранее в работах [15–19]. В работах [15, 18] изложены и обоснованы варианты массопереноса, когда необходимо учитывать зависимость коэффициента диффузии от пространственных координат, в частности экспоненциальную зависимость. Так, например, имеются модели, в которых учитывается влияние на коэффициент диффузии напряжений в листовых материалах, вызванных набуханием полимеров при сорбции [18].
В этом случае математическая задача формулируется следующим образом. Уравнение для коэффициента диффузии имеет вид
где $k{\text{*}}$ и b – положительные постоянные. В таком случае классическое уравнение диффузии запишется следующим образом:(2)
$\frac{\partial }{{\partial x}}\left[ {k{\text{*}}\exp \left( { - bx} \right)\frac{{\partial C}}{{\partial x}}} \right] = \frac{{\partial C}}{{\partial t}},$Они выражают соответственно то, что в начальный момент времени концентрация ЦК в теле постоянная всюду величина, равная С0, и что на входе системы х = 0 поддерживается в любой момент времени постоянная, нулевая концентрация ЦК, т.е. на входе нет ЦК. Последнее условие (3) выражает ограниченность решения при стремлении координаты х к бесконечности. Фактически на границе х = 0 из тела ЦК быстро удаляется. Принято приближение полубесконечного тела, т.е. переменная х ∈ [0, ∞). То же самое имеет место и для диапазона изменения времени t ∈ [0, ∞). Следует заметить, что фактически реальное тело можно считать бесконечно протяженным только в ограниченном диапазоне времени, пока влияние второй границы незначительно.
РАЗМЕРНОСТНЫЙ АНАЛИЗ
Теория размерности позволяет делать определенные выводы о структуре решения задачи без проведения расчетов. Основной интерес для практики представляет величина $Q{\text{*}}$ – выход вещества на единицу площади через границу области х = 0. В нашей задаче размерности параметров и переменных следующие:
(4)
$\begin{gathered} \left[ {Q{\text{*}}} \right] = {{{\text{кг}}} \mathord{\left/ {\vphantom {{{\text{кг}}} {{{{\text{м}}}^{{\text{2}}}}}}} \right. \kern-0em} {{{{\text{м}}}^{{\text{2}}}}}},\,\,\,\,\left[ {{{C}_{0}}} \right]~\,\, = {{{\text{кг}}} \mathord{\left/ {\vphantom {{{\text{кг}}} {{{{\text{м}}}^{{\text{3}}}},}}} \right. \kern-0em} {{{{\text{м}}}^{{\text{3}}}},}} \\ \left[ {k{\text{*}}} \right]~\,\, = {{{{{\text{м}}}^{{\text{2}}}}} \mathord{\left/ {\vphantom {{{{{\text{м}}}^{{\text{2}}}}} {\text{c}}}} \right. \kern-0em} {\text{c}}},\,\,\,\,[t] = c,\,\,\,\,[b] = {{{\text{м}}}^{{ - 1}}}, \\ \end{gathered} $определение которой является главной целью работы. На данном этапе ясно, что задача описывается одной кривой [графиком функции Q(τ)]. На основании соотношений (4) и (5) можно делать определенные выводы об извлечении вещества. Например, при увеличении величины k* в два раза можно получить тот же выход вещества за вдвое меньшее время, поскольку параметр τ остается тем же. Более интересно, когда меняется только параметр b. В этом случае, если, например, его уменьшить в два раза, то можно получить выход вещества $Q{\text{*}}$ в два раза больше за время в четыре раза большее, чем при прежнем значении b.
В частном случае b = 0 приходим к хорошо известной задаче с постоянным значением коэффициента диффузии. В таком варианте из набора переменных (4) следует исключить параметр b. Соображения теории размерности приводят к одному безразмерному комплексу $Q{\text{*}}$/[C0($k{\text{*}}$t)1/2], которой должен равняться некоторой постоянной величине А:
Постоянную А нельзя найти из соображений размерности. Аналитическое решение задачи, приведенное в [15, 20], приводит к равенству А = 2/π1/2. Данное обстоятельство позволяет определить поведение функции Q при малых значениях времени. Поскольку при b → 0 будет τ → 0, а значит, выражение
позволяет согласовать формулы (5) и (6) при τ → 0. В общем случае, поскольку из физических соображений очевидно, что уменьшение коэффициента диффузии k* (а значит, и параметра τ = b2$k{\text{*}}$t) приводит к уменьшению выхода вещества, приходим к выводу, что функция Q(τ) – монотонно возрастающая.
ТЕОРЕТИЧЕСКИЙ АНАЛИЗ
Целесообразно переписать задачу (2), (3) в безразмерных координатах:
(7)
$z = bx,~\,\,\,\,\tau {\text{ }} = {{b}^{2}}k{\text{*}}t,\,\,\,\,~G = 1 - {C \mathord{\left/ {\vphantom {C {{{C}_{0}}}}} \right. \kern-0em} {{{C}_{0}}}},$что приводит к следующим выражениям:
(8)
$\frac{\partial }{{\partial z}}\left[ {\exp \left( { - z} \right)\frac{{\partial G}}{{\partial z}}} \right] = \frac{{\partial G}}{{\partial {\tau }}},$Применим к уравнению (8) преобразование Лапласа по переменной τ:
В результате получим следующую задачу для преобразованной функции G*:
(11)
$G{\text{*}}(0,p) = {1 \mathord{\left/ {\vphantom {1 p}} \right. \kern-0em} p},~\,\,\,\,G{\text{*}}(\infty ,p) < \infty .$Заменим независимую переменную в уравнении (10) как ξ = exp(z). Это преобразует задачу (10), (11) к виду
(12)
${\xi }\frac{{{{d}^{2}}G{\text{*}}}}{{d{{{\xi }}^{{\text{2}}}}}} = pG{\text{*}},\,\,\,G{\text{*}}(1,p) = {1 \mathord{\left/ {\vphantom {1 p}} \right. \kern-0em} p},\,\,\,\,G{\text{*}}(\infty ,p) < \infty .$Уравнение (12) сводится к уравнению Бесселя [21]. Следует заметить, что структура зависимости концентрации G от координаты z в уравнении (8) была представлена в работе [18]. Однако там не были конкретизированы дополнительные условия по переменной z.
Нашим дополнительным условиям удовлетворяет следующее решение уравнения (12):
(13)
$G{\text{*}}(p,{\xi }) = \sqrt {\xi } \frac{{{{\operatorname{K} }_{1}}\left( {2\sqrt {p{\xi }} } \right)}}{{p{{\operatorname{K} }_{1}}\left( {2\sqrt p } \right)}},$Основной интерес представляет выход целевого компонента к моменту времени τ через сечение z = 0. Имеем зависимость
(14)
$\begin{gathered} Q = - \int\limits_0^{\tau } {{{{\left. {\frac{{\partial G}}{{\partial z}}} \right|}}_{{z = 0}}}d{\tau }} = - \int\limits_0^{\tau } {{{{\left. {\frac{{\partial G}}{{\partial {\xi }}}} \right|}}_{{{\xi } = 1}}}d{\tau }} \Rightarrow Q{\text{*}}(p) = \\ = - \frac{1}{p}{{\left. {\frac{{dG*}}{{d{\xi }}}} \right|}_{{{\xi } = 1}}} = \frac{1}{{p\sqrt p }}\frac{{{{\operatorname{K} }_{0}}\left( {2\sqrt p } \right)}}{{{{\operatorname{K} }_{1}}\left( {2\sqrt p } \right)}}. \\ \end{gathered} $При выводе соотношения (14) использовалось выражение (13) и рекуррентные формулы [22], связывающие функции Макдональда различных порядков. Здесь K0(z) – функция Макдональда нулевого порядка.
Используя разложения функций Макдональда при малых и больших значениях аргумента [22] для зависимости (14), и далее выполнив почленное обращение полученных выражений, найдем
(15)
$\begin{gathered} Q\left( {\tau } \right) \cong {{Q}_{\infty }}\left( {\tau } \right) = \\ = \,\,\ln \left( {{\gamma \tau }} \right) + \frac{2}{{\tau }}\left( {\ln {\tau } - 1 + {\gamma }} \right),\,\,\,\tau \to \infty , \\ \end{gathered} $В работах [13, 14] приведены (в размерном виде) слагаемые до порядка τ2 для разложения (16) и, фактически, только одно слагаемое lnτ для разложения (15). Приведенная там формула Q(t) ≈ ln(t/4) обязана своим происхождением упрощению разложения K0(z) ≈ ln(2/z), z → 0 [19] до величины −ln(z).
Обращение выражения (14) выполним при помощи формулы Римана–Меллина [23], предварительно уменьшив сингулярность подынтегральной функции в окрестности точки p = 0. В окрестности этой точки главная часть разложения (14) следующая:
(17)
$\frac{1}{{p\sqrt p }}\frac{{{{\operatorname{K} }_{0}}\left( {2\sqrt p } \right)}}{{{{\operatorname{K} }_{1}}\left( {2\sqrt p } \right)}} \cong - \frac{{\ln p}}{p},\,\,\,\,p \to 0.$Известно [24, с. 404], что обращение выражения −ln(p)/p дает формула ln(γτ). Фактически мы это использовали при определении первого слагаемого разложения (15). Представим зависимость (14) следующим образом:
(18)
$\begin{gathered} \frac{1}{{p\sqrt p }}\frac{{{{\operatorname{K} }_{0}}\left( {2\sqrt p } \right)}}{{{{\operatorname{K} }_{1}}\left( {2\sqrt p } \right)}} = - \frac{{\ln p}}{p} + \left[ {\frac{1}{{p\sqrt p }}\frac{{{{\operatorname{K} }_{0}}\left( {2\sqrt p } \right)}}{{{{\operatorname{K} }_{1}}\left( {2\sqrt p } \right)}} + \frac{{\ln p}}{p}} \right], \\ F\left( p \right) = \frac{1}{{p\sqrt p }}\frac{{{{\operatorname{K} }_{0}}\left( {2\sqrt p } \right)}}{{{{\operatorname{K} }_{1}}\left( {2\sqrt p } \right)}} + \frac{{\ln p}}{p} \\ \end{gathered} $Обратить осталось только выражение для функции F(p). Для этого используем стандартный контур (рис. 1) в плоскости комплексной переменной p. Имеем
(19)
$Q\left( {\tau } \right) = \ln \left( {{\gamma \tau }} \right) + \frac{1}{{2{\pi }i}}\int\limits_L {\exp \left( {p{\tau }} \right)F\left( p \right)dp} ,$Свойства функций Макдональда при больших и малых значениях параметра p [22, 23] позволяют, используя лемму Жордана [23] (R → ∞, см. рис. 1) и теорему Коши, свести интегрирование по прямой L к интегрированию по двум берегам разреза по отрицательной части оси плоскости p (см. рис. 1). В результате получим следующее выражение для функции Q(τ):
(20)
$\begin{gathered} Q\left( {\tau } \right) = \ln \left( {{\gamma \tau }} \right) + \\ + \,\,2\int\limits_0^\infty {\frac{{\exp ( - {{{\xi }}^{2}}{\tau })}}{{\xi }}\left\{ {1 - \frac{1}{{{\pi \xi }}}\operatorname{Re} \left[ {\frac{{{{\operatorname{K} }_{0}}(2{\xi }i)}}{{{{\operatorname{K} }_{1}}(2{\xi }i)}}} \right]} \right\}} d{\xi }, \\ \end{gathered} $Интеграл в соотношении (20) легко численно вычисляется, поскольку современные программы, в частности Mathcad, могут оперировать с комплексными числами, и в Mathcad представлены необходимые нам функции Макдональда. Результат расчета представлен на рис. 2. Формулу (20) представляет кривая 5. Кривую 1 на этом рисунке представляет формула [13, 14]
(21)
${{Q}_{1}} = \int\limits_0^{\tau } {{{{\left. {\frac{{\partial G}}{{\partial z}}} \right|}}_{{z = 0}}}d{\tau } = 2\sqrt {\frac{{\tau }}{{\pi }}} - \frac{{\tau }}{4} + \frac{{\tau }}{8}\sqrt {\frac{{\tau }}{{\pi }}} - \frac{{3{{{\tau }}^{2}}}}{{128}}} .$Видим, что выписанные в (21) слагаемые разложения могут хорошо представлять решение до τ = 1.
Предложенная в работах [13, 14] формула ln(τ/4) в реализованном на рис. 2 диапазоне времени (кривая 6) не может аппроксимировать решение. Оценки показывают, что с точностью 10% и менее эта формула может быть использована при τ > 8000, что свидетельствует о ее слабой практической пригодности. В то же время зависимость (15) (кривая 4) пригодна для расчета с достаточной для практики точностью с момента времени τ = 10. Отметим, что рассмотренная задача может быть проиллюстрирована всего одной кривой 5. При этом при τ < 1 в области резкого изменения функции Q можно использовать разложение (16), а за пределом графика (τ > 30) следует использовать зависимость (15).
Видим, что увеличение числа слагаемых в разложении типа (16) приводит к улучшению аппроксимации точного решения, но это улучшение (продление области пригодности ряда при малых значениях времени) происходит слишком “медленно”, т.е. подобная форма расчета выглядит бесперспективно.
Поскольку представленные на рис. 2 кривые в широком диапазоне времени значительно различаются, для большей наглядности целесообразно область вблизи начала координат проиллюстрировать еще раз. Это проделано на рис. 3. На этом рисунке более точно выражены различия графиков представленных функций.
Использование асимптотических формул для функций Макдональда [22] для соотношения (14), приводящее к ряду по степеням p–1 –n/2, n = = 1, 2, … , и последующее обращение результата, что использовано для построения отрезка ряда (16), показывает, что, в принципе, можно получить бесконечный ряд, представляющий функцию Q(τ) по дробным степеням τn/2. Известно [22, 25], что функция K1(z) не имеет корней в области |argz| < π/2, а функция K0(z) регулярна вне любой окрестности |z| > r > 0 начала координат плоскости z. Таким образом, функция (14) не имеет особенностей вне начала координат p = 0 (|argp| < π), а значит, ее ряд по степеням p–1 –n/2 будет сходиться вне любой окрестности начала координат. При переходе к оригиналам коэффициенты при степенях τn/2 будут, по сравнению с коэффициентами разложения изображения по Лапласу, разделены на Γ(1 + n/2), где Γ – гамма-функция. Учитывая асимптотическое поведение гамма-функции при n → ∞ [22, 23, 25], можно легко убедиться, что разложение функции Q(τ) по степеням τn/2 будет сходиться абсолютно.
Для инженерных целей можно предложить аппроксимационную формулу
(22)
$\begin{gathered} {{Q}_{a}}\left( {\tau } \right) = {{Q}_{0}}\left( {\tau } \right){{H}_{0}}\left( {\tau } \right)\Phi \left( {6.45 - {\tau }} \right) + \\ + \,\,{{Q}_{\infty }}\left( {\tau } \right){{H}_{\infty }}\left( {\tau } \right)\Phi \left( {{\tau } - 6.45} \right), \\ \end{gathered} $Функции H0(τ) и H∞(τ), при определенном влиянии зависимости (23), подобраны так, чтобы при малых значениях τ было H0(τ)Φ(τ) ≈ 1 (при этом H∞(τ)Φ(τ) ≈ 0). В таком случае при представлении функции Qa(τ) главное значение имеет функция Q0(τ). В то же время при больших значениях τ (τ → +∞) должна “работать” функция Q∞(τ). В этом случае H0(τ)Φ(τ) ≈ 0, а H∞(τ)Φ(τ) ≈ 1. Формула (22) осуществляет плавный переход от функции Q0(τ) к функции Q∞(τ) с ростом безразмерного времени τ. Предложены следующие выражения для функций H0(τ) и H∞(τ):
Расчет по зависимости (22) показывает, что график функции Qa(τ) на экране монитора практически не отличается от графика точного решения (кривая 5 на рис. 2). Поэтому она не представлена на рис. 2.
ЗАКЛЮЧЕНИЕ
Результаты работы показывают, что получение асимптотических формул без установления их области пригодности для практических расчетов имеет ограниченное значение. Подобные формулы следует понимать как предельные при выполнении математических выкладок. Однако они могут служить основой для построения удобных аппроксимационных зависимостей. Повышение качества программного обеспечения (типа Mathcad) позволяет в расчетах опираться на сингулярные интегралы, что ранее часто было проблематично.
ОБОЗНАЧЕНИЯ
А | постоянная в зависимости (6) |
b | положительная постоянная, м–1 |
C | концентрация ЦК, кг/м3 |
C0 | концентрация ЦК в начальный момент времени t = 0, кг/м3 |
G | величина, равная 1 – C/C0 |
H0(τ), H∞(τ), Φ(τ) | аппроксимационные функции в выражении (22) |
K0(z), K1(z) | функции Макдональда нулевого и первого порядка соответственно |
$k{\text{*}}$ | положительная постоянная, м2/c |
L | прямая линия на рис. 1 |
p | переменная преобразования Лапласа |
$Q{\text{*}}$ | выход вещества на единицу площади через границу области, кг/м2 |
Q | безразмерное значение ${{Q}_{*}}$ |
R | радиус дуги на рис. 1 |
t | время, c |
х | пространственная координата, м |
z | величина, равная bx |
ξ = exp(z) | независимая переменная в уравнении (10) |
τ | безразмерное время |
ИНДЕКСЫ
Список литературы
Аксельруд Г.А., Лысянский В.М. Экстрагирование (система твердое тело–жидкость). М.: Химия, 1974.
Аксельруд Г.А., Альтшуллер М.А. Введение в капиллярно-пористую технологию. М.: Химия, 1983.
Романков П.Г., Курочкина М.И. Экстрагирование из твердых материалов. Л.: Химия, 1983.
Бабенко Ю.И., Иванов Е.В. Диффузионно-конвективное экстрагирование в замкнутый объем // Теор. осн. хим. технол. 2008. Т. 42. № 1. С. 48.
Бабенко Ю.И., Иванов Е.В. Экстрагирование в движущуюся жидкость с градиентом скорости // Теор. осн. хим. технол. 2008. Т. 42. № 5. С. 504.
Мошинский А.И. Математическая модель пропитки и экстрагирования в случае бидисперсного пористого материала // Теор. осн. хим. технол. 2009. Т. 43. № 4. С. 401.
Бабенко Ю.И., Иванов Е.В. Экстрагирование из тела с бидисперсной пористой структурой // Теор. осн. хим. технол. 2009. Т. 43. № 4. С. 408.
Мошинский А.И. Одномерные дискретные математические модели экстрагирования из пористого материала // Теор. осн. хим. технол. 2010. Т. 44. № 1. С. 45.
Мошинский А.И. Описание начальной стадии процессов пропитки и экстрагирования на основе бидисперсной модели // Теор. осн. хим. технол. 2012. Т. 46. № 1. С. 94.
Мошинский А.И. Начальная стадия пропитки пористого материала при учете конвективного массопереноса // Инж.-физ. ж. 2013. Т. 86. № 1. С. 14.
Мошинский А.И. Начальная стадия массопереноса в пористом материале с двумя типами пор // Мат. модел. 2013. Т. 25. № 5. С. 29.
Мошинский А.И. Извлечение вещества из пористого тела в движущуюся жидкость // Прикл. мех. тех. физ. 2019. № 4. С. 100.
Бабенко Ю.И., Иванов Е.В. Экстрагирование из капилляра с переменным коэффициентом диффузии // Теор. осн. хим. технол. 2008. Т. 42. № 4. С. 355.
Бабенко Ю.И., Иванов Е.В. Экстрагирование. Теория и практические приложения. Санкт-Петербург: Профессионал, 2009.
Лыков А.В. Теория теплопроводности. М.: Высшая школа, 1967.
Мошинский А.И. Анализ проблем теплопроводности при экспоненциальной зависимости коэффициента температуропроводности от координаты // Изв. Акад. наук. Энерг. 1999. № 1. С. 160.
Мошинский А.И. Решение задач теплопроводности при экспоненциальной зависимости коэффициента температуропроводности от координаты // Изв. Акад. наук. Энерг. 2000. № 3. С. 62.
Рудобашта С.П., Карташов Э.М. Химическая технология: диффузионные процессы. Часть 1. М.: Юрайт, 2018.
Рудобашта С.П., Карташов Э.М. Химическая технология: диффузионные процессы. Часть 2. М.: Юрайт, 2018.
Фролов В.Ф. Лекции по курсу “Процессы и аппараты химической технологии”. СПб.: Химиздат, 2003.
Мошинский А.И. Теория размерности в проблемах химической технологии. Saarbrücken: Lambert Academic, 2017.
Лебедев Н.Н. Специальные функции и их приложения. М.: Физматгиз, 1963.
Лаврентьев М.А., Шабат Б.В. Методы теории функций комплексного переменного. М.: Наука, 1973.
Диткин В.А., Прудников А.П. Интегральные преобразования и операционное исчисление. М.: Наука, 1974.
Справочник по специальным функциям / Под ред. Абрамовица М., Стиган И. М.: Наука, 1979.
Дополнительные материалы отсутствуют.
Инструменты
Теоретические основы химической технологии