# The Form of Waiting Time Distributions of Continuous Time Random Walk in Dead-End Pores

Authors: Yusong Hou; Jianguo Jiang; J. Wu
Journal: Geofluids
 (2018)

---

## Abstract

Anomalous dispersion of solute in porous media can be explained by the power-law distribution of waiting time of solute particles. In this paper, we simulate the diffusion of nonreactive tracer in dead-end pores to explore the waiting time distributions. The distributions of waiting time in different dead-end pores show similar power-law decline at early time and transit to an exponential decline in the end. The transition time between these two decline modes increases with the lengths of dead-end pores. It is well known that power-law distributions of waiting time may lead to anomalous (non-Fickian) dispersion. Therefore, anomalous dispersion is highly dependent on the sizes of immobile zones. According to the power-law decline, we can directly get the power index from the structure of dead-end pores, which can be used to judge the anomalous degree of solute transport in advance.


---

## Body

## 1. Introduction

Anomalous dispersion widely exists in solute transport in ground water or in other porous media [1–7], which makes it difficult to simulate solute transport by traditional advection-dispersion equation. In recent years, scholars have developed a few models to describe the anomalous transport such as continuous time random walk (CTRW) [8–12] and fractional advection-dispersion equation [13–16]. The simulations by CTRW can agree well with the experimental data by fitting the transfer probability density function φx,t [2, 17]. Assume φx,t has the form φx,t=ϕxψt, where ψt denotes the waiting time distribution between two successive jumps. Corresponding to anomalous transport, ψt should follow a power-law decline, that is, ψt~t-1-α with 0<α<2 [17]. Many researchers have investigated the form of ψt by simulating the transport of solute in porous media [18–21], which may show power-law decline. However, fitting parameters from transport data in complex fluid fields are not suitable to uncover the physical mechanism of the power-law decline. An alternative way for the investigation of the form of ψt is to directly simulate the motion of solute particles in immobile zones. In this paper, we construct the dead-end pores as immobile zones and simulate the waiting time distributions to explore the modes of solute transport.ψ
t is more likely to show an exponential decline than a power-law decline because particles in stagnant zones are determined by diffusion equation. As many experiments show, an anomalous transport eventually converts to a Fickian process [2, 19]. Consequently, in the mathematical models, the power-law decline is often truncated arbitrarily by an exponential decline in the assumption of ψt [4, 9]. However, why does there exist a power-law decline of waiting time distribution? When does a power-law decline begin to convert to an exponential decline? In addition, it is not clear which factors determine the value of the power index α. To explore the general form of waiting time distribution function, we directly simulate ψt in dead-end pores. In this paper, we select three types of dead-end pores. By these dead-end pores with different sizes and shapes, we aim to get a universal form of ψt, which may show both the power-law part (corresponding to anomalous transport) and the exponential decline part (Fickian transport). In addition, from our simulations, we expect to get the index α which describes the anomalous degree.


## 2. Method of Simulation

In the framework of CTRW model for solute transport, the distribution of waiting time between two successive jumps determines whether the dispersion is anomalous or not. For simplicity, skipping over the displacement distributionϕx of jumps, we only investigate the waiting time distribution ψt. In our physical picture, waiting events happen when nonreactive particles are trapped in immobile zones. Dead-end pores are typical immobile zones. The flow velocity in a dead-end pore is zero and solute particles escape only by diffusion. The waiting time is the period from entering a dead-end pore to moving out of it. We consider three types of dead-end pores: straight tube, bent tube, and hemisphere as illustrated in Figure 1. These three dead-end pores lie above the flow tube which represents the mobile zone. Tracer particles on the interface between the flow tube and dead-end pores may enter the dead-end pore by diffusion. We begin to count waiting time when the particles pass through interfaces to dead-end pores and end counting when particles pass through interfaces to the flow tube.Figure 1
The schematic diagram of mobile zone and immobile zone. Three types of dead-end pores, that is, straight tube, bent tube, and hemisphere, represent the immobile zones.For simplicity, we use the surviving probabilityPt in the immobile zone to describe the waiting time distribution ψ(t); the relation between them can be written as(1)Pt=1-∫0tψτdτ=∫t∞ψτdτ,or Ψt=-dP(t)/dt. Compared to Pt, Ψt fluctuates more violently as the derivative enlarges the fluctuations of simulated values of Pt.P
t can be directly calculated by counting the number of particles in the dead-end pores and the particles will not be tracked when they move out of immobile zones. Then, ψt is obtained from the derivative of Pt. However, the probability Pt is very small when the time is large, especially if Pt decreases in an exponential form. By conventional particle tracing method with each particle of the same weight, Pt is proportional to the number of particles staying in the immobile zone. The number corresponding to small Pt is zero or fluctuates violently. To simulate Pt in a wide range, variable weight for tracing particles is adopted in this paper. In recent years the pruned-enriched method has been developed to simulate the random walk with variable weight [22–24]. We use the pruned-enriched method to adjust the weight according to Pt. At the beginning of simulation, each particle is released at the interface between the flow tube and dead-end pores with the same weight W0=1. Consider the case that the simulated particle is still in the dead-end pore at the time t=nΔt. Define r=Wt/Pt. If r>1, this particle is split in Intr (the integer part of r) copies and each copy continues to diffuse with a new weight Wnew=Wt/Intr. If r≤1, we continue to simulate the diffusion of the particle with probability r′=1/Int1/r and the new weight Wnew is adjusted to Wt/r′. In other words, the simulation of the particle is ceased with the probability 1-r′. Pt is updated in every step. By the scheme of variable weight, the smallest probability we could simulate can be reduced by many orders of magnitude. In the whole range of simulated time, the number of samples for each t has a nearly uniform distribution. Without the scheme of variable weights, the number of visited samples is proportional to Pt and decreases sharply with t, which leads to violent oscillation of simulated Pt.


## 3. Results and Discussion

For most heavy metal ions, the order of magnitude of diffusion coefficient is about 10−10 m2/s. In this paper, the diffusion coefficient of solute is set to be 4×10-10 m2/s. When tracer particles enter the dead-end pores through the interface, we start the time. When particles escape from the dead-end pores, the tracing process of this particle will be terminated. For simplicity, particles after escaping from a dead-end pore cannot be allowed to enter the same dead-end pore again because they drift in the flow tube (mobile zone). Firstly, we consider the surviving probability Pt in large dead-end pores. The lengths of straight tube and bent tube are both taken 10 mm and the radius of hemisphere is taken 15 mm. Figure 2 shows the power-law decrease of simulated Pt and ψ(t). The Pt and ψ(t) in bent tube and in hemisphere are reduced by the factor 1/2 and 1/4, respectively, which makes the comparison in Figure 2 clear without altering the power law. We take the values of power from the data of Pt because ψ(t) has more violent fluctuations. During the time interval 10 s–104 s, Pt decays in a power law, that is, Pt~t-α. Corresponding to straight tube, bent tube, and hemisphere, the value of α is 0.50, 0.49, and 0.52, respectively. They are close to 0.5. Thus, it is reasonable to suppose that the power α is approximately equal to 0.5 in these three types of dead-end pores.Figure 2
The power-law form of the surviving probabilityP(t) and the waiting time distribution Ψ(t) in dead-end pores at early time. The power α approximates to 0.5 for each type of dead-end pore.Let us discuss the physical meaning of the simulated results. This power-law form ofPt with α=0.5 results in anomalous dispersion, which can also be described by time-fractional advection-dispersion equation (tFADE). If the waiting time distribution has the form ψt~t-1-α, the concentration Cx,t can be written as [25](2)∂Cx,t∂t=Dt1-α0LCx,t,where Dt1-α0 is the fractional Riemann-Liouville operator; that is,(3)Dt1-α0ft=1Γα∂∂t∫0tdτfτt-τ1-α.Here, LCx,t is the operator representing dispersion. Therefore, the fractional order of time derivation is about 0.5 if anomalous dispersion is caused by solute trapping in long dead-end pores.The sizes of above dead-end pores are set large in the simulations. Let us consider the dead-end pores with short length. Tracer particles can jump out of them easily. Which form does surviving probability densityPt show? To investigate it, we set the length of straight tubes to be 1 mm, 2 mm, 3 mm, and 4 mm, respectively. The diffusion coefficient is the same as above. The simulated Pt and ψt are illustrated in Figure 3. At the early time, Pt and ψt in each straight tube decrease with the power law. As time increases, Pt and ψt deviate from the power-law decline. The shorter the length is, the earlier the deviation happens. Therefore, anomalous dispersion strongly depends on the sizes of dead-end pores. The duration of anomalous dispersion becomes longer as the length of dead-end pore increases.Figure 3
The deviation from the power-law decline ofP(t) and Ψ(t) in straight tubes of different length L as t increases. The shorter the length is, the faster the deviation happens.Now let us consider the form ofPt at large time in small dead-end pores. We also investigate three types of dead-end pores, that is, straight tube, bent tube, and hemisphere, but their sizes are shortened. The length of straight tube and bent tube is set to be 2 mm and 3 mm, respectively. The radius of hemisphere is set to be 3 mm. Using the particle tracing method with viable weight, we can simulate the surviving probability Pt with a high resolution. Figure 4 shows that log⁡Pt and log⁡ψt are linear with time t in every type of dead-end pore, which indicates the exponential decline when t is large. The exponential decline can be easily explained by immobile-mobile model [26] by neglecting the variability of concentration in immobile zones, which reads(4)∂Cimx,t∂t=-βCim-Cm,where Cim and Cm represent the concentration in dead-end pores (immobile zones) and in flow tube (mobile zones), respectively. Cimx,t is proportional to Pt. If Cm=0 at large t because of drift with water flow, Cim decreases as Cim~e-βt, that is, an exponential decline.Figure 4
The exponential decrease ofP(t) and Ψ(t) at large time for each type of dead-end pore.From simulations above, we can see that surviving probabilityPt decreases with time in a power law at early time and in an exponential form at large time. The power law with α≈0.5 corresponds to an anomalous dispersion, while the exponential decline corresponds to a Fickian dispersion. Therefore, it is important to determine the transition time Tt at which the decline transits from a power law to an exponential form. To get the transition time Tt, we first calculate the exponentβin (4). Tt has a close relation with the reciprocal of β. For simplicity, we use a straight tube as the dead-end pore. The diffusion equation about the concentration Cx,t in the tube can be written as(5)∂Cx,t∂t=Dm∂2Cx,t∂x20≤x≤L,where L is the length of the straight tube and Pt=∫0LCx,tdx. The interface between the straight tube and the mobile zone is located at x=0. As tracer particles move out of the straight tube, they are not allowed to enter again because of drift with flow. Consequently, the boundary condition at x=0 can be set as C0,t=0. The condition of reflecting boundary at x=L is written as ∂Cx,t/∂xx=L=0. By the method of variable separation, Cx,t has the form Cx,t=Xxe-βt and we can get the minimum of β; that is,(6)βmin=π2Dm4L2.Other eigenvalues of β are much larger than βmin and the modes with large β decrease more rapidly than the mode with βmin. Neglecting the modes with large β, β in Cx,t=Xxe-βt approximates to βmin. We fit the β of straight tubes of different length according to Pt~e-βt and compare them with the values of βmin calculated from (6). Table 1 shows the values of β and βmin. The two sets of values are very close. In hemisphere dead-end pores, β can also be approximated by βmin as the tube despite replacing L with R and adding a coefficient.Table 1
The comparison between theβ fitted by Pt~e-βt and the βmin of (6).
L (mm)
β (s−1)
β
m
i
n (s−1)
β
/
β
m
i
n
1
9.72
×
1
0
-
4
9.87
×
1
0
-
4
0.98
2
2.44
×
1
0
-
4
2.47
×
1
0
-
4
0.99
3
1.08
×
1
0
-
4
1.10
×
1
0
-
4
0.98
4
6.10
×
1
0
-
4
6.20
×
1
0
-
4
0.99Now let us discuss the transition timeTt at which the decline of Pt transits from a power-law form to an exponential decline form. In Figure 3, Tt is the time when the curve deviates from the straight line. Comparing the transition time Tt to 1/β in Table 1, we can find that Ttis located in the range of time 1/2β,1/β. Since β is very close to βmin (illustrated in Table 1), Tt is approximately proportional to the square of length L and inversely proportional to Dm as shown in (6). As Tt represents the transition time from power-law decline to exponential decline, the enduring time of anomalous dispersion has the same relation with the length L and Dm as that in (6).


## 4. Conclusions

In this paper, we investigate the waiting time distributionψt of nonreactive solute particles in immobile zones, which are represented by dead-end pores. Three types of dead-end pores, straight tube, bent tube, and hemisphere, are simulated. When a dead-end pore is large enough, ψt follows a power-law distribution, that is, ψt~t-1-α. All the simulated dead-end pores in this paper show that the values of α are very close to 0.5. According to CTRW theory, this power-law distribution results in anomalous dispersion [17, 25]. If dead-end pore is finite, ψt follows an exponential decrease finally. Therefore, the transport is anomalous at early time and is Fickian at later time. The transition time Tt is proportional to the square of the length of dead-end pores and inversely proportional to the diffusion coefficient of solute. From CTRW theory, we could use Tt to sign the transition time from anomalous transport to Fickian transport. Thus, we can draw the conclusion that anomalous transport phenomena will be more easily observed in the porous media with large immobile zones.In practice, dead-end pores may be very long, such as the fractures in fracture porous media. In this case, the waiting time distribution can seem as a power-law distribution during the whole observation period. Therefore, anomalous dispersion lasts for the whole observation period. The power-law distribution ofψt has a close relation with fractional dispersion equation [25]. According to our simulations, the fractional order of time derivative approximates to 0.5 when using tFADE (i.e., (2)) to describe this class of anomalous dispersion.



---

*Source: 8329406-2018-01-08.xml*