Boundary conditions at a thin membrane for normal diffusion equation which generate subdiffusion
aa r X i v : . [ c ond - m a t . s t a t - m ec h ] J a n Boundary conditions at a thin membrane for normal diffusion equation whichgenerate subdiffusion
Tadeusz Koszto lowicz ∗ Institute of Physics, Jan Kochanowski University,Uniwersytecka 7, 25-406 Kielce, Poland
Aldona Dutkiewicz † Faculty of Mathematics and Computer Science,Adam Mickiewicz University,Uniwersytetu Pozna´nskiego 4, 61-614 Pozna´n, Poland (Dated: January 19, 2021)We consider a particle transport process in a one-dimensional system with a thin membrane,described by a normal diffusion equation. We consider two boundary conditions at the membranethat are linear combinations of integral operators, with time dependent kernels, which act on thefunctions and their spatial derivatives define on both membrane surfaces. We show how boundaryconditions at the membrane change the temporal evolution of the first and second moments of parti-cle position distribution (the Green’s function) which is a solution to normal diffusion equation. Asthese moments define the kind of diffusion, an appropriate choice of boundary conditions generatesthe moments characteristic for subdiffusion. The interpretation of the process is based on a particlerandom walk model in which the subdiffusion effect is caused by anomalously long stays of theparticle in the membrane.
PACS numbers:
I. INTRODUCTION
Anomalous diffusion in a one-dimensional system isusually characterized by the following relation defined inthe long time limit [1–4] (cid:10) (∆ x ) ( t ) (cid:11) ∼ t α , (1)where (cid:10) (∆ x ) ( t ) (cid:11) is the mean square displacement of dif-fusing particle, 0 < α < α = 1 is fornormal diffusion, and α > < α ≤
1. Eq. (1) characterizes a kind of diffusionwhen the parameter α is uniquely defined. When thereis a probability distribution of α [5], the particle meansquare displacement is described by a more complicatedequation. In the following we assume that α is unique.Different models of subdiffusion lead to Eq. (1) in thelong time limit [1–3]. We mention here diffusion in a sys-tem having comb–like structure and diffusion on fractals.We focus our attention on models based on differentialequations. Subdiffusion can be described by a differentialequation with a fractional time derivative [2–4, 6] ∂P ( x, t | x ) ∂t = D α ∂ − α ∂t − α ∂ P ( x, t | x ) ∂x , (2)where P ( x, t | x ) is the Green’s function which is inter-preted as probability density that a diffusing particle is at ∗ Electronic address: [email protected] † Electronic address: [email protected] a point x at time t , D α is a subdiffusion coefficient mea-sured in the units of m /second α , and x is the initialposition of the particle. The initial condition is P ( x, | x ) = δ ( x − x ) , (3) δ is the Dirac delta function. The Riemann-Liouvillefractional derivative is defined for 0 < γ < d γ f ( t ) dt γ = 1Γ(1 − γ ) ddt Z t dt ′ f ( t ′ )( t − t ′ ) γ . (4)The physical interpretation of subdiffusion within theContinuous Time Random Walk model that leads to Eq.(1) is that a diffusing particle waits an anomalously longtime for its next jump. The probability density of thewaiting time ψ α has a heavy tail, ψ α ( t ) ∼ /t α [2–4].The other example is the subdiffusion differential equa-tion with derivatives of natural orders [7, 8] ∂P µ ( x, t ) ∂t = ∂∂x D ( x, t ) ∂P ν ( x, t ) ∂x , (5) µ, ν >
0. When D ( x, t ) = const. the solution P pro-vides Eq. (1) with α = 2 µ/ ( µ + ν ); when µ < ν we havesubdiffusion. The physical interpretation of this processis based on the non-additive Sharma–Mittal entropy [7].When D ( t ) ∼ t α − and µ = ν = 1 one gets P which leadsto Eq. (1) [9]. For diffusion in a box bounded by impen-etrable walls assuming D ( x, t ) = D | x | − Θ , Θ >
0, onegets the Green’s function which provides (cid:10) (∆ x ) ( t ) (cid:11) ∼ ( Dt ) Θ / (2+Θ) [10].The Continuous Time Random Walk model of subdif-fusion assumes that particle jumps are significantly hin-dered at each point of the system. However, in someprocesses particle diffusion can be very hindered at amembrane only. Considering diffusion of a particle alongthe x -axis, we have diffusion in a one-dimensional systemdisturbed at a single point at which the perpendicular tothe x axis membrane is placed. Obstruction of a particlepassage through the membrane may affect the nature ofdiffusion. An example is breaking the Markov propertyfor normal diffusion due to specific boundary conditionsat the membrane [11]. The change of the character of dif-fusion can also be caused by the presence of an adsorbingwall in a system in which the process is described by thenormal diffusion equation. A boundary condition at thewall involves an integral operator with a time dependentkernel [12].The mechanisms of a particle transport through themembrane may be very complicated. Some of them leadto great difficulties in particle transport inside the mem-brane, which affect the process in the outer regions. Froma mathematical point of view, these mechanisms providespecific boundary conditions at the membrane [13, 14],see also the discussion in Ref. [11] and the referencescited therein, the list of references regarding this issuecan be significantly extended. In particular, the bound-ary conditions may contain fractional derivatives [15–17].The diffusing particle can stay in the membrane for a longtime, which can happen, among others, in a lipid bilayermembrane [18].The question considered in this paper is whether thereare boundary conditions at the membrane that changethe nature of the diffusion process described by the nor-mal diffusion equation in such a way that the processhas subdiffusion properties. In our considerations we arebased on the Laplace transforms of the Green’s functions.We consider the boundary conditions for which Laplacetransforms are linear combination of probabilities andfluxes defined on both membrane surfaces with coeffi-cients depending on the Laplace transform parameter.As it is argued in Ref. [11], such boundary conditions of-ten occur in models of diffusion in a membrane system. Inthe time domain the boundary conditions are expressedby integral operators with time–dependent kernels. Weshow that appropriately chosen boundary conditions atthe membrane lead to Green’s functions for the normaldiffusion equation providing Eq. (1) with 0 < α < II. METHOD
In this section we consider how boundary conditionsat the membrane are related to the first and second mo-ments of distribution of particle location. This distribu-tion (Green’s function) is a solution to normal diffusionequation with the initial condition Eq. (3).
A. Boundary conditions at a membrane
The normal diffusion equation with constant diffusioncoefficient D is ∂P ( x, t | x ) ∂t = D ∂ P ( x, t | x ) ∂x . (6)In the following we use the Laplace transform L [ f ( t )] =ˆ f ( s ) = R ∞ e − st f ( t ) dt . In terms of the Laplace transformEq. (6) is s ˆ P ( x, s | x ) − P ( x, | x ) = D ∂ ˆ P ( x, s | x ) ∂x . (7)We assume that a thin membrane is located at x = 0. Athin membrane means that the particle can stop insidethe membrane, but its diffusive motion is not possible init. We additionally assume that x <
0. The regionsbounded by the membrane are denoted as A = ( −∞ , B = (0 , ∞ ). In the following the function P and adiffusive flux J are marked by the indexes A and B whichindicate the location of the point x . In the time domainthe flux is defined as J i ( x, t | x ) = − D ∂P i ( x, t | x ) ∂x , (8)its Laplace transform isˆ J i ( x, s | x ) = − D ∂ ˆ P i ( x, s | x ) ∂x , (9) i ∈ { A, B } .We consider boundary conditions at a thin membranewhich in terms of the Laplace transform areˆ P B (0 + , s | x ) = ˆΦ( s ) ˆ P A (0 − , s | x ) , (10)ˆ J B (0 + , s | x ) = ˆΞ( s ) ˆ J A (0 − , s | x ) . (11)Assuming that the system is unbounded, the aboveboundary conditions are supplemented byˆ P A ( −∞ , s | x ) = ˆ P B ( ∞ , s | x ) = 0 . (12)In the time domain the boundary conditions (10)–(12)are P B (0 + , t | x ) = Z t dt ′ Φ( t − t ′ ) P A (0 − , t ′ | x ) , (13) J B (0 + , t | x ) = Z t dt ′ Ξ( t − t ′ ) J A (0 − , t ′ | x ) , (14) P A ( −∞ , t | x ) = P B ( ∞ , t | x ) = 0 . (15)The question arises whether Eqs. (10) and (11) do notconstitute too narrow set of linear boundary conditions ata thin membrane. Let us consider the following boundaryconditions γ ( s ) ˆ P A (0 − , s | x ) + γ ( s ) ˆ J A (0 − , s | x ) (16)= γ ( s ) ˆ P B (0 + , s | x ) + γ ( s ) ˆ J B (0 + , s | x ) ,λ ( s ) ˆ P A (0 − , s | x ) + λ ( s ) ˆ J A (0 − , s | x ) (17)= λ ( s ) ˆ P B (0 + , s | x ) + λ ( s ) ˆ J B (0 + , s | x ) . Eqs. (16) and (17) are more general that Eqs. (10) and(11). However, as it is shown in Appendix I, the bound-ary conditions (16) and (17) and the ones (10) and (11)provide the same Green’s functions whenˆΦ( s ) = 2 √ DsW B ( s ) W ( s ) + 2 √ DsW A ( s ) , (18)ˆΞ( s ) = 2 √ DsW B ( s ) W ( s ) − √ DsW A ( s ) , (19)where W ( s ) = ( λ ( s ) − √ Dsλ ( s ))( γ ( s ) + √ Dsγ ( s )) (20) − ( λ ( s ) + √ Dsλ ( s ))( γ ( s ) − √ Dsγ ( s )) ,W A ( s ) = 12 (cid:20)(cid:18) γ ( s ) √ Ds + γ ( s ) (cid:19)(cid:18) λ ( s ) + √ Dsλ ( s ) (cid:19) (21) − (cid:18) λ ( s ) √ Ds + λ ( s ) (cid:19)(cid:18) γ ( s ) + √ Dsγ ( s ) (cid:19)(cid:21) ,W B ( s ) = 12 (cid:20)(cid:18) γ ( s ) √ Ds + γ ( s ) (cid:19)(cid:18) λ ( s ) − √ Dsλ ( s ) (cid:19) (22) − (cid:18) λ ( s ) √ Ds + λ ( s ) (cid:19)(cid:18) γ ( s ) − √ Dsγ ( s ) (cid:19)(cid:21) , under conditions W ( s ) = 0 and W A ( s ) = ± W ( s ) / √ Ds .Since the boundary conditions determine the solutions tothe diffusion equation uniquely, the boundary conditionsEqs. (16) and (17) can be written as Eqs. (10) and (11)under the above mentioned conditions which interpreta-tion is given in Appendix I. In general, the boundaryconditions (16) and (17) depend on eight functions γ i and λ i , i ∈ { , , , } , while the boundary conditionsEqs. (10) and (11) are generated by two functions ˆΦ andˆΞ only. Thus, due to Eqs. (18) and (19), the boundaryconditions Eqs. (10) and (11) are uniquely determinedby Eqs. (16) and (17) but the opposite is not true.For example, one of the most used boundary condi-tions at the membrane is J A (0 , t | x ) = λ P A (0 − , t | x ) − λ P B (0 + , t | x ), λ , λ >
0, supplemented by thecondition that the flux is continuous J A (0 − , t | x ) = J B (0 + , t | x ). These boundary conditions can be writ-ten in the form of Eqs. (13) and (14) with Φ( t ) = J A (0 - ,t|x ) x (cid:1)(cid:2) t (cid:3)(cid:4)(cid:2) t (cid:3) P A (0 - ,t|x ) P B (0 + ,t|x )J B (0 + ,t|x ) x (cid:5) (cid:6) FIG. 1: Illustration of the boundary conditions at a thin mem-brane. The operator Φ changes the probabilities that theparticle is located at the membrane surface, the operator Ξchanges the flux flowing through the membrane. λ √ D (cid:20) √ Dt − λ √ D e λ tD erfc (cid:16) λ √ t √ D (cid:17)(cid:21) and Ξ( t ) = δ ( t ), whereerfc( u ) = (2 / √ π ) R ∞ u e − τ dτ is the complementary errorfunction [11]. For this case we have ˆΦ( s ) = λ / ( λ + √ Ds ) and ˆΞ( s ) = 1.The Laplace transform of Green’s functions for normaldiffusion equation obtained for the boundary conditions(10)–(12) are [11]ˆ P A ( x, s | x ) = 12 √ Ds e −| x − x | √ sD (23) − ˆΦ( s ) − ˆΞ( s )ˆΦ( s ) + ˆΞ( s ) ! √ Ds e ( x + x ) √ sD , ˆ P B ( x, s | x ) = ˆΦ( s )ˆΞ( s )ˆΦ( s ) + ˆΞ( s ) ! √ Ds e − ( x − x ) √ sD . (24)In the following we use the function P M defined as P M ( t | x ) = 1 − Z −∞ P A ( x, t | x ) dx (25) − Z ∞ P B ( x, t | x ) dx. Eqs. (23), (24), and the Laplace transform of Eq. (25)provideˆ P M ( s | x ) = e x √ sD s ˆΦ( s ) (cid:16) − ˆΞ( s ) (cid:17) ˆΦ( s ) + ˆΞ( s ) . (26)The function P M is the probability of not finding the par-ticle in the regions A or B at time t . The Green’s func-tions Eqs. (23) and (24) are normalized when P M ( t | x ) ≡
0. Thus, the normalization condition is met when theflux through the membrane is continuous, ˆΞ( s ) ≡
1, orwhen ˆΦ( s ) ≡ + is still zero with a non-zero fluxflowing from the region A to B .In Sec.II B we consider a model of a random walk ofa particle as it passes through a membrane. This modelgives a stochastic interpretation of the boundary condi-tions. It also imposes a certain condition on the functionsˆΦ and ˆΞ. B. Random walk model of particle passing throughthe membrane
We consider a model in which a diffusing particle canbe inside a thin membrane for a very long time. x (cid:1) a (cid:2) t (cid:3) a b x - x + P b (x + ,t;x )P a (x - ,t;x )J(x,t;x ) (cid:1) b (cid:2) t (cid:3) FIG. 2: Illustration of the transport process described by Eq.(27). The diffusive flux J at the point x depends on thedistribution of waiting times ψ a and ψ b for the particle tojump between the neighbouring points x − and x + located inthe media a and b , respectively. xx (cid:1) (cid:2) - + (cid:3)(cid:4) t (cid:5) (cid:3) (cid:6) (cid:4) t (cid:5) J B J A (cid:3)(cid:4) t (cid:5) P A P B P M FIG. 3: Transport of a particle through the membrane. Point0 represents the inside of the membrane where the particlecan stay even for a long time, points 0 − and 0 + mark the po-sitions of the particle on membrane surfaces, a more detaileddescription is in the text. We define the Laplace transform of diffusive flux thatflows through the boundary between two media a and b located at x asˆ J ( x, s | x ) = ǫs ˆ ψ a ( s )2(1 − ˆ ψ a ( s )) ˆ P a ( x − , s | x ) (27) − ǫs ˆ ψ b ( s )2(1 − ˆ ψ b ( s )) ˆ P b ( x + , s | x ) , where ˆ ψ i ( s ) is the Laplace transform of probability den-sity of time which is needed to take a particle next stepin the medium i , i ∈ { a, b } , ǫ = x + − x − is a length ofparticle step, see Fig. 2, the derivation of Eq. (27) is inAppendix II. The function ˆ ψ is expressed by the formula[15] ˆ ψ ( s ) = 11 + ǫ η ( s ) , (28)where the function η , which in practice determines a kindof diffusion, fulfils the condition η ( s ) → s → ǫ we have ˆ ψ ( s ) = 1 − ǫ η ( s ). Weassume that the particle can stay inside the membraneat the point 0. Let the points 0 − and 0 + represent pointslocated on the membrane surfaces. Applying Eq. (27) tothe system presented in Fig. 3 we getˆ J A (0 − , s | x ) = s ǫη ( s ) ˆ P A (0 − , s | x ) (29) − s ǫη M ( s ) ˆ P M ( s | x ) , ˆ J B (0 + , s | x ) = s ǫη M ( s ) ˆ P M ( s | x ) (30) − s ǫη ( s ) ˆ P B (0 + , s | x ) , where ˆ ψ M ( s ) = 11 + ǫ η M ( s ) . (31)For normal diffusion the distribution of time to take theparticle next step is given by Eq. (28) with η ( s ) = s D . (32)We are going to find the function η M which togetherwith Eqs. (29), (30) provide Eq. (11). The probabilitythat the particle is inside the membrane, represented bythe point 0, is P M ( t | x ). From Eqs. (23) and (24) we getˆ P A (0 − , s | x ) = ˆΞ( s )ˆΦ( s ) + ˆΞ( s ) ! e x √ sD √ Ds , (33)ˆ P B (0 + , s | x ) = ˆΦ( s )ˆΞ( s )ˆΦ( s ) + ˆΞ( s ) ! e x √ sD √ Ds . (34)Combining Eqs. (11), (26), and (29)–(34) we obtain η M ( s ) = ˆΦ( s )(1 − ˆΞ ( s ))2ˆΞ( s )( ˆΦ( s ) + ˆΞ( s )) r sD . (35)The boundary conditions at the membrane Eqs. (10) and(11) are generated by the residence time of the particlein the membrane with distribution Eq. (31) in which η M is expressed by Eq. (35). However, due to the normal-ization condition ˆ ψ M (0) = 1, there is η M ( s ) → s →
0. This condition and Eq. (35) provide the followingcondition for the functions ˆΦ and ˆΞ √ s ˆΦ( s )(1 − ˆΞ ( s ))ˆΞ( s )( ˆΦ( s ) + ˆΞ( s )) → s → C. First and second moments of P ( x, t | x ) We derive the relations between the moments of parti-cle locations at time t , generated by Green’s functions P A and P B , and the functions Φ and Ξ that define boundaryconditions at the membrane. The moments are calcu-lated by means of the formula (cid:10) x i ( t ) (cid:11) = Z −∞ x i P A ( x, t | x ) dx (37)+ Z ∞ x i P B ( x, t | x ) dx. From Eqs. (23), (24), and the Laplace transform of Eq.(37) we get L [ h x ( t ) i ] = x s + e x √ sD ˆ v ( s ) , (38) L (cid:2)(cid:10) x ( t ) (cid:11)(cid:3) = x s + 2 Ds + e x √ sD ˆ w ( s ) , (39)where ˆ v ( s ) = √ Ds / (cid:16) ˆΦ( s ) − (cid:17) ˆΞ( s )ˆΦ( s ) + ˆΞ( s ) , (40)ˆ w ( s ) = 2 Ds (cid:16) ˆΞ( s ) − (cid:17) ˆΦ( s )ˆΦ( s ) + ˆΞ( s ) . (41)We consider the first and second moments in the limitof long time which corresponds to the limit of smallparameter s . If s ≪ D/ | x | , which corresponds to t ≫ | x | /D , we can use the approximation e x √ s/D ≈ z ( s ) = ˆ w ( s ) + 2 Ds . (42)Then, Eqs. (38) and (39) read L [ h x ( t ) i ] = x s + ˆ v ( s ) , (43) L (cid:2)(cid:10) x ( t ) (cid:11)(cid:3) = x s + ˆ z ( s ) . (44) From Eqs. (41) and (42) we getˆ z ( s ) = 2 Ds (cid:16) ˆΞ( s ) + 1 (cid:17) ˆΞ( s )ˆΦ( s ) + ˆΞ( s ) . (45)From Eqs. (40) and (45) we obtainˆΦ( s ) = ˆ z ( s ) + 2 q Ds ˆ v ( s )ˆ z ( s ) − q Ds ˆ v ( s ) , (46)ˆΞ( s ) = ˆ z ( s ) + 2 q Ds ˆ v ( s ) Ds − ˆ z ( s ) + 2 q Ds ˆ v ( s ) . (47)Thus, knowing the boundary conditions at the mem-brane we can determine the time evolution of the firstand second moments of the particle position distributionin the long time limit putting Eqs. (40) and (45) to Eqs.(43) and (44), respectively, and then calculating the in-verse Laplace transforms of the obtained functions. Con-versely, the temporal evolution of these moments definesthe boundary conditions at the membrane by Eqs. (46)and (47). D. Boundary conditions at the membranegenerated by the first and second moments
The boundary conditions at the membrane generatedby Eqs. (10), (11), (46), and (47) read (cid:18) s ˆ z ( s )2 D − s / ˆ v ( s ) √ D (cid:19) ˆ P B (0 + , s | x ) (48)= (cid:18) s ˆ z ( s )2 D + s / ˆ v ( s ) √ D (cid:19) ˆ P A (0 − , s | x ) , (cid:18) − s ˆ z ( s )4 D + s / ˆ v ( s )2 √ D (cid:19) ˆ J B (0 + , s | x ) (49)= (cid:18) s ˆ z ( s )4 D + s / ˆ v ( s )2 √ D (cid:19) ˆ J A (0 − , s | x ) . Due to the formula L − h ˆ g ( s )ˆ h ( s ) i = Z t g ( t ′ ) h ( t − t ′ ) dt ′ , (50)in the time domain the boundary conditions Eqs. (48)and (49) take the forms of integral operators with thekernels depending on the functions v ( t ) and z ( t ). E. Green’s functions generated by the first andsecond moments
From Eqs. (23), (24), (26), (46), and (47) we getˆ P A ( x, s | x ) = e −| x − x | √ sD √ Ds (51) − (cid:18) − s ˆ z ( s )2 D + s / ˆ v ( s ) √ D (cid:19) e ( x + x ) √ sD √ Ds , ˆ P B ( x, s | x ) = (cid:18) s ˆ z ( s )4 D + s / ˆ v ( s )2 √ D (cid:19) e − ( x − x ) √ sD √ Ds , (52)we also obtainˆ P M ( s | x ) = (cid:18) − s ˆ z ( s )2 D (cid:19) e x √ sD s . (53) III. BOUNDARY CONDITIONS AT A THINMEMBRANE WHICH GENERATESUBDIFFUSION
We consider how the temporal evolution of the firstand second moments that are power functions of timeaffects the boundary conditions and Green’s functions.These moments lead to the relation Eq. (1).
A. Moments as power functions of time
We consider time evolution of the first and second mo-ments, and consequently the mean square displacement,as power functions of time. We use Eqs. (43) and (44)assuming ˆ v ( s ) = Bs β , (54)ˆ z ( s ) = As α , (55)where α, β, A >
0. In the time domain we have h x ( t ) i = x + B ′ t β , (56) (cid:10) x ( t ) (cid:11) = x + A ′ t α , (57)where A ′ = A/ Γ(1 + α ) and B ′ = B/ Γ(1 + β ). Using theequation (cid:10) (∆ x ) ( t ) (cid:11) = (cid:10) x ( t ) (cid:11) − h x ( t ) i , (58)we get (cid:10) (∆ x ) ( t ) (cid:11) = A ′ t α − B ′ t β − x B ′ t β . Since (cid:10) (∆ x ) ( t ) (cid:11) >
0, we suppose α ≥ β , but if α = 2 β we assume that A ′ > B ′ . Under these conditions forsufficiently long times this relation can be approximatedas (cid:10) (∆ x ) ( t ) (cid:11) = ˜ At α , (59)where ˜ A = A ′ when α > β and ˜ A = A ′ − B ′ when α = 2 β . B. Boundary conditions at the membrane
Combining Eqs. (48), (49), (54), (55), and using thefollowing formula valid for bounded function g L − [ s γ ˆ g ( s )] = d γ g ( t ) dt γ , < γ < , (60)we get the boundary conditions at the membrane withRiemann–Liouville fractional time derivatives (cid:18) A D ∂ − α ∂t − α − B √ D ∂ / − β ∂t / − β (cid:19) P B (0 + , t | x ) (61)= (cid:18) A D ∂ − α ∂t − α + B √ D ∂ / − β ∂t / − β (cid:19) P A (0 − , t | x ) , (cid:18) − A D ∂ − α ∂t − α + B √ D ∂ / − β ∂t / − β (cid:19) J B (0 + , t | x ) (62)= (cid:18) A D ∂ − α ∂t − α + B √ D ∂ / − β ∂t / − β (cid:19) J A (0 − , t | x ) . The discussion in Sec.III A shows that 0 < α ≤ ≤ β ≤ /
2. Thus, all fractional derivatives in the aboveboundary conditions are of non-negative orders which arenot greater than one.
C. Solutions to diffusion equation
From Eqs. (51)–(55) we getˆ P A ( x, s | x ) = 12 √ Ds h e −| x − x | √ sD − e ( x + x ) √ sD i (63)+ (cid:18) As − α +1 / D / − Bs − β D (cid:19) e ( x + x ) √ sD , ˆ P B ( x, s | x ) = (cid:18) As − α +1 / D / + Bs − β D (cid:19) e − ( x − x ) √ sD , (64)ˆ P M ( s | x ) = (cid:18) − As − α D (cid:19) e x √ sD s . (65)We calculate the inverse Laplace transforms of Eqs.(63)–(65) using the formulas L − [e − x √ s/D / √ Ds ] =e − x / Dt / √ πDt , L − [e − x √ s/D /s ] = erfc( x/ √ Dt ), x >
0, and [19] L − h s ν e − as β i ≡ f ν,β ( t ; a ) (66)= 1 t ν +1 ∞ X k =0 k !Γ( − kβ − ν ) (cid:16) − at β (cid:17) k , a, β >
0. In this way we obtain the following solutions tothe diffusion equation Eq. (6) with the boundary condi-tions Eqs. (61) and (62) P A ( x, t | x ) = 12 √ πDt (cid:20) e − ( x − x Dt − e − ( x + x Dt (cid:21) (67)+ A D / f − α +1 / , / (cid:18) t ; − ( x + x ) √ D (cid:19) − B D f − β, / (cid:18) t ; − ( x + x ) √ D (cid:19) ,P B ( x, t | x ) = A D / f − α +1 / , / (cid:18) t ; x − x √ D (cid:19) (68)+ B D f − β, / (cid:18) t ; x − x √ D (cid:19) . The inverse Laplace transform of Eq. (65) reads P M ( t | x ) = erfc (cid:18) − x √ Dt (cid:19) − A D f − α, / (cid:18) t ; − x √ D (cid:19) . (69) -30 -20 -10 0 10 20 300.000.020.040.060.08 P ( x ,t | x ) x =0.6, t=10 =0.6, t=8 =0.6, t=6 =0.6, t=10 =0.6, t=8 =0.6, t=6 FIG. 4: Plots of the Green’s functions Eqs. (70) and (71)which are solutions to the normal diffusion equation withfractional boundary conditions Eqs. (61) and (62) (lines withopen symbols) and the Green’s function Eq. (73) for the subd-iffusion equation (lines with filled symbols), for times given inthe legend, the other parameters are α = 0 . D = D α = 10,and x = −
1, the values of parameters are given in arbitrarilychosen units.
D. Comparison of two models
We compare the Green’s functions for the diffusionequation (6) and for the fractional subdiffusion equation -30 -20 -10 0 10 20 300.000.020.040.060.08 P ( x ,t | x ) x =0.9, t=10 =0.9, t=8 =0.9, t=6 =0.9, t=10 =0.9, t=8 =0.9, t=6 FIG. 5: The description is similar to the one in Fig. 4, buthere α = 0 . P M ( t | x ) t =0.6=0.7=0.8=0.9=0.95=0.99 FIG. 6: Plots of P M ( t | x ) Eq. (72) for different α , the otherparameters are D = D α = 10 and x = − (2). In both cases we assume the boundary conditionsthat the functions are continuous at the membrane, butthe flux is continuous for the solutions to Eq. (2) only.The discontinuity of the flux at the membrane in thefirst case generates a subdiffusion effect. We also assumethat the Green’s functions for both equations generatethe same relation (cid:10) (∆ x ) ( t ) (cid:11) = 2 D α t α Γ(1 + α ) . Thus, we solve the normal diffusion equation with theboundary conditions (61) and (62) with A = 2 D α / Γ(1 + α ) and B = 0. We obtain P A ( x, t | x ) = 12 √ πDt (cid:18) e − ( x − x Dt − e − ( x + x Dt (cid:19) (70)+ D α D / Γ(1 + α ) f / − α, / (cid:18) t ; | x + x |√ D (cid:19) ,P B ( x, t | x ) = D α D / Γ(1 + α ) (71) × f / − α, / (cid:18) t ; x − x √ D (cid:19) , the function P M is P M ( t | x ) = erfc (cid:18) − x √ Dt (cid:19) (72) − D α D Γ(1 + α ) f − α, / (cid:18) t ; − x √ D (cid:19) , The solution to fractional diffusion equation in termsof the Laplace transform isˆ P ( x, s | x ) = s − α/ √ D α e −| x − x | q sαDα . In the time domain we get P ( x, t | x ) = 12 √ D α f − α/ ,α/ (cid:18) t ; | x − x |√ D α (cid:19) . (73)The plots of the Green’s functions Eqs. (70), (71) forthe model considered in this paper and for the ones Eq.(73) being solutions to the fractional subdiffusion equa-tion are shown in Figs. 4 and 5. The Green’s functionsare assumed to be continuous at the membrane. How-ever, as opposed to Eq. (73), the flux is assumed tobe discontinuous at the membrane for the functions Eqs.(70) and (71). Then, the particle can stay inside themembrane as it passes through it. The plots show thatthe subdiffusion effect is achieved by anomalous long res-idence times within the membrane. The effect is strongerfor less α . In Fig. 6 we can see that the probability offinding a particle inside the membrane strongly dependson α . If α is greater, the mobility of the particle is greaterand it is less likely to remain in the membrane. From Eqs.(35), (46), (47), (54), and (55) we obtain η M ( s ) = 2 √ DA s α − / (cid:18) − A D s − α (cid:19) (74) × − B √ D s − β +1 / B √ DA s α − β − / ! , In the limit of small s we get η M ( s ) ≈ √ Ds α − / . Usingthe approximation ˆ ψ M ( s ) ≈ − ǫ η M ( s ) ≈ e − ǫ η M ( s ) andEq. (66) with ν = 0 we find that ψ M has the heavy tail ψ M ( t ) ≈ κt α +1 / , t → ∞ , (75) where κ = 2 ǫ √ D ( α − / /A Γ(3 / − α ). This tail is”heavier” than the one ψ α ( t ) ∼ /t α , t → ∞ , for themodel provides the fractional subdiffusion equation Eq.(2) [2, 4]. IV. FINAL REMARKS
We have shown how boundary conditions at a thinmembrane affect the first and second moments of proba-bility density P ( x, t | x ) of a particle position at x at time t . This probability is a solution to the normal diffusionequation for the initial condition P ( x, | x ) = δ ( x − x ).We also considered the inverse problem, how knowing thetime evolution of these moments we can find the bound-ary conditions and the Green’s functions. The first andsecond moments, considered in the long time limit, alsodetermine the temporal evolution of (cid:10) (∆ x ) ( t ) (cid:11) which isusually considered as the definition of the kind of diffu-sion. We have shown that assuming appropriate bound-ary conditions we can change the kind of diffusion in themembrane system despite the fact that outside the mem-brane the process is described by the normal diffusionequation. The other remarks are as follows.(1) Whether the relation (1) defines a kind of diffusionalone has been treated by some authors rather as an openproblem. It has been shown in Ref. [20] that an appropri-ate combination of subdiffusion and superdiffusion leadsto Green’s functions that generate Eq. (1) with α = 1which is characteristic for normal diffusion, although theprocess is non–Gaussian and non–Markovian. The con-clusion is that, in addition to the relation (1), the char-acteristics of the diffusion process should be based on itsstochastic interpretation. We have presented a stochasticrandom walk model in which, if the particle enters themembrane, the waiting time for its jump has a heavy tail ψ M ( t ) ∼ /t α +1 / when t → ∞ , the waiting time for aparticle jump in the regions external to the membraneis the same as for normal diffusion. This tail is heavierthan the tail of distribution of waiting time for the par-ticle to jump ψ α ( t ) ∼ /t α +1 in a model providing thefractional subdiffusion equation Eq. (2). The function ψ M affects diffusion of a particle at only one point cor-responding to the position of the membrane, while thefunction ψ α affects particle diffusion at each point in thesystem. However, both determine the relation Eq. (1)with the same α in the long time limit. Thus, in thepresented model subdiffusion is generated by the effectof the long retention of the diffusing particle inside themembrane.(2) Possible application of the particle random walkmodel in a system with a subdiffusive thin membranecould be diffusion of antibiotic through a thin layer ofbacterial biofilm. The bacteria in the biofilm have manydefense mechanisms against the action of the antibiotic.One of them is the thickening of the biofilm which causesthat antibiotic particles can be trapped in the biofilm fora long time [21].(3) As an example, we have considered first and secondmoments that are power functions of time. However, theresults obtained in this paper can be applied to otherforms of the temporal evolution of the moments. Forexample, assuming that the functions ˆ v and ˆ z are slowlyvarying, we obtain the temporal evolution of the meansquare of the particle displacement which is characteristicfor slow subdiffusion (ultraslow diffusion), see [15, 16, 22].(4) The relations between the moments and the bound-ary conditions at the membrane has the following prop-erties. (a) When the Green’s function is continuousat the membrane, ˆΦ( s ) ≡
1, then ˆ v ( s ) ≡
0, see Eq.(40). Due to Eq. (43) there is h x ( t ) i = x . The sec-ond moment evolves over time according to the formula (cid:10) x ( t ) (cid:11) = L − [( x + 2 D ˆΞ) /s ]. (b) When the flux iscontinuous at the membrane, ˆΞ( s ) ≡
1, then Eq. (47)provides ˆ z = 2 D/s . Thus, the flux is continuous atthe membrane only if (cid:10) x ( t ) (cid:11) = x + 2 Dt . Due to Eq.(26), the probability of a particle becoming trapped inthe membrane is zero. Eq. (35) shows that η M ( s ) ≡ ψ M ( s ) ≡ ψ M ( t ) = δ ( t ). This means that evenwhen a particle enters the membrane, it will immediatelyleave it. In this case the first moment evolves in timeas long as the Green’s function is not continuous at themembrane, ˆΦ( s ) = 1. (c) When the probability density P and flux J are continuous at the membrane, ˆΦ( s ) ≡ s ) ≡
1, then in time domain we have h x ( t ) i = x and (cid:10) x ( t ) (cid:11) = x + 2 Dt . In this case we get the standardrelation for normal diffusion (cid:10) (∆ x ) ( t ) (cid:11) = 2 Dt . Thisresult is obvious as the continuity of the Green’s func-tion and flux means that there is no membrane effect onparticle diffusion. Acknowledgments
This paper was partially supported bythe Jan Kochanowski University under grantSMGR.RN.20.222.628.
Appendix I
The Laplace transforms of solutions to the diffusionequation with boundary conditions Eq. (12) readˆ P A ( x, s | x ) = 12 √ Ds e −| x − x | √ sD (76)+ A e ( x + x ) √ sD , ˆ P B ( x, s | x ) = B e − ( x − x ) √ sD . (77)From Eqs. (9), (16), (17), (76), and (77) we get thefollowing system of linear equations with respect to A and BA (cid:18) γ ( s ) − √ Dsγ ( s ) (cid:19) − B (cid:18) γ ( s ) + √ Dsγ ( s ) (cid:19) (78)= − (cid:18) γ ( s ) √ Ds + γ ( s ) (cid:19) ,A (cid:18) λ ( s ) − √ Dsλ ( s ) (cid:19) − B (cid:18) λ ( s ) + √ Dsλ ( s ) (cid:19) (79)= − (cid:18) λ ( s ) √ Ds + λ ( s ) (cid:19) . The determinants W ( s ), W A ( s ), and W B ( s ) for the sys-tem of equations (78) and (79) are given by Eqs. (20),(21), and (22), respectively. Solutions to Eqs. (78) and(79) A = W A ( s ) /W ( s ) and B = W B ( s ) /W ( s ) are uniqueonly if W ( s ) = 0. Under this condition the solutionsto diffusion equation are determined by the membraneboundary conditions uniquely. Comparing Eqs. (23) and(24) with (76) and (77), respectively, we get Eqs. (18)and (19) if A = ± / √ Ds . Since boundary conditionsdetermine the solution to diffusion equation uniquely, theequivalence of solutions (23), (24) and (76), (77) meansthe equivalence of the boundary conditions (10), (11) and(16), (17). If A = ± / √ Ds , from Eq. (76) we getˆ P A ( x, s | x ) = 12 √ Ds e −| x − x | √ sD (80) ± √ Ds e ( x + x ) √ sD . The + sign before the second term on the right–hand sideof Eq. (80) gives the Green’s function for a system withfully reflecting wall, in this case the boundary conditionat the membrane is J A (0 − , t | x ) = 0. The sign - givesthe Green’s function for a system with fully absorbingwall, the boundary condition is P A (0 − , t | x ) = 0. Inboth cases the diffusion is considered in region A only. Appendix II
We present how to get Eq. (27), here we use the no-tation as shown in Fig 2. Within the Continuous TimeRandom Walk model the Laplace transform of diffusionflux reads [16]ˆ J ( x, s | x ) = − ǫ s ˆ ψ − ˆ ψ ( s )) ∂ ˆ P ( x, s | x ) ∂x . (81)The mean number of particle jumps in the time inter-val [0 , t ] is h n ( t ) i = P ∞ n =1 nQ n ( t ), where Q n is theprobability that the particle jumps n times in the timeinterval. In terms of the Laplace transform we haveˆ Q n ( s ) = ˆ ψ n ( s )(1 − ˆ ψ ( s )) /s , then L [ h n ( t ) i ] = ˆ ψ ( s ) /s (1 − ˆ ψ ( s )). The frequency of particle jumps ν is defined as ν ( t ) = d h n ( t ) i /dt . Since h n (0) i = 0 we get ˆ ν ( s ) =0ˆ ψ ( s ) / (1 − ˆ ψ ( s )). Using the above formula and approxi-mating the derivative as ∂ ˆ P ( x, s | x ) /∂x = [ ˆ P ( x + , s | x ) − ˆ P ( x − , s | x )] /ǫ we define the probability flux by the unidi-rectional fluxes. The unidirectional flux J x − → x + controlsthe probability that a particle jumps from x − to x + ina time unit, similar interpretation is of J x + → x − whichcontrols a particle jump in the opposite direction. Fromthe above equations we obtainˆ J ( x, s | x ) = ˆ J x − → x + ( x − , s | x ) − ˆ J x + → x − ( x − , s | x ) , (82)where J x − → x + ( x − , s | x ) = ǫs ˆ ν ( s )2 ˆ P ( x − , s | x ) , (83) J x + → x − ( x + , s | x ) = ǫs ˆ ν ( s )2 ˆ P ( x + , s | x ) . (84) By adapting the above equations to the system presentedin Fig. 2, we change the particle jump frequency intofrequencies defined in the media a and b . We get J x − → x + ( x − , s | x ) = ǫs ˆ ν a ( s )2 ˆ P a ( x − , s | x ) , (85) J x + → x − ( x + , s | x ) = ǫs ˆ ν b ( s )2 ˆ P b ( x + , s | x ) , (86)where ˆ ν i ( s ) = ˆ ψ i ( s ) / (1 − ˆ ψ i ( s )), i ∈ { a, b } . From Eqs.(82), (85), and (86) we obtain Eq. (27). [1] J.P. Bouchaud and A. Georgies, Phys. Rep. , 127(1990).[2] R. Metzler and J. Klafter, Phys. Rep. , 1 (2000).[3] R. Metzler and J. Klafter, J. Phys. A , R161 (2004).[4] J. Klafter and I.M. Sokolov, First step in random walks.From tools to applications , (Oxford UP, NY, 2011).[5] T. Sandev, R. Metzler, and A. Chechkin, Fract. Calc.Appl. Analys. , 10 (2018); A. Chechkin, R. Gorenflo,and I.M. Sokolov, Phys. Rev. E , 046129 (2002); Frac.Calc. Appl. Anal. , 259 (2003); A. Chechkin, J. Klafter,and I.M. Sokolov, In: Fractional Dynamics: Recent Ad-vances , World Scientific, Singapore (2011); T. Sandev,I.M. Sokolov, R. Metzler, and A. Chechkin, Chaos Solit.Fract. , 210 (2017); C.H. Eab and S.C. Lim, Phys.Rev. E , 031136 (2011);[6] A. Compte, Phys. Rev. E , 4191 (1996).[7] T. D. Frank, Nonlinear Fokker-Planck Equations. Fun-damental and Applications , (Springer, Berlin, 2005); T.Koszto lowicz and K.D. Lewandowska, Phys. Rev. E ,021108 (2012).[8] E.K. Lenzi, R.S. Mendes, and C. Tsallis, Phys. Rev. E , 031104 (2003).[9] S.C. Lim and S.V. Muniandy, Phys. Rev. E , 021114(2002).[10] K.S. Fa and E.K. Lenzi, Phys. Rev. E , 012101 (2005).[11] T. Koszto lowicz, Phys. Rev. E , 022123 (2020).[12] V.G. Guimar˜aes, H.V. Ribeiro, Q. Li, L.R. Evangelista,E.K. Lenzi, and R.S. Zola, Soft Matter , 1658 (2015).[13] T. Zhang, B. Shi, Z. Guo, Z. Chai, and J. Lu, Phys.Rev. E , 016701 (2012); T. Koszto lowicz, K. Dworecki,and K.D. Lewandowska, ibid. , 021123 (2012); T.Koszto lowicz, Physica A , 285 (2001); D.K. Singhand A.R. Ray, J. Membr. Sci. , 107 (1999); Y.D. Kim,J. Y. Kim, H. K. Lee, and S. C. Kim, ibid.
69 (2001);R. Ash, ibid. , 9 (2004); S.M. Huang, M. Yang, W.-F.Zhong, and Y. Xu, ibid. , 8 (2013); A. Adrover, M.Giona, M. Grassi, R. Lapasin, and S. Pricl, ibid. ,7 (1996); M.J. Abdekhodaie, ibid. , 81 (2000); P.Taveira, A. Mendes, and C. Costa, ibid. , 123 (2003);M.I. Cabrera, J.A. Luna, and R.J.A. Grau, ibid. , 693 (2006); T. Koszto lowicz, ibid. , 492 (2008); I. Goy-chuk and P. H¨anggi, Phys. Rev. E , 051915 (2004);N. Korabel and E. Barkai, ibid. , 051113 (2011); Phys.Rev. Lett. , 170603 (2010); M.A. Lomholt, I.M. Zaid,and R. Metzler, ibid. , 200603 (2007); I.M. Zaid, M.A.Lomholt, and R. Metzler, Biophys. J. , 710 (2009);D.S. Grebenkov, J. Chem. Phys. , 104108 (2019); ibid. , 034104 (2010).[14] A. Bobrowski, Convergence of one–parameter operatorsemigroups in models of mathematical biology and else-where (Cambridge UP, 2016).[15] T. Koszto lowicz and A. Dutkiewicz, Math. Meth. Appl.Sci. , 10500 (2020).[16] T. Koszto lowicz, Phys. Rev. E , 022127 (2019).[17] T. Koszto lowicz, S. W¸asik, and K.D. Lewandowska,Phys. Rev. E , 010101(R) (2017); T. Koszto lowicz, ibid. , 022102 (2015); Int. J. Heat Mass Transf. ,1322 (2017).[18] E. Awoonor–Williams and Ch.N. Rowley, Biochim. Bio-phys. Acta , 1672 (2016); W. Shinoda ibid. ,2254 (2016).[19] T. Koszto lowicz, J. Phys. A , 10779 (2004).[20] B. Dybiec and E. Gudowska–Nowak, Phys. Rev. E ,061122 (2009).[21] T. Koszto lowicz and R. Metzler, Phys. Rev. E ,032408 (2020); T. Koszto lowicz, R. Metzler, S. W¸asik,and M. Arabski, PLoS One , e0243003 (2020).[22] A.V. Chechkin, J. Klafter, and I.M. Sokolov, Europhys.Lett. , 326 (2003); S.I. Denisov and H. Kantz, Phys.Rev. E , 041132 (2011); S.I, Denisov, S.B. Yuste, Yu.S.Bystrik, H. Kantz, and K. Lindenberg, ibid. , 061143(2011); R. Metzler, J.H. Jeon, A.G. Cherstvy, and E.Barkai, Phys. Chem. Chem. Phys. , 24128 (2014); L.P.Sanders, M.A. Lomholt, L. Lizana, K. Fogelmark, R.Metzler, and T. Abj¨ornsson, New J. Phys. , 113050(2014); A.S. Bodrova, A.V. Chechkin, A.G. Cherstvy,and R. Metzler, ibid. , 063038 (2015); A.V. Chechkin,H. Kantz, and R. Metzler, Eur. Phys. J. B90