# Permeability impairment by hydrodynamic pore bridging: Probabilistic pore-network modeling and microfluidic experiments

Cyprien Soulaine, Walid Okaybi, Laurez Fogouang Maya, Emmanuel Le Trong, Sophie Roman

<figure><img src="image_1.png" /></figure>

## To cite this version:

Cyprien Soulaine, Walid Okaybi, Laurez Fogouang Maya, Emmanuel Le Trong, Sophie Roman. Permeability impairment by hydrodynamic pore bridging: Probabilistic pore-network modeling and microfluidic experiments. International Journal of Rock Mechanics and Mining Sciences, 2026, 200, pp.106423. 10.1016/j.ijrmms.2026.106423. insu-05512557

HAL Id: insu-05512557

https://insu.hal.science/insu-05512557v1

Submitted on 16 Feb 2026

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
<figure><img src="image_3.png" /></figure>

# Permeability impairment by hydrodynamic pore bridging: Probabilistic pore-network modeling and microfluidic experiments☆

Cyprien Soulaine<sup>1 ID*</sup>, Walid Okaybi<sup>1 ID</sup>, Laurez Fogouang Maya, Emmanuel Le Trong<sup>1 ID</sup>, Sophie Roman<sup>1 ID</sup>

Institut des Sciences de la Terre d'Orléans, Université d'Orléans, CNRS, OSUC, UMR7327, Orléans, France

## ARTICLE INFO

### Keywords:

- Pore-network modeling

- Hydrodynamic bridging

- Clogging

- Microfluidics

- Permeability

## ABSTRACT

Permeability impairment caused by the migration and retention of suspended particles is a critical issue in numerous industrial and environmental processes. While pore-network models (PNMs) have successfully described clogging by sieving and particle aggregation, they have failed to capture hydrodynamic bridging – a mechanism where particle arches form and block pore throats. This study introduces a novel probabilistic PNM that incorporates a stochastic law for arch formation, accounting for the particle-to-throat size ratio, particle concentration, and pore geometry. The probability law is calibrated using high-fidelity CFD-DEM simulations of single-pore bridging. Microfluidic experiments in heterogeneous micromodels representative of the rock microstructure are carried out to investigate the effect of particle size and concentration, and flow rate on permeability reduction. The proposed probabilistic framework successfully reproduces experimental trends in clogging dynamics and permeability decline, thereby extending the capability of PNM to capture all three pore-clogging mechanisms.

## 1. Introduction

Permeability impairment caused by the flow of suspended particles is a critical mechanism in various industrial processes, including filtration, groundwater flow, and subsurface engineering. For instance, the detachment and mobilization of fines — small sand or clay particles loosely adhered to pore walls — by hydrodynamic forces and change in fluid salinity<sup>1–3</sup> can significantly reduce permeability near wells. This phenomenon poses a serious threat to the productivity of geothermal and oil reservoirs, often requiring costly remediation treatments or even leading to the abandonment of affected wells. Three primary mechanisms are responsible for pore clogging during the flow of particle suspension: (i) sieving, (ii) bridging, and (iii) particle aggregation.<sup>4</sup> Sieving occurs when particles are physically excluded from entering a pore due to their size.<sup>5</sup> Bridging refers to the formation of arch-like structures by a few particles arriving simultaneously at a pore entrance where the flow converges.<sup>6</sup> Aggregation involves the gradual accumulation of colloidal particles, driven by long-range electrochemical interactions between the fluid, the particles, and the solid matrix.<sup>7</sup> A secondary process associated with aggregation involves the progressive growth of dendritic structures that ultimately result in pore plugging.<sup>8</sup> During operation, these clogging mechanisms may co-occur, and their combined effects can exacerbate system performance degradation.<sup>9</sup>

This study focuses on clogging by bridging – a phenomenon observed in micromodel experiments since the late 1970s, though it was not well understood at the time.<sup>10</sup> Later experiments using straight and constricted microtubes demonstrated that stable particle arches can form and obstruct flow, even at low particle concentrations.<sup>11</sup> Bridging occurs when hydrodynamic forces overcome repulsive interactions and balance contact forces between particles and between particles and the pore walls. Numerous studies have investigated the influence of particle size and concentration, flow rate, and pore geometry on the onset of arch formation. A breakthrough was achieved by Ref. 12, who found that clogging by arch formation typically occurs after a certain amount of volume of fluid is extracted. Based on that, Ref. 9 defines a critical number of particles that need to pass through a pore throat before clogging. Notably, this number was found to be largely independent of flow rate but highly sensitive to pore geometry—specifically the particle-to-throat size ratio, the pore-to-throat radius ratio, and the smoothness of the throat.<sup>13–15</sup> These findings have been confirmed by numerous experimental and numerical studies on clogging at the single-pore scale.<sup>16–18</sup> They laid the foundation for probabilistic laws linking arch formation to particle-to-throat size ratio and particle concentration.<sup>19,20</sup> More recently, Ref. 21 highlighted the role of throat angle in facilitating bridging, which

<footnote>☆ This article is part of a Special issue entitled: 'CouFrac 2024' published in International Journal of Rock Mechanics and Mining Sciences.</footnote>

<footnote>* Corresponding author.
E-mail address: cyprien.soulaine@cnrs-orleans.fr (C. Soulaine).</footnote>
Nomenclature

<table><tr><td>i, j</td><td>subscript of pores</td></tr><tr><td>k, l</td><td>subscript of particles</td></tr><tr><td>(·)<sup>T</sup></td><td>Transpose operator</td></tr><tr><td>α<sub>ij</sub></td><td>Geometric conductance factor of throat ij</td></tr><tr><td>ω<sub>k</sub></td><td>Particle rotational velocity</td></tr><tr><td>F<sup>c</sup></td><td>Particle contact force</td></tr><tr><td>F<sup>f</sup></td><td>Particle-fluid force</td></tr><tr><td>F<sup>g</sup></td><td>Particle gravitational force</td></tr><tr><td>M<sup>c</sup></td><td>Particle contact torque</td></tr><tr><td>M<sup>hyd</sup></td><td>Particle hydrodynamic torque</td></tr><tr><td>r<sub>k</sub></td><td>Position vector</td></tr><tr><td>U<sup>p</sup><sub>k</sub></td><td>Particle translational velocity</td></tr><tr><td>v<sup>f</sup></td><td>Fluid velocity vector</td></tr><tr><td>v<sup>p</sup><sub>k</sub></td><td>Particle velocity vector</td></tr><tr><td>δP<sub>ij</sub></td><td>Pressure drop in throat ij</td></tr><tr><td>[·]</td><td>Floor operator</td></tr><tr><td>D<sub>ij</sub></td><td>Effective diffusion in throat ij</td></tr><tr><td>D<sub>m</sub></td><td>Particle diffusion</td></tr><tr><td>μ<sub>f</sub></td><td>Fluid viscosity</td></tr><tr><td>ν<sub>p</sub></td><td>Particle-averaged velocity in a cell</td></tr><tr><td>ρ<sub>f</sub></td><td>Fluid density</td></tr><tr><td>ρ<sub>p</sub></td><td>Particle density</td></tr><tr><td>θ</td><td>Throat angle</td></tr><tr><td>ε</td><td>Local-cell porosity</td></tr><tr><td>C</td><td>Particle concentration</td></tr><tr><td>C<sub>m</sub></td><td>Maximal packing concentration</td></tr><tr><td>C<sub>jam</sub></td><td>Jamming transition concentration</td></tr><tr><td>D</td><td>Throat diameter</td></tr><tr><td>d<sub>p</sub></td><td>Particle diameter</td></tr><tr><td>D<sub>ch</sub></td><td>Channel diameter</td></tr><tr><td>H(·)</td><td>Heaviside operator</td></tr><tr><td>h</td><td>Micromodel thickness</td></tr><tr><td>I<sub>k</sub></td><td>Particle moment of inertia</td></tr><tr><td>k</td><td>Fitting parameter of the probability law</td></tr><tr><td>K<sub>0</sub></td><td>Microstructure permeability</td></tr><tr><td>k<sub>B</sub></td><td>Boltzmann constant</td></tr><tr><td>L</td><td>Length of the throat</td></tr><tr><td>L<sub>ij</sub></td><td>Length of throat ij</td></tr><tr><td>m<sub>k</sub></td><td>Particle mass</td></tr><tr><td>n</td><td>Number of particles of the arch</td></tr><tr><td>N<sup>*</sup></td><td>Number of particles escaping the throat</td></tr><tr><td>N<sub>l</sub></td><td>Number of pores connected to pore i</td></tr><tr><td>N<sup>c</sup><sub>k</sub></td><td>Number of objects in contact with particle k</td></tr><tr><td>p</td><td>Fluid pressure</td></tr><tr><td>P<sub>clog</sub></td><td>Probability of clogging the throat</td></tr><tr><td>Q<sub>ij</sub></td><td>Volumetric flow rate in throat ij</td></tr><tr><td>R<sub>i</sub></td><td>Radius of pore i</td></tr><tr><td>R<sub>ij</sub></td><td>Radius of throat ij</td></tr><tr><td>Re<sub>p</sub></td><td>Particle Reynolds number</td></tr><tr><td>T</td><td>Temperature</td></tr><tr><td>V<sub>i</sub></td><td>Volume of pore i</td></tr></table>

led Ref. 22 to propose a new stochastic model that predicts pore clogging probabilities based on suspension properties, throat angle, and the

particle-to-throat size ratio. The model incorporates two fitting parameters that capture the influence of pore geometry and has demonstrated strong predictive capability across a wide range of simulation and experimental conditions. These probabilistic laws of clogging have been developed for individual pore throats, but a significant gap remains in using them to predict permeability impairment at the scale of a rock, where the microstructure consists of a heterogeneous network of interconnected pores.

The prediction of permeability impairment by particle suspensions through numerical simulation remains in its infancy. While high-fidelity pore-scale simulations are now feasible, they are still highly computationally demanding. Among the most advanced approaches in this field are four-way coupled Computational Fluid Dynamics–Discrete Element Methods (CFD-DEM), which combine the grid-based – or lattice-based – solution of the Navier–Stokes equations in the pore space with Lagrangian tracking of solid particle transport and retention.<sup>23,24</sup> Modern CFD-DEM frameworks account for the feedback of particle presence on the fluid flow, and incorporate the mechanical interactions between particles and pore walls<sup>25–27</sup> as well as colloidal interaction forces.<sup>28–31</sup> Therefore, these frameworks are capable of simulating a wide range of pore-clogging mechanisms, including sieving, bridging, and particle aggregation. However, because they explicitly resolve all contact and long-range interactions, such simulations are computationally intensive and, at present, are not feasible for modeling particle migration and retention across large domains representative of the three-dimensional structure of a core sample.<sup>32</sup>

An alternative approach widely adopted in porous media research and engineering to compute permeability from digital rock images is the use of pore-network models (PNMs).<sup>33</sup> In this approach, the rock microstructure is approximated as a network of interconnected nodes – corresponding to pore bodies – and cylindrical throats.<sup>34</sup> A wide variety of PNMs have been proposed to simulate particulate transport and retention in porous media. Eventually, they should be able to predict not only the permeability response and the effluent concentration during deep bed filtration, but also the mapping of pores plugged by particles under varying particle size distributions and concentrations, flow rates, fluid pH, and salinity. The first PNM to simulate formation damage caused by injection of particles has been proposed by Ref. 35 using a simplified approach to particle transport and retention. Since this pioneering work, PNMs for particulate transport have evolved significantly and become increasingly sophisticated. Two main modeling approaches have emerged. On the one hand, particle tracking methods explicitly simulate the motion of individual particles through the network.<sup>36–38</sup> On the other hand, Eulerian approaches compute the particle concentration by solving an advection–diffusion-retention equation within the network.<sup>39–41</sup>

The particle tracking approach in PNM provides a discrete, Lagrangian framework in which particle trajectories and retention behavior are governed by local flow conditions and probabilistic rules.<sup>36</sup> A key challenge in early implementations was the redistribution of particle paths at network bifurcations. Ref. 37 addressed this using a flow-biased probabilistic algorithm, where particles choose exit throats based on the relative flow distribution, favoring pathways with higher flow rates. While their model accounted for a range of particle sizes, it was limited to simulating clogging by sieving (size exclusion). This framework was later extended by Refs. 42, 43, who introduced both straining of particles larger than pore throats and probabilistic deposition at solid walls for smaller particles, thereby capturing a broader spectrum of retention mechanisms. Follow-up studies have built on this approach, refining deposition models to assess how particle and network characteristics influence permeability damage.<sup>38,44</sup> More recent advancements, such as the model by Ref. 45, incorporate surface roughness and detailed particle–surface interactions to improve predictions of particle capture. Collectively, these efforts have significantly enhanced our understanding of how microscale transport phenomena contribute to permeability decline. However, while particle-tracking PNMs have
proven effective at simulating pore clogging by sieving and deposition, modeling clogging by bridging remains a persistent challenge.

Eulerian approaches in PNMs solve advection-dispersion equations for particle concentrations within pore bodies and throats, coupled with mechanistic descriptions of particle retention processes. Originally introduced to simulate reactive transport at the pore scale – such as surface reaction kinetics and molecular adsorption onto pore walls<sup>46,47</sup> – this framework has since been extended to account for particulate transport and clogging phenomena. A key advantage of the Eulerian formulation is its ability to capture hydrodynamic dispersion within throats, which is critical for accurately predicting the spreading and mixing of particles. Ref. 48 demonstrated that Eulerian PNM approaches can reliably reproduce results from direct numerical simulations of passive scalar transport, validating their suitability for transport modeling in heterogeneous porous structures. Building on this foundation, Ref. 40 developed a PNM framework that integrates advective-diffusive transport with size-dependent retention mechanisms, enabling the simulation of aggregation-driven clogging. Ref. 49 employed a similar PNM to investigate the impact of colloid deposition on breakthrough curves and retention profiles across pore networks. In addition to deposition kinetics, Ref. 39 introduced straining and a hydrodynamic bridging criterion: clogging is predicted if the cumulative projected area of retained colloids exceeds the cross-sectional area of the pore throat. With their framework, they were able to investigate, by simulation, the impact of particle size and concentration, flow rate, and fluid chemistry on permeability reduction. Ref. 41 further advanced this modeling framework by considering the number of particles passing through a throat before forming an arch that plugs the pores, as proposed by Ref. 9. Ref. 50 proposed a probabilistic PNM in which the bridging probability is controlled by particle concentration, pore-to-particle size ratio, and friction, and is calibrated against CFD-DEM simulations. While these studies have significantly enhanced our understanding of pore-scale clogging mechanisms, the probability law governing this particle threshold does not explicitly account for pore geometry. Moreover, their predictive capabilities have yet to be validated against experimental observations, leaving a gap between modeling and empirical evidence.

In this paper, we investigate permeability impairment due to particle bridging using both microfluidic experiments and pore-scale modeling. The use of microfluidic devices that replicate the complexity of natural porous media,<sup>51</sup> combined with state-of-the-art imaging and post-processing methods,<sup>52</sup> allows for the assessment of permeability impairment caused by particle clogging. We present a novel pore-network model that incorporates the probabilistic nature of hydrodynamic pore bridging. It accounts for key parameters influencing bridging, including particle size and concentration, the particle-to-throat size ratio, and the pore-to-throat radius ratio. The structure of the paper is as follows. In Section 2, we introduce the developed pore-network model along with the CFD-DEM-based approach used to calibrate the bridging probability law. We also describe the microfluidic experiments conducted to support and validate the numerical results. In Section 3, we compare experimental and numerical results on the injection of suspensions into porous media, highlighting the emergence of pore bridging and the resulting permeability impairment. Finally, we conclude with a summary of the main findings and an outlook on future research directions.

## 2. Material and methods

In this section, we present the numerical and experimental framework developed to investigate hydrodynamic bridging in porous media. We begin by describing the geometry of the porous medium used in both simulations and experiments. Next, we introduce our novel pore-network model, which incorporates a probabilistic law for arch formation at pore entrances. We then detail the computational approach used to calibrate this probability law using CFD-DEM simulations. Finally, we describe the microfluidic experiments designed to assess the robustness and predictive capability of the model.

### 2.1. Heterogeneous porous geometry

We design a porous geometry that mimics the complexity and heterogeneity of subsurface environments, which is critical for studying clogging and particle transport. This structure is generated using a MATLAB-based image quilting algorithm that enables the creation of large, non-repetitive porous patterns by stitching together smaller texture patches.<sup>53</sup> The porous domain created with this method spans 12 mm in length and 6 mm in width, with an approximate porosity of 47%, and is shown at the center of Fig. 1a, along with an inset highlighting the detailed microstructure of grains (black) and pore space (white). To further characterize the domain, we use built-in MATLAB scripts to extract both the two-dimensional throat size distribution (Fig. 1b) and pore size distribution (Fig. 1c). These analyses reveal an arithmetic mean throat size of approximately 20 µm and a mean pore diameter of 73 µm. The distributions follow log-normal trends, with most probable (modal) values around 5 µm for throat size and 50 µm for pore diameter, as indicated by the red fits. It is important to note that the pore and throat size distributions are two-dimensional, as they are extracted from a planar geometry with a fixed channel depth of 20 µm in the z-direction. To ensure uniform flow distribution across the domain, we include intricate inlet and outlet networks upstream and downstream of the porous structure. The complete micromodel, including the central porous region and the boundary networks (Fig. 1a), is designed and assembled in KLayout, a CAD software for photomask design, with a total length of 36 mm.

### 2.2. Pore-network model with probability of bridging

**Flow model.** In a pore-network model, the microstructure is represented as a collection of pore bodies interconnected by cylindrical throats. The volumetric flow rate, $Q_{ij}$, through the throat connecting pore i to pore j is approximated using a Hagen-Poiseuille-type relation for rectangular sections:

$$Q_{ij} = \frac{\alpha_{ij} (D_{ij}, h)}{\mu L_{ij}} (P_i - P_j), \quad (1)$$

where $\alpha_{ij}$ is the geometric hydraulic conductance factor, h is the micro-model thickness, $D_{ij}$ is the two-dimensional diameter of the throat, $L_{ij}$ is its length, $\mu$ is the suspension viscosity, and $P_i$ and $P_j$ are the fluid pressure in pore i and j, respectively. For a throat with a rectangular cross-section, as commonly encountered in microfluidic devices, the geometric hydraulic conductance factor is $\alpha_{ij} = \frac{h D_{ij}^3}{12}$ if $D_{ij} < h$ and $\alpha_{ij} = \frac{h^3 D_{ij}}{12}$ if $D_{ij} > h$ corresponding to the parallel-plate (wide-slit) limit in each orientation. For solid volume fraction below the jamming transition, the suspension viscosity reads,<sup>54</sup>

$$\mu(C) = \mu_f \frac{\exp(-2.34C)}{(1 - C/C_m)^3}, \quad (2)$$

where $\mu_f$ is the fluid viscosity, C is the particle concentration and $C_m \approx 0.62$ is the maximum packing density.<sup>54</sup>

Due to mass conservation at each pore, the net flow into pore i must be zero, leading to the condition:

$$\sum_{j}^{N_i} Q_{ij} = 0, \quad (3)$$

where $N_i$ denotes the number of neighboring pores connected to pore i.

Constant pressures are applied at the inlet and outlet boundaries. These boundaries may correspond either to the inlet/outlet of the porous domain itself or to those of the intricate injection channels. In the latter case, the injection channels are represented using pore bodies and throats. All other boundaries parallel to the main flow direction are treated as no-flow boundaries.

The flow equations, combined with mass conservation at each pore and the imposed boundary conditions, result in a linear system that determines the pressure in all pore bodies. This system is solved numerically using MATLAB and is updated at every time step.
<figure><img src="image_4.png" /></figure>

<figure><img src="image_5.png" /><figcaption>Fig. 1. (a) An image of the micromodel showing intricate networks upstream and downstream of the central porous domain. The inset presents a magnified view of the porous structure, where solid grains appear in black and pore spaces in white. (b) Throat size distribution within the porous domain, with a red log-normal fit. The distribution is right-skewed, with a most probable (modal) throat size near 5 µm and an arithmetic mean of approximately 20 µm. (c) Pore diameter distribution of the porous domain, fitted with a red log-normal curve. The most probable pore diameter is around 50 µm, while the arithmetic mean is approximately 73 µm. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)</figcaption></figure>

**Transport model.** We use the Mixed-Cell Method to solve the advection-dispersion equation governing the evolution of the particle volume fraction throughout the network. In this method, pores are seen as control volumes over which mass balance is imposed. We have,

$$V_i \frac{dC_i}{dt} = \sum_{j=1}^{N_i} Q_{ij}^c, \qquad (4)$$

where $Q_{ij}^c$ is the mass flux through the pore throat $ij$. This flux includes both advective and diffusive contributions. Using a first-order upwind scheme for the advective term, the total flux is expressed as:

$$Q_{ij}^c = H(\delta P_{ij}) Q_{ij} C_j + H(-\delta P_{ij}) Q_{ij} C_i + D_{ij} S_{ij} \frac{C_j - C_i}{L_{ij}}, \qquad (5)$$

where $H(\cdot)$ is the Heaviside function used to implement upwinding based on the sign of pressure difference, $\delta P_{ij}$, in the throat,⁵⁵ $S_{ij} = hD_{ij}$ is the throat cross-sectional area, and $D_{ij}$ is the effective diffusion coefficient. The latter accounts for Taylor-Aris dispersion in the throat⁵⁶⁻⁵⁸ using,

$$D_{ij} = D_m + \frac{\min(h, D_{ij})^2}{4D_m} U_{ij}^2 \left( \frac{2}{105} - \sum_{n=1}^{\infty} \frac{18}{(n\pi)^6} \exp \left( -4 \left( \frac{n\pi}{\min(h, D_{ij})} \right)^2 D_m t \right) \right),$$

(6)

where $U_{ij} = \frac{Q_{ij}}{S_{ij}}$ is the average velocity in the throat and the particle diffusion coefficient, $D_m$, is given by Stokes-Einstein relationship, $D_m = \frac{k_B T}{3\pi\mu_f d_p}$, $k_B$ and $T$ being the Boltzmann constant and temperature, respectively.

The jamming transition, $C_{jam}$, bounds the particle concentration and is defined as the concentration at which the suspension becomes immobile. Its value is always lower than the packing density, $C_m$, and depends on the shear stress, and therefore varies with pore geometry, particle size and concentration, and flow rate.⁵⁹⁻⁶¹

**Hydrodynamic bridging criterion.** Bridging occurs when a stable arch of particles forms across a pore throat. The number of particles, $n$, required to form such an arch depends on the particle size and the throat diameter. It is estimated as:

$$n = \left\lfloor \frac{D}{d_p} \right\rfloor + 1, \qquad (7)$$

where $\lfloor \cdot \rfloor$ is the floor function (i.e., the greatest integer less than or equal to the argument) and $D/d_p$ is the throat-to-particle diameter ratio.
The probability of forming a stable arch that clogs a pore is governed by $N^*$, the number of particles that pass through a throat before clogging occurs. Ref. 9 has demonstrated that this number is independent of the flow rate. This number depends on the probability of clogging according to Ref. 62,

$$N^* = \frac{(1 - P_{clog})}{P_{clog}}. \qquad (8)$$

Building on our recent work,²² the probability of forming a stable arch of n particles reads,

$$P_{clog}(n) = k \frac{[(n-1)C]^{n-1}}{(n-1)! \times C_m^{n-1} \times \exp\left(\frac{(n-1)C}{C_m}\right)} \left(\frac{d_p}{D}\right)^n \left(\sin \theta\right)^a \exp\left(\left(\frac{C}{C_m}\right)^n\right). \qquad (9)$$

This discrete stochastic model for predicting pore-clogging probabilities was established under the assumption of a creeping-flow regime (low Reynolds number) for a suspension of monodisperse spherical particles. We also assume smooth particle surfaces and smooth pore walls (no roughness effects). The model was established and validated in Ref. 22 using more than 700 CFD-DEM simulations in which particle size, particle concentration, particle-to-throat diameter ratio, and constriction angle were systematically varied. We verified that the predictions of this law are consistent with the probability density function of particle escaping before clogging obtained from both the CFD-DEM simulations and the experimental datasets of Refs. 14, 20, 21. In contrast with the probability laws proposed in the literature,¹⁴,¹⁹,²⁰ Eq. (9) depends on the particle concentration C, the maximal packing concentration C<sub>m</sub>, and geometric parameters including the throat angle θ, and the particle-to-throat diameter ratio D/d<sub>p</sub>. The PNM computes the particle concentration within the throat at each time step. The throat angle is determined from the pore geometry and is given by tan θ = R<sub>i</sub>/R<sub>ij</sub>/L<sub>ij</sub>/2, where R<sub>i</sub> and R<sub>ij</sub> are the pore and throat radius, respectively. Eq. (9) involves two fitting parameters, k and α, that reflect the packing morphology, the shape of the channel, and the particle friction. We have demonstrated that these parameters are independent of flow rate, particle size, throat angle, and particle concentration, in the limit of the semidilute suspensions.²² They will be calibrated in Section 3.1 using CFD-DEM simulations of suspension injection through a single pore.

In the PNM, we begin by computing the number of particles required to form a clogging arch in each throat ij, using Eq. (7). This calculation is performed only once at the beginning of the simulation, as the suspension is monodisperse. At each time step, we evaluate the clogging probability for every throat using Eq. (9), and determine the number of particles that must pass the constriction before clogging, $N^*$, using Eq. (8). If a throat is smaller than the particle diameter, we set $N^* = 1$. Simultaneously, we track the cumulative number of particles that have passed through each pore. Once this number exceeds $N^*$ the likelihood of clogging is considered high, and the throat conductance factor is replaced by $\alpha_{ij} = D_{ij}h \frac{d_p^2(1-C_{jam})^3}{180C_{jam}^2}$ which corresponds to an estimate of the permeability of a throat clogged by packed bed of beads, based on the Kozeny-Carman relationship.

The time step is selected as the minimum value required to update the mass flux in each pore, subject to a Courant-Friedrichs-Lewy (CFL) stability condition with $CFL \le 0.5$.

## 2.3. CFD-DEM for particulate transport

The fitting parameters in the clogging probability law (Eq. (9)) are calibrated numerically using four-way coupling CFD-DEM simulations of particle bridging at the entrance of a single pore throat. The simulations are carried out using an in-house CFD-DEM coupling framework developed within the open-source finite-volume library OpenFOAM

(version 9, https://www.openfoam.org). The key features of the implementation are summarized below, and a detailed description of the numerical algorithms can be found in Ref. 26.

The fluid motion is modeled in two dimensions by solving the depth-integrated incompressible Volume-Averaged Navier-Stokes (VANS) equations,⁶³ which consist of the mass conservation equation for the fluid phase:

$$\frac{\partial(\epsilon\rho_f)}{\partial t} + \nabla \cdot (\epsilon\rho_f \mathbf{v}_f) = 0, \qquad (10)$$

and the volume-averaged momentum conservation equation:

$$\frac{\partial (\epsilon \rho_f \mathbf{v}_f)}{\partial t} + \nabla \cdot (\epsilon \rho_f \mathbf{v}_f \mathbf{v}_f) = -\epsilon \nabla p + \nabla \cdot (\epsilon \mu_f (\nabla \mathbf{v}_f + (\nabla \mathbf{v}_f)^T)) - \frac{12\mu_f}{h^2} \epsilon \mathbf{v}_f - \frac{\epsilon \mu_f}{K_0 (1 - \epsilon)^2} (\mathbf{v}_f - \bar{\mathbf{v}}_p), \qquad (11)$$

where ε is the volume fraction of fluid in a grid block, ρ<sub>f</sub> is the fluid density, v<sub>f</sub> is the cell-averaged fluid velocity, p is the fluid pressure, μ<sub>f</sub> is the fluid dynamic viscosity, v̄<sub>p</sub> is the mean particle velocity on the Eulerian grid, and K<sub>0</sub> is a predefined local cell-permeability. The Hele-Shaw correction term, 12μ<sub>f</sub>v<sub>f</sub>/h<sup>2</sup>, accounts for the hydrodynamic effects in the thickness, h, of the micromodel.⁶⁴ The last term represents the drag force exerted by the fluid on suspended particles. This force vanishes in pure fluid regions and increases sharply as the local solid volume fraction approaches unity, effectively penalizing fluid motion in fully occupied cells.⁶⁵

Particle motion is governed by Newton's second law, solved in a Lagrangian framework where trajectories evolve under the combined influence of fluid forces, particle-particle contacts, and interactions with solid boundaries. The total velocity, v<sub>k</sub><sup>p</sup>, of particle k reads,⁶⁶

$$\left\{ \begin{aligned} \mathbf{v}_k^p &= U_k^p + \omega_k \times \mathbf{r}_k, \\ m_k \frac{dU_k^p}{dt} &= \sum_l \mathbf{F}_{kl}^c + \mathbf{F}_k^f, \\ I_k \frac{d\omega_k}{dt} &= \sum_l \mathbf{M}_{kl}^c + \mathbf{M}_k^{hyd} \end{aligned} \right. \qquad (12)$$

where U<sub>k</sub><sup>p</sup> and ω<sub>k</sub> are the translational and the angular velocity of particle k, respectively, and r<sub>k</sub> is the position vector. m<sub>k</sub> and I<sub>k</sub> are the particle mass and moment of inertia, respectively. The drag force, F<sub>k</sub><sup>f</sup>, is computed from the fluid-particle shear stress acting on the particle k surfaces following the immersed boundary method of Ref. 67. F<sub>kl</sub><sup>c</sup> and M<sub>kl</sub><sup>c</sup> represent the contact forces and torques resulting from particle-particle and particle-wall interactions, while M<sub>k</sub><sup>hyd</sup> denotes the hydrodynamic torque. The summations run over N<sub>k</sub><sup>c</sup> objects in contact with particle k, including other particles and the confining geometry. The specific force and torque models used are described in detail in Ref. 26.

The time-stepping strategy accounts for three numerical stability criteria: (i) a criterion for the pressure-velocity coupling algorithm, (ii) a constraint due to the calculation of colliding particles, and (iii) a stability criterion related to the fluid-particle interaction. A full description is found in Ref. 26.

## 2.4. Microfluidic design and experiments

**Microfabrication of the chip.** A master wafer of the micromodel geometry is fabricated by standard photolithography (micro and nanotechnologies platform, LAAS-CNRS). Our microfluidic devices are then fabricated using PDMS (polydimethylsiloxane), a widely used elastomer in soft lithography techniques.⁶⁸ The PDMS is prepared by mixing the base and curing agent in a 10:1 weight ratio. After thorough mixing, the mixture is degassed under vacuum to remove air bubbles, then poured over the master wafer and cured in an oven at 65 °C for 1.5 h. Once cured, the PDMS replica is trimmed and punched to create the inlet and outlet holes. It is then cleaned in an ultrasonic bath for 30 min. For
<figure><img src="image_6.png" /><figcaption>**Fig. 2.** Schematic of the experimental setup used for flow experiments. A pressure controller drives the injection of particle suspension from a reservoir through tubing equipped with a one-way valve, which facilitates fluid switching without introducing air bubbles. The flow enters the microfluidic device mounted on an inverted microscope platform, where image sequences are recorded over time using a high-speed camera. A flow rate sensor positioned downstream continuously monitors the decline in flow rate due to pore clogging. The fluid then exits into a disposal vial.</figcaption></figure>

final assembly, the PDMS replica and a glass slide are placed in an air plasma chamber for 35 s to activate their surfaces and enable strong, irreversible bonding. To ensure stable and consistent hydrophobicity during the experiments, the assembled device is placed in an oven at 120 °C for at least 3 h prior to use.¹⁸

**Particle suspension.** The particle suspensions consist of 78% Millipore Milli-Q water and 22% glycerol by weight, forming a density-matched fluid with a density of approximately 1.05 g/cm³, equal to that of the polystyrene particles. This prevents sedimentation during the experiments. The particles in use are monodisperse Polybead® Carboxylate Microspheres (polystyrene particles) of either 4.5 µm or 6 µm in diameter, manufactured by Polysciences, Inc., and are prepared at volume fractions of 0.2% and 0.3%, depending on the experimental conditions. The viscosity of the final mixture is approximately 1.45 mPa s.

**Flow experiments.** The experiment begins by saturating the microfluidic device with a particle-free, density-matched fluid to eliminate trapped air and ensure consistent initial conditions (see Fig. 2 for a schematic of the experimental setup). Flow is driven using a pressure controller (FLUIGENT EZ FLOW) connected to a reservoir containing the particle suspension, with experiments conducted at two different pressure drops: 100 mbar and 150 mbar. Since deformation of the PDMS structure could influence particle deposition and clogging behavior, ultimately altering the hydraulic resistance of the microfluidic device by dilating the channel walls, we conducted a series of calibration experiments to identify the deformation threshold. Based on these observations, all experiments were conducted at inlet pressures

below this threshold, identified to be around 250 mbar. The suspension is delivered through 0.5 mm inner diameter Tygon tubing, which passes through a one-way valve. This valve is essential when switching between particle-free fluid and particle suspension, as it prevents air from entering the flow path. Once the particle suspension reaches the entrance of the porous domain, the experiment is initiated. The flow rate is continuously monitored over time using a Coriolis flow sensor (Bronkhorst) positioned downstream of the microfluidic device. The flow sensor records data at 1-second intervals, but due to sensor-induced fluctuations, we apply an additional smoothing step using a Python script. A moving average with a 50-second window is used to consistently smooth the flow rate data across all experiments.⁶⁹ The experiment continues until full clogging occurs, which is defined as the point at which the flow rate drops to 85% of its initial value.

**Image analysis.** To monitor the evolution of clogging over time, image sequences are recorded at selected time points throughout the experiment using an inverted microscope equipped with a ×2 objective and a FASTEC high-speed camera. The images are captured at a frame rate of 300 fps, allowing us to track the progression of particle deposition and clogging dynamics. Since the field of view of the camera covers less than the full length of the porous domain, three separate sequences are acquired at consecutive locations along the flow direction (x-axis). Each sequence is temporally averaged to highlight only the immobile (deposited) particles while removing those in motion. The resulting averaged images are then stitched together using a standard image stitching algorithm⁷⁰ to reconstruct a complete 2D view of the porous domain. To quantify particle deposition along the x-axis of the porous domain, we develop a Python code that processes the final stitched images of the domain obtained under different experimental conditions. For consistency, we apply the same gray-scale contrast settings to all stitched images, regardless of the experimental condition, once clogging is reached — defined here as a decrease in flow rate to 85% of its initial value. In the analysis, all pixels with intensity values below 30 are classified as particles. The stitched image of the porous domain is then divided into ten equal segments along the x-direction, and for each segment, the number of particle pixels is calculated and normalized by the total number of pixels in the entire stitched image. The resulting distribution is plotted, allowing direct comparison of particle deposition patterns along the x-direction for different experimental conditions.

### 3. Results and discussion

In this section, we investigate permeability impairment in a heterogeneous medium during the injection of a monodisperse suspension. We begin by calibrating the clogging probability law using CFD-DEM simulations. Next, we conduct microfluidic experiments to investigate the impact of injection pressure, particle size, and concentration on clogging dynamics. Finally, this high-resolution experimental dataset is used to assess the predictive robustness of the proposed pore network model.

#### 3.1. Calibration of the arch formation probability law

The probability law used to predict pore plugging (Eq. (9)) includes two fitting parameters, namely *k* and *a*. In previous work, we demonstrated that they are relatively independent of flow rate and particle concentration, but strongly dependent on the pore geometry.²² To calibrate these parameters, we perform CFD-DEM simulations of particle bridging in a single pore for a wide range of throat angles, $\theta$, and particle-to-throat diameter ratio, $D/d_p$.

The computational domain consists of a two-dimensional channel made of a converging and diverging section with a central constriction (see Fig. 3a). The walls of the throats follow the trajectory $y(x) = (D_{ch} + D)/4 + (D_{ch} - D)\cos(2\pi x/L)/4$, where $D_{ch}$ is the channel diameter, *D* is the throat diameter, and *L* is the length of the throat.
<figure><img src="image_7.png" /><figcaption>**Fig. 3.** CFD-DEM simulations of hydrodynamic bridging in a single pore. (a) Geometry of the depth-integrated pore used to investigate particle bridging. The suspension is injected from left to right. (b) Influence of the throat angle, $\theta$, on the formation of particle arches. (c) Effect of the particle-to-throat size ratio, $D/d_p$, on bridging dynamics. (d) Average number of particles, $N^*$, that escape before clogging as a function of the throat angle. (e) Variation of $N^*$ with particle-to-throat size ratio, $D/d_p$. The CFD-DEM simulations enable calibration of the fitting parameters in the stochastic bridging model, yielding $k = 0.3$ and $\alpha = 0.6$.</figcaption></figure>

The channel geometry is discretized using a structured hexahedral grid with 200 × 40 cells. The grid size is at least 3 times finer than the particle diameter to ensure numerical convergence, following the guidelines established in Ref. 26. The suspension is injected from the left (inlet) to the right (outlet) under a constant pressure difference, $\Delta p$. Particles are introduced randomly at the inlet, allowing for statistical analysis through multiple simulation replicas. For each configuration of $(\theta, D/d_p)$, we perform 15 simulation replicates. The randomized particle positions at the inlet lead to variability between runs, allowing us to compute ensemble-averaged clogging probabilities along with their standard deviations. The simulation parameters are set to $C = 0.125$, and $Re_p = 0.05$.

We first investigate the influence of the throat angle, $\theta$, on bridging by varying it within the range [3°:50°] while the particle-to-throat diameter ratio fixed at $D/d_p = 1.72$. To maintain both this ratio and a constant channel length, $L$, across different angles, we adjust the channel diameter, $D_{ch}$, while keeping the throat diameter, $D$, constant. Fig. 3b shows pore geometries for selected values of $\theta = 5.7°$, 17°, and 31°. In all cases, arches are principally made of $n = 2$ particles. However, no clogging is observed for $\theta = 3°$ even after 800 particles have escaped through the throat. As shown in Fig. 3d, the average number of particles that escape before clogging, $N^*$, decreases nonlinearly with increasing throat angle until it reaches a plateau for $\theta > 40°$, in agreement with previous experimental observations.²¹ Notably, this plateau value is nearly four times lower than that observed at $\theta \approx 5.7°$, indicating that larger throat angles promote earlier bridging. The stochastic model for $N^* = (1 - P_{clog}/P_{clog})$ in which $P_{clog}$ is provided by Eq. (9) shows good agreement with the CFD-DEM results when the fitting parameters are set to $k = 0.3$ and $\alpha = 0.6$.

The role of the particle-to-throat size ratio, $D/d_p$, on bridging is analyzed by varying the particle diameter $d_p$ while keeping the throat

angle fixed at $\theta = 17°$. We consider values of $D/d_p$ ranging in [1.25; 1.92] which corresponds to 260 µm $\le$ $d_p$ $\le$ 400 µm, with most clogging events involving primarily two-particle arches. Snapshots of particle bridging for $D/d_p = 1.93$, 1.67, and 1.43 are presented in Fig. 3c. Although two-particle arches dominate across this range, we occasionally observe three-particle arches for configurations approaching $D/d_p \approx 2$. Fig. 3e presents $N^*$ as a function of $D/d_p$. It follows a staircase-like increase along with an increase in $D/d_p$, consistent with the experimental observations reported by Ref. 20. Our stochastic model captures this behavior when using $k = 0.3$ and $\alpha = 0.6$.

These CFD-DEM simulations, performed on a single depth-integrated two-dimensional pore geometry, allow for calibration of the stochastic model with $k = 0.3$ and $\alpha = 0.6$. The calibrated probability laws will be used in the PNM simulations to model permeability impairment in porous media representative of our microfluidic chip. This calibration procedure should be repeated for other pore configurations, including fully three-dimensional geometries and non-spherical particles.

## 3.2. Permeability impairment of hydrodynamic pore-bridging via microfluidics

Microfluidic experiments are conducted to investigate how injection pressure, particle size, and particle concentration affect permeability impairment and particle deposition patterns resulting from hydrodynamic pore bridging. A total of five experiments are performed: a base case with a replicate to assess reproducibility, and three additional cases in which flow conditions and suspension properties are systematically varied. Starting from the base case, each parameter is varied by a factor of 1.5, except for particle size, which is varied by a factor of 1.33, to assess the influence of each parameter on permeability impairment fairly. The parameters for each case is summarized in Table 1. Each
<table><caption>**Table 1**

Cases investigated both using microfluidic experiment and pore-network modeling.</caption><thead><tr><th>Plot color</th><th>Pressure drop<br>Δp (in mbar)</th><th>Concentration<br>C (in %)</th><th>Particle size<br>d<sub>p</sub> (in μm)</th><th>Jamming concentration<br>C<sub>jam</sub> (in %)</th></tr></thead><tbody><tr><td>Blue</td><td>150</td><td>0.3</td><td>6</td><td>0.6115</td></tr><tr><td>Dashed blue</td><td>150</td><td>0.3</td><td>6</td><td></td></tr><tr><td>Red</td><td>150</td><td>0.3</td><td>4.5</td><td>0.6135</td></tr><tr><td>Green</td><td>150</td><td>0.2</td><td>6</td><td>0.6152</td></tr><tr><td>Purple</td><td>100</td><td>0.3</td><td>6</td><td>0.6125</td></tr></tbody></table>

experiment is initiated by saturating the porous micromodel with a particle-free, density-matched fluid to ensure consistent initial conditions. The particle-laden suspensions are then injected at controlled pressure drops using a pressure-driven system.

**Base configuration.** A suspension containing 0.3% of 6 µm particles is injected into the heterogeneous micromodel by applying a pressure difference of ΔP = 150 mbar, and the injection is maintained until complete clogging occurs. Full clogging is observed approximately 10 min after the start of the injection. Time-lapse images of the experiment are presented in Fig. 4. Yellow arrows in the magnified insets of the experimental images indicate the last constriction to become clogged before complete blockage of the micromodel occurs.

At early times, particles flow freely through the porous domain. Some particles are deposited at stagnation points around solid grains, while a few others become trapped by sieving in pore throats narrower than 6 µm, as shown in the magnified inset. After 131 s, several black regions appear, indicating throats that have become blocked, followed by the progressive filling of the associated pore spaces with particles depositing on both sides of the constriction marked by yellow.

As time progresses, the accumulation intensifies, and larger regions become densely packed with particles. Clogging becomes increasingly pronounced near the inlet, where a filter cake begins to form after 394 s. Regarding the constriction of interest, we observe particles depositing and spanning the 2D cross-section of the constriction; however, a light black streamline is still visible downstream in the averaged images at 394 s and 525 s, indicating residual particle flow.

Multiple filter cakes subsequently develop and expand until, by 660 s, the inlet is completely obstructed. In addition, the constriction of interest shows a clear and complete blockage of particle flow. The spatial profile of deposited particle concentration along the longitudinal axis at the time of clogging is shown by the blue curve in Fig. 5, confirming that most clogging occurs near the inlet, within the first quarter of the micromodel's length, which corresponds to the field of view shown in the time-lapse images.

This experiment was repeated once to assess the reproducibility of the results. The two experiments show excellent agreement. As shown in Fig. 6, the permeability reduction follows nearly identical trends in both cases, with complete clogging occurring after approximately 10 min. A comparison of the particle deposition maps at full clogging (see Fig. 5) reveals a surprisingly high level of similarity between the two replicates: the clogged regions and blocked pores are nearly identical. This consistency is likely due to the disordered nature of the porous medium, which tends to yield more reproducible flow patterns than ordered or highly structured systems.⁷¹ This replication reinforces our confidence in the reproducibility of results obtained for other experimental configurations.

**Effect of injection pressure, particle size, and concentration on deposition patterns.** In addition to the base case, three further experiments are carried out, each varying a single parameter — flow rate, particle size, or particle concentration — by a factor of approximately 1.5. The results are summarized in Fig. 5. For each experimental condition, the progression of particle deposition over time is shown for the entire porous domain. We also show the mean deposition profiles along the x-axis after complete clogging. They correspond to the dotted boxes.

The experiment performed at a lower pressure difference, shown in violet, exhibits only a minimal filter cake near the inlet, with a greater amount of deposition occurring in the second quarter of the domain. The reduced pressure difference slows the clogging process, enabling particles to penetrate deeper into the porous medium and accumulate within internal pores before eventual blockage at the inlet. This observation aligns with the concept that a constant number of particles, N*, must traverse a pore before an arch forms to clog a throat.⁹ Since N* remains constant, a lower flow rate results in a longer time required for clogging to occur.

Among all experiments, the one carried out with a particle concentration reduced by a factor of 1.5 (see green plots) exhibits the highest overall deposition in the porous domain. It has the most pronounced accumulation at the inlet, forming a large filter cake. This behavior can be explained by the fact that a high pressure difference combined with 6 µm particles promotes the formation of a filter cake, while the reduced concentration delays the onset of clogging in certain constrictions. As a result, particles continue to deposit on the previously formed filter cakes near the inlet while simultaneously migrating deeper into the domain, progressively filling pores, particularly within the first half of the porous structure, until complete clogging of the flow is eventually reached.

The experiment shown in red uses particles with a diameter of 4.5 µm instead of 6 µm in the base case. The deposition profile follows a similar overall trend with a more pronounced filter cake near the inlet and a significant deposition in the second quarter of the domain. Since the particles are smaller by a factor of 1.33, a larger number of them, n, are required to form an arch. It follows that the likelihood, $P_{clog}(n)$, to form a stable arch decreases, resulting in a longer delay to achieve complete clogging.

From these experiments, a general trend emerges: higher injection pressures lead to greater particle deposition near the inlet,⁷²,⁷³ resulting in the formation of filter cakes whose intensity varies depending on the other parameters specific to each experiment. Regarding particle size and particle concentration, both parameters strongly influence the intensity of deposition in the porous domain before complete clogging occurs. A higher amount of deposition is observed in the case of a lower suspension concentration, which is consistent with the understanding that a larger number of particles is required to form a stable bridging arch across a single pore throat.²⁰ This delays the onset of complete clogging, allowing more particles to migrate and accumulate behind already clogged pores. In contrast, reducing the particle size relative to the base configuration results in a higher overall deposition than the base case, yet still lower than the deposition observed in the experiment with a reduced particle concentration. This is attributed to the smaller particle diameter, which increases the ratio of pore width to particle size and thus requires more time for clogging to develop.⁹,²⁰,⁷⁴

**Permeability impairment.** Real-time flow rate measurements for the different experimental conditions are acquired downstream and processed using a moving-average smoothing algorithm to reduce sensor noise. The resulting evolution of permeability impairment is presented in Fig. 6 for the 4 experimental conditions.

The base configuration depicted in blue shows the fastest permeability decline, with complete clogging occurring after approximately 750 seconds. The replication of this experiment yields an identical
<figure><img src="image_8.png" /><figcaption>Fig. 4. Time-lapse images of the base configuration in the first section of the porous domain, showing how particle deposition progressively develops over time until the domain becomes fully clogged. The insets show a constant magnified view of a selected region highlighted by a yellow arrow, indicating a constriction where particles gradually accumulate over time until the pore is completely blocked after 660 s. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)</figcaption></figure>

permeability impairment curve. Reducing the injection pressure by a factor of 1.5 (violet curve) results in a decline that is approximately twice as slow. A 1.5-fold reduction in particle concentration (green curve) leads to an even slower decline, about seven times slower than in the base case. The slowest permeability decline is observed for the experiment using particles that are 1.33-fold smaller (red curve), with complete clogging occurring approximately seventeen times later than in the base case.

The stepwise sharp declines observed in the permeability curves are likely associated with the progressive clogging of major flow paths within the porous domain. Each steep drop corresponds to the complete blockage of a main flow line, with the final steep drop marking the clogging of the last remaining primary flow path. This pattern appears in almost all experiments, except in the baseline case, where the permeability shows a consistently steep decline from the start.

These results indicate that, among the parameters investigated, particle size has the strongest influence on permeability impairment,

followed by particle concentration, and finally injection pressure. This experimental dataset will serve as a reference to assess the ability of our probabilistic PNM to capture permeability impairment caused by particle bridging.

### 3.3. Permeability impairment in PNM with probabilistic bridging law

Our novel pore-network model is employed to simulate pore clogging induced by particle bridging. The two-dimensional network is extracted from the micromodel geometry defined in Fig. 1. The resulting network spans 12 mm × 6 mm and is made of 5835 pores and 8952 throats. Consistent with the microfluidic experiments, four cases are considered: a base configuration and three additional cases in which a single parameter is varied relative to the base case (see Table 1).

The base case involves the injection of a suspension containing 0.3% particles with a diameter of 6 µm under a constant pressure difference of 150 mbar applied between the inlet and outlet. The time step is constraint by CFL ≤ 0.5. The jamming transition, $C_{jam}$, is obtained
<figure><img src="image_9.png" /><figcaption>**Fig. 5.** Time-lapse averaged stitched experimental images of the porous domain under different experimental conditions, with each condition distinguished by a separate color. The corresponding plot shows the normalized particle deposition fraction from the final stitched images, highlighted by dotted boxes, as a function of position along the x-direction for each experimental condition. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>**Fig. 6.** Evolution of permeability impairment over time for the different experimental conditions, each represented by a different color. The light black horizontal line indicates the 85% decline in initial permeability, while the vertical lines in corresponding colors mark the time at which each experiment reaches this threshold, indicating the onset of complete particle clogging. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)</figcaption></figure>

by trial and error. Snapshots of the baseline simulation are shown in Fig. 7. At the beginning of the simulation, the suspension is transported through the network. The probability of clogging and the expected number of particles, $N^*$, passing through a throat before clogging are computed at every time step for every throat. Once the cumulative number of particles per throat goes above $N^*$, the corresponding throat is clogged, and the conductance of the throat decreases significantly. After a throat is clogged, the upstream adjacent pore accumulates particles until the packing value is reached. Over time, clusters of clogged pores and throats develop within the network. The evolution

of permeability is presented in plain blue line in Fig. 8. We observe an excellent agreement with the microfluidic measurements for $C_{jam} = 0.6115$: the curves almost overlap, and both predict an 80% reduction in permeability after approximately 600 s.

Three additional PNM simulations are performed, each varying one parameter — particle size, particle concentration, or the pressure difference between the inlet and outlet. First, we use particle of 4.5µm instead of 6µm. The resulting permeability impairment develops more slowly than in the base case, consistent with experimental observations (see red curve in Fig. 8). This behavior is expected because a greater
<figure><img src="image_11.png" /><figcaption>**Fig. 7.** Snapshots of pore clogging obtained by probabilistic pore-network modeling for the baseline case. The throats in red are clogged. $t_1$ = 0 s, $t_2$ = 50 s, $t_3$ = 100 s, and $t_4$ = 700 s.. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)</figcaption></figure>

<figure><img src="image_12.png" /><figcaption>**Fig. 8.** Evolution of permeability impairment for the different conditions, each represented by a different color, obtained with PNM simulations (dashed line) and microfluidic experiments (plain line). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)</figcaption></figure>

number of smaller particles, $n$, is required to form a stable bridging arch, thereby reducing the overall probability of clogging. Similar to microfluidic results, the curve displays a stepwise decline, reflecting the fact that when a primary flow pathway is blocked, the fluid is forced to reroute through secondary channels with higher hydraulic resistance. Then, we inject a suspension with 0.2% of particles instead of 0.3%. The resulting permeability decline curve is plotted in green in Fig. 8. Finally, the pressure difference between the network inlet and outlet was reduced by a factor 1.5 (see violet curve in Fig. 8). In agreement with microfluidic experiments, the permeability reduction is delayed compared with the baseline, and the effect is less pronounced than for the small-particle case.

In all the cases, the permeability exhibits a stair-wise decline reflecting the blockage of preferential flow paths and the subsequent diversion of the suspension toward more resistive pore throats. Some discrepancies are nevertheless observed in the sequence of local clogging events, as evidenced by abrupt permeability drops that do not perfectly

coincide with those measured experimentally. These differences can be attributed to several factors. First, the heterogeneous microstructure – in particular the pore geometry – is only approximately captured during the network extraction, and the transport laws used in the 2D model are themselves subject to modeling uncertainties. Second, additional physical mechanisms, such as shear-induced migration, may influence particle transport at bifurcations between adjacent throats but are not accounted for in the PNM simulations. Moreover, the jamming criterion in the PNM is imposed globally over the entire network, whereas it should ideally be defined locally at the throat scale, depending on constriction geometry and local shear conditions. Despite these limitations, the predicted clogging timescale remains in the same order of magnitude as that observed experimentally.

Overall, the simulations reproduce permeability impairment trends consistent with those observed in the microfluidic experiments. This set of results demonstrates the capability of the PNM approach, when coupled with an appropriate probability law for arch formation, to
robustly capture the clogging dynamics by hydrodynamic bridging. Most previous PNMs were restricted to simulating clogging through sieving and the accumulation of particle layers driven by aggregation. By introducing a probabilistic treatment of arch formation, our approach now enables PNMs to capture the full spectrum of pore-clogging processes including bridging.

## 4. Conclusion

Permeability impairment in heterogeneous porous media remains a long-standing challenge with significant implications for numerous industrial and environmental processes. In this study, we used both microfluidic experiments and pore-scale modeling to investigate the mechanisms of clogging induced by particle injection.

We introduce a novel probabilistic pore-network model that quantifies the likelihood of stable particle arch formation within pore throats. The law of probability accounts for the particle-to-throat diameter ratio, the particle concentration and the pore constriction angle. Its fitting parameters are calibrated in a data-driven manner using high-fidelity CFD-DEM simulations. This constitutes a key advantage over previous probabilistic models, as it eliminates ad-hoc or case-specific tuning and yields a framework that is more predictive, generalizable, and physically interpretable.

Both the microfluidic experiments and PNM simulations exhibit similar trends in permeability impairment. First, the permeability reduction follows a stepwise pattern, reflecting the sudden blockage of flow paths caused by hydrodynamic bridging. Second, the effects of flow rate, concentration and particle size are well captured by the model and align with the experimental observations: smaller particles significantly delay permeability reduction, and this effect is more pronounced than that of a reduced pressure difference. Overall, the simulations show very good agreement with the microfluidic results, confirming the validity of the probabilistic PNM developed in this study to model clogging driven by hydrodynamic bridging.

State-of-the-art PNMs have so far been able to simulate pore clogging through sieving and the progressive aggregation of particles at pore walls, but they have not captured hydrodynamic bridging properly. This study demonstrates that probabilistic PNMs can successfully represent this mechanism, thereby extending the capability of PNMs to describe all three pore-clogging processes during the injection of particle-laden flows in porous media. The simulations, however, were conducted in two dimensions. While extending the PNM framework to three dimensions is straightforward for solving flow and transport, the validity of the probability law in 3D remains uncertain. Further work is needed to clarify this aspect. In this context, CFD-DEM simulations can be employed to calibrate the probability law in three dimensions, and in-situ imaging techniques, such as time-lapse X-ray microtomography, will provide high-resolution datasets to validate PNM predictions.

### CRediT authorship contribution statement

**Cyprien Soulaine:** Writing – original draft, Validation, Supervision, Software, Project administration, Methodology, Investigation, Fund-ing acquisition, Formal analysis, Conceptualization. **Walid Okaybi:** Writing – original draft, Visualization, Validation, Methodology, Investi-gation, Formal analysis, Conceptualization. **Laurez Fogouang Maya:** Writing – original draft, Visualization, Validation, Software, Method-ology, Investigation, Formal analysis, Conceptualization. **Emmanuel Le Trong:** Writing – review & editing, Visualization, Validation, Soft-ware, Methodology, Investigation, Formal analysis, Conceptualization. **Sophie Roman:** Writing – review & editing, Validation, Supervision, Methodology, Conceptualization.

### Declaration of Generative AI and AI-assisted technologies in the writing process

During the preparation of this manuscript, the author(s) used Chat-GPT in order to improve language and readability. After using this tool/service, the author(s) reviewed and edited the content as needed and take(s) full responsibility for the content of the publication.

### Declaration of competing interest

The authors declare the following financial interests/personal relationships which may be considered as potential competing interests: Cyprien Soulaine reports financial support was provided by European Research Council. If there are other authors, they declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

### Acknowledgments

This research was funded by the European Research Council (ERC) through the ERC COCONUT grant (101043288), awarded to Cyprien Soulaine. Experiments have been performed at NanoU Lab, a nano-microfluidic laboratory in Orléans (France) that is part of the ECCSEL consortium and of the project StoRIES (STorage Research Infrastructure Eco-System)-European Research Infrastructures capacities and services to address European Green Deal challenges. This work was also supported by LAAS-CNRS micro and nanotechnologies platform, a member of the Renatech French national network. The authors benefited from the Centre de Calcul Scientifique en région Centre-Val de Loire's (CaSciMoDOT) high-performance computational resources.

### Data availability

The pore-network model the MATLAB script used in this article are available in Zenodo at REF<sup>75</sup>. The microfluidic dataset is available on the Digital Rock Portal REF<sup>76</sup>.

### References

1. Khilar KC, Fogler HS. The existence of a critical salt concentration for particle release. *J Colloid Interface Sci.* 1984;101(1):214-224. http://dx.doi.org/10.1016/0021-9797(84)90021-3, URL https://ui.adsabs.harvard.edu/abs/1984JCIS.101.214K.

2. Cerda CM. Mobilization of kaolinite fines in porous media. *Colloids Surf.* 1987;27(4):219-241.

3. Khilar KC, Vaidya RN, Fogler HS. Colloidally-induced fines release in porous media. *J Pet Sci Eng.* 1990;4(3):213-221. http://dx.doi.org/10.1016/0920-4105(90)90011-Q, URL https://ui.adsabs.harvard.edu/abs/1990JPSE...4..213K.

4. Dressaire E, Sauret A. Clogging of microfluidic systems. *Soft Matter.* 2017;13:37-48.

5. Sauret A, Barney EC, Perro A, Villeraux E, Stone HA, Dressaire E. Clogging by sieving in microchannels: Application to the detection of contaminants in colloidal suspensions. *Appl Phys Lett.* 2014;105(7). http://dx.doi.org/10.1063/1.4893459.

6. Ramachandran V, Fogler HS. Plugging by hydrodynamic bridging during flow of stable colloidal particles within cylindrical pores. *J Fluid Mech.* 1999;385(1):129-156. http://dx.doi.org/10.1017/S0022112098004121, URL https://ui.adsabs.harvard.edu/abs/1999JFM...385.129R.

7. Dersoir B, Schofield A, Robert de Saint Vincent M, Tabuteau H. Dynamics of pore fouling by colloidal particles at the particle level. *J Membr Sci.* 2019;573:411-424. http://dx.doi.org/10.1016/j.memsci.2018.12.025.

8. Okaybi W, Roman S, Soulaine C. Progressive colloidal clogging mechanism by dendritic build-up in porous media. *Soft Matter.* 2025. http://dx.doi.org/10.1039/d5sm00285k.

9. Wyss HM, Blair DL, Morris JF, Stone HA, Weitz DA. Mechanism for clogging of microchannels. *Phys Rev E.* 2006;74. http://dx.doi.org/10.1103/physreve.74.061402.

10. Muecke T. Formation fines and factors controlling their movement in porous media. *J Pet Technol.* 1979;31(02):144-150. http://dx.doi.org/10.2118/7007-PA, arXiv:https://onepetro.org/JPT/article-pdf/31/02/144/2226523/spe-7007-pa.pdf.

11. Sharp KV, Adrian RJ. On flow-blocking particle structures in microtubes. *Microfluid Nanofluidics.* 2005;1:376-380. http://dx.doi.org/10.1007/s10404-005-0043-x.
12. Goldsztein GH, Santamarina JC. Suspension extraction through an opening before clogging. *Appl Phys Lett.* 2004;85(19):4535-4537. http://dx.doi.org/10.1063/1.1818342.

13. Souzy M, Zuriguel I, Marin A. Transition from clogging to continuous flow in constricted particle suspensions. *Phys Rev E.* 2020;101(6):060901. http://dx.doi.org/10.1103/physrev.e.1010.060901.

14. Marin A, Lhuissier H, Rossi M, Kähler CJ. Clogging in constricted suspension flows. *Phys Rev E.* 2018;97(2):0211102. http://dx.doi.org/10.1103/physrev.e.97.0211102.

15. Janda A, Zuriguel I, Garcimartín A, Pugnaloni LA, Maza D. Jamming and critical outlet size in the discharge of a two-dimensional silo. *EPL (Europhysics Letters).* 2008;84(4):44002. http://dx.doi.org/10.1209/0295-5075/84/44002.

16. Agbangla GC, Climent É, Bacchin P. Experimental investigation of pore clogging by microparticles: Evidence for a critical flux density of particle yielding arches and deposits. *Sep Purif Technol.* 2012;101:42-48. http://dx.doi.org/10.1016/j.seppur.2012.09.011.

17. Agbangla GC, Bacchin P, Climent E. Collective dynamics of flowing colloids during pore clogging. *Soft Matter.* 2014;10:6303-6315. http://dx.doi.org/10.1039/c4sm00869c.

18. Dersoir B, Vincent MRd, Abkarian M, Tabuteau H. Clogging of a single pore by colloidal particles. *Microfluid Nanofluid.* 2015;19:953-961. http://dx.doi.org/10.1007/s10404-015-1624-y.

19. Ortega-Roano E, Souzy M, Weinhart T, van der Meer D, Marin A. Clogging of noncohesive suspensions through constrictions using an efficient discrete particle solver with unresolved fluid flow. *Phys Rev E.* 2023;108(064905). http://dx.doi.org/10.1103/physrev.e.108.064905.

20. Vani N, Escudier S, Sauret A. Influence of the solid fraction on the clogging by bridging of suspensions in constricted channels. *Soft Matter.* 2022;18(36):6987-6997. http://dx.doi.org/10.1039/d2sm00962e.

21. Vani N, Escudier S, Jeong D-H, Sauret A. Role of the constriction angle on the clogging by bridging of suspensions of particles. *Phys Rev Res.* 2024;6(3):1032060. http://dx.doi.org/10.1103/physrevresearch.6.1032060.

22. Fogouang Maya L, André L, Soulaine C. Numerical investigation and stochastic modeling of particle bridging undervarious flow, pore, and particle properties. *Adv Water Resour.* 2025;206(105158). http://dx.doi.org/10.1016/j.adwatrems.2025.105158.

23. Tsuji Y, Kawaguchi T, Tanaka T. Discrete particle simulation of two-dimensional fluidized bed. 1993:79-87. 77, http://dx.doi.org/10.1016/0032-5910(93)85010-7, 76.

24. Ladd AJC. Numerical simulations of particulate suspensions via a discretized Boltzmann equation. I: Theoretical foundation. II: Numerical results. *J Fluid Mech.* 1994;271:285-339.

25. Sun D. Hydrodynamics and effect of velocity on particle filtration due to bridging in water-saturated porous media using CFD-DEM simulation. *Environ Sci Pollut Res.* 2023;30:1001-124-100136. http://dx.doi.org/10.1007/s11356-023-29398-5.

26. Maya Fogouang L, André L, Soulaine C. Particulate transport in porous media at pore-scale. Part 1: Unresolved-resolved four-way coupling CFD-DEM. *J Comput Phys.* 2025;521(113540). http://dx.doi.org/10.1016/j.jcp.2024.113540.

27. Mondal S, Wu C-H, Sharma MM. Coupled CFD-DEM simulation of hydrodynamic bridging at constrictions. *Int J Multiphase Flow.* 2016;84:245-263. http://dx.doi.org/10.1016/j.ijmultiphaseflow.2016.05.001. URL https://ui.adsabs.harvard.edu/abs/20161JMF...84.245M.

28. Samari-Kermani M, Jafari S, Rahnama M, Raoof A. Ionic strength and zeta potential effects on colloid transport and retention processes. *Colloid Interface Sci Commun.* 2021;42:100389. http://dx.doi.org/10.1016/j.colcom.2021.100389.

29. Samari-Kermani M, Jafari S, Rahnama M, Raoof A. Direct pore scale numerical simulation of colloid transport and retention. Part I: Fluid flow velocity, colloid size, and pore structure effects. *Adv Water Resour.* 2020;144:103694. http://dx.doi.org/10.1016/j.adwatrems.2020.103694. URL https://ui.adsabs.harvard.edu/abs/2020AdWR...14403694S.

30. Maya Fogouang L, André L, Leroy P, Soulaine C. Particulate transport in porous media at pore-scale. Part 2: CFD-DEM and colloidal forces. *J Comput Phys.* 2024;519:113439. http://dx.doi.org/10.1016/j.jcp.2024.113439. URL https://www.sciencedirect.com/science/article/pii/S0021999124006879.

31. Li H, Wang S, Chen X, Xie L, Shao B, Ma Y. CFD-DEM simulation of aggregation and growth behaviors of fluid-flow-driven migrating particle in porous media. 2023;212343. 231, http://dx.doi.org/10.1016/j.geen.2023.212343.

32. Mirabolghasemi M, Prodanović M, DiCarlo D, Ji H. Prediction of empirical properties using direct pore-scale simulation of straining through 3D microtomography images of porous media. *J Hydrol.* 2015;529:768-778. http://dx.doi.org/10.1016/j.jhydrol.2015.08.016. URL https://ui.adsabs.harvard.edu/abs/2015JHyd...529.768M.

33. Blunt MJ, Bijeljic B, Dong H, Gharbi O, et al. Pore-scale imaging and modelling. *Adv Water Resour.* 2013;51:197-216.

34. Fatt I. The network model of porous media. *Trans AIME.* 1956;207(SPE-574-G):144-181. http://dx.doi.org/10.2118/574-G.

35. Todd AC, Somerville JE, Scott G. The application of depth of formation damage measurements in predicting water injectivity decline. 1984:568-MS. In: SPE international conference and exhibition on formation damage control; SPE Formation Damage Control Symposium, http://dx.doi.org/10.2118/12498-MS, arXiv:https://onepetro.org/SPEFD/proceedings-pdf/84FD/84FD/SPE-12498-MS/4134028/spe-12498-ms.pdf.

36. Meng X, Yang D. Pore-network modeling of particle dispersion in porous media. *Colloids Surfaces A.* 2019;580:123768. http://dx.doi.org/10.1016/j.colsurfa.2019.123768.

37. Rege SD, Fogler HS. Network model for straining dominated particle entrapment in porous media. *Chem Eng Sci.* 1987;42(7):1553-1564. http://dx.doi.org/10.1016/0009-2509(87)80160-4, URL https://ui.adsabs.harvard.edu/abs/1987ChEnS.42.1553R3.

38. Feng Q, Li S, Han X, Wang S. Network simulation for formation impairment due to suspended particles in injected water. *J Pet Sci Eng.* 2015;133:384-391. http://dx.doi.org/10.1016/j petrol.2015.06.027.

39. Lin D, Hu L, Bradford SA, Zhang X, Lo IMC. Pore-network modeling of colloid transport and retention considering surface deposition, hydrodynamic bridging, and straining. *J Hydrol.* 2021;603:127020. http://dx.doi.org/10.1016/j.jhydrol.2021.127020, URL https://ui.adsabs.harvard.edu/abs/2021JHyd.60327020L.

40. Seetha N, Raoof A, Kumar MSM, Hassanizadeh SM. Upscaling of nanoparticle transport in porous media under unfavorable conditions: Pore scale to Darcy scale. *J Contam Hydrol.* 2017;200:1-14. http://dx.doi.org/10.1016/j.jconhyd.2017.03.002.

41. Tang Q, Raoof A. Exploring particle clogging mechanisms through deterministic pore network modeling. (5047338). 2024 http://dx.doi.org/10.2139/ssrn.5047338, SSRN.

42. Rege SD, Fogler HS. A network model for deep bed filtration of solid particles and emulsion drops. *AIChE J.* 1988;34:1761-1772. http://dx.doi.org/10.1002/aic.690341102.

43. Rege SD, Fogler HS. Development of radial models for formation damage in porous media. *Chem Eng Commun.* 1991;108:67-83. http://dx.doi.org/10.1080/00986449108910951.

44. Siqueira AG, Bonet EJ, Scheeca FS. A 3D network model of rock permeability impairment due to suspended particles in injection water. 2003:SPE-82232-MS. In: Vol. SPE European Formation Damage Conference of SPE European Formation Damage Conference and Exhibition; http://dx.doi.org/10.2118/82232-MS, arXiv:https://onepetro.org/SPEEFDC/proceedings-pdf/03EFDC/03EFDC/SPE-82232-MS/2891666/spe-82232-MS.pdf.

45. Yang H, Balhoff MT. Pore-network modeling of particle retention in porous media. *AIChE J.* 2017;63(7):3118-3131.

46. Raoof A, Hassanizadeh SM, Leijnse A. Upscaling transport of adsorbing solutes in porous media: Pore-network modeling. *Vadose Zone J.* 2010;9(3):624-636. http://dx.doi.org/10.2136/vzj2010.0026, URL https://ui.adsabs.harvard.edu/abs/2010VZJ...9.624R.

47. Raoof A, Nick HM, Hassanizadeh SM, Spiers CJ. PoreFlow: A complex pore-network model for simulation of reactive transport in variably saturated porous media. *Comput Geosci.* 2013;61:160-174. http://dx.doi.org/10.1016/j.cageo.2013.08.005, URL https://ui.adsabs.harvard.edu/abs/2013CG...61.160R.

48. Mehmani Y, Tchelepi HA. Minimum requirements for predictive pore-network modeling of solute transport in micromodels. *Adv Water Resour.* 2017;108:83-98. http://dx.doi.org/10.1016/j.adwatrems.2017.07.014, URL https://ui.adsabs.harvard.edu/abs/2017AdWR.108...83M.

49. Lin D, TANG M, ZHANG X, BRADFORD SA, Hu L. Exploring coupled effects of attachment, blocking, straining, and clogging in a pore network model. (4896182). 2024 http://dx.doi.org/10.2139/ssrn.4896182, SSRN.

50. Li Z, Yang H, Sun Z, Espinoza DN, Balhoff MT. A probability-based pore network model of particle jamming in porous media. *Transp Porous Media.* 2021;139:419-445. http://dx.doi.org/10.1007/s11242-021-01673-4.

51. Roman S, Soulaine C, AbuAlSaud M, Kovsek A, Tchelepi H. Particle velocimetry analysis of immiscible two-phase flow in micromodels. *Adv Water Resour.* 2016;95:199-211. http://dx.doi.org/10.1016/j.adwatrems.2015.08.015, URL http://www.sciencedirect.com/science/article/pii/S0309170815002018.

52. Roman S, Rembert F, Kovsek A, Poonoosamy J. Microfluidics for geosciences: metrological developments and future challenges. *Lab Chip.* 2025;25:4273-4289. http://dx.doi.org/10.1039/d5lc00108K.

53. Efros AA, Freeman WT. Image quilting for texture synthesis and transfer. In: *ACM SIGGRAPH conference on computer graphics, los angeles*. 2001.

54. Zarraga IE, Hill DA, Leighton DT. The characterization of the total stress of concentrated suspensions of noncolloidal spheres in Newtonian fluids. *J Rheol.* 2000;44:185-220. http://dx.doi.org/10.1122/1.551083.

55. Mehmani Y, Balhoff MT. Eulerian network modeling of longitudinal dispersion. *Water Resour Res.* 2015;51(10):8586-8606. http://dx.doi.org/10.1002/2015WR017543, URL https://ui.adsabs.harvard.edu/abs/2015WRR...51.8586M.

56. Aris R. On the dispersion of a solute in a fluid flowing through a tube. *Proc Roy. Soc A.* 1956;235:67-77.

57. Taylor. Dispersion of soluble matter in solvent flowing slowly through a tube. *Proc Roy. Soc A.* 1953;219:186-203.

58. Dentz M, Carrera J. Mixing and spreading in stratified flow. *Phys Fluids.* 2007;19(1). http://dx.doi.org/10.1063/1.2427089.

59. Blanc F, D'Ambrosio E, Lobry L, Peters F, Lemaire E. Universal scaling law in frictional non-Brownian suspensions. *Phys Rev Fluids.* 2018;3(11):114303. http://dx.doi.org/10.1103/physrevfluids.3.114303.

60. Lobry L, Lemaire E, Blanc F, Gallier S, Peters F. Shear thinning in non-Brownian suspensions explained by variable friction between particles. *J Fluid Mech.* 2018;860:682-710. http://dx.doi.org/10.1017/jfm.2018.881.
61. d'Ambrosio E, Gilbert D, Blanc F, Cohen C, Lemaire E. Rheology of suspensions of cubic particles. *J Rheol.* 2025;69(6):807-828. http://dx.doi.org/10.1122/8.0001045.

62. Zuriguel I, Garcimartín A, Maza D, Pugnaloni LA, Pastor JM. Jamming during the discharge of granular matter from a silo. *Phys Rev E.* 2005;71(5):051303. http://dx.doi.org/10.1103/physreve.71.051303.

63. Whitaker S. The forchheimer equation: A theoretical development. *Transp Porous Media*. 1996;25(1):27-61. http://dx.doi.org/10.1007/bf00141261.

64. Roman CSAJMAS. Computational microfluidics for geosciences. *Front Water*. 2021;3:1-11. http://dx.doi.org/10.3389/frwa.2021.643714.

65. Angot P, Bruneau C-H, Fabrie P. A penalization method to take into account obstacles in incompressible viscous flows. *Numer Math.* 1999;81(4):497-520. http://dx.doi.org/10.1007/s002110050401.

66. Cundall PA, Strack ODL. A discrete numerical model for granular assemblies. *Géotechnique*. 1979;29(1):47-65. http://dx.doi.org/10.1680/geot.1979.29.1.47.

67. Uhlmann M. An immersed boundary method with direct forcing for the simulation of particulate flows. *J Comput Phys.* 2005;209(2):448-476. http://dx.doi.org/10.1016/j.jcp.2005.03.017.

68. Whitesides GM, Ostuni E, Takayama S, Jiang X, Ingber DE. Soft lithography in biology and biochemistry. *Annu Rev Biomed Eng.* 2001;3(1):335-373.

69. Dincau B, Tang C, Dressaire E, Sauret A. Clog mitigation in a microfluidic array via pulsatile flows. *Soft Matter*. 2022;18(9):1767-1778.

70. Preibisch S, Saalfeld S, Tomancak P. Globally optimal stitching of tiled 3D microscopic image acquisitions. *Bioinformatics*. 2009;25(11):1463-1465.

71. Zami-Pierre F, de Loubens R, Quintard M, Davit Y. Effect of disorder in the pore-scale structure on the flow of shear-thinning fluids through porous media. *J Non-Newton Fluid Mech.* 2018;261:99-110. http://dx.doi.org/10.1016/j.jnnfm.2018.08.004, URL https://ui.adsabs.harvard.edu/abs/2018JNNFM.261...99Z.

72. Elmøe TD, Tricoli A, Grunwaldt J-D, Pratsinis SE. Filtration of nanoparticles: evolution of cake structure and pressure-drop. *J Aerosol Sci.* 2009;40(11):965-981.

73. Chen Y, Liu H, Gao Y, Ye X, et al. Visualization of the microscopic mechanism of slurry infiltration and filter cake formation for slurry shield tunneling in saturated sand: A microfluidic chip experiment. *Transp Geotech.* 2025;51:101529.

74. Dersoir B, de Saint Vincent MR, Abkarian M, Tabuteau H. Clogging of a single pore by colloidal particles. *Microfluid Nanofluidics*. 2015;19:953-961.

75. Maya Fogouang L, Soulaine C, Okaybi W, LE TRONG E, Roman S. Dataset and model for "permeability impairment by hydrodynamic pore bridging: probabilistic pore-network modelling and microfluidic experiments" [data set]. *Zenodo*. 2026. http://dx.doi.org/10.5281/zenodo.18300225.

76. Okaybi W, Roman S, Maya L, Le Trong E, Soulaine C. Hydrodynamic pore-bridging and permeability impairment in heterogeneous porous media: microfluidic data-set [dataset]. *Digital Porous Media Portal*. 2026. https://www.doi.org/10.17612/d1fv-1z91.