Role of the Number of Adsorption Sites and Adsorption Dynamics of Diffusing Particles in a Conﬁned Liquid with Langmuir Kinetics

: In this work, we investigate the effect of the number of available adsorption sites for diffusing particles in a liquid conﬁned between walls where the adsorption (desorption) phenomena occur. We formulate and numerically solve a model for particles governed by Fickian’s law of diffusion, where the dynamics at the surfaces obey the Langmuir kinetic equation. The ratio between the available number of adsorption sites and the number of total particles are used as a control parameter. The investigation is carried out in terms of characteristic times of the system for different initial conﬁgurations, as well as the cases of identical or non-identical surfaces. We calculate the bulk and surface densities dynamics, as well as the variance of the system, and demonstrate that the number of sites affects the bulk, surface distributions, and diffusive regimes.


Introduction
The kinetics of diffusing particles in confined space with adsorption (desorption) by solid substrates represent an important class of problems that is readily applied from basic sciences to industrial separation processes [1][2][3][4][5][6][7]. One basic approach to model adsorptiondesorption phenomena is to use the simple, yet very powerful, Langmuir adsorption [8], often referred to as Langmuir kinetics [9,10]. At the same time, the diffusion process is described by Fick's law [11]. Indeed, in any practical application where any transport takes place, adsorption-desorption is likely to occur, even in very simple geometries [12,13]. However, from an analytical point of view, it is not easy to solve the coupled equations representing diffusion and Langmuir's kinetic altogether. In many cases, due to the nature of adsorbents, first-and second-order kinetic equations are employed to analyze experimental data [14][15][16][17]. Such equations account, for example, for diffusion within the adsorbing wall [18]. In other cases, a linearization process is performed to the Langmuir kinetic equation, in which the number of available sites for adsorption is much larger than the number of particles in the system [3,13,19,20]. Such an approach is convenient for studying the relation between diffusing particles in bulk samples adsorbed by confining walls.
Several systems of interest have been studied within the scope of the linearized Langmuir kinetic equation coupled to the diffusion equation. For example, it was used to study general aspects of bulk and surface dynamics [3], characterize diffusion in different geometries [12,13], to probe time-dependent diffusion coefficients [20], to insert memory effects and distinguish adsorption kind [21,22], to study adsorption and diffusion in systems with non-identical surfaces [23] and in systems with augmented surfaces [24], analyze adsorption in systems with space-dependent diffusion coefficient [25], and to study adsorption effects in electrolyte cells in the context of impedance spectroscopy [3,26]. However, in the above cases, the importance of the number of particles to the number of adsorbing sites ratio is neglected.
This article aims to provide a simple yet general model to study the diffusion of neutral particles in an isotropic liquid confined by adsorbing (desorbing) walls that obey Langmuir's kinetic. In this sense, a broad spectrum of phenomena arises by combining diffusion and a limited number of adsorbing sites, including how the bulk and surface dynamics take place and the effect of Langmuir's kinetic on diffusive regimes on the bulk. Furthermore, we demonstrate how memory effects can be included within Langmuir's kinetic equation to account for more general scenarios related to the number of particles, number of adsorbing sites, and the dependence of the previous state of the particle on the next one.

Model
Since we want to study diffusing particles diluted in a liquid in contact with adsorbing (desorbing) walls, we assume a sample in the shape of a slab where the only relevant direction is the z−axis. Such geometry is similar to experimental situations such as those found in liquid crystal displays, for example, where the adsorption and desorption of particles is known to play critical roles in the organization of the material [3]. The cell has thickness L, and the substrates are located at z = −L/2 (left side) and z = +L/2 (right side), as shown in Figure 1. Particles in bulk obey the diffusion equation, here assumed to be diluted so the Fickian approach can be used, that is ∂ρ ∂t where D is the diffusion coefficient, here assumed to be constant in space, and ρ(z, t) is the time-and space-dependent bulk density of particles. At the walls limiting the sample, we assume that for t = 0, all the particles are in the bulk (all adsorption sites are vacant), and that particles may be adsorbed-desorbed according to the well-known Langmuir's kinetic equation: In Equation (2), the subscript i = ±L/2 represents the set of parameters characterizing either the surface on the left side or the surface on the right side of the sample. Moreover, σ i (t) is the density of adsorbed particles and κ ai is a parameter representing the rate of adsorption, while κ di represents the rate (or time) of desorption. Furthermore, σ 0i is the number of adsorption sites available at surface i. The Langmuir kinetics assumes that the adsorption rate (dσ i /dt) is proportional to the difference between adsorption and desorption rates and that both rates are the same at equilibrium. Figure 1 shows the system studied here, including the occupied and free sites at the walls. Figure 1. System studied in this work. It consists of a liquid with diluted neutral particles (diffusing particles) that diffuse in the z−direction and may be adsorbed (desorbed) following Langmuir's equation. Notice that the left surface is located at z = −L/2, while the right surface is at Z = L/2. Gray spheres represent occupied sites on the surfaces, while green spheres represent free adsorption sites.
Notice that if we consider an infinity number of sites available for adsorption (σ 0i → ∞), we recover the linearized version previously used in different works [21]. Moreover, if the total density of particles is ρ 0 , and are initially in the bulk, i.e., ρ(z, t = 0) = ρ 0 , we can introduce the normalized quantities ρ R = ρ/ρ 0 and σ Ri = σ i /σ 0i , and, at equilibrium, arrive at the Langmuir isotherm: where α i = κ ai ρ 0 /(κ di σ 0i ) is a parameter governing the steady-state equation for each substrate. The set of Equations (1) and (2) are solved with the aid of the current density at the walls, that is, which implies that the number of particles must be conserved, that is To account for memory effects, we modify Equation (2). We introduce a kernel dependence as follows: where K(τ) is a kernel introduced to represent distinct scenarios related to the adsorptiondesorption phenomena. An equation such as (6) has been widely used to describe memory effects and non-Debye relaxations [21,22,27,28]. It has been used to describe chemisorption, physisorption, or a combination of both depending on the choice of K(τ) [21,22]. Although we introduced the kernel in the desorption term of Equation (2), it is important to stress that a memory effect represents that the previous state of the particle is important for the next state, so the kernel modifies the whole dynamics of the surface. Unfortunately, the process of solving Equations (2), (5) and (6) is difficult and does not have an analytical solution, so we employ a numerical method to solve it (see, for example, refs. [29,30]). We first introduce reduced quantities in such a way that Equation (1) becomes where ρ R = ρ/ρ 0 , Z = 2z/L, t * = 4t/τ D , and τ D = L 2 /D is the diffusion time. Regarding Equation (2), we first notice that the parameter 1/κ di has the dimension of time; thus, we write it as τ i (now i = ±1, which we call either l or r, meaning left or right side, respectively) from now on. Second, we have to choose the form of the kernel K(t) to proceed with the solution procedure. As proposed in [21], the form K(t) = 1/(τ ai )e −t/τ ai , which is a nonlocalized function of time, is an excellent choice to represent memory effects, such as may occur during multiple collisions of particles with the surfaces, in which energy is lost after each collision, and therefore the previous state of the particle is important in determining the next state. This memory time, represented by τ ai , reproduces the adsorption-desorption phenomena often found in physisorption or mixed processes [21]. If τ ai → 0, we recover the original kinetic equation, Equation (2), which is better suited to describe chemisorption processes. Numerically, it is more convenient to work with differential equations than with integral equations. We therefore take the time derivative on both sides of Equation (6) and apply the kernel K(t) = 1/(τ ai )e −t/τ ai to arrive at: where τ κi = L/2κ ai is the adsorption time and β i = ρ 0 L/σ 0i is a parameter that relates the number of particles available in bulk to be adsorbed to the number of sites available at the surfaces. Thus, if β i < 1, the specific surface has more sites than particles to be adsorbed, whereas if β i > 1, there are not enough sites at the surface for all the particles in bulk. The conservation of the number of particles, Equation (5), now becomes which also implies: 2 A summary of the main parameters and a brief description is given in Appendix A. Now, in order to solve the set of Equations (7)-(9) (or (10)), we employ a numerical method based on finite differences [31]. We use a mesh with n z points separated by a fixed distance δz = L/(n z − 1). All derivatives were taken using a central difference with a secondorder approximation, using ghost points at the borders. Thus, including the two walls, it results in a system of ordinary equations in t * , with n z + 4 equations. The time integration was performed with the 8th-order Dormand-Prince method implemented by the GSL library [32]. To ensure that the simulations were performed without instabilities, we monitored the conservation of the number of particles during each small increment in time by checking if Equation (9) is still satisfied. Through all simulations, we used n z = 500, and set the absolute error to be smaller than 10 −9 , allowing a free time step size, and retrieved data at every ∆t * = 10 −6 . Moreover, we assumed the initial surface density for both walls was zero (σ Ri (t * = 0) = 0). For the initial distribution in bulk, we used two different initial configurations: 1-particles are uniformly distributed across the cell (ρ R (Z, t * = 0) = 1), or 2-particles are initially concentrated in a plane set in the center of the sample, that is, the initial bulk distribution obeys a Dirac delta configuration. Finally, the characteristic times used in this work come from experimental data of similarly confined systems [3]. The adsorption parameter, κ ai , is usually in the order of 10 −6 ms −1 [3,33]. At the same time, the desorption time was estimated to be nearly 0.01 s for liquid crystal samples and close to 1 s for other isotropic samples [3,33]. A typical slab sample such as the one studied here is around 10 µm thick, while the diffusion coefficient is in the order of D ≈ 10 −11 m 2 s −1 , so τ D ≈ 10 s. Last, the memory time was estimated from experimental results to be τ a ≈ 1 s [21]. It is important to notice that, as long as the system obeys Fick's law and Langmuir's kinetic, it in principle can be investigated within the scope of this model.

Results and Discussion
We start by analyzing how bulk and surface distributions change over time as the characteristic parameters of the system change. In particular, we want to understand how the parameter β, which gives the ratio between particles available to be adsorbed and the number of sites available for adsorption, affects the bulk distribution of particles. Indeed, previous attempts to model adsorption-desorption with diffusing particles in a limited system [21] have used a kinetic equation that considers the number of sites to be infinity, so we can check here the importance of having a limited number of sites on the dynamics of the system. Furthermore, in our model, we can treat both surfaces as being non-identical; that is, each surface has its dynamics, so the model becomes more closely related to experimental situations (where it is difficult to assure both surfaces are completely identical) and to other cases where adsorption is present and substrates are non-identical, such as in fuel cells (hybrid microfluidic fuel cells, for example) [34], hybrid aligned liquid crystal cells [35], wetting layers [36], liquid crystals doped with dyes [37], and polymer adsorption in confined regions [38].
We start by showing a simple example where the role played by the parameter β can be clearly understood. For this first case, we consider both surfaces to be identical (so σ Rl = σ Rr , β l = β r and so on) and that the initial bulk distribution is uniform, i.e., ρ R (Z, t * = 0) = 1. Figure 2a reproduces the surface density (in Z = −1) vs. t * for τ D /τ l = 100, τ κl /τ l = 10 and for two values of τ al /τ l , that is, 0.1 or 20, meaning in the first case very short memory time and in the second case a long memory time, which results in the oscillations seen during the adsorption phenomena [21]. As discussed elsewhere [21], memory effects occur when the adsorption process depends on the previous state of the particle being adsorbed. For example, if a particle falls into an adsorbing well, it may be desorbed and eventually adsorbed again. Since the energy landscape changes after each process, the next phenomenon depends on the previous state, which in such a dynamic model it is represented by the parameter τ a [21,22]. We used two different values of β, 0.2 and 2, meaning in the first case there are five times more sites on the surfaces to adsorb than particles to be adsorbed, while in the second case, there are two times more particles in bulk to be adsorbed than available sites. Since τ κi < τ i , the desorption rate is larger than adsorption, so the equilibrium density is fairly low. However, in the case of β > 1, the surfaces reach a larger value, as expected (see Equation (3)). It is also interesting to note that the adsorption peak position (t * ), related to memory effects, also changes depending on β, which is indicative that memory effects are also sensitive to the ratio between the number of particles to the available number of adsorbing sites. To the best of our knowledge, this is the first time memory effects are considered in Langmuir kinetics. The inset shows the bulk density for t * = 0.6, indicating that for larger β there is a tendency of larger bulk density at all times. Figure 2b shows the left-side surface density for τ D /τ l = 100, τ κl /τ l = 1 and τ al /τ l = 0.1 and several values of β. Since the adsorption and desorption rates are the same, the surface coverage is more extensive and determined by the parameter β. For example, if the number of particles is ten times larger than the number of sites, the surfaces are nearly filled to maximum once the equilibrium is achieved. Clearly, the surface density for β = 10.0 reaches saturation quicker when compared with β = 0.2. We now start investigating the cases in which the surfaces are not identical. This is particularly interesting because one can see how one surface affects the other and the bulk distribution. Moreover, it is helpful to understand the effect of β on such distributions if the substrates are seen individually. For our study, we fixed the left substrate (Z = −1) with the following parameters: β l = 1.0, τ κl /τ l = 1 and τ al /τ l = 0.01. Furthermore, we kept the diffusion time the same for all the analyses, that is, τ D /τ l = 10.0. Figure 3a shows a case in which all the parameters are the same for both substrates, except the parameter β. We chose τ D /τ i = 10.0, τ κ i /τ i = 1.0, τ ai /τ i = 1.0, β l = 1.0, and varied β r . The substrate on the left side (parameters are always the same) is shown in the main figure, while the inset shows the substrate in Z = 1. As β r gets larger, the coverage in Z = 1 grows, which is expected since larger β represents more particles to be adsorbed compared with the number of sites on that specific surface. The left-side surface also saturates at higher concentration values when the parameter β of the rightside surface increases. Since the right side has fewer sites for adsorption (as β increases), the left-side surface has more particles to adsorb, since β l = 1.0. This trend can also be understood by Equation (9). If β r = 0.2, σ Rl = 1 − 5σ Rr − bulk/2, whereas for β r = 10.0, σ Rl = 1 − 0.1σ Rr − bulk/2, where "bulk" stands for the final concentration in the bulk. Since the number of particles in bulk is never larger than 1, σ Rl saturates at higher values as β r increases. In Figure 3b, we keep the surface in Z = −1 with the same parameters as in Figure 3a but used the parameters τ κr /τ r = 10, τ ar /τ r = 20 and β r = 0.2 and β r = 2.0. Notice that the adsorption time is ten times larger than the desorption time, so the overall coverage is low for the right surface. In this case, changing β r does not affect the left surface. Thus, for β r = 0.2, the surface has five times more sites than particles available, and the resulting coverage is small. If β r = 2.0, the right surface has room for only half of the particles, so the coverage is higher. Moreover, since the memory time is longer, the right surface displays oscillations in the adsorption-desorption process. The inset shows the bulk density for three different times, where the black curves represent β r = 0.2 and the red curves β r = 2.0. In the initial moments, both cases are identical. As time passes, for t * = 0.5, the bulk density near the right surface for β = 2.0 is slightly higher because the right surface has fewer sites to adsorb, so there is an excess of particles when compared with the case where β r = 0.2. As time increases, the bulk densities for the two values of β become very similar, being only slightly higher for β r = 0.2 due to the limited capability of the right surface to adsorb the particles.  To demonstrate that a larger β represents a lower overall density of particles compared with lower β cases, in Figure 4, we show the density of adsorbed particles for both surfaces normalized by their respective β. For this figure, we used the same parameters for the left surface as in Figure 3. For the right surface, we used τ κr /τ r = 0.01 and τ ar /τ r = 0.01, so the right surface has a high adsorption rate and neglectable memory effects. As expected, for both β r = 0.2 and β r = 2.0, the right surface reaches higher coverage values when compared with the left side, which is a consequence of the high adsorption rate. However, when β r = 2.0, both surfaces reach the equilibrium (the left coverage is still slowly decreasing for t * = 3.5) much more quickly when compared with β r = 0.2, which is due to the quick filling of the available sites when β r = 2.0. On the other hand, for β r = 0.2, the system is still evolving (filling the adsorption sites) when t * = 3.5. Figure 4 clearly shows that the equilibrium point depends on the number of available sites rather than just adsorption and desorption rates. The inset of Figure 4 shows the bulk dynamics three different times for both values of β. Initially, when most particles are still in bulk, the distribution for both cases is the same. For t * = 1.0, the distribution near the left side is the same, whereas, near the right surface, the bulk distribution is lower for β r = 0.2 when compared with β r = 2.0. As the system continues evolving, the case for β r = 2.0 quickly reaches an equilibrium value, which is larger when compared with β r = 0.2, where even for t * = 5.0, the bulk distribution is still not at equilibrium (which is also seen from the surface curves). Now, we use the mean square displacement (MSD) to probe how the diffusive regimes affect the number of sites at the surfaces. This is an important parameter because it is related to the spreading of particles across the cell, which characterizes the time-dependent diffusion coefficient, D(t). Thus, we are now interested in understanding how the number of sites, hence the parameter β, changes the way particles diffuse in bulk. Indeed, it is known that limiting surfaces may affect how particles diffuse in the bulk [24], which is of particular importance mainly for transport in living cells [39,40], but in general, the number of sites is neglected during the modeling. The MSD, which takes into account the conservation of particles (bulk distribution does not remain normalized at all times as particles leave the volume to be adsorbed at the walls), as in Equation (9), is given by In other words, the MSD is calculated as: (11)  In Figure 5, shows the MSD vs. time for different values of the parameter β for the situations in which the surfaces are identical and the case where the surfaces are nonidentical. For all cases, we used τ D /τ l = 100. For the case of identical surfaces, depicted by the solid lines, we used τ k /τ = 1 and τ a /τ = 20. Notice that when the MSD increases with time, the bulk distribution is spreading up, while MSD decreasing with time means that the distribution is shrinking; that is, particles are returning to the bulk. We first notice that for all the curves, the initial behavior of the MSD is equal, which corresponds to the initial spreading of particles from the center of the cell toward the surfaces. This initial diffusing process is not affected by the surfaces and corresponds to the usual diffusion that occurs for boundless samples. The black solid line represents β i = 0.2, while the solid red line shows β i = 5.0, and the solid blue curve depicts the case β i = 10.0. For β i = 0.2, after the initial free diffusion, the MSD decreases until t * ≈ 0.6. This is a consequence of the desorption process that takes place together with the adsorption process (τ κ = τ), which is favored by the small number of particles remaining in bulk. After this minimum occurs for the MSD, the adsorption process continues, favored by the larger number of sites compared with particles in bulk, until an equilibrium is reached. For β i = 5.0, this minimum happens much faster and is much less pronounced. At the same time, for β i = 10.0, it does not appear at all, which is a consequence of a large number of particles compared with the available sites, resulting in a much weaker desorption process of the high particle concentration in bulk at all times. The dashed lines in Figure 5a represent the same values of the parameter β but for non-identical walls. Here, the right surface is the same as used to produce the solid curves, only β r is changed, but the left surface is fixed at β l = 1.0, τ κl /τ l = 1 and τ al /τ l = 0.01. Interestingly, the bulk dynamic is also affected by the non-identical surface, not only by changing how the spreading and shrinking of the distribution occurs but also by changing the overall amount of particles to reach equilibrium in bulk, a direct consequence of the ratio between the number of particles and number of sites. Figure 5b shows the cases in which τ κ /τ = 0.01, when the surfaces are identical, and for the right surface (Z = 1) in the case of non-identical surfaces. The left surface for non-identical surfaces remains the same as Figure 5a. Here, since the rate of adsorption is much larger than the rate of desorption, it is expected that the surface quickly adsorbs most of the particles, leaving the bulk with a low concentration. This is true for β = 0.2, where the MSD is nearly zero for t * ≈ 2, meaning that the distribution continuously decreases where most of the particles are trapped at the surfaces. However, as the number of particles increases when compared with the number of available sites, this is no longer true, which significantly affects how diffusion takes place in the bulk, indicating that indeed the parameter β is essential, rather than just looking at the rates of adsorption and desorption as usually performed. This is true for both identical and non-identical cases. Figure 5b shows the cases in which τ κ /τ = 0.01, when the surfaces are identical, and for the right surface (Z = 1) in the case of non-identical surfaces. The left surface for non-identical surfaces remains the same as Figure 5a. Here, since the rate of adsorption is much larger than the rate of desorption, it is expected that the surface quickly adsorbs most of the particles, leaving the bulk with a low concentration. This is true for β = 0.2, where the MSD is nearly zero for t * ≈ 2, meaning that the distribution continuously decreases, where most of the particles are trapped at the surfaces. However, as the number of particles increases when compared with the number of available sites, this is no longer true, which significantly affects how diffusion takes place in the bulk, indicating that indeed the parameter β is essential, rather than just looking at the rates of adsorption and desorption as usually performed. This is true for both identical and non-identical cases. Notice how the dynamics of the MSD closely resemble the data observed for the diffusive characteristics observed in some systems that are essentially anomalous and with limited amount of adsorption sites to diffusing particles, such as in living cells [41] and in the diffusion of gold-labeled dioleoylPE in the plasma membrane of fetal rat skin keratinocyte cells [39]. Finally, noticing how the diffusive regimes change with the parameter β is interesting. Accordingly, (∆Z) 2 ∼ t * a , where the power a is related to how particles diffuse, that is, if a < 1, the diffusive regime is said to be subdiffusive, if a = 1, it is called usual, whereas for a > 1, the diffusive regime is called superdiffusive. Usually, confined systems are subjected to adsorption-desorption phenomena. However, the ratio β is not considered when studying the diffusive regime. We observe that in the initial moments before particles reach the surfaces, the diffusion is usual, that is, a = 1. Nonetheless, after interacting with the surfaces, the diffusion becomes essentially subdiffusive and is heavily affected by the parameter β. In Figure 5a, we show, as dotted lines, a few examples of exponents a for certain time intervals in which the distribution spreads after a quick desorption process. For example, for β = 0.2, a = 0.52, while for β = 5.0, a = 0.80, indicating that β changes the diffusion regime. As it is well-known, molecular crowding is one source of anomalous diffusion, especially in biological fluids [42], which, within this model, may be altered by changing the parameter β. This fact shall be further explored in future studies.
In conclusion, we modeled a system composed of an isotropic liquid limited by adsorbing surfaces where particles diffuse and may be adsorbed/desorbed following the complete Langmuir's kinetic equation. In the modeling process, we incorporated, in addition to the natural parameters arising from the Langmuir kinetic equation such as the ratios of adsorption, desorption, and number of available sites, memory effects that may occur depending on the inherent nature of adsorbing walls, such as occurring during chemisorption, physisorption, or a mixed process. We scaled all the variables in terms of characteristic times of the system and kept as main parameter the ratio between the available sites to the amount of particles in the bulk. It turns out that not only the dynamics but the diffusive regimes are heavily affected by such ratio, which is often neglected during adsorption dynamics. We hope our results may be helpful in describing separation processes and other systems, such as in living mater, where the limited amount of adsorption sites plays a crucial role.

Conflicts of Interest:
The authors declare no conflict of interest.

Appendix A. Table of Symbols
Here we present a table summarizing the main quantities used in this article (Table A1). Memory time t * = 4t/τ D Dimensionless time (∆Z) 2 Mean Square Displacement β = ρ 0 L/σ 0 Ratio of particles in the bulk to the number of available sites Subscript i May be l (left) or r (right). Indicates the substrate considered.