Propagation and adsorption of nanoparticles in porous medium as traveling waves Confinement effects -->

Generally, one attempts to globally predict or interpret, from numerical simulations, the history of the effluents exiting porous media (i.e., the breakthrough curve), without a clear view of the detailed evolutions of deposition inside the medium. We developed a simple physical frame of description of the colloidal particle transport and adsorption, which allows to predict the main characteristics of transport and deposition in porous media from a set of directly measurable (macroscopic) physical parameters. More precisely, we show that the deposition distribution is basically a traveling wave propagating in the medium with a shape (frontal or extended) and velocity depending on the flow rate and the availability of particles with regards to the adsorption capacity. This in particular makes it possible to predict or interpret the breakthrough curve shape from a physical approach. We also show that additional effects may be included, such as a multiporosity leading to confinement effects (delayed deposition in less accessible regions). The validity of the model is checked from original direct visualizations by confocal microscopy of particle adsorption in time and space for nanoparticle suspensions flowing through a bead packing. This makes it possible to measure the evolution of the deposition profiles in time distinguishing the deposition in confined regions. The model appears to successfully predict the different trends: traveling wave, global deposition profile shape, profiles of deposition in confined regions.


I. INTRODUCTION
Colloidal particles from industrial or natural sources propagate and may deposit and alter the environment they flow through. Common problems include the accumulation of particles impacting industrial [1][2][3][4], medical and biological [5][6][7][8] processes (injection, filtration, storage, cleaning, sorting, etc.), or the leaking of contaminants from industrial or hydrologic processes transported in ground water and accumulating in soils [9][10][11][12]. Predicting particle transport and deposition in these porous media is the key to solve these problems.
All these applications involve suspensions of submicronic particles (chemicals, heavy metals, microorganisms, engineered remediation agents) flowing through saturated porous matrices with a diversity of geometries and porosities, with some possibility of adsorption thanks to van der Waals or electrostatic attraction. Obviously, the physical problem is complex since particle transport and retention on the walls of the porous media are a function of the specific interactions and of the way the particles access the surfaces. Many studies described in details the possible interactions relying on the classical DLVO theory and/or taking into account specific hydrodynamic processes around the solid walls of porous systems [13][14][15][16][17][18][19][20][21]. Geometrical aspects were also sometimes considered as important (theoretical impact of pore size polydispersity [22,23], shadow effect [24], local confinement models [25]) but the evaluation of their impact on macroscopic transport modeling is somewhat indirect.
Deposition models generally propose empirical corrections of the common advection/diffusion system of equations [9,[26][27][28][29] to take into account the different types of interactions between the particles and the solid walls of the porous structure, generally through a modification of the adsorption kinetic constant. However, the validity of these approaches is often only indirectly tested through the analysis of the breakthrough curve [20,21,[30][31][32][33], corresponding to the exited concentration as a function of time. The impact of a multiporosity, which is a crucial point for transport in soils, has been considered within the same frame [20,29,30]. We are thus faced with the paradoxical situation that the validity of a detailed description of the adsorption processes inside the medium, for a macroscopic model of transport and deposition, is tested through the observation of the transported particles outside the porous medium.
On the other side, there has been several direct observations of some processes, in particular clogging effects, in simplified systems (2D and/or at pore scale) [23,[34][35][36][37][38][39], and qualitative static imaging or indirect imaging in realistic systems [40][41][42]. Besides, detailed visualizations at the pore scale inside 3D porous media provided a full set of information on the local mechanisms of bacterial motion and trapping [43,44]. At last, it was possible from MRI (magnetic resonance imaging) to get a view of the adsorbed and suspended particles during the flow through a bead packing, in a limited range of flow rates and for fast adsorption: a regular progression of the front of adsorption was observed [45]. These different works helped clarifying our understanding of the local phenomena but were in general not directly used for transport and deposition modeling at a macroscopic scale. Nevertheless, recent quantitative internal visualizations, with NMR (nuclear magnetic resonance) [46] and confocal microscopy [47], appeared useful to develop probabilistic models of deposits and/or clogging at a macroscopic scale and check their agreement with internal measurements.
Here, our objective is to establish a straightforward relation between the events at a local scale and the macroscopic behavior of the system. We focus on deposition processes leaving the porous medium unchanged, i.e., clogging or even thick deposit as compared to the pore scale, which have been considered in Refs. [46][47][48], are negligible. We show that a minimal number of ingredients may be used to describe the principles of the propagation and deposition of colloids in porous medium, and that such an approach leads a straightforward quantification of the process in simple physical terms. This relies on a description of the local processes which overcomes the complexity of the porous structure and the hydrodynamics, by using a probability of adsorption for each particle. This probability can be considered as constant on average precisely thanks to the complexity of the flow which rapidly induces a homogeneous dispersion of the particles.
We show that, within the frame of a simple model of colloidal particle transport and adsorption through porous medium, the deposition distribution is a traveling wave propagating in the medium with a shape and velocity depending on the flow rate and the availability of particles with regards to the adsorption capacity. The breakthrough curve is then predicted as a straightforward consequence of this shape. This model is validated from direct visualizations through confocal microscopy of particle adsorption in time and space for nanoparticle suspension flowing through a bead packing. This model can serve as a reference to which various more complex effects may be added. As an example, we show that it can be extended to media with multiporosities, leading to "confinement effects."

A. Homogeneous porous medium
We consider the transport and deposition of colloidal particles initially suspended in a liquid with negligible sedimentation, the suspension (liquid + particles) being injected at constant flow rate in porous medium. This process is usually described by an advection-dispersion equation containing a blocking function which can express in more or less details the physical effects leading to deposits [48][49][50][51]. Here we built the description from more global and straightforward physical assumptions. We assume that the particles are essentially transported by the fluid but each particle entering a pore has some given probability to be adsorbed to the walls of the pore. We consider the case for which only one layer of particles can be adsorbed to the walls. Another fundamental assumption (in agreement with our experimental observations under various conditions: see below) is that a particle which has been adsorbed will never start again to move. It seems in fact rather reasonable to consider that a deposited particle will not easily be moved again if the flow conditions, which allowed the particle to fix onto the solid surface, do not change, as expected under steady injection. Remark that, on the contrary, some particle removal can result from a significant evolution of the pore structure due to deposits, which has led to a significant increase of the local velocity [47].
We then define the deposited amount s(x, t ) as the number of adsorbed particles per unit volume at time t in the pores situated at the distance x from the entrance. This variable must be understood as a mean value over a thin cross-sectional layer at the distance x. We also define n(x, t ) as the number of particles in suspension per unit fluid volume.
Finally we assume that over a short distance, typically a few pores, the suspension is sufficiently mixed so that the particle distribution can be considered as homogeneous and the suspension characteristics are fully described through the value of n. Thus we consider that the complex fluid flow tends to rapidly disperse again the particles still flowing at all possible distances from the pore surfaces. It follows that, for a given velocity of the fluid, we can describe the process in terms of a constant probability of adsorption for the particles moving through some elementary volume of the porous medium (possibly of several pore lengths). Actually, we have also to take into account that this probability of adsorption is decreased if the solid surface is already partly covered by particles. Note that here we neglect the impact of adsorption on the porous structure, which amounts to consider that the particle size is much smaller than the pore size.
During some elementary flow step, i.e., over an elementary time period dt, during which the suspension flows through the elementary volume of porous medium (see above), the gain in deposition ds is proportional to (1) The number of particles entering the pores, i.e., nV dt , in which V is the mean flow velocity through the pores; V is the effective mean velocity of the fluid and particles, which must not be confused with the mean velocity through the porous medium, equal to εV , where ε is the porosity of the medium.
(2) The probability of adsorption of a particle entering a pore, i.e., k, which depends on the characteristics of the system (particle size and interaction with wall, structure of the porous medium, flow characteristics).
(3) The fraction of surface sites available at that time: 1 − s/s 0 , in which s 0 is the number of surface adsorption sites per unit volume (i.e., the maximum number of particles adsorbed at saturation). Thus we assume that only one layer of adsorbed particles can form on the bead.
We thus get ds = kn(1 − s/s 0 )V dt. Note that this expression is close to the classical "Langmuirian" blocking function [51] except that here we give a precise definition of the adsorption probability which will lead to take into account the flux of particles entering the pore.
Besides, the relationship between the evolution of n and s is described through the standard advection equation, which expresses that the variation of concentration of suspended particles in the moving fluid (with the velocity V ) results from the loss of particles due to adsorption. Here we assume 033074-2 negligible dispersion, as this effect essentially tends to smoothen or expand the profiles. Moreover, it was recently shown [41] that dispersion of colloidal suspension flowing through porous medium is much weaker than generally assumed. The full set of equations then writes This set of equations is similar to that proposed by Leij et al. [52] under the same assumptions. The only difference concerns the definition of k, but this does not affect the solution. Usually, such a system, including a dispersion term, is solved numerically. Leij et al. [52] provided the full analytical solution of (1) for an injected pulse of suspension (see Appendix A). They were then able to predict the breakthrough curve and the setback distance (where a particular concentration is reached) as a function of time, and compared these results to experimental data. However, from this solution, general trends can hardly be inferred as the predictions mainly concern specific parameter variations. In the following, we will use a different approach leading to more straightforward descriptions of the various characteristics of the deposits, as a function of a reduced number of parameters. This approximate dimensionless approach can then be used for a physical analysis of the processes. We will nevertheless refer to the general solution (see Appendix A) to confirm the range of validity of our solutions. Finally, our approach will also allow to find an approximate solution for the case of a biporous medium. We consider the case of a continuous injection of the suspension, i.e., "a train" of suspension. Under these conditions, our fundamental assumption is that after some distance from the entrance we get a stationary distribution of deposited particles progressing through the medium. As a consequence, we look at a simple solution to these equations in the form of traveling waves of constant velocity (V 1 ) for s and n. In such a case, some time after the entrance of the suspension in the medium, the profiles of adsorbed particles [i.e., s(x, t )] collapse onto a single master curve thanks to a shift proportional to the time. We thus introduce X = x − V 1 t, and now look at s(x, t ) and n(x, t ) as respectively s(X ) and n(X ). For this theoretical description, it is useful to define these variables over the whole domain of X , i.e., ] − ∞; ∞[, with the boundary conditions s(−∞) = s 0 ("adsorption capacity," i.e., for saturated medium) and n(−∞) = n 0 (initial concentration of particles in suspension at entrance), and s(+∞) = n(+∞) = 0 (medium not yet reached by any particle).
The mass conservation imposes that the number of particles entering the domain per unit time and unit open section area (i.e., sample cross-section times porosity) be equal to the rate of variation of the total number of particles inside the domain. We thus have Using this relation between V and V 1 and the rescaled concentrations S = s/s 0 and N = n/n 0 , the set of equations (1) then becomes with a = (1 + n 0 /s 0 ). This parameter expresses the availability of particles with regards to the adsorption capacity. The solution of Eq. (3a) taking into account the boundary condition at −∞ is This result, which is intimately related to our description of the process as traveling waves, in particular expresses the fact that the adsorption increases with the number of suspended particles, and cannot exist otherwise, i.e., when there is no particle in suspension. Since S continuously varies between 1 and 0 along the X axis it reaches 0.5 in some point, that we use as the position of reference (X = 0). Then the solution of Eq. (2b) is with λ = 1/ak. This deposition distribution corresponds to a sigmoidal shape with a characteristic length λ. Note that this description is obviously an approximation of reality since, to be strictly exact under the above boundary conditions, the expression (5) would require the train of suspension to have covered an infinite distance along the porous medium. The conditions of validity of this approximation are analyzed in the next section. However, we can already remark one particular aspect of this approximation: the real distribution is in fact limited by the distance (say, x 1 ) reached by the suspension at the time considered, so that S(x > x 1 ) = 0. If we note x 0 = V 1 t 0 , the center of the wave reached at time t 0 , then x 1 = V t 0 , and S abruptly drops to zero beyond a distance x 1 − x 0 = x 0 s 0 /n 0 from the wave center. We can neglect this effect if this value if much larger than λ, which implies that V t 0 1/k, i.e., if the probability of adsorption is such that most of the particles have been adsorbed over the distance V t 0 .
Finally, the deposited distribution during the train motion takes the form of a wave front propagating at the velocity V 1 along the medium axis, and with a more vertical shape for smaller λ (see Fig. 1), i.e., for larger probability of adsorption and larger availability of suspended particles at entrance.

B. Validity of the traveling wave assumption
For the traveling wave assumption to be valid, a sufficient volume of suspension must have penetrated the medium. A criterion of validity is that the deposited amount at the porous medium entrance, as predicted by the model (1), be closed to that predicted by the traveling wave model. At the entrance (i.e., x = 0), the exact model [i.e., Eq. (1b)] with the boundary condition n(0, t ) = n 0 tells us that the deposited amount evolves as S(x = 0) = 1 − exp(−kV (a − 1)t ). On the other side the approximate solution [Eq. the two expressions differ by less than 10% when kn 0 V t > s 0 . This means that, for the adsorption to approach saturation around the entrance, the length of train already injected (V t) must be larger than a characteristic length depending on the adsorption capacity and suspension concentration (s 0 /kn 0 ). A similar conclusion is reached by approaching the traveling wave solution from the general solution of Eq. (1) (see Appendix A).

C. Distribution of adsorbed particles after a pulse passing
The above description, through a traveling wave propagation, fundamentally relies on the assumption that the suspension is continuously injected upstream in the medium from some initial time, in the form of a train without end. However, in practice, we generally inject a pulse, i.e., a train of suspension of finite length, or a series of pulses. This may in particular allow to measure the distribution of adsorbed particles after a pulse.
Since in our description there is no diffusion effect and no possibility for a deposited particle to start again moving, a succession of pulses will be equivalent, in terms of deposition characteristics, to a train of total length equal to the sum of pulse lengths. However, the deposited distribution after the passage of the last pulse differs from that associated with the traveling wave, as the end of particle injection at the back of the pulse induces a transient deposit which alters the stationary distribution (during flow).
Actually, we can deduce the distribution after the passage of the pulse, from the stationary distribution (thus assuming that the traveling wave is established). It may indeed be remarked that the deposited concentration reached at the point of the back of the pulse will not change any more during the further advance of the pulse. Let us consider a pulse of time duration t (i.e., injection of suspension during a time t). After the complete passing of the pulse, the deposited fraction at the point x is that reached in the train (full wave) at the time t for which the back of the pulse is situated in that point. This implies that the deposit is S * = S(X * ), as given by (5) Thus we find that deposited distributions following a pulse of suspension are similar to the distributions associated with a moving train [i.e., Eq. (5)], with a characteristic length 1/k instead of λ = 1/ak. Since the factor a is larger than 1, the deposited distribution S * is a sigmoid of larger characteristic length than S. This solution is effectively the general solution of the set of equations (1) under the condition that the duration of the pulse be sufficient large, i.e., kn 0 V t > s 0 (see Appendix A).
Since the medium remains intact after the passage of a pulse, it starts evolving again as inside a train as soon as the next pulse arrives. Thus, if we consider a series of pulses (successively numbered as p = 1, 2, …) of identical duration t we will get similar successive deposited distributions given by Eq. (6) with now a shift factor along the x-axis equal to the length p(n 0 V t/s 0 ). Finally, the deposits resulting from this series of pulses will appear as traveling waves at a velocity V 2 = n 0 V /s 0 . Note that this value is larger than the traveling waves of the deposited distribution during flow (i.e., V 1 ) because in the latter case there is still a fraction of particles in suspension.

D. Breakthrough curve
We can also predict the breakthrough curve under the condition of a sufficiently long pulse for the traveling wave to be established. We consider a pulse of length l injected from the initial time t = 0 in the medium at the velocity V and a medium exit located at the distance L. In that case the breakthrough curve is given by the suspension concentration at the exit of the medium: This may be expressed with dimensionless variables as in which K = kL is a factor of probability of adsorption in the system, and T = V t/L is a dimensionless time, i.e., the time rescaled by the time needed for the pulse to reach the medium end. An illustration of different shapes of breakthrough curves obtained for different values of a is provided in Fig. 2. The interest of this model clearly appears here. Within the frame of this approach, it is possible to describe or predict the trends with the help of a limited number of parameters with a physical meaning, namely the availability of particles for adsorption and the probability of adsorption. This might also constitute a first approach of any deposition problem, from which further aspects could be considered as alterations of the predictions from this simple model.

E. Multiporous medium
Let us now consider a more complex porous medium (than the above homogeneous one) made of q different regions at the pore scale, with different probabilities k i , with 1 i q, for the particles to be adsorbed. As a consequence, we now follow the deposited fractions s i in each region. We can extend the above approach for a homogeneous medium by assuming that the complex flow allows to fully mix the suspension over a short distance, so that at any time, in each region, the suspended particle concentration is n in the different region types. The characteristic velocity of the flow through each region (i.e., V i ) may differ from V and we write it as β i V . Note that for a laminar (Stokes) flow the coefficients β i are independent of V . The set of equations describing the process is now: where s = p 1 s i and s 0,i is the maximum number of deposited particles in the region i.
This approach appears particularly useful in the case of a medium with a double porosity, i.e., with locally (at the pore scale) two types of regions with respective probabilities of adsorption k 1 and k 2 . We consider the solution in the form of traveling waves of same constant velocity (V 1 ) for s, n, and every s i . Due to mass conservation the relation (2) between V 1 and V is still valid. This leads to the following dimensionless set of equations: where the same definitions for S and N hold, with now s 0 = s 0,1 + s 0,2 , and with S i = s i /s 0,i and ξ i = s 0 /s 0,i for i = 1, 2. It is important to keep the different factors in the right-hand side of Eqs. (9b) and (9c) as they express different major physical effects: ξ i represents the inverse of the fraction of adsorbed particles in region of type i, β i the relative velocity in this region, and k i the specific probability of adsorption in this region.
Coupling equations (9b) and (9c) in particular provides the following relationship between S 1 and S 2 : where u = β 2 ξ 2 k 2 /β 1 ξ 1 k 1 . The solution for S 1 as a function of X may be obtained by inserting the expression (10) for S 2 as a function of S 1 in Eq. (9b) then solving for S 1 .
It is also possible to get an approximate solution in the case of a small number of particles in regions of type 2 compared to the number of particles in regions of type 1 (i.e., s 0,2 s 0,1 ) and a lower probability of adsorption in regions 2 leading to a delay for the adsorption in these regions, i.e., S 2 < S 1 . In that case, we have s 2 /s 1 = (S 2 /S 1 )(s 0,2 /s 0,1 ) 1, so that S 1 ≈ S, and the solution of the system (9a) and (9b) is similar to (5). We finally get The distribution of particle deposits in regions of type 2 is a transformation of the sigmoid (11a) into a curve of shape close to that of a sigmoid and positioned either behind that described by (11a), i.e., S 2 (X ) < S 1 (X ), for u < 1, or ahead of that described by (11a), i.e., S 2 (X ) < S 1 (X ), for u > 1. In the case u 1, S 2 approaches 1 only for −X values much larger than the characteristic length of the sigmoid S 1 , i.e., 1/β 1 ξ 1 k 1 , which means that the adsorption in regions of type 2 is considerably delayed with regards to the adsorption in regions of type 1. The shape of the profile for region of type 2 is shown in Fig. 1 for some values of λ and u. It is worth noting that the profiles for both regions are very sensitive to the values of these two parameters, which means that a comparison between theory and experiment might allow a good evaluation of the physical processes.
The deposited distributions after a pulse may be deduced, in the same way as for a homogeneous medium [see Eq. (6)], from the expressions (11a) and (11b). The deposited distributions for a series of pulses then appear as traveling waves of velocity V 2 = n 0 V /s 0 and of shapes given by (11a) and (11b) after dropping the factor a.

A. Porous media and working fluid
The porous medium is made of monodisperse, spherical, borosilicate beads, with a diameter d g = 63 μm (see Fig. 12 in Appendix B), and a strongly negative (SiO − groups) surface charge. The beads are packed in disorder in a 1 × 1 × 20 mm 3 square capillary at a volume concentration of 62%. Hence, we have typically 15 bead diameters per sample width. Although the solid phase at a local scale (bead) has a simple structure, the porous structure at a larger scale, i.e., including several particles, is complex. One may nevertheless estimate an order of magnitude of the mean pore size (d p ) from the diameter of the largest sphere that can be inserted in the hole at the center of four neighboring spheres whose centers form a tetrahedron: d p ≈ d g /6 ≈ 10 μm. We checked that the bead packing structure did not change over the duration of our tests, so that we can rely on direct imaging of the local evolution of the porous media to measure particle deposition.
This porous medium is saturated with a fluid made of 90% of DMSO (Dimethyl sulfoxide) in distilled water, leading to a fluid of viscosity μ ≈ 3.3 mPa s and density ρ = 1090 kg m −3 . This fluid matches the refractive index of the beads and allows for internal imaging of the 3D pack with the help of confocal microscopy, a technique which proved to be relevant for the study of fluorescent particles in 3D porous media [47,53].

B. Charged nanoparticles
The colloidal particles dispersed in the fluid are nanoparticles (Molday ION EverGreen CL-50Q02-6A-51 from BioPAL Inc.) made of a magnetite core and a thick Dextran shell. These particles have a density ρ S = 1250 kg m −3 and a high positive surface charge (C6 amine surface functionalization) so that we expect strong electrostatic attractive forces between them and the surface of the beads (see Fig. 3 for a view of the adsorbed particles). On the contrary, since their charges are of the same sign, two particles strongly repel each other when they approach one another. As a consequence, we only expect the formation of a monolayer of particles adsorbed onto the bead surface, but with particles somewhat dispersed on the surface due to electrostatic repulsion between them. The fluorescence of these particles is achieved thanks to additional surface groups (excitation/emission wavelengths 505/515 nm). To estimate the mean particle size, the nanoparticles are imaged under SEM (scanning electron microscopy) after dispersion then drying on a clean surface (see Ap- pendix B). The particle diameter appears to range from 35 to 120 nm with a mean size d = 60 nm (higher than given by the supplier, i.e., 35 nm). Such a particle size is more than two orders of magnitude smaller than the typical pore size, which ensures negligible clogging due to geometrical confinement (see Ref. [46]).
We prepared the suspension at a particle volume fraction n 0 = 7.4 × 10 −4 , which is sufficiently large for a significant deposition to occur over the typical duration of our tests, i.e., a few hours, and sufficiently small for hydrodynamic interactions between them to be negligible [54].

C. Confocal imaging
A single argon-ion laser source is used to excite the fluorophores of the nanoparticles, in a single plane 400 µm above the capillary bottom. We image a large portion of the sample (1 mm × 10 mm) by stitching together successive 1 × 1 mm 2 square acquisitions (see Fig. 4). The time required to perform a complete scan over the medium being a few minutes, significant evolution of the porous structure due to particle adsorption may occur during this time. In order to fix the structure during imaging we injected pulses of suspension, i.e., injection over a limited time, followed by pure fluid to rinse the sample and subsequent imaging (see Fig. 4). We checked independently that the rinsing period does not modify the deposition: the fluorescence signal during the injection of pure fluid following a pulse does not evolve, meaning that the flushing step with pure fluid does not detach particles already deposited. Note that this approach allows to directly and strictly image the deposited particles, instead of having to distinguish between suspended and deposited particles in the images. This also proves the fundamental point that when they are adsorbed somewhere the particles do not move anymore or do not detach at any time during the flow.
The particles being smaller than our imaging resolution, we do not follow the number of adsorbed particles but the 033074-6 fluorescence (summing the intensity over all pixels) in crosssectional layers which can be considered as proportional to the particle concentration. We thus obtain a distribution (profile) of deposited particle concentration along the porous medium axis.

D. Flow characteristics
The volume of the pulses is adjusted to have a significant change of the signal between two successive pulses, and is constant for all our tests, i.e., 20 µL. The volume used for rinsing is chosen so as to apparently remove all suspended particles, i.e., 80 µL. The flow rate (Q) was varied in the range 1.3−208 μl min −1 , giving a mean flow velocity through the pores V = Q/S (with S the sample cross-section area) in the range 0.055−8.8 mm s −1 . With the successive injection cycles the total duration of an experiment ranged from about one to nine hours.
The Reynolds number of the flow, i.e., ρV d p /μ, is in the range 0.00017-0.027, which means that the flow is laminar. The particles have a density larger than that of the fluid but their characteristic velocity of sedimentation through the fluid at rest [i.e., 2(ρ S − ρ)g(d 2 /4)/9μ] is equal to 10 −10 m s −1 . This means that it would take about 10 5 s for a particle to fall over a distance equal to the pore size under gravity. Such a value is much larger than the average duration of its motion through the porous medium (L/V , in which L is the sample length), i.e., between 2 and 360 s. Thus, even over the largest possible distance inside the medium sedimentation is negligible. Moreover, the ratio of the particle inertia to the typical work of viscous forces over a distance equal to the particle size, is equal to about 1/10 of the Reynolds number. This means that this ratio is far much smaller than 1 and thus particle inertia plays a negligible role on the particle motion. At last, the diffusion coefficient of the particles, estimated from the Stokes-Einstein expression, i.e., D = k B T /3πμd, where k B = 1.38 × 10 −23 J K −1 is the Boltzmann constant and T ≈ 293 K is the room temperature, is equal to 2.2 × 10 −12 m 2 s −1 . The characteristic time for a particle to move through a pore thanks to fluid motion is d p /V , while the characteristic time for a particle to move over a distance of the order of the pore size thanks to thermal agitation is d 2 p /D. The ratio of these two times is the Peclet number, i.e., V d p /D, which is here in the range 250-40 000. This means that particle diffusion due to thermal agitation is fully negligible in comparison to advection processes.
Turbulence, inertia, sedimentation, and diffusion due to thermal agitation being negligible, the particles are essentially carried by the liquid and follow the paths imposed by the fluid flow. There is nevertheless one additional effect which can make them deviate from these paths, namely the electrostatic interactions between the particles and the bead surfaces.

E. Adsorption
Let us consider the interaction between a particle and a charged bead surface. The electrostatic field induced by a single bead can be determined from the Gauss theorem, it is equal to E = −q/4πε 0 (R + x) 2 , in which ε 0 = 8.9 × 10 −12 F m −1 , R = d/2, −q the total charge on the bead, and x the distance from the bead surface. The resulting (attractive) force on the charged particle is then eE (with e the particle charge). It follows that a particle situated at some distance from the bead will be progressively dragged towards this surface. A particle is in fact surrounded by beads which finally induce a more complex electrostatic field resulting from the sum of the elementary fields. It follows that the force acting on the particle whose position is r, is equal to F e = −(eq/4πε 0 ) i |r − r i | −3 (r − r i ), where the sum is taken on all beads in the sample, whose center positions are given by r i . On the other side, a particle moving relatively to the fluid at rest experiences a drag force equal to F d = −6πμRṙ. Finally, a particle initially at rest in the fluid and situated at r 0 , will have a motion described by the fundamental law of mechanics: mr = F d + F e . The solution of this equation taking into account the initial condition provides the time t (r 0 ) for the particle to reach the solid surface of one of the closest beads, i.e., those forming the wall of the pore containing the particle. This is more generally the characteristic time for a particle to be adsorbed onto a solid surface.
On the other side, the characteristic time spent in a pore by a particle moving with the fluid is t v = d p /V . We can then determine the domain (a surface) of the different positions r 0 at the entrance of the pore for which t v < t (r 0 ). When the particle appears in this domain it will ultimately reach a solid surface within the duration of the flow through the pore. The ratio of this surface to the total surface of entrance in the pore, finally represents the probability to reach a solid surface during the particle motion through the pore.
Obviously, we expect that this probability decreases when the flow velocity is increased. For example, let us consider the simplistic case of a particle moving parallel to a bead surface (i.e., its trajectory remains tangential to the bead surface) as a result of liquid flow, and initially situated at a distance x 0 . Neglecting inertia and the other beads we find that the approach to the bead is described from the equation: 0 = −6π Rμẋ − eq/4πε 0 (R + x) 2 , in which x is the current distance from the bead surface. Integrating this equation and assuming at first order that all distances remain small compared to R we get the time for the particle to reach the bead surface: t (x 0 ) = Kx 0 , with K = 24π 2 R 3 με 0 /eq. Since the time spent by the particle in the pore is approximately d p /V , all particles initially situated at a distance smaller than x c = (d p /K )V −1 should reach and stop on the pore surface before they can exit the pore, while the others will exit the pore. However, when x c is larger than the pore size, i.e., for sufficiently low velocity, the probability to reach the pore wall is equal to 1. In this asymptotic case, we should observe a straight front of deposition. Otherwise, i.e., in the general case, the probability of adsorption of a particle entering the pore is proportional to During the flow of a suspension through the porous medium the particles thus progressively adsorb to the bead surface. However, as soon as a particle encounters another particle, they strong repel each other. As a consequence, we expect only the formation of a monolayer of adsorbed particles. This implies that the concentration of particles increases until no new particle reaching the bead can find a sufficient room to be submitted to the attraction of the bead surface and at the same time sufficiently far from other particles. This situation (i.e., maximum concentration) corresponds to what we will call the saturation, even if we ignore the exact value of the concentration of particles in such a case.

A. Global flow characteristics
Let us look at a typical series of confocal images taken after successive cycles (see Fig. 4). Each image is composed of dark areas (no fluorescence, i.e., no adsorbed particles) and bright areas (containing adsorbed particles). At the inlet (far left, dotted white line), a sharp transition from dark to bright highlights the transition between the inlet tubing (no beads, particles cannot adsorb) and the porous system (particles meet solid surfaces and start adsorbing). The image is the brightest around the porous medium entrance and decreases towards the right (flow direction) up to a completely black region, where a priori no particles deposited (see Fig. 5). Globally, after each new cycle (pulse + rinse), the image become brighter at larger distance, indicating that the particles progressively deposit farther in the sample. This is confirmed by the corresponding distributions of fluorescence intensity, which extend towards larger distance for successive cycles (see Fig. 6). Interestingly, the intensity near the entrance does not evolve with further injections: it seems to saturate (see Figs. 5 and 6). This observation is consistent with the assumption of a mono-layer formed by adsorption on the bead surface. Here we can consider that the signal at saturation corresponds to the maximum coverage that may be reached under such circumstances (i.e., suspension flow).
We then represent the fluorescence distribution rescaling the signal by the mean value of the plateau level progressively formed from the entrance (see Fig. 6). One can now observe more precisely the progressive spreading towards larger distance along with the extension of a plateau from the medium entrance. It is interesting to note that, in the plateau region, the profiles converge towards the same curve exhibiting (fixed) fluctuations due to slight density heterogeneities of the bead packing: this proves that the number of deposited particles in the saturated region vary negligibly from one test to another, and finally shows the reliability of our measurements. Moreover, the shape of the front of the profile remains approximately similar for the different cycles (see Fig. 6), in agreement with our assumption of traveling waves propagating through the sample. At last, the profiles can be rather well represented by the expression (6) corresponding to the model for a homogeneous sample (see Fig. 6). Data for different flow rates were collected, again after different cycles of injection. Similar results were obtained for the different flow rates: the profiles of fluorescence distribution extend at larger depth with successive cycles; these different profiles can be shifted and superimposed along a master curve for each flow rate (see Fig. 7). The fundamental impact of velocity thus clearly appears in this representation: the master profile extends over larger distance when the flow rate increases. Moreover, our model well represents these profiles for the different velocities in this range covering more than two orders of magnitude, which provides a good support to the theory. At last, the parameter 1/k fitted to data appears to vary approximately with V 1/2 , which means that 1/k would be proportional to V −1/2 . This result differs from that obtained in the case of a flow around a single bead (i.e., V −1 ) (see above) but is qualitatively consistent with the expectations in the general case, i.e., decrease with velocity increase.

B. Flow in confined regions
So far, the confocal images have been treated as a whole to extract average deposition profiles over time and distance along the porous medium. Actually, the same data may be used to study some local specificities of adsorption. A closer look [see Fig. 8(b)] to the adsorption distribution from an image taken in a region of significant deposition for a typical experiment [see Fig. 8(a)] shows, as expected, deposited particles outlining the surface of the silica beads in the form of circles (i.e., along the intersection between the imaging plane and the beads). The interior of these circles appears naturally darker, since the signal here comes from close particles not situated in the main optical focus plane.
We can nevertheless remark on the same picture some small regions much darker than the rest of the frame [see Fig. 8 between neighboring beads (see Fig. 9). These regions exhibit different shapes depending on the relative orientation of the plane tangent to the two beads in contact and the plane of observation (confocal image). Large (circular) areas appear when these two planes are parallel, while elliptic regions are observed when they are not parallel, these ellipses tending to a line when the two planes are perpendicular. Remarkably, this shows that such dark regions may be obtained along any orientation of the contact between neighboring particles, i.e., they are mostly independent of the flow direction. Nevertheless note that we can observe a slight effect of the flow: these regions often exhibit a thin tail in the flow direction (see Fig. 8); but we will not discuss further this effect here.
From images taken at successive times in the same window one can see that the apparent size of these regions changes as deposition progresses [see Fig. 10(a)], as if the adsorbed particles invaded progressively these regions. A further look at the local distribution of the signal around the point of contact between the beads provides a deeper insight in the process [see Fig. 10(b)]. At a sufficient distance (say, beyond ≈8 μm, i.e., in region A) from the point of contact the deposition is fast. At shorter distance (i.e., in the region B) the deposition is much slower. These data thus show that the region of the bead surface beyond some distance (i.e., ≈8 μm) from the point of contact is much more rapidly saturated than the inner region, which suggests that the coverage of the dark areas has a different dynamics than the rest of the bead surface. Finally, looking at the deposition profile we see that it progressively spreads towards the point of contact, in the form of a saturated front (maximum adsorption) [see Fig. 10(b)]. This in particular means that the description of the adsorption process in these regions, in the form of a progressive reduction of the size of the empty domain, is a good approximation.
It is clear that the dark regions correspond to the regions of the surface of the beads for which the surrounding liquid is the most confined, i.e., situated at the smallest distance from the solid surfaces. This obviously suggests that the effect observed, i.e., the slower dynamics of adsorption in these regions, is due to this confinement. Actually, several effects resulting from this confinement may tend to slow down the adsorption dynamics. The first effect of confinement is a ratio of liquid volume to solid surface lower than in the rest of the sample. As a consequence, for a similar particle concentration in the liquid volume, the covering of a unit bead surface in the confined region will be slower than elsewhere. Another effect of confinement is that particles necessarily need to move close to already adsorbed particles to reach deeper regions (closer to the contact point). Thus they can be somewhat repelled by these particles, which would slow down their access to deeper regions. Another effect directly reduces the concentration of adsorbed particles in these regions: the particles are stopped at some distance from the contact point because the distance between bead surfaces is too small. At last, in such a confined region the flow rate is on average smaller than elsewhere because of the proximity of the solid walls (bead surfaces). This has two opposite consequences: the rate of entrance of suspended particles is lower, which reduces the number of particles adsorbed per unit time; the time spent by the particles along unit surface is larger, which increases the probability of adsorption per particle.
To sum up, we have regions in which the probability of adsorption is lower than in the rest of the material. We can then test our theoretical predictions for a double porosity material (see Sec. II) on these data. In that aim we look again at the distribution of deposited particles after each injection cycle, now distinguishing the deposit in the confined regions. We determine the deposition profiles in the confined regions by directly measuring the apparent diameter of each dark area and translating it as an effective non-covered area over the bead (i.e., in 3D). The current saturation is then estimated from the relative difference between the current area and the initial area of each of these regions. Such an approach is justified by the observation that the coverage essentially develops in the form of a front progressing towards the point of contact (see above). This provides a discrete distribution of saturation (see Fig. 13 in Appendic C) which is then somewhat smoothened through running average to get the profiles shown in Fig. 11.
The deposited distribution obtained in this way decreases and tends to zero at a distance which increases with the injection cycle number (see Fig. 11). Then these profiles appear to fall along the same master when shifted of the same values as for the global distributions corresponding to injected volumes (see Fig. 11). This means that the progression of the deposited distribution in the confined regions also takes the form of a traveling wave progressing through the sample and with the same velocity. Moreover, as expected from the above observation, the deposition in the confined regions occurs with some delay with regards to the global deposition: the master curve for confined regions is situated behind the master curve for global deposition (see Fig. 11).

033074-10
Let us now consider the number of particles expected to be adsorbed in these regions. The maximum area of each of these regions typically corresponds to a spherical cap of apparent radius r c ≈ 8 μm, and the number of contacts of one bead with other beads is around 6 in a disordered packing. As a consequence, the fraction of bead surface associated with such regions expresses as 3(r c /R) 2 / 1 − (r c /R) 2 ≈ 0.2. This might represent the fraction of coverage associated with these regions, i.e., s 0,2 /s 0 , but we can expect that the effective value is lower, due to a more difficult access to some adsorption sites at the approach of the contact point, as a result of electrostatic interactions and geometrical confinement. Assuming for example that the deposit of particles is hindered in a spherical cap of radius 4 μm [see Fig. 10(b)], we get a new fraction of coverage of these regions equal to 0.15. In addition, we clearly observe experimentally a delay between the deposit in these regions and on the rest of the surface, thus S 2 < S 1 . Under these conditions, considering that s ≈ s 1 constitutes a reasonable approximation.
In this context it is justified to use Eq. (5) to represent the global deposition profile as we did so far, and then transform it through (10) by fitting the parameter u. We can see that this allows to well represent the stationary deposition distribution (see Fig. 11). In this way we get u = 0.35. The physical origin of this value may now be analyzed. Let us recall that the theoretical expression for u is β 2 ξ 2 k 2 /β 1 ξ 2 k 1 . We have ξ 2 /ξ 1 = s 0,1 /s 0,2 , i.e., the ratio of maximum coverage of the two region types, which according to the above analysis should be of the order of 5. k 2 /k 1 is the ratio of adsorption probabilities on the two surface types. As we excluded the region the closest to the contact point in our estimation of s 0,2 we kept a domain in which the impact of confinement on the probability of adsorption under given flow velocity should be rather low with regards to that in the rest of the material, the only effect being a smaller ratio of volume to surface. Thus, it is considered here that the difference between k 2 and k 1 essentially comes from the difference of velocities in the different regions, so that we expect k 2 /k 1 ≈ (V 1 /V 2 ) 1/2 . At last we have β 2 /β 1 = V 2 /V 1 , so that u ≈ 5(V 2 /V 1 ) 1/2 . The point is now to estimate the characteristic velocity in each region V 1 and V 2 . At as first approximation, from the Poiseuille law, we can consider that under similar pressure drop along one bead the mean velocity will vary with the square size of the pore. A simple result is then obtained for a "standard" pore, whose size was estimated to be 1/6 of the bead diameter. But a similar estimation is not straightforward for the confined regions, which makes us realize the limitations of the model: in these regions the distance between walls continuously decreases towards zero at the approach of the point contact. This suggests that for a more precise approach we should consider such regions as a continuum of media of decreasing pore size. Such a sophisticated approach nevertheless does not seem relevant with regards to the precision of our data. Coming back to the estimation of a single pore size for the confined regions we can consider that it corresponds to the gap between two beads at a distance of 5 microns of their point of contact. This gives a characteristic pore size of d p2 = 0.8 μm, and finally u ≈ 5(V 2 /V 1 ) 1/2 = 5d p2 /d p = 0.4, a value almost exactly equal to that found from the data (0.35), which confirms the validity of our approach.
Finally, we can now evaluate the validity of the model to predict the data concerning confined regions for other flow velocities. For V 5.5 m s −1 , the delay to fill the confined regions was too large consistently with the above model (see Fig. 1) so that we could not observe it clearly within our test duration. For V = 0.05 m s −1 we would observe again that the profiles for the confined regions advance as traveling waves with the same velocity as the global profiles. Moreover, at our scale of observation the profiles are rather close to the global profile and their shape is well represented by the model keeping approximately the same value (0.4) for u, which confirms the consistency of the above model (see expected profiles in Fig. 1), but the scattering on data does not allow to give a precise estimation of this parameter in that case.

V. CONCLUSION
We proposed a simple probabilistic model for the transport and desposition of colloidal particles in porous medium. This model predicts that after some distance of injection the deposition distribution is a traveling wave propagating in the medium with a shape (frontal or extended) and velocity depending on the flow rate and the availability of particles with regards to the adsorption capacity. This model can serve as a reference model in which various additional effects may be included, such as a multiporosity leading to confinement effects. Besides we carried out original direct visualizations through confocal microscopy of particle adsorption in time and space for nanoparticle suspensions flowing through a bead packing. This made it possible to measure the evolution of the deposition profiles in time distinguishing the deposition in confined regions. The model appears to successfully predict the different trends: traveling wave, deposition profile shape, profiles of deposition in confined regions.
The model relies on simplifying assumptions, such as no particle removal after adsorption, negligible diffusion, inertia and sedimentation, sufficient length of the initial pulse, etc, which have been verified in the experimental setup. However, the interest of the model is that it describes the main physical effects through a reduced set of parameters which encompasses several physicochemical effects. There is indeed one fundamental parameter, i.e., the characteristic length of the deposition profile, which appears to encompass the suspension concentration, the adsorption probability and the number of adsorption sites. A second parameter appears when confinement effects are taken into account, which describes the relative capacities of adsorption of the different media. In this context it is possible to directly test the impact of varying one physicochemical aspect, by looking at the variation it induces on each of the two above parameters.
The general objective of this approach is to have a global approach of colloid transport, leaving apart as a first step the details of the physicochemistry, and instead taking into account their effects globally, in the spirit of the approach of Ref. [55], but in a different way, by identifying the main physical mechanisms. In particular, we suggest that the shape of breakthrough curves could be studied with the help of such a model.