Stochastic modeling of non-linear adsorption with Gaussian kernel density estimators

Adsorption is a relevant process in many fields, such as product manufacturing or pollution remediation in porous materials. Adsorption takes place at the molecular scale, amenable to be modeled by Lagrangian numerical methods. We have proposed a chemical diffusion-reaction model for the simulation of adsorption, based on the combination of a random walk particle tracking method involving the use of Gaussian Kernel Density Estimators. The main feature of the proposed model is that it can effectively reproduce the nonlinear behavior characteristic of the Langmuir and Freundlich isotherms. In the former, it is enough to add a finite number of sorption sites of homogeneous sorption properties, and to set the process as the combination of the forward and the backward reactions, each one of them with a prespecified reaction rate. To model the Freundlich isotherm instead, typical of low to intermediate range of solute concentrations, there is a need to assign a different equilibrium constant to each specific sorption site, provided they are all drawn from a truncated power-law distribution. Both nonlinear models can be combined in a single framework to obtain a typical observed behavior for a wide range of concentration values.



There are no comments yet.


page 1

page 2

page 3

page 4


An Analytical Model for Molecular Communication over a Non-linear Reaction-Diffusion Medium

One of the main challenges in diffusion-based molecular communication is...

Highly scalable numerical simulation of coupled reaction-diffusion systems with moving interfaces

A combination of reaction-diffusion models with moving-boundary problems...

Modeling Firn Density through Spatially Varying Smoothed Arrhenius Regression

Scientists use firn (compacted snow) density models as a function of dep...

Deep Adversarial Koopman Model for Reaction-Diffusion systems

Reaction-diffusion systems are ubiquitous in nature and in engineering a...

Analysis of Flux Corrected Transport Schemes for Evolutionary Convection-Diffusion-Reaction Equations

We report in this paper the analysis for the linear and nonlinear versio...

Correlated pseudo-marginal schemes for time-discretised stochastic kinetic models

Performing fully Bayesian inference for the reaction rate constants gove...

Unsupervised Machine Learning Based on Non-Negative Tensor Factorization for Analyzing Reactive-Mixing

Analysis of reactive-diffusion simulations requires a large number of in...
This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.

1 Introduction

Adsorption is defined as the binding of atoms or molecules from a gas or liquid to a surface. It is a phenomenon well-described in many physical, biological, and chemical systems and processes, and that has been widely employed in industrial applications such as pharmaceutical industry, chillers and air conditioning systems, water purification, coatings, and resins, to name a few. To design and optimize an adsorption-based process, it is necessary to characterize accurately the adsorption equilibria and their dependence on the experimental conditions.

Equilibrium relations are described by adsorption isotherms, relating the equilibrium concentration of a solute on the surface of an adsorbent to the concentration of the solute in the liquid/gas being in contact. In 1916, Langmuir introduced the first scientifically based nonlinear isotherm by assuming a homogeneous surface with a specific number of sites where the solute molecules could be adsorbed (Langmuir, 1918). Furthermore, he assumed that the adsorption involves the attachment of only one layer of molecules to the surface, i.e. mono-layer adsorption. However, the Langmuir model deviates from the experimental observations in the presence of a rough inhomogeneous surface where multiple site-types or layers are available for adsorption and some parameters vary from site to site. This problem was tackled by Freundlich who proposed the first mathematical fit to a nonlinear isotherm, leading to a purely empirical formula for adsorption on heterogeneous surfaces (Freundlich, 1906). The Freundlich isotherm was established by assuming that adsorption varies directly with pressure without reaching saturation, while experimentally the rate of adsorption saturates by applying very high pressures. Therefore, the use of the Freundlich isotherm is appropriate when dealing with low/medium pressures/concentrations.

The Langmuir and Freundlich models are the two most commonly used isotherms due to their simplicity and their ability to properly fit a variety of adsorption experimental data. Several isotherms have been also proposed to combine the features of both the Langmuir and Freundlich isotherms, including the BET (Brunauer et al., 1938), Sips(Sips, 1948) and Redlich-Peterson (Redlich and Peterson, 1959) models. Apart from these classical isotherms, other thermodynamically consistent adsorption models can be derived from a fundamental integral equation relating the experimental isotherm, the adsorption energy distribution, and the local isotherm, see Quinonesa and Guiochon (1996) for a review of these models.

A major challenge to control and predict an adsorption-based process is the heterogeneity and complexity of interactions between the adsorbate and adsorbent surfaces in a dynamic solid/fluid or solid/gas system. This complexity is the main reason behind proposing various isotherm models to predict an adsorption process based on macroscopic experimental data. To tackle this complexity, computational models have been developed by incorporating the interaction of the reactants into the dynamic system at the molecular level. A number of these models employ Density Functional Theory (DFT) methods, force-field techniques and Molecular Dynamics (MD) methods to simulate adsorption at the molecular level, see Costa et al. (2015) for a recent review. These models are ideal to study adsorption close to the surface; however, due to a high computational cost, they cannot be employed for large simulation domains, i.e., at the pore scale or above.

To span over a wide range of scales, continuum models based on a diffusion-reaction equation have been developed to study reactive solute transport driven by diffusion. This classical problem can be addressed using both Eulerian and Lagrangian approaches, taking into account their limitations and advantages. In an Eulerian approach, the problem is defined in terms of reactants concentration which is used to describe the reaction rate by macroscopic laws. A macroscopic mass balance equation such as the diffusion-reaction equation (DRE) is then obtained in terms of the concentration of a solute that is in instantaneous equilibrium at each point in space with the adsorbed concentration [S], both expressed as mass per unit volume of solute. The system is then formulated in terms of two coupled equations:


where is the diffusion coefficient. Alternative expressions for Eq. (1) involve defined as mass of sorbed species per unit mass of solid. These two equations can be combined to write down a single one


Except in the particular case of Eq. (2

) being a linear relationship, the system results in a nonlinear partial differential equation, thus limiting the usefulness of Eulerian approaches to tackle the mathematical problem. An approach to simulate nonlinear sorption with particle tracking based on Eq. (

3) has been presented by Tompson (1993). This author estimated by reconstructing concentrations from particles at each time step without controlling the statistical errors involved. In this context, the method presented here does not require the computation of solute concentrations during the course of the simulation, therefore the statistical errors cannot propagate overtime.

In the presence of chemical heterogeneity, small enough volumes are required to capture the reactions while a macroscopic quantity such as concentration (mass per unit volume) loses its integrity in this limit. An alternative approach is to address this problem at the molecular level within a Lagrangian framework, where the movement of each individual molecule is tracked (Gillespie, 1976, 1977, 2000). However, since an extremely large number of molecules can exist even at very small concentrations, this approach is computationally unfeasible for simulating reactive transport in porous media. Moreover, at the molecular level, the continuum assumption is no longer valid and chemical kinetics should be derived by molecular collision.

Particle Tracking Methods (PTM) offer a practical and efficient alternative to overcome these problems by combining the key features of both approaches. In this case, the solute plume is divided into manageable number of particles, typically ranging between . This way, each particle does not represent a molecule itself but a certain fraction of mass containing numerous molecules. This advantage comes with the need to translate reaction-rate equations into particle relationships. In non-linear reactive transport problems, the latter includes the chemical interaction between particles, which is normally defined based on the particle area of influence. A concept that has different interpretations among researchers.

Kernel Density Estimators (KDE) are typically used in statistics as a non-parametric approach to estimate probability distribution function from a finite data sample. Its application in Lagrangian reactive transport problems has several advantages. Since it is non-parametric, KDE allows the identification of complex solute plume distributions (multimodal or non-Fickian behaviors). It provides not only adequate estimates of concentrations but also of their functionals (e.g. mixing, Human-health risk). And more importantly, it provides an adequate mathematical framework to select an optimal choice of the particle area of influence

(Fernàndez and Sanchez, 2011). This parameter will essentially dictate if particle of appropriate kinds are mutually in contact for the occurrence of chemical reactions. In this context, Rahbaralam et al. (2015) has demonstrated that this approach avoids the segregation of particles and the resulting incomplete mixing which is common in classical diffusion-based models. A variety of kernel functions can be used, among them, Gaussian kernels are usually preferred for mathematical advantages (Pedretti and Fernàndez, 2013).

Motivated by this potential, here we extend the Gaussian KDE model to simulate nonlinear adsorption. We show that the model proposed is able to reproduce the results of the Langmuir and Freundlich isotherms and to combine the features of these two classical models. This approach opens up a new way to predict and control an adsorption-based process using a particle-based method with a finite number of particles. The paper is structured as follows. First, Section 2 sets out the background, the problem and the numerical approach. Section 3 introduces the proposed adsorption model using PTM and Gaussian KDEs. Simulation results and discussion are presented in Section 4. Finally, Section 5 provides a summary of the main contributions of this paper.

2 Background and statement of the problem

2.1 The Langmuir isotherm

The Langmuir isotherm explains adsorption by assuming that the adsorbent is an ideal solid surface composed of series of identical sites capable of binding the adsorbate. This binding is treated as a chemical reaction between the adsorbate molecule and an empty site, . This reaction yields an adsorbed complex . This process can be reversed through desorption whereby the adsorbed molecule is released from the surface and the complex is transformed to and . This dynamic equilibrium existing between the adsorbate and the adsorbent can be expressed as


This model assumes adsorption and desorption as being elementary processes, where the rate of forward adsorption and the rate of backward desorption are given by:


where is the forward adsorption reaction constant, is the backward desorption reaction constant, and denotes the concentration of species (, , or ). At equilibrium, the rate of adsorption equals the rate of desorption, i.e. , then by rearranging the terms we obtain


where is the equilibrium constant. By adding up the concentration of free sites and of occupied sites , the concentration of all sites , assumed constant in time, is obtained as . Combining this relation and Eq. (7) yields the Langmuir adsorption isotherm:


In the presence of a high concentration of the adsorbate , Eq. (8) leads to the saturation of surface sites, . In other words, the surface reaches a saturation point where the maximum adsorption capacity of the surface will be achieved.

2.2 The Freundlich isotherm

The Freundlich isotherm is an empirical model, which is commonly used to describe the adsorption performance of heterogeneous surfaces. This isotherm is mathematically expressed as


where and are called the Freundlich constants. The constant is an adsorption coefficient, while is a measure of the deviation from the linearity of the adsorption. Unlike in the Langmuir model, in this one there is no adsorption maximum or saturation. Equation (9) can be linearized as


where represents the slope of the line in a plot.

While the Freundlich model has been found to fit most existing adsorption experimental data (Frankenburg, 1944; Davis, 1946; Thomas, 1957; Urano et al., 1981), its theoretical basis is under scrutiny. In principle, the Freundlich isotherm can be derived from the fundamental integral equation for the overall adsorption isotherm stated as


where is the number of sites having adsorption energy and is the local coverage of individual sites. The integration region is over all possible adsorption energies. It has been shown that Eq. (11

) leads to the Freundlich isotherm when an exponential distribution of adsorption energies is assumed

(Sips, 1948; Sheindorf et al., 1981):


where and are constants, and is the temperature (assumed constant as implied by the word isotherm). Furthermore, it is assumed that the local coverage follows the Langmuir isotherm as


with the equilibrium constant depending on the adsorption energy as


Plugging Eqs. (12) - (14) into Eq. (11) and assuming that yields


This integral has a solution in the form of the Freundlich isotherm given by


with the Freundlich constant . Therefore, Eqs. (11)-(16) establish a theoretical basis for the Freundlich isotherm. In Section 3, we demonstrate that this approach, which is based on an exponential distribution of adsorption energy, is mathematically equivalent to consider

, in the Langmuir model, as a random variable that follows a truncated power-law distribution.

2.3 The KDE-based diffusion-reaction model

This section presents a summary of our recently proposed chemical reaction model using KDEs (Rahbaralam et al., 2015). This model has been developed based on a single forward bimolecular irreversible reaction , where and represent two species in the dissolved phase, reacting kinetically with unitary stoichiometric reaction coefficients. A diffusion-reaction equation governs the transport of the reactants:


where represents , the concentration of species is indicated by , is the reaction rate constant, and is the diffusion coefficient.

To simulate this problem, the random walk particle tracking method is used whereby an equal fraction of the total mass is carried by each particle. For simplicity and without loss of generality, the 1-D form of the reactive-advective-dispersive equation is considered. Given as the particle number, the locations of the A and B particles, and

, are initialised using a statistical uniform distribution. Then, a Brownian random walk motion governs the diffusion of the species during a time interval between

and +, formally written as


where is the flow velocity,

is the variance, and

is a normally distributed random number drawn at each time step.

As stated in the introduction, Kernels provide an influential area around each particle. This area is controlled by a parameter whose optimal value is obtained based on minimising the variance error while maintaining smoothness. The optimum bandwidth is found as (Park and Marron, 1990)


where the value of is time dependent and can be formally derived from the second derivative of the concentration spatial function. So, from Eq. (19), the optimal bandwidth size increases inversely proportional to the number of particles. This feature is particularly beneficial to avoid incomplete mixing due to the segregation of particles during chemical reactions.

In a Reactive Particle Tracking (RPT) method, the reaction of particles is simulated through probabilistic rules,


where is the particle mass and is the (forward) probability that two particles, and , separated by a distance react within the time interval

. The co-location probability density function (pdf) is represented by

which defines the probability that two particles separated by a distance occupy the same position after a time interval

. By attributing a Gaussian density with a standard deviation

to each particle, the convolution of two particles density results in the co-location pdf, representing the reaction zone of the two particles. This probabilistic rule was demonstrated by Benson and Meerschaert (2008). The fundamental difference in the KDE model is that the area of influence, , around a particle is not only attributed to diffusion but also depends on the distribution of particles (shape and number of particles). This dynamic area of influence changes the probability of reaction at each time step, avoiding the formation of segregated areas of particles. Based on the principles of the law of mass action, the probability of forward reaction for the KDE-based model is obtained as (Rahbaralam et al., 2015)


This equation implies that by increasing the area of influential of each particle, , the probability of reaction decreases while the unitary area under the co-location pdf is preserved. Such expansion, in turn, would increase the number of potential reactive pairs.

3 The KDE-based adsorption model

The probability of forward reaction provided in Eq. (21) was originally proposed to model a simple chemical model, i.e., the single forward bimolecular irreversible reaction. To extend this approach for modeling adsorption, is considered as the probability of forward adsorption reaction in Eq. (4) with the particle as the product of the reaction. In contrast to the bimolecular irreversible reaction, an adsorption reaction is reversible, so that a particle can transform back to a particle and another one . The probability of backward desorption is given by


To numerically implement the probability of forward and backward reactions in Eqs. (21) and  (22), we use the following approach. In the forward case, we first computed from Eq. (21), and then we compared this value with a random number generated from a uniform distribution . Then, if , the reaction was supposed to have occurred, and both and particles were removed from the system and substituted by a particle; otherwise, the reaction was not supposed to have taken place at that particular time step and both particles were kept. This procedure was repeated for every pair of particles (A,B) in the system.

A reversed procedure was repeated for the backward reaction by comparing a random number with the probability of backward reaction obtained from Eq. (22). If , the reaction was supposed to have occurred, the particle was removed and substituted by an and a particles located at the same point was originally considered; otherwise, no action was taken, representing that the reaction had not occurred. After the loop for all particles were performed, the simulation continued to the next time step. We also assumed that the locations of the adsorbent particles and the reaction product particles were fixed (the location of the sorption sites did not change with time), so that Eq. (18) was only applied to the displacement of the adsorbate particles . At each time step, the number of remaining particles multiplied by and divided by the volume resulted in the concentrations and . The simulation was carried out until equilibrium was clearly achieved, as indicated by the stabilization in the ratio of concentrations .

As explained later, we set up two simulation runs to assess the capacity of the model to simulate adsorption based on both the Langmuir and the Freundlich isotherms. In order to account for the former, it is just enough to set a finite number of sites and then fix the and values to honor the relationship set in Eq. (7). In the case of the Freundlich isotherm, it is necessary to incorporate an exponential distribution of adsorption energies, as given in Eq. (12), and to set the individual value of the equilibrium constant at each site. The generation of the equilibrium constant requires several considerations. As indicated in Section 2.2, the assumption of in Eq. (15) is necessary for the validation of the Freundlich model in Eq. (16). However, several authors (Hill, 1949; DaRocha et al., 1949) have indicated that the integration limits over all possible adsorption energies in Eq. (11) should have a minimum cut-off to have a physical significance. To include this in our model, we consider the site-energy distribution given by Eq. (12

) as a truncated exponential distribution with cumulative distribution function given as


where is the minimum energy bound. Considering Eq. (14), the cumulative distribution is obtained as a function of the equilibrium constant ,


where . Notice that Eq. (24) is equivalent to a truncated power-law function, given as


with the corresponding probability density function,


Therefore, the generation of values can be determined directly from the cumulative distribution function of . Given to be a random value of a uniform distribution between 0 and 1, the corresponding value is obtained as


Rewriting Eq. (11) in terms of the probability density function of equilibrium constant and the concentration of the adsorbed chemical compound , we have


Introducing Eq. (13) and Eq. (26) into Eq. (28) results in


This integral has two limiting solutions depending on the concentration of . When is small, the solution is obtained by applying a change of variable in the integration as , so that


Now, we can see that the lower limit of integration approaches to zero and the solution is equivalent to


Thus, when , the concentration of the adsorbed chemical compound follows the Freundlich isotherm written as




On the contrary, when is large, Eq. (30) simplifies to


which gives


Thus, at large concentration of , the solution approaches the solution limit as given by the Langmuir isotherm. In short, a truncated power law distribution of equilibrium constant gives a sorption isotherm that follows the Freundlich model for small concentrations and changes to the Langmuir saturation limit for large concentrations. To understand the transition between the two limiting models, one can analyze the deviation of Eq. (29) from the Freundlich model. The relative deviation from the Freundlich model can be written as


To find an approximate solution of the relative deviation we expand the integral around =0 in an ascending series in


This is a relatively good approximation given that is small and the integral limits are very close to each other. Truncating the series expansion at the first term results in


where is the concentration of at which the isotherm starts deviating from the Freundlich model.

The random walk model of sorption proposed here has three parameters , , and . The parameter represents directly the exponent of the Freundlich isotherm. The parameter determines the maximum concentration of the chemical compound in the liquid phase for which sorption follows the Freundlich isotherm. The parameter is the concentration of sorption sites in the system which relates to the Freundlich coefficient according to Eq. (33).

4 Numerical simulations

4.1 The Langmuir Model

The first test was design to show the performance of the proposed model to reproduce adsorption based on the Langmuir model. The constants of forward adsorption and background desorption were set to and , respectively. A 1-D domain of size , a particle mass of , and initial concentrations of and are considered. Simulations are performed for 2000 time steps with a time interval of and a diffusion coefficient of . All values throughout the paper are reported in normalized units.

At each time step, all the adsorbate particles are moved following Eq. (18) to account for the effect of diffusion. To effectively search for potentially reactive pairs of particles and , the 1-D domain is divided by elements of size , then we find at each element and its two neighbors, all pairs of particles and . The distance between each pair is first obtained, then the procedure mentioned below Eq. (22) is followed to implement the forward reaction. After repeating this procedure for every pair of particles in all elements, the backward reaction for each particle is implemented following a similar procedure mentioned in Section 3.

Figure 1 shows the concentration ratio as a function of time for different initial concentration values . It can be observed that all the simulations approach equilibrium at late times. Fluctuations around the equilibrium value decrease as increases. For a concentration of , which corresponds to an initial 40,000 particles, the actual equilibrium constant is closely approximated by the model at late times while the number of remaining adsorbate particles in the system is very low, in the order of 1000 particles. This reasonable approximation shows the potential of our proposed model to simulate adsorption with a low number of particles. We note that the classical diffusion-based model is unable to reach this approximation due to the segregation of particles at late times (Rahbaralam et al., 2015).

Figure 1: Concentration ratio as a function of time for different initial concentration values . At equilibrium this ratio should be equal to , which is plotted for comparison purposes.

To further examine the performance of the proposed model, we performed twenty simulations considering different initial concentration values = 40, 50, …, 250. For each simulation, we obtained the concentrations and from the average of their values at the last 100 time steps (so that equilibrium was achieved). Figure 2 shows these concentrations together with the results of the analytical model in Eq. (8) for comparison purposes. It is clear that the simulation results have a good agreement with the analytical Langmuir isotherm even for a concentration of , where the initial number of particles is 8000 and the number of remaining adsorbate particles is below 100. By increasing the concentration , the product concentration tends to saturate towards which is the maximum capacity of free sites for adsorption. This saturation level implies the mono-layer adsorption in the Langmuir isotherm.

Figure 2: Concentration of the adsorbed species as a function of adsorbate concentration at equilibrium. Analytical solution of the Langmuir isotherm in Eq. (8) is also plotted for comparison purposes.

4.2 The Freundlich Model

This is a challenging model due to the need to specify different equilibrium constants to each individual sorption site. We first set up all the parameters, except , identical to that of the previous section, to be consistent with the simulations for the Langmuir model. Then, according to the theoretical developments for deriving Eq. (27), the constant of forward adsorption can be obtained in terms of the Freundlich constant for each pairs of and particles by , where in the minimum equilibrium constant given by Eq. (38) considering a relative deviation = and is a random number generated from a standard uniform distribution. Thus, different values were used in the simulations. Having these parameters, we followed the same procedures explained in the previous section to implement the forward and backward reactions.

To obtain a good agreement between the model and the theory, the sorption sites should follow the exponential distribution of energy in Eq. (12). However, in a range of low concentrations, where the number of particles is small, there is a risk of undersampling. For this reason, we performed a larger number of simulations in this range. Figure 3 shows the adsorbate concentration and the product concentration at equilibrium. It is clear that, in a range of low concentrations, the points in the graph follow closely a linear logarithmic slope whose value is controlled by the Freundlich constant . By increasing the value of , this agreement holds for a smaller range of concentrations .

Figure 3: Concentration as a function of concentration at equilibrium. Logarithmic slope of the solid lines correspond to the Freundlich constant values used in the simulations. The dashed line shows the saturation concentration .

As mentioned in the introduction, a drawback of the Freundlich isotherm is that it cannot be employed for large adsorbate concentrations since this model does not taken into account the saturation limit of free adsorption sites. However, Fig. 3 indicates that the proposed model is also able to incorporate this physical limitation without the need to use a complex model that transitions to a Langmuir-type behavior. Again in Figure 3 it is observed that as the adsorbate concentration increases, the number of available adsorbent sites decreases and the product concentration tends to saturate to the maximum adsorbent concentration , in agreement with Eq. (35). This saturation occurs at lower adsorbate concentration by increasing the value of . Therefore, these results show that our model combines the features of both the Langmuir and Freundlich isotherms in all range of adsorbate concentrations.

5 Conclusions

Nonlinear adsorption in diffusion-reaction problems is a challenging problem. Eulerian methods lead to nonlinear partial differential equations, and so there is a need to develop efficient Lagrangian methods to tackle such a problem. While some methods have been proposed in the literature to address the Langmuir model, none so far is capable of addressing the Freundlich model. In this work we proposed a numerical method that combines a simple Particle Tracking Methods (PTM) with some predefined rules for particle interaction that properly reproduces nonlinear adsorption. The method uses Gaussian Kernel Density Estimators (KDEs), enabling it to avoid the effects of incomplete mixing due to the segregation of particles which is common in classical diffusion-based models with finite (and small) number of particles.

For the Langmuir model, the method involves writing the adsorption process as the combination of a forward and a backward reaction. The former relies on the concept of particle distance, where one of the particles corresponds to the adsorbate (mobile) and the second one indicates a free site for adsorption (immobile). The probability of reaction is based on a function that involves the distance between both particles, and at each time step a random function is drawn to see if the reaction had taken place. The backward reaction is first-order. It was found that the PTM-KDE method provides a good reproduction of nonlinear adsorption following the Langmuir model by simply adding a finite number of sorption sites of homogeneous sorption properties. The advantage of using KDEs is that it is possible to obtain very good results in terms of adsorbed versus adsorbate concentrations with quite a small number of particles, while other PTM methods would rely on a very large number of tracked particles to obtain good approximations. The system quickly gets to equilibrium and it is very stable in terms of the very low number of particles that can be used. The effect of the nonlinearity resulted in a different equilibrium time depending on the initial concentration of the dissolved species (adsorbate).

The approach was then extended to address nonlinear sorption modeled with a Freundlich isotherm. This involves some theoretical considerations regarding heterogeneity of the sorption properties at each specific sorption site. When such properties follow a truncated exponential statistical distribution, adsorption follows a power law in the ratio of sorbed and dissolved concentrations. This is valid for low to intermediate concentration values.

The saturation of adsorbent sites has been also observed for high concentrations, pointing out the ability of our model to combine the features of the classical Langmuir and Freundlich isotherms. Our proposed approach opens up a new way to predict and control an adsorption-based process using a particle-based method with a finite (and actually quite low) number of particles.


  • Benson and Meerschaert (2008) Benson, D. A., Meerschaert, M. M. (2008). Simulation of chemical reaction via particle tracking: Diffusion-limited versus thermodynamic rate-limited regimes. Water Resour. Res. 44, W12201
  • Benson and Meerschaert (2009) Benson, D. A., Meerschaert, M. M. (2009). A simple and efficient random walk solution of multi-rate mobile/immobile mass transport equations. Adv. Water Resour. 32, 532–539.
  • Brunauer et al. (1938) Brunauer, S., Emmett, P. H., Teller, E. (1938). Adsorption of gases in multimolecular layers. J. Am. Chem. Soc. 60, 309–319.
  • Bolster et al. (2012) Bolster, D. de Anna P., Benson, D.  A., Tartakovsky, A. M. (2012). Incomplete mixing and reactions with fractional dispersion. Adv. Water Resour. 37, 86–93.
  • Boso et al. (2013) Boso, F., Bellin, A., Dumbser, M. (2013). Numerical simulations of solute transport in highly heterogeneous formations: A comparison of alternative numerical schemes. Adv. Water Resour. 52, 178–189.
  • Costa et al. (2015) Costa, D., Pradier, C. M., Tielens, F., Savio, L. (2015). Adsorption and self-assembly of bio-organic molecules at model surfaces: A route towards increased complexity. Surf. Sci. Rep. 70, 449–553.
  • Davis (1946) Davis, R. T. (1946). The activated adsorption of Nitrogen on a finely divided Tungsten powder. J. Am. Chem. Soc. 68, 1395.
  • Dentz et al. (2011) Dentz, M., Le Borgne, T., Englert, A., Bijeljic, B. (2011). Mixing, spreading and reaction in heterogeneous media: a brief review. J. Contam. Hydrol. 120-121,1–17
  • Fernàndez and Sanchez (2011) Fernàndez-Garcia, D., Sanchez-Vila, X. (2011). Optimal reconstruction of concentrations, gradients and reaction rates from particle distributions. J. Contam. Hydrol. 120-121, 99–114.
  • Frankenburg (1944) Frankenburg, W. G. (1944). The adsorption of hydrogen on Tungsten. J. Am. Chem. Soc. 66, 1827–1838.
  • Freundlich (1906) Freundlich, U. (1906). Die adsorption in lösungen, 385 – 470.
  • Gillespie (1976) Gillespie, D. T. (1976). A general method for numerically simulating the stochastic time evolution of coupled chemical reactions, J. Comput. Phys. 22, 403–434.
  • Gillespie (1977) Gillespie, D. T. (1977) Exact stochastic simulation of coupled chemical reactions. J. Phys. Chem. 81, 2340–2361.
  • Gillespie (2000) Gillespie, D. T. (2000). The chemical Langevin equation. J. Chem. Phys. 113(1),297–306.
  • Gramling et al. (2002) Gramling, C. M., Harvey, C. F., Meigs, L. C. (2002). Reactive transport in porous media: a comparison of model prediction with laboratory visualization. Environ. Sci. Technol. 36, 2508–2514
  • Henri et al. (2014) Henri, C. V., Fernàndez-Garcia, D. (2014). Toward efficiency in heterogeneous multispecies reactive transport modeling: A particle-tracking solution for first-order network reactions, Water Resour. Res. 50, 7206–7230
  • Jones et al. (1996) Jones, M. C., Marron, J. S., Sheather, S. J. (1996). A brief survey of bandwidth selection for density estimation. J. Am. Stat. Assoc. 91, 401–407.
  • Kang et al. (1985) Kang, K., Redner, S. (1985). Fluctuation-dominated kinetics in diffusion-controlled reactions. Phys. Rev. A. 32, 435–447.
  • Langmuir (1918) Langmuir, I. (1918). The adsorption of gases on plane surfaces of glass, mica, and platinum. J. Am. Chem. Soc. 40, 1361.
  • Moroni et al. (2007) Moroni, M., Kleinfelter, N,. Cushman, J. H. (2007). Analysis of dispersion in porous media via matched-index particle tracking velocimetry experiments. Adv. Water Resour. 30, 1–15.
  • Ovchinnikov et al. (1978) Ovchinnikov, A. A., Zeldovich, Y. B. (1978). Role of density fluctuations in bimolecular reaction kinetics. Chem. Phys. 28, 215–218.
  • Park and Marron (1990) Park, B. U., Marron, J. S. (1990). Comparison of data-driven bandwidth selectors. J. Am. Stat. Assoc. 85, 66–72.
  • Paster et al. (2014) Paster, A., Bolster, D., Benson, D. A. (2014). Connecting the dots: Semi-analytical and random walk numerical solutions of the diffusion-reaction equation with stochastic initial conditions. J. Comput. Phys. 263, 91–112.
  • Pedretti and Fernàndez (2013) Pedretti, D., Fernàndez-Garcia, D. (2013). An automatic locally-adaptive method to estimate heavily-tailed breakthrough curves from particle distributions. Adv. Water Resour. 59, 52–65.
  • Pollock (1998) Pollock, D. (1998). Semianalytical computation of path lines for finite-difference models. Ground Water 26, 743–750.
  • Quinonesa and Guiochon (1996) Quinonesa, I., Guiochon, G. (1996). Derivation and application of a Jovanovic-Freundlich isotherm model for single-component adsorption on heterogeneous surfaces. J. Colloid Inter. Sci. 183, 57–67.
  • Rahbaralam et al. (2015) Rahbaralam, M., Fernandez-Garcia, D., Sanchez-Vila, X. (2015). Do we really need a large number of particles to simulate bimolecular reactive transport with random walk methods? A kernel density estimation approach. J. Comput. Phys. 303, 95 –104.
  • Raje and Kapoor (2000) Raje, D. S., Kapoor, V. (2000). Experimental Study of Bimolecular Reaction Kinetics in Porous Media. Environ. Sci. Technol. 34, 1234–1239.
  • Redlich and Peterson (1959) Redlich, O. (1959). Peterson, D. L. A useful adsorption isotherm. J. Chem. Phys. 63, 1024 – 1024.
  • Salamon et al. (2006a) Salamon, P., Fernàndez-Garcia D, Gomez-Hernandez J. J. (2006). Modeling mass transfer processes using random walk particle tracking. Water Resour. Res. 42, W11417.
  • Salamon et al. (2006b) Salamon, P., Fernàndez-Garcia D, Gomez-Hernandez J. J. (2006). A review and numerical assessment of the random walk particle tracking method. J. Contam. Hydrol. 87, 277–305.
  • Sanchez et al. (2010) Sanchez-Vila, X., Fernàndez-Garcia, D., Guadagnini, A. (2010). Interpretation of column experiments of transport of solutes undergoing an irreversible bimolecular reaction using a continuum approximation. Water Resour. Res. 46, 1–7.
  • Sheindorf et al. (1981) Sheindorf, C., Rebhum, M., Sheintuch, M. A. (1981). Freundlich-type multicomponent isotherm. J. Colloid Inter. Sci. 79,136–42.
  • Sips (1948) Sips R. (1948). On the structure of a catalyst surface. J. Chem. Phys. 16, 490–495.
  • Silverman (1986) Silverman, B. W. (1986). Density Estimation for Statistics and Data Analysis, London:Chapman and HallLondon.
  • Tartakovsky and Meakin (2005) Tartakovsky, A. M., Meakin, P. (2005). A smoothed particle hydrodynamics model for miscible flow in three-dimensional fractures and the two-dimensional Rayleigh-Taylor instability. J. Comput. Phys. 207, 610–624.
  • Thomas (1957) Thomas, W.  J. (1957). Chemisorption of ethane at iron and nickel oxide. Trans. Faraday Soc. 53, 1124–1131
  • Tompson (1993) Tompson A. F. B. (1993). Numerical Simulation of Chemical Migration in Physically and Chemically Heterogeneous Porous Media. Water Resour. Res. 29, 3709–3726.
  • Tompson and Dougherty (1992) Tompson, A. F. B., Dougherty D. E. (1992). Particle-grid methods for reacting flows in porous media with application to Fisher’s equation. Appl. Math. Model. 16, 374–383.
  • Urano et al. (1981) Urano, K., Koichi, Y., Nakazawa,Y. (1981). Equilibria for adsorption of organic compounds on activated carbons in aqueous solutions I. Modified Freundlich isotherm equation and adsorption potentials of organic compounds. J. Colloid Interface Sci. 81, 477–485.
  • Willmann et al. (2010) Willmann, M., Carrera, J., Sanchez-Vila, X., Silva, O., Dentz, M. (2010). Coupling of mass transfer and reactive transport for nonlinear reactions in heterogeneous media. Water Resour. Res. 46, W07512.
  • Hill (1949) Hill, T. L. (1949). Statistical Mechanics of Adsorption. VI. Localized Unimolecular Adsorption on a Heterogeneous Surface. J. Chem. Phys. 17, 762–771.
  • DaRocha et al. (1949) DaRocha, M. S., Iha, K., Faleiros, A. C., Corat, E. J., Suarez-Iha, M. E. V. (1997). Freundlich’s Isotherm Extended by Statistical Mechanics. J. Coll. Inter. Sci. 185, 493–496.
  • Garabedian et al. (1988) Garabedian S. P., Gelhar L. W., Celia, M. A.(1988). Large-scale dispersive transport in aquifers: Field experiments and reactive transport theory. Rep. 315, W10103-17.
  • Fernàndez-Garcia et al. (2004) Fernàndez-Garcia, D., Illangasekare, T. H., Rajaram, H. (2004). Conservative and sorptive forced-gradient and uniform flow tracer tests in a three-dimensional laboratory test aquifer. Water Resour. Res. 40, W10103.
  • Engel, et al. (1994) Engel, J., Herrmann, E., Gasser, T., (1994). An iterative bandwidth selector for kernel estimation of densities and their derivatives. Nonparametric Stat. 4, 21-34.
  • Sheather (1992) Sheather, S. J. (1992). The performance of six popular bandwidth selection methods on some real data sets (with discussion). Comp. Stat., 7, 225-50, 271-281.
  • Fernàndez et al. (2012) Fernàndez-Garcia, D., Bolster, D., Sanchez-Vila, X., Tartakovskyc, D. M. A Bayesian approach to integrate temporal data into probabilistic risk analysis of monitored NAPL remediation. Adv. Water Resour. 36, 108-120.