Cybernetics and Systems Analysis, Vol. 47, No. 5, September, 2011
SYSTEMS ANALYSIS IDENTIFYING KINETIC PARAMETERS OF MASS TRANSFER IN COMPONENTS OF MULTICOMPONENT HETEROGENEOUS NANOPOROUS MEDIA OF A COMPETITIVE DIFFUSION SYSTEM V. S. Deineka,a M. R. Petryk,b and J. Fraissardc
UDC 519.6
Abstract. The identification of diffusion parameters of a two-component solution in heterogeneous nanoporous materials is analyzed. The gradient of the residual functional is obtained based on optimal control theory. The results of numerical experiments are presented. Keywords: mathematical model, competitive diffusion, initial–boundary-value problems, modified Crank–Nicholson scheme, identification of parameters, gradient method, heterogeneous nanoporous medium. INTRODUCTION New theoretical developments in systems analysis and mathematical modeling constitute the basis for information technologies of the control of research experiment and the analysis of the state of complex physical objects. The latter include multicomponent systems of competitive mass transfer in heterogeneous nanoporous media; studying their kinetics is an important problem of the modern nanophysics and nanodiffusion. Heterogeneous nanoporous media widely used in various branches of industry (medicine, petrochemistry, catalysis, partition of liquids and gases) consist of thin layers of particles of ramified porous structure with different physical and chemical (including diffusion) properties. Each layer is a multilevel system of pores with two most important subsystems (spaces): system of micropores and nanopores with high adsorption capacity and low diffusion penetration rate (intraparticle space) and system of macropores and cavities among particles with low capacity and high penetration rate (interparticle space) [1–4]. As a rule, numerous studies in this domain concerned molecular transport of isolated substances in a porous medium, where mass transfer was mainly considered on a macrolevel without significant influence of micro- and nanotransfer in particles [1–8], which is a limiting and governing factor of the general kinetics. The major problems of intermolecular interaction, based on the Langmuir–Hinshelwood principle [4], which take place in real systems of diffusion “competition” (competitive diffusion of two and more substances) are not investigated. As recent nanophysical experiments show [9–11], the physical state of a competitive diffusion system is incomparable to similar results for monodiffusion systems. Identifying kinetic transfer parameters that determine the rate of the process at macro- and microlevels and the conditions of their equilibrium is an important scientific problem, which arises along with determining the concentration and gradient fields for each diffusing substance. In the present paper, we will employ mathematical models of competitive diffusion in heterogeneous nanoporous media and methods of deriving analytical solutions with the use of integral transforms [12–19] and experimental results a
V. M. Glushkov Institute of Cybernetics, National Academy of Sciences of Ukraine, Kyiv, Ukraine,
[email protected]. bI. Pulyui National Technical University, Ternopil, Ukraine,
[email protected]. cPierre and Marie Curie University, Paris, France,
[email protected]. Translated from Kibernetika i Sistemnyi Analiz, No. 5, pp. 45–64, September–October 2011. Original article submitted September 23, 2010. 1060-0396/11/4705-0705
©
2011 Springer Science+Business Media, Inc.
705
[9, 11] to obtain an analytical and numerical solution of a model that describes a nanoporous system. Based on the optimal control theory developed for multicomponent distributed systems [20–22] and following [13, 23–26], we will substantiate the formulations of the direct and conjugate boundary-value coefficient identification problems, implement the gradient procedure of parameter identification of internal transfer kinetics, and obtain the distributions of coefficients of competitive diffusion for intraparticle transfer in a heterogeneous nanoporous medium. MATHEMATICAL MODEL OF A COMPETITIVE DIFFUSION SYSTEM The model under study is similar to the biporous model [2–10, 13, 15]. This system of complex competitive mass transfer of two components diffusing with each other in a heterogeneous medium (catalytic bed) of particles of nanoporous structure is characterized by two forms of mass transfer: diffusion in the space of macropores (due to interparticle space) and diffusion in the system of micro- and nanopores of intraparticle space. The heterogeneous catalytic domain of transfer consists of a large number of ( n + 1) thin layers, nanoporous spherical particles arranged perpendicularly to input flows (Fig. 1). It plays a defining role for thin nanoporous heterogeneous samples, especially in the case of gas diffusions before the state of absorption equilibrium in view of the system of multi-interface interactions. The competitive two-component transfer is carried out through the permeable surface of the bed of particles of nanoporous structure in two directions: axial, i.e., in the space of macropores (oz is the direction along the bed height perpendicular to the layer surfaces) and radial, i.e., in the space of micro- and nanopores. The system evolves toward the equilibrium by the gradients of concentration in macro-, micro-, and nanopores of particles (from surface to center). The molecular transport of two mutually diffusing components in the system of nanoporous particles (crystallites) is characterized by two forms of joint mass transfer: intercrystallite competitive diffusion in macropores and intracrystallite competitive diffusion in the system of micro- and nanopores. To determine the influence of each form of diffusion on the system of general joint competitive mass transfer (including their interaction at the macro- and microlevel), it is necessary to know the parameters of the system that limit the molecular transport and adsorption equilibrium conditions at the macro- and microlevels. Such key parameters are primarily the matrices of diffusions coefficients obtained based on the system representation of the kinetics of competitive diffusion. Constructing the model of competitive transfer involves systemic approach, which implies the decomposition of the complex competitive two-component diffusion D into more simple and investigated models (monodiffusion) followed by their composition to obtain the general characteristic of competitive diffusion. According to the proposed scheme, the monodiffusion of one component regardless of other components (for example, the first diffusing component (benzene)) is considered first. In the matrices of the coefficients of diffusion in micropores of particles [ Dinter ij ] and [ Dintra ij ] , i, j = 1, 2, this assumption is taken into account by components Dinter 11 and Dintra 11 . The situation for the macrolevel is similar. Given the first component, the diffusion in an element of the porous medium of the second diffusing component (hexane) is determined by components of the matrix of diffusion coefficients Dintra 12 . Given the first component, the situation for the diffusion of the second component is similar. Thus, for an individual element of the porous medium, the competitive éDintra 11 Dintra 12 ù diffusion is determined by the matrix of four components of diffusion coefficient: Dintra = ê ú. ë Dintra 21 Dintra 22 û The mathematical model of the system of competitive mass transfer in a heterogeneous nanoporous medium in view of the above physical factors can be described as an initial–boundary-value problem as follows. With allowance for [5, 7], concentrations U1m ( t , z ) and U 2 m ( t , z ) on domains W mT = ( 0, T ) ´ W m (W m = ( lm -1 , lm ) , m = 1, n + 1, l0 = 0 < l1 0 , system (14) is strictly parabolic according to Petrovskii [29]. After reducing system (14) to the diagonal form and obtaining the fourth-order differential equation ( Dintra 11 Dintra 22 - Dintra 12 Dintra 21 )
d 4 u* dr
4
- ( Dintra 11 + Dintra 22 ) p
d 2 u* dr
2
+ p 2 u * = 0,
(16)
whose general solution is the function u * ( p, r, z ) = C1 ( p, z )ch w1 pr + C 2 ( p, z )sh w1 pr + C 3 ( p , z )ch w 2
pr + C 4 ( p , z )sh w 2
pr,
we find the general solution of system (14) as Qi* ( p, r ) = ( Dintra 2s - Dintra 1s )
d 2 u*
+ E iintra [C 3 ch w 2 2
dr 2
- p u * = p( E iintra [C1 ch w1 pr + C 2 sh w1 pr] 1
pr + C 4 sh w 2
ì 1, s = 2, pr]), i = í î 2, s = 1,
ì 1, s = 2, where E ijintra = ( Dintra 2s - Dintra 1s )w 2j - 1, i = í j = 1, 2, w1,2 = b1,2 / î 2, s = 1, b1,2 = ±
(17)
p, and
( Dintra 11 + Dintra 22 ) ± [( Dintra 11 - Dintra 22 ) 2 + 4Dintra 12 Dintra 21 ]1/ 2 2( Dintra 11 Dintra 22 - Dintra 12 Dintra 21 )
p
are the roots of the characteristic equation ( Dintra 11 Dintra 22 - Dintra 12 Dintra 21 ) b 4 - ( Dintra 11 + Dintra 22 ) p b 2 + p 2 = 0.
(18)
The conditions at the boundary r = 0 (15) generate the algebraic system of equations intra intra intra intra intra intra B11 × C1 + B12 × C 3 º ( Dintra 11 E11 + Dintra 12 E 21 ) C1 + ( Dintra 11 E12 + Dintra 12 E 22 ) C 2 = 0, intra
intra intra B 21 × C1 + B 22 × C 3 º ( Dintra 21 E11
(19)
intra intra intra + Dintra 22 E 21 )C1 + ( Dintra 21 E12 + Dintra 22 E 22 )C 2 = 0,
where B ijintra = Dintra i 1 × E1intra . + Dintra i 2 × E 2intra j j LEMMA 1 (resolvability condition for system (19) and determination of the coefficients C1 and C 3 of the solution of the boundary-value problem (14), (15)). System (19) has a unique solution, its components are C1 = C 3 = 0. Proof. We shall establish the resolvability condition for system (19) by finding its determinant intra intra intra intra + Dintra 12 E 22 ) D = ( Dintra 11 E11 + Dintra 12 E 21 )( Dintra 21 E12 intra intra intra intra - ( Dintra 21 E12 + Dintra 22 E 21 )( Dintra 11 E12 + Dintra 12 E 22 ) intra intra intra intra E 22 = [Dintra 11 Dintra 21 × E11 E12 + Dintra 11 Dintra 22 × E11 intra intra intra intra + Dintra 12 Dintra 21 × E12 E 21 + Dintra 12 Dintra 22 × E 21 E 22 ] intra intra intra intra E 21 - [Dintra 11 Dintra 21 × E11 E12 + Dintra 11 Dintra 22 × E12 intra intra intra intra + Dintra 12 Dintra 21 × E11 E 22 + Dintra 12 Dintra 22 × E 21 E 22 ] intra intra intra intra E 21 ), = ( Dintra 11 Dintra 22 - Dintra 12 Dintra 21 )( E11 E 22 - E12
or by obtaining from E ijintra 710
intra intra D = ( Dintra 11 E11 - Dintra 12 E 21 )[ ([ Dintra 22 - Dintra 12 ] w12 - 1)([ Dintra 11 - Dintra 21 ]
´ w 22 - 1) - ([ Dintra 22 - Dintra 12 ] w 22 - 1)([ Dintra 11 - Dintra 21 ] w12 - 1)] .
(20)
Substituting the values of w1,2 and, respectively, the values of the roots b1,2 of the characteristic equation (16) into (20) and making transformations yield D = (( Dintra 11 - Dintra 22 ) + ( Dintra 12 - Dintra 21 ))[( Dintra 11 - Dintra 22 ) 2 + 4Dintra 12 Dintra 21 ] -1 ¹ 0 .
(21)
Hence, system (19) has a unique solution, where C1 = C 3 = 0, as was to be shown. As a result, the structure of the solution of the boundary-value problem (14), (15) becomes intra intra Q1* ( p, r ) = ( pE11 sh w1 pr ) C 2 + ( pE12 sh w 2 intra intra Q2* ( p, r ) = ( pE 21 sh w1 pr ) C 2 + ( pE 22 sh w 2
(22)
pr ) Ñ 4 ,
pr ) Ñ 4 , i = 1, 2 .
The conditions at the boundary r = R for the unknown integration constants C 2 and C 4 yield the algebraic system of equations k1 * ì intra intra ï ( E11 sh w1 pR ) C 2 + ( E12 sh w 2 pR ) C 4 = p f1 ( p, z ), ï (23) í ï ( E intra sh w1 pR ) C 2 + ( E intra sh w 2 pR ) C 4 = k 2 f 2* ( p, z ). 22 ïî 21 p The unique solution of the algebraic system (23) (nonzero determinant) and substitution of the calculated values of C 2 and C 4 into Eq. (22) yield the unique solution of the boundary-value problem (14), (15) é E intra E intra sh w pr E intra E intra sh w pr ù 2 1 21 22 ú k1 f1* Q1* ( p, r, z ) = ê 11 - 12 D D sh w1 pR sh w 2 pR ú ê ë û -
Q2* ( p, r, z ) =
intra intra E11 E 22 é sh w1 pr sh w 2 pr ù * ê ú k 2 f2 , D êësh w1 pR sh w 2 pR úû
(24)
intra intra é E intra E intra sh w pr E intra E intra sh w pr ù E 21 E 22 é sh w1 pr sh w 2 pr ù 2 1 * ê 12 21 22 ú k 2 f 2* . - 11 k f ê ú 1 1 D D D sh w1 pR sh w 2 pR ú ê êësh w1 pR sh w 2 pR úû ë û
Recovering the Original Functions. To recover the original functions, let us consider a function sh w j pr , j = 1, 2. According to the definition [28] of original function, F *j ( p, r ) = sh w j pR s +i×¥ s +i×¥ é sh w j pr ù sh w j pr pt sin i pw j r pt 1 0 1 0 e dp = e dp. F j ( t , r ) º L-1 ê ú= ò ò 2pi sin i pw j R êësh w j pR úû 2pi s 0 - i × ¥ sh w j pR s 0-i × ¥ Denote i p = l, then the transcendental equation sin l w j R = 0 has the roots l jn = points p jn = -l2j being simple poles for the n
¥
F j (t , r ) = å
( -i )sin
n =1 æç d sh w j ç dp è
pn , j = 1, 2, n = 0, ¥ , n = 1, ¥ , w jR
function F *j . Based on the Heaviside expansion theorem [19, 28],
pn r -p t 2p R e jn = ö w 2j R 2 pR ÷÷ ø p = p jn
¥
å
n =1
-
( -1) n ne
p 2n 2 w 2j R 2
t
rö æ sin ç pn ÷ . Rø è
711
According to (24), we obtain t éq1 ( t , r, z )ù = êq ( t , r, z )ú ò û ë 2
0
æ é K11 ( t - t , r ) - K12 ( t - t , r )ù é k1 f1 ù R ö ÷ çê ç K ( t - t , r ) - K ( t - t , r )ú êk f ú r ÷ dt, 22 ûë 2 2û ø è ë 21
(25)
intra intra é E intra E intra ù E12 E 21 11 22 F1 ( t , r ) F 2 ( t , r )ú , K11 ( t , r ) = ê D D ê ú ë û
K12 ( t , r ) =
intra intra E12 E 21
D
[ F1 ( t , r ) - F 2 ( t , r )] ,
K 21 ( t , r ) =
intra intra E12 E 22
D
[ F1 ( t , r ) - F 2 ( t , r )] ,
intra intra é E intra E intra ù E11 E 22 12 21 ê F1 ( t , r ) F 2 ( t , r )ú . K 22 ( t , r ) = D D ê ú ë û
Integrating and considering that the functions f j , j = 1, 2, do not depend on integration variable, we obtain intra intra æ é E intra E intra öR ù E intra E 21 E intra E12 22 E 1 ( t , r ) - 12 E 2 ( t , r )ú k1 f1 - 11 [E 1 ( t , r ) - E 2 ( t , r )] k 2 f 2 ÷ , q1 ( t , r, z ) = ç ê 11 çç ê ÷÷ r D D D ú û ø èë
(26)
ö æ intra intra intra intra é E intra E intra ù E11 E 22 ÷R ç E12 E 22 12 21 E 1 (t , r ) E 2 ( t , r )ú k 2 f 2 ÷ , q 2 ( t , r, z ) = ç [E 1 ( t , r ) - E 2 ( t , r )] k1 f1 - ê D D D ê ú ÷ r ç ë û ø è ¥
t
2 1é where E j ( t , r ) º ò F j ( t - t , r ) dt = å ( -1) n ê1- e p n =1 në 0
-
p 2n 2 w 2j R 2
t
rö ù æ ú sin ç pn R ÷ , j = 1, 2. è ø û
In this case, the residual functional has the form J ( Dintra sp ) =
T
1 2ò
2
å ( q s ( t, R / 2, Dintra sp s ) - g s )2 dt.
(27)
0 s=1
Formulation of the Initial–Boundary-Value Problem for Increments. In view of the increments of diffusion n éDD n ù DDintra intra 11 n n n 12 ú ê coefficients Dintra + DDintra , where DDintra = is the matrix of increments, and based on problem n êDD n ú DDintra 22 û ë intra 21 (8)–(11) we obtain the corresponding increments q for concentrations q. Following [25] and neglecting the terms of the second order of smallness, we obtain the following initial–boundary-value problem for the increment q of state of the system at a fixed point z ÎW: ¶ ¶ ö 1 ¶æ 2 1 ¶æ 2 n ¶ ö n (28) q ( t , r, z ) = q÷ + q ÷ , r Î ( 0, R ), t Î ( 0, T ). ç r Dintra ç r DDintra 2 ¶r è 2 ¶r è ¶t ¶ r ¶ r ø ø r r The initial conditions are (29) q( t , r, z ) t= 0 = 0, r Î ( 0, R ), z Î W . The boundary conditions with respect to the variable r are defined by the equalities n Dintra
¶ q ¶r
r= 0
n = -DDintra
¶ q ¶r
r= 0
, q
r= R
= 0, t Î ( 0, T ).
(30)
Let us introduce the generalized solution for the initial–boundary-value problem (28)–(30). Definition 1. The generalized solution of the initial–boundary-value problem (28)–(30) is a function q( t , r ) Î L2 ( 0, T ; V ) that satisfies the following equalities "w( r ) ÎV0 : 712
R
òr
2
0
¶q w dr + a( q , w ) = lq ( w ), t Î ( 0, T ), ¶t
q
t=0
= 0, r Î ( 0, R ),
q
r= R
= 0, t Î ( 0, T ),
(31)
where V = {u( t , r ) : u Î (W21 (( 0, R ))) 2 , u( t , R ) = 0, t Î ( 0, T )}, V0 = {u( r ) Î (W21 (( 0, R ))) 2 : u( R ) = 0}, R
n a( q , w ) = ò r 2 Dintra 0
R
¶q ö ¶q ¶w ¶æ 2 n dr, lq ( w ) = ò ç r DDintra ÷ w dr. ¶r ø ¶r ¶r ¶r è 0
Conjugate Problem. Following [25, 30], let us introduce the notation p ( u, u ) = (Y ( u ) - Y ( u n ), Y ( u ) - Y ( u n )) r ,
(32)
L( u ) = ( g - Y ( u n ), Y ( u ) - Y ( u n )) r , where
"u Î U = R +4
we
have
u n Î U,
Y ( u ) = q( u; t , R / 2) ,
n , u n = Dintra
T
( j, y ) r = ò j × y dt ,
j = {j i ( t )}i2=1 ,
0
2
y = {y i ( t )}2i =1 , j × y = å j i y i , and g = {g s ( t )}s2=1 . i =1
Considering (32), we obtain "u Î U 2J ( u ) = p ( u, u ) - 2L( u ) + || g - Y ( 0)|| 2r .
(33)
Let è ï + Dè ï Î U . Then "l Î ( 0, 1) we have u n + lDu n Î U. In view of problems (9)–(11) and (28)–(30), y( u n + lDu n ; t , r ) » ( u n ; t , r ) + lq ( t , r ),
(34)
where y( u n ; t , r ) and q ( t , r ) are the solutions of problems (9)–(11) and (28)–(30), respectively. Neglecting the terms of the second order of smallness and considering (33) and (34), we obtain J ( u n + lDu n ) - J ( u n ) ~ » (Y ( u n ) - g , Y ( u n + 1 ) -Y ( u n )) r , l® 0 l
J u¢ n , Du n = lim
(35)
~ where Y ( u n + 1 ) = ~y( u n + 1 ; t , R / 2) = y( u n ; t , R / 2) + q ( t , R / 2) and J ¢un is the gradient of the functional J ( u ) at the point n Î U. u = u n = Dintra
Following [25] and considering (35), we obtain the conjugate problem as ¶y 1 ¶ æ 2 n T ¶y ö + ç r ( Dintra ) ÷ = 0, r Î W d , t Î ( 0, T ) , 2 ¶t r ¶r è ¶r ø [ y] é n T ¶y ù êë( Dintra ) ¶r úû
=r= R /2
n - ( Dintra )T
¶y ¶r y
where W d = ( 0, R / 2) È ( R / 2, R ) , [ j]
r= R /2
r= R /2
1 r2
r= 0 t =T
= 0, t Î ( 0, T ), ( q( u n ; t , R / 2) - g ( t )), t Î ( 0, T ),
= 0, y
r= R
(36)
= 0, t Î ( 0, T ),
= 0, r Î W d ,
n . = j( t , R / 2 + 0) - j( t , R / 2 - 0) , u n = Dintra
713
Instead of the classical solution of the initial–boundary-value problem (36), we will use its generalized solution. Definition 2. The generalized solution of the initial–boundary-value problem (36) is a function y( t , r ) Î L2 ( 0, T ; Vd ) that satisfies the following equalities "w( r ) ÎVd 0 : T
- ò r2 0
¶y n w dt + a( y , w ) = ( q( Dintra ; t , R / 2) - g ( t ))w( R / 2), t Î ( 0, T ), ¶t y
r= R
= 0, t Î ( 0, T ), y ( t , r )
t =T
= 0, r Î W d .
(37) (38)
Taking the difference ~y( u n + 1 ; t , r ) - y( u n ; t , r ) instead of the function w in identity (37), neglecting the terms of the second order of smallness, and considering (35), (31), and (38) we obtain T
n n n+ 1 n » ò ( y( Dintra J u¢ n , DDintra ; t , R / 2) - g ( t ))( ~y( Dintra ; t , R / 2) - y( Dintra ; t , R / 2)) dt 0
T
T
TR
0
00
¶q y dt + ò a( q , y ) dt = ò ¶t
=ò
0
Hence,
ò
¶q ö ¶æ 2 n ç r DDintra ÷y drdt . ¶r è ¶r ø
~ J ¢un » y n , TR
~ ~ ~ where y n = {y nij }2i , j =1 , y nii = ò
ò
00
(39) (40)
TR ¶ æ 2 ¶q i ö ¶ æ 2 ¶q j ö ~n ÷ y i drdt for i ¹ j , i, j = 1, 2 . çr ÷ y i drdt , i = 1, 2 ; and y ij = ò ò çç r ¶r è ¶r ø ¶r è ¶r ÷ø 00
Solving the Conjugate Problem. We will solve the conjugate problem (36) by analytical methods. To this end, we rearrange problem (36) in the equivalent form ¶ 1 ¶æ 2 n ö T ¶ y ( t , r, z ) + y÷ ç r ( Dintra ) 2 ¶t ¶r ø r ¶r è =
(q
r= R /2 2
-g )
r
× d( r - R / 2), r Î ( 0, R ), z Î W , t Î ( 0, T ),
y( t , r, z ) n - ( Dintra )T
DTintra
¶ y ¶r
t =T
r= 0
(41)
= 0, r Î ( 0, R ), z Î W,
= 0; y
r= R
= 0, t Î ( 0, T ), z Î W.
Note that in contrast to the direct problem, the conjugate problem employs the transposed coefficient matrix éDintra 11 Dintra 21 ù =ê ú. ëDintra 12 Dintra 22 û
We obtain the solution of the conjugate problem (41) (as well as the solution of the direct problem) analytically with the use of the Heaviside’s operational calculus after passing to the new time variable t = T - t . In the Laplace transform [28] ¥
for functions y *ik ( p , r, z ) º L[ q i ] = ò y i ( t , r, z )e - p t d t, i = 1, 2, problem (41) has the form 0
éy * ù æ ¶ 2 2 ¶ ö éDintra 11 ÷ + pê 1 ú +ç * ç ¶ r 2 r ¶ r ÷ êDintra y 12 ë 2û è øë where g = ( q
r = R /2
(42)
- g ) × d( r - R / 2), with the boundary conditions with respect to r éD n intra - ê n 11 êD ë intra 12
714
Dintra 21 ù éy * ( p, r, z )ù ég * ù 1 = 1 , Dintra 22 úû êy *2 ( p, r, z )ú êg * ú û ë 2û ë
n ù é * Dintra ¶ y1 ( p , r, z )ù 21 ú ú ê n Dintra ú ¶ r ëy *2 ( p , r, z )û 22 û
éy * ( p , r, z )ù = 0, ê 1 ú * ëy 2 ( p , r, z )û r= 0
= 0, z ÎW . r =R
(43)
Applying the replacement y i = R × r -1Yi , i = 1, 2 , to the boundary-value problem (42), (43) yields the problem é ù d2 Dintra 21 êDintra 11 2 + p úé *ù é *ù dr ê ú êY1 ú = êg 1 ú , 2 ê ú ëY2* û ëg 2* û d D D + p intra intra ê ú 12 22 dr 2 ë û é1 Dintra s 2 ê 2 êë r
é 1 æ dY* æ d Y1* öù öù 2 çr - Y1* ÷ú - Y2* ÷ú - Dintra s 1 ê ç r = 0, 2 ç dr ÷ú ÷ú êë r çè dr è øû r = 0 øû r = 0 éY1* ( p , r, z )ù = 0, z Î W . ê * ú ëY2 ( p , r, z )û r = R
(44)
Reducing the heterogeneous system of problem (44) to the diagonal form, consisting of the fourth-order differential equation d 4 w* d 2 w* (45) ( Dintra 11 Dintra 22 - Dintra 12 Dintra 21 ) + ( Dintra 11 + Dintra 22 ) p + p 2 w * = Gi* ( p, r ), i = 1, 2, 4 2 dr dr yields the solution éR ù (46) Yi* ( p, r ) = L*i [ w * ( p, r )] = L*i ê ò H * ( p, r, r )Gi* ( p, r ) drú , êë 0 úû where L*i = ( Dintra ss - Dintra is )
¶2 ¶r2
+ p,
æ ö ¶2 ¶2 * ì 2, i = 1, Gi* ( p, r ) = ç Dintra ss + p ÷ g i* - Dintra si g , i = 1, 2; s = í 2 2 s ç ÷ ¶r ¶r î 1, i = 2, è ø H * ( p, r, r ) is the Cauchy fundamental function of the form [19] ì H - * º C 12 sh w1 pr + C 14 sin w 2 pr, 0 < r < r < R , ïï H * ( p, r, r ) = í H + * º C12 ch w1 pr + C 22 sh w1 pr ï 2 2 0 < r < r < R. ïî + C 3 cos w 2 pr + C 4 sin w 2 pr,
(47)
The unknown coefficients C ks , k = 1, 4, s = 1, 2 , are determined from the conditions, which take into account the behavior of the Cauchy function inside the domain H * ( p, r, r )
r = r- 0
¶ * H ( p, r, r ) ¶r
r = r- 0
¶2 ¶r
2
¶2 ¶r
2
H * ( p, r, r ) H * ( p, r, r )
r = r- 0
r = r- 0
- H * ( p, r, r ) -
r = r+ 0
¶ * H ( p, r, r ) ¶r ¶2 ¶r
-
2
¶3 ¶r
4
= 0,
r = r+ 0
H 0* ( p, r, r ) H * ( p, r, r )
= 0,
r = r+ 0
r = r+ 0
(48)
= 0, =1
and on the boundaries of domain [19] é1 æ d ö ù Dintra s 1 ê ç r - 1÷H * ú = 0, H * 2 è dr ø û r= 0 ër
r= R
= 0,
¶2 ¶r
2
= 0, r Î ( 0, R ), z Î W .
H*
(49)
r= R
715
Substituting expression (47) into the first boundary condition, we obtain C11 = 0 and C 13 = 0. Solving the system of fourth-order equations, we obtain from the four conditions (48) 1 1 sin ( w 2 C12 = sh ( w1 p r ), C 32 = 2 2 3/ 2 2 w1 ( w1 + w 2 ) p w 2 ( w1 - w 22 ) p 3 / 2 (C 22 - C 12 ) = (C 42 - C 14 ) = -
1 w1 ( w12
- w 22 ) p 3 / 2 1
w 2 ( w12
+ w 22 ) p 3 / 2
p r ),
(51)
ch ( w1 p r ), cos ( w 2
(50)
p r) .
The second boundary conditions (49) (for r = R) with allowance for (50) make it possible to determine C 22 and C 42 : sh ( w1 pR )C 22 + sin ( w 2 w12
pR )C 22
sin( w1
+ w 22
sin( w 2
2
1
pR ) C 42 =
p
pR )C 42
3/ 2
2 k =1 w k ( w1
p
3/ 2
- ( -1)
k -1
w 22
)
M *k ,
w 2k ( -1) k -1
2
1
=
1
å
å
2 k =1 w k ( w1
- ( -1)
k -1
w 22
)
M *k ,
ì sh ( w1 p r )ch ( w1 pR ), k = 1, where M *k = í î sin( w 2 p r )cos( w 2 pR ), k = 2. Respectively, we obtain k -1 k -1 2 1 2 ( -1) (1- ( -1) w k ) sh ( w1 p r )ch ( w1 pR ) , C 22 = å p k =1 w k ( w 2 - ( -1) k -1 w 2 ) sh ( w1 pR ) 1
C 42
=-
1
2
å p
2 k =1 w k ( w1
(53)
2
( -1) k -1 (1- ( -1) k -1 w k2 ) sin ( w 2 - ( -1)
k -1
w 22
p r )cos( w 2
sin ( w 2
)
pR )
(52)
.
pR )
As a result of the unique solution of system (52) and due to substituting constants C 22 and C 42 calculated according to (53) and, respectively, C12 and C 32 into (47), since the Cauchy function is symmetric, we obtain ì -* sh ( w1 p ( R - r ))sh ( w1 p r ) + ï H º j11 sh ( w1 pR ) ï ï sin ( w 2 p ( R - r ))sin( w 2 p r ) 0 < r < r < R, , ï + j12 sin( w 2 pR ) 1 ï * H ( p, r, r ) = í p pï + s h ( w1 p ( R - r )) sh ( w1 pr ) * + º j 21 ïH sh ( w1 pR ) ï sin( w 2 p ( R - r ))sin( w 2 pr ) ï , 0 < r < r < R, ï + j 22 sin( w 2 pR ) î ( -1) k -1 (1- ( -1) s-1 w 2k )
sh ( w1 p ( R - r ))sh ( w1 pr )
- ( -1)
sh ( w1 pR )
where
j ks =
sin( w 2
p ( R - r ))sin( w 2 sin( w 2
w k ( w12 pR )
k -1
pr )
w 22
ì 1, k = 2, The components , k = 1, 2, s = í ) î 2, k = 1.
(54)
, and
of the Cauchy function H * ( p, r, r ) (54) satisfy the existence conditions for Laplace
original functions. Hence, it is possible to pass to the original function with respect to the variable ð [28]: 1 H ( t , r, r ) º L [H ( p, r, r )] = 2pi -1
716
*
s 0+ i × ¥
ò
H * ( p, r, r ) e pt dp .
s 0-i × ¥
(55)
The generalized Heaviside expansion theorem [28] can be applied to the first component of the Cauchy function. The original function of the second component is defined by the expression [33]: ésin( w 2 p ( R - r ))sin( w 2 L-1 ê p × sin( w 2 pR ) p êë
pr ) ù ( R - r)r . ú = w2 R úû
Some transformations yield the complete expression of the original of the Cauchy function ì 2w1 R ¥ ( -1) n + 1 - l2kn t æ rö R -rö æ e sinç np ÷ sinç np ÷ ï j11 2 å 2 R R è ø è ø p n= 0 n ï ï w (R - r )r , 0 < r < r < R, ï+ j12 2 ï R H( t , r, r ) = í 2w1 R ¥ ( -1) n + 1 - l2kn t æ R - rö æ rö ï e sinç np j ÷ sinç np ÷ 21 ï 2 2 å R ø è Rø è p n= 0 n ï w 2 ( R - r )r ï , 0 < r < r < R. ïî + j 22 R
(56)
The solution of the conjugate problem will respectively have the form y i (t , r ) =
r R ù R Ré Yi ( t , r ) = êò H + (T - t , r, r )Gi ( r ) dr + ò H - (T - t , r, r )Gi ( r ) drú , i = 1, 2. r rê úû r ë0
ì 1, s = 2, In this problem, E ijintra = ( Dintra s 2 - Dintra s 1 )w 2j - 1, i = í j = 1, 2, w1,2 = b1,2 / î 2, s = 1, b1,2 = ± i
- (Dintra 11 + Dintra 22 ) ± [( Dintra 11 - Dintra 22 ) 2 + 4Dintra 12 Dintra 21 ]1/ 2 2( Dintra 11 Dintra 22 - Dintra 12 Dintra 21 )
(57)
p, where
p
are the roots of the characteristic equation (45).
ALGORITHM OF IMPLEMENTING THE GRADIENT IDENTIFICATION METHOD FOR INTRAPARTICLE DIFFUSION COEFFICIENTS OF A COMPETITIVE TRANSFER SYSTEM Let us describe the general procedure of implementing the gradient method to identify intraparticle diffusion coefficients for intraparticle space (Dintra sp , s, p = 1, 2), which is decisive and renders limiting influence on the general q transfer in the system. The matrix of the state of the system M s ( t k , z i , Dintra ) corresponds to the total accumulation of the sp
mass of the jth component in the layers of nanoporous particles for the interparticle space and intraparticle space (Fig. 2) [9]. To identify this distribution (vector) Dintra sp , one of gradient methods is used. See [25, 31] for the mathematical substantiation of applying such methods to parameter identification of multicomponent distributed systems. Based on the specific features of the problem, the minimum-error method is best suited. Taking this into account, to determine the ( q + 1)th approximation of the diffusion coefficient in the intraparticle space Dintra sp , we will apply the following identification gradient procedure defined in the matrix form || e q || 2 q+ 1 q . (58) Dintra = Dintra -Ñ J ( Dintra sp ) q || ÑJ ( Dintra )|| 2 Figure 3 presents the framework of the identification algorithm for intraparticle diffusion coefficients Dintra sp m , m = 1, n + 1, s, p = 1, 2 . 717
Benzene + Hexane
Hexane + Benzene
5.88 4.58
15.30
2.54
4.36
0.84 0.56
2.41 0.84
0.31
0.39
0.09
0.27 0.17
0.04 0.02 Sample
Sample
0.11 0.06
Fig. 2. Experimental data of the analysis of competitive mass transfer in a heterogeneous nanoporous catalytic medium.
Numerical Simulation and Identification of Kinetic Parameters of the System. The identification process here is to determine components of the diffusion coefficients Dintra 12 , m in micropores of particles with the use of the values of the already identified components of the matrix of diffusion coefficients D inter 11 , m and D intra 11 , m obtained during the monodiffusion of benzene [11]. As the experimental data, the curves of the adsorption of a two-component (benzene and hexane) mixture in heterogeneous porous zeolite catalysts (Fig. 4) [9, 10] of the gradient method described above were used, and the procedure of identifying the diffusion coefficient Dintra 12 , m (as the solution of problem (8)–(11)), which determines the influence of the mutual diffusion (diffusion component 1 (benzene) in the presence of component 2 (hexane)) in particles was implemented. Elements of the matrix of experimental data [ M exp k ] i =1, M are the values of the distribution of the total absorbed i
k =1, N
mass along the coordinate z for different time intervals of the process of two-component competitive diffusion. The results obtained in the identification of the kinetics of intraparticle competitive diffusion with the use of the identification technique described above are presented in Figs. 5–7 along the coordinate of the nanoporous medium layer thickness z for different times t = 0.02 hour (72 sec), t = 0. 31 hour (1116 sec), t = 0.84 hour (3024 sec), t = 2. 54 hour (9144 sec), and t = 5.88 hour (21170 sec) according to the program of the conducted physical experiments. The other parameters were assumed l = 0.1m, R = 0.001 m , and e = 0.8 . For all the considered times, the curves in Figs. 5a–7a represent the numerical results of the solution of the inverse problem, i.e., identified distributions of competitive components Dintra 12 of the matrices of diffusion coefficients Dintra for the intraparticle space. As was mentioned above, the competitive components Dintra 12 determine the influence of the diffusion of the first component on the total transfer in the presence of the second component. The values of the distributions of the monodiffusion components Dintra 11 of the matrices of diffusion coefficients Dintra (diffusion of the first component regardless of the others) were determined from the results of modeling of the monodiffusion for the heterogeneous nanoporous media obtained in [10, 11].
718
0 [ Dintra
Íà÷àëüíîå ïðèáëèæåíèå Initial approximation êîýôôèöèåíòîâ äëÿ of the coefficients for each mth layer êàæäîãî m-ãî ñëîÿ
n +1 ] sp , m m =1
q=0
Íîìåð èòåðàöèè Number of iteration Îïðåäåëåíèå çíà÷åíèé qsmq,smysmand Determining the values of for èòåðàöèè the qth iteration y sm q-é äëÿ
q q qsm , y sm
q ÑJ sm
R
=
ò
0
hqsm =
q +1 Dintra
sp , m
× y qsm dx 2
¶x
|| ensm ||2
-
Finding the value of the residual ensm and coefficient h
q 2 || ÑJ sm ||
q = Dintra
sp , m
q +1 d = max Dintra
q=q+1
Âû÷èñëåíèåthe ãðàäèåíòà Calculating gradient îò of the ôóíêöèîíàëà íåâÿçêè residual functional (8) (8) withñ èñïîëüthe use çîâàíèåì ðåøåíèé ïðÿìîé qsmand of solutions of the direct qsm conjugate y sm problems è ñîïðÿæåííîé ysm çàäà÷
q ¶ qsm
sp
d
